3930bfdec7a0f99d7ae2eb605a266541ca0405a4
[ghc-hetmet.git] / ghc / lib / std / PrelNum.lhs
1 %
2 % (c) The AQUA Project, Glasgow University, 1994-1996
3 %
4
5 \section[PrelNum]{Module @PrelNum@}
6
7 Numeric part of the prelude.
8
9 It's rather big!
10
11 \begin{code}
12 {-# OPTIONS -fno-implicit-prelude -#include "cbits/floatExtreme.h" #-}
13 {-# OPTIONS -H20m #-}
14
15 #include "../includes/ieee-flpt.h"
16
17 \end{code}
18
19 \begin{code}
20 module PrelNum where
21
22 import PrelBase
23 import PrelGHC
24 import {-# SOURCE #-} PrelErr ( error )
25 import PrelList
26 import PrelMaybe
27
28 import PrelArr          ( Array, array, (!) )
29 import PrelIOBase       ( unsafePerformIO )
30 import Ix               ( Ix(..) )
31 import PrelCCall        ()      -- we need the definitions of CCallable and 
32                                 -- CReturnable for the _ccall_s herein.
33                 
34
35 infixr 8  ^, ^^, **
36 infixl 7  /, %, `quot`, `rem`, `div`, `mod`
37 \end{code}
38
39
40 %*********************************************************
41 %*                                                      *
42 \subsection{Standard numeric classes}
43 %*                                                      *
44 %*********************************************************
45
46 \begin{code}
47 class  (Num a, Ord a) => Real a  where
48     toRational          ::  a -> Rational
49
50 class  (Real a, Enum a) => Integral a  where
51     quot, rem, div, mod :: a -> a -> a
52     quotRem, divMod     :: a -> a -> (a,a)
53     toInteger           :: a -> Integer
54     toInt               :: a -> Int -- partain: Glasgow extension
55
56     n `quot` d          =  q  where (q,r) = quotRem n d
57     n `rem` d           =  r  where (q,r) = quotRem n d
58     n `div` d           =  q  where (q,r) = divMod n d
59     n `mod` d           =  r  where (q,r) = divMod n d
60     divMod n d          =  if signum r == negate (signum d) then (q-1, r+d) else qr
61                            where qr@(q,r) = quotRem n d
62
63 class  (Num a) => Fractional a  where
64     (/)                 :: a -> a -> a
65     recip               :: a -> a
66     fromRational        :: Rational -> a
67
68     recip x             =  1 / x
69
70 class  (Fractional a) => Floating a  where
71     pi                  :: a
72     exp, log, sqrt      :: a -> a
73     (**), logBase       :: a -> a -> a
74     sin, cos, tan       :: a -> a
75     asin, acos, atan    :: a -> a
76     sinh, cosh, tanh    :: a -> a
77     asinh, acosh, atanh :: a -> a
78
79     x ** y              =  exp (log x * y)
80     logBase x y         =  log y / log x
81     sqrt x              =  x ** 0.5
82     tan  x              =  sin  x / cos  x
83     tanh x              =  sinh x / cosh x
84
85 class  (Real a, Fractional a) => RealFrac a  where
86     properFraction      :: (Integral b) => a -> (b,a)
87     truncate, round     :: (Integral b) => a -> b
88     ceiling, floor      :: (Integral b) => a -> b
89
90     truncate x          =  m  where (m,_) = properFraction x
91     
92     round x             =  let (n,r) = properFraction x
93                                m     = if r < 0 then n - 1 else n + 1
94                            in case signum (abs r - 0.5) of
95                                 -1 -> n
96                                 0  -> if even n then n else m
97                                 1  -> m
98     
99     ceiling x           =  if r > 0 then n + 1 else n
100                            where (n,r) = properFraction x
101     
102     floor x             =  if r < 0 then n - 1 else n
103                            where (n,r) = properFraction x
104
105 class  (RealFrac a, Floating a) => RealFloat a  where
106     floatRadix          :: a -> Integer
107     floatDigits         :: a -> Int
108     floatRange          :: a -> (Int,Int)
109     decodeFloat         :: a -> (Integer,Int)
110     encodeFloat         :: Integer -> Int -> a
111     exponent            :: a -> Int
112     significand         :: a -> a
113     scaleFloat          :: Int -> a -> a
114     isNaN, isInfinite, isDenormalized, isNegativeZero, isIEEE
115                         :: a -> Bool
116
117     exponent x          =  if m == 0 then 0 else n + floatDigits x
118                            where (m,n) = decodeFloat x
119
120     significand x       =  encodeFloat m (negate (floatDigits x))
121                            where (m,_) = decodeFloat x
122
123     scaleFloat k x      =  encodeFloat m (n+k)
124                            where (m,n) = decodeFloat x
125 \end{code}
126
127 %*********************************************************
128 %*                                                      *
129 \subsection{Overloaded numeric functions}
130 %*                                                      *
131 %*********************************************************
132
133 \begin{code}
134 even, odd       :: (Integral a) => a -> Bool
135 even n          =  n `rem` 2 == 0
136 odd             =  not . even
137
138 {-# SPECIALISE gcd ::
139         Int -> Int -> Int,
140         Integer -> Integer -> Integer #-}
141 gcd             :: (Integral a) => a -> a -> a
142 gcd 0 0         =  error "Prelude.gcd: gcd 0 0 is undefined"
143 gcd x y         =  gcd' (abs x) (abs y)
144                    where gcd' x 0  =  x
145                          gcd' x y  =  gcd' y (x `rem` y)
146
147 {-# SPECIALISE lcm ::
148         Int -> Int -> Int,
149         Integer -> Integer -> Integer #-}
150 lcm             :: (Integral a) => a -> a -> a
151 lcm _ 0         =  0
152 lcm 0 _         =  0
153 lcm x y         =  abs ((x `quot` (gcd x y)) * y)
154
155 {-# SPECIALISE (^) ::
156         Integer -> Integer -> Integer,
157         Integer -> Int -> Integer,
158         Int -> Int -> Int #-}
159 (^)             :: (Num a, Integral b) => a -> b -> a
160 x ^ 0           =  1
161 x ^ n | n > 0   =  f x (n-1) x
162                    where f _ 0 y = y
163                          f x n y = g x n  where
164                                    g x n | even n  = g (x*x) (n `quot` 2)
165                                          | otherwise = f x (n-1) (x*y)
166 _ ^ _           = error "Prelude.^: negative exponent"
167
168 {-# SPECIALISE (^^) ::
169         Double -> Int -> Double,
170         Rational -> Int -> Rational #-}
171 (^^)            :: (Fractional a, Integral b) => a -> b -> a
172 x ^^ n          =  if n >= 0 then x^n else recip (x^(negate n))
173
174 {-# SPECIALIZE fromIntegral ::
175     Int         -> Rational,
176     Integer     -> Rational,
177     Int         -> Int,
178     Int         -> Integer,
179     Int         -> Float,
180     Int         -> Double,
181     Integer     -> Int,
182     Integer     -> Integer,
183     Integer     -> Float,
184     Integer     -> Double #-}
185 fromIntegral    :: (Integral a, Num b) => a -> b
186 fromIntegral    =  fromInteger . toInteger
187
188 {-# SPECIALIZE fromRealFrac ::
189     Double      -> Rational, 
190     Rational    -> Double,
191     Float       -> Rational,
192     Rational    -> Float,
193     Rational    -> Rational,
194     Double      -> Double,
195     Double      -> Float,
196     Float       -> Float,
197     Float       -> Double #-}
198 fromRealFrac    :: (RealFrac a, Fractional b) => a -> b
199 fromRealFrac    =  fromRational . toRational
200
201 atan2           :: (RealFloat a) => a -> a -> a
202 atan2 y x       =  case (signum y, signum x) of
203                         ( 0, 1) ->  0
204                         ( 1, 0) ->  pi/2
205                         ( 0,-1) ->  pi
206                         (-1, 0) ->  (negate pi)/2
207                         ( _, 1) ->  atan (y/x)
208                         ( _,-1) ->  atan (y/x) + pi
209                         ( 0, 0) ->  error "Prelude.atan2: atan2 of origin"
210 \end{code}
211
212
213 %*********************************************************
214 %*                                                      *
215 \subsection{Instances for @Int@}
216 %*                                                      *
217 %*********************************************************
218
219 \begin{code}
220 instance  Real Int  where
221     toRational x        =  toInteger x % 1
222
223 instance  Integral Int  where
224     a@(I# _) `quotRem` b@(I# _) = (a `quotInt` b, a `remInt` b)
225     -- OK, so I made it a little stricter.  Shoot me.  (WDP 94/10)
226
227     -- Following chks for zero divisor are non-standard (WDP)
228     a `quot` b  =  if b /= 0
229                    then a `quotInt` b
230                    else error "Integral.Int.quot{PreludeCore}: divide by 0\n"
231     a `rem` b   =  if b /= 0
232                    then a `remInt` b
233                    else error "Integral.Int.rem{PreludeCore}: divide by 0\n"
234
235     n `div` d
236      | n > 0 && d < 0 = mk_neg (quotInt (n-d-1) d)
237      | n < 0 && d > 0 = mk_neg (quotInt (n-d+1) d)
238      | otherwise      = quotInt n d
239       where
240        {-
241          - the result of (integral) division is
242            defined as being truncated towards
243            negative infinity. (see Sec 6.3.2 of
244            the Haskell 1.4 report.)
245
246          - in the case of Int, if either nominator or
247            denominator is negative, we adjust the nominator
248            to account for the above property before
249            computing the quotient.
250
251          - in the case of Int, the adjustment of the
252            nominator runs the risk of overflowing. If
253            we make the assumption that arithmetic is
254            modulo word size, and adjust the final result
255            to account for this.
256        -}
257
258        mk_neg r 
259         | r <= 0    = r
260         | otherwise = -(r+1)
261
262     x `mod` y = if x > 0 && y < 0 || x < 0 && y > 0 then
263                     if r/=0 then r+y else 0
264                 else
265                     r
266               where r = remInt x y
267
268     divMod x@(I# _) y@(I# _) = (x `div` y, x `mod` y)
269     -- Stricter.  Sorry if you don't like it.  (WDP 94/10)
270
271 --OLD:   even x = eqInt (x `mod` 2) 0
272 --OLD:   odd x  = neInt (x `mod` 2) 0
273
274     toInteger (I# n#) = int2Integer# n#  -- give back a full-blown Integer
275     toInt x           = x
276
277 \end{code}
278
279
280 %*********************************************************
281 %*                                                      *
282 \subsection{Type @Integer@}
283 %*                                                      *
284 %*********************************************************
285
286 These types are used to return from integer primops
287
288 \begin{code}
289 data Return2GMPs     = Return2GMPs     Int# Int# ByteArray# Int# Int# ByteArray#
290 data ReturnIntAndGMP = ReturnIntAndGMP Int# Int# Int# ByteArray#
291 \end{code}
292
293 Instances
294
295 \begin{code}
296 instance  Eq Integer  where
297     (J# a1 s1 d1) == (J# a2 s2 d2)
298       = (cmpInteger# a1 s1 d1 a2 s2 d2) ==# 0#
299
300     (J# a1 s1 d1) /= (J# a2 s2 d2)
301       = (cmpInteger# a1 s1 d1 a2 s2 d2) /=# 0#
302
303 instance  Ord Integer  where
304     (J# a1 s1 d1) <= (J# a2 s2 d2)
305       = (cmpInteger# a1 s1 d1 a2 s2 d2) <=# 0#
306
307     (J# a1 s1 d1) <  (J# a2 s2 d2)
308       = (cmpInteger# a1 s1 d1 a2 s2 d2) <# 0#
309
310     (J# a1 s1 d1) >= (J# a2 s2 d2)
311       = (cmpInteger# a1 s1 d1 a2 s2 d2) >=# 0#
312
313     (J# a1 s1 d1) >  (J# a2 s2 d2)
314       = (cmpInteger# a1 s1 d1 a2 s2 d2) ># 0#
315
316     x@(J# a1 s1 d1) `max` y@(J# a2 s2 d2)
317       = if ((cmpInteger# a1 s1 d1 a2 s2 d2) ># 0#) then x else y
318
319     x@(J# a1 s1 d1) `min` y@(J# a2 s2 d2)
320       = if ((cmpInteger# a1 s1 d1 a2 s2 d2) <# 0#) then x else y
321
322     compare (J# a1 s1 d1) (J# a2 s2 d2)
323        = case cmpInteger# a1 s1 d1 a2 s2 d2 of { res# ->
324          if res# <# 0# then LT else 
325          if res# ># 0# then GT else EQ
326          }
327
328 instance  Num Integer  where
329     (+) (J# a1 s1 d1) (J# a2 s2 d2)
330       = plusInteger# a1 s1 d1 a2 s2 d2
331
332     (-) (J# a1 s1 d1) (J# a2 s2 d2)
333       = minusInteger# a1 s1 d1 a2 s2 d2
334
335     negate (J# a s d) = negateInteger# a s d
336
337     (*) (J# a1 s1 d1) (J# a2 s2 d2)
338       = timesInteger# a1 s1 d1 a2 s2 d2
339
340     -- ORIG: abs n = if n >= 0 then n else -n
341
342     abs n@(J# a1 s1 d1)
343       = case 0 of { J# a2 s2 d2 ->
344         if (cmpInteger# a1 s1 d1 a2 s2 d2) >=# 0#
345         then n
346         else negateInteger# a1 s1 d1
347         }
348
349     signum n@(J# a1 s1 d1)
350       = case 0 of { J# a2 s2 d2 ->
351         let
352             cmp = cmpInteger# a1 s1 d1 a2 s2 d2
353         in
354         if      cmp >#  0# then 1
355         else if cmp ==# 0# then 0
356         else                    (negate 1)
357         }
358
359     fromInteger x       =  x
360
361     fromInt (I# n#)     =  int2Integer# n# -- gives back a full-blown Integer
362
363 instance  Real Integer  where
364     toRational x        =  x % 1
365
366 instance  Integral Integer where
367     quotRem (J# a1 s1 d1) (J# a2 s2 d2)
368       = case (quotRemInteger# a1 s1 d1 a2 s2 d2) of
369           Return2GMPs a3 s3 d3 a4 s4 d4
370             -> (J# a3 s3 d3, J# a4 s4 d4)
371
372 {- USING THE UNDERLYING "GMP" CODE IS DUBIOUS FOR NOW:
373
374     divMod (J# a1 s1 d1) (J# a2 s2 d2)
375       = case (divModInteger# a1 s1 d1 a2 s2 d2) of
376           Return2GMPs a3 s3 d3 a4 s4 d4
377             -> (J# a3 s3 d3, J# a4 s4 d4)
378 -}
379     toInteger n      = n
380     toInt (J# a s d) = case (integer2Int# a s d) of { n# -> I# n# }
381
382     -- the rest are identical to the report default methods;
383     -- you get slightly better code if you let the compiler
384     -- see them right here:
385     n `quot` d  =  if d /= 0 then q else 
386                      error "Integral.Integer.quot{PreludeCore}: divide by 0\n"  
387                    where (q,r) = quotRem n d
388     n `rem` d   =  if d /= 0 then r else 
389                      error "Integral.Integer.quot{PreludeCore}: divide by 0\n"  
390                    where (q,r) = quotRem n d
391     n `div` d   =  q  where (q,r) = divMod n d
392     n `mod` d   =  r  where (q,r) = divMod n d
393
394     divMod n d  =  case (quotRem n d) of { qr@(q,r) ->
395                    if signum r == negate (signum d) then (q - 1, r+d) else qr }
396                    -- Case-ified by WDP 94/10
397
398 instance  Enum Integer  where
399     toEnum n             =  toInteger n
400     fromEnum n           =  toInt n
401     enumFrom n           =  n : enumFrom (n + 1)
402     enumFromThen m n     =  en' m (n - m)
403                             where en' m n = m : en' (m + n) n
404     enumFromTo n m       =  takeWhile (<= m) (enumFrom n)
405     enumFromThenTo n m p =  takeWhile (if m >= n then (<= p) else (>= p))
406                                       (enumFromThen n m)
407
408 instance  Show Integer  where
409     showsPrec   x = showSignedInteger x
410     showList = showList__ (showsPrec 0) 
411
412 instance  Ix Integer  where
413     range (m,n)         =  [m..n]
414     index b@(m,n) i
415         | inRange b i   =  fromInteger (i - m)
416         | otherwise     =  error "Integer.index: Index out of range."
417     inRange (m,n) i     =  m <= i && i <= n
418
419 integer_0, integer_1, integer_2, integer_m1 :: Integer
420 integer_0  = int2Integer# 0#
421 integer_1  = int2Integer# 1#
422 integer_2  = int2Integer# 2#
423 integer_m1 = int2Integer# (negateInt# 1#)
424 \end{code}
425
426
427 %*********************************************************
428 %*                                                      *
429 \subsection{Type @Float@}
430 %*                                                      *
431 %*********************************************************
432
433 \begin{code}
434 instance Eq Float where
435     (F# x) == (F# y) = x `eqFloat#` y
436
437 instance Ord Float where
438     (F# x) `compare` (F# y) | x `ltFloat#` y = LT
439                             | x `eqFloat#` y = EQ
440                             | otherwise      = GT
441
442     (F# x) <  (F# y) = x `ltFloat#`  y
443     (F# x) <= (F# y) = x `leFloat#`  y
444     (F# x) >= (F# y) = x `geFloat#`  y
445     (F# x) >  (F# y) = x `gtFloat#`  y
446
447 instance  Num Float  where
448     (+)         x y     =  plusFloat x y
449     (-)         x y     =  minusFloat x y
450     negate      x       =  negateFloat x
451     (*)         x y     =  timesFloat x y
452     abs x | x >= 0.0    =  x
453           | otherwise   =  negateFloat x
454     signum x | x == 0.0  = 0
455              | x > 0.0   = 1
456              | otherwise = negate 1
457     fromInteger n       =  encodeFloat n 0
458     fromInt i           =  int2Float i
459
460 instance  Real Float  where
461     toRational x        =  (m%1)*(b%1)^^n
462                            where (m,n) = decodeFloat x
463                                  b     = floatRadix  x
464
465 instance  Fractional Float  where
466     (/) x y             =  divideFloat x y
467     fromRational x      =  fromRat x
468     recip x             =  1.0 / x
469
470 instance  Floating Float  where
471     pi                  =  3.141592653589793238
472     exp x               =  expFloat x
473     log x               =  logFloat x
474     sqrt x              =  sqrtFloat x
475     sin x               =  sinFloat x
476     cos x               =  cosFloat x
477     tan x               =  tanFloat x
478     asin x              =  asinFloat x
479     acos x              =  acosFloat x
480     atan x              =  atanFloat x
481     sinh x              =  sinhFloat x
482     cosh x              =  coshFloat x
483     tanh x              =  tanhFloat x
484     (**) x y            =  powerFloat x y
485     logBase x y         =  log y / log x
486
487     asinh x = log (x + sqrt (1.0+x*x))
488     acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))
489     atanh x = log ((x+1.0) / sqrt (1.0-x*x))
490
491 instance  RealFrac Float  where
492
493     {-# SPECIALIZE properFraction :: Float -> (Int, Float) #-}
494     {-# SPECIALIZE truncate :: Float -> Int #-}
495     {-# SPECIALIZE round    :: Float -> Int #-}
496     {-# SPECIALIZE ceiling  :: Float -> Int #-}
497     {-# SPECIALIZE floor    :: Float -> Int #-}
498
499     {-# SPECIALIZE properFraction :: Float -> (Integer, Float) #-}
500     {-# SPECIALIZE truncate :: Float -> Integer #-}
501     {-# SPECIALIZE round    :: Float -> Integer #-}
502     {-# SPECIALIZE ceiling  :: Float -> Integer #-}
503     {-# SPECIALIZE floor    :: Float -> Integer #-}
504
505     properFraction x
506       = case (decodeFloat x)      of { (m,n) ->
507         let  b = floatRadix x     in
508         if n >= 0 then
509             (fromInteger m * fromInteger b ^ n, 0.0)
510         else
511             case (quotRem m (b^(negate n))) of { (w,r) ->
512             (fromInteger w, encodeFloat r n)
513             }
514         }
515
516     truncate x  = case properFraction x of
517                      (n,_) -> n
518
519     round x     = case properFraction x of
520                      (n,r) -> let
521                                 m         = if r < 0.0 then n - 1 else n + 1
522                                 half_down = abs r - 0.5
523                               in
524                               case (compare half_down 0.0) of
525                                 LT -> n
526                                 EQ -> if even n then n else m
527                                 GT -> m
528
529     ceiling x   = case properFraction x of
530                     (n,r) -> if r > 0.0 then n + 1 else n
531
532     floor x     = case properFraction x of
533                     (n,r) -> if r < 0.0 then n - 1 else n
534
535 instance  RealFloat Float  where
536     floatRadix _        =  FLT_RADIX        -- from float.h
537     floatDigits _       =  FLT_MANT_DIG     -- ditto
538     floatRange _        =  (FLT_MIN_EXP, FLT_MAX_EXP) -- ditto
539
540     decodeFloat (F# f#)
541       = case decodeFloat# f#    of
542           ReturnIntAndGMP exp# a# s# d# ->
543             (J# a# s# d#, I# exp#)
544
545     encodeFloat (J# a# s# d#) (I# e#)
546       = case encodeFloat# a# s# d# e# of { flt# -> F# flt# }
547
548     exponent x          = case decodeFloat x of
549                             (m,n) -> if m == 0 then 0 else n + floatDigits x
550
551     significand x       = case decodeFloat x of
552                             (m,_) -> encodeFloat m (negate (floatDigits x))
553
554     scaleFloat k x      = case decodeFloat x of
555                             (m,n) -> encodeFloat m (n+k)
556     isNaN x = 
557       (0::Int) /= unsafePerformIO (_ccall_ isFloatNaN x) {- a _pure_function! -}
558     isInfinite x =
559       (0::Int) /= unsafePerformIO (_ccall_ isFloatInfinite x) {- ditto! -}
560     isDenormalized x =
561       (0::Int) /= unsafePerformIO (_ccall_ isFloatDenormalized x) -- ..
562     isNegativeZero x =
563       (0::Int) /= unsafePerformIO (_ccall_ isFloatNegativeZero x) -- ...
564     isIEEE x    = True
565
566 instance  Show Float  where
567     showsPrec   x = showSigned showFloat x
568     showList = showList__ (showsPrec 0) 
569 \end{code}
570
571 %*********************************************************
572 %*                                                      *
573 \subsection{Type @Double@}
574 %*                                                      *
575 %*********************************************************
576
577 \begin{code}
578 instance Eq Double where
579     (D# x) == (D# y) = x ==## y
580
581 instance Ord Double where
582     (D# x) `compare` (D# y) | x <## y   = LT
583                             | x ==## y  = EQ
584                             | otherwise = GT
585
586     (D# x) <  (D# y) = x <##  y
587     (D# x) <= (D# y) = x <=## y
588     (D# x) >= (D# y) = x >=## y
589     (D# x) >  (D# y) = x >##  y
590
591 instance  Num Double  where
592     (+)         x y     =  plusDouble x y
593     (-)         x y     =  minusDouble x y
594     negate      x       =  negateDouble x
595     (*)         x y     =  timesDouble x y
596     abs x | x >= 0.0    =  x
597           | otherwise   =  negateDouble x
598     signum x | x == 0.0  = 0
599              | x > 0.0   = 1
600              | otherwise = negate 1
601     fromInteger n       =  encodeFloat n 0
602     fromInt (I# n#)     =  case (int2Double# n#) of { d# -> D# d# }
603
604 instance  Real Double  where
605     toRational x        =  (m%1)*(b%1)^^n
606                            where (m,n) = decodeFloat x
607                                  b     = floatRadix  x
608
609 instance  Fractional Double  where
610     (/) x y             =  divideDouble x y
611     fromRational x      =  fromRat x
612     recip x             =  1.0 / x
613
614 instance  Floating Double  where
615     pi                  =  3.141592653589793238
616     exp x               =  expDouble x
617     log x               =  logDouble x
618     sqrt x              =  sqrtDouble x
619     sin  x              =  sinDouble x
620     cos  x              =  cosDouble x
621     tan  x              =  tanDouble x
622     asin x              =  asinDouble x
623     acos x              =  acosDouble x
624     atan x              =  atanDouble x
625     sinh x              =  sinhDouble x
626     cosh x              =  coshDouble x
627     tanh x              =  tanhDouble x
628     (**) x y            =  powerDouble x y
629     logBase x y         =  log y / log x
630
631     asinh x = log (x + sqrt (1.0+x*x))
632     acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))
633     atanh x = log ((x+1.0) / sqrt (1.0-x*x))
634
635 instance  RealFrac Double  where
636
637     {-# SPECIALIZE properFraction :: Double -> (Int, Double) #-}
638     {-# SPECIALIZE truncate :: Double -> Int #-}
639     {-# SPECIALIZE round    :: Double -> Int #-}
640     {-# SPECIALIZE ceiling  :: Double -> Int #-}
641     {-# SPECIALIZE floor    :: Double -> Int #-}
642
643     {-# SPECIALIZE properFraction :: Double -> (Integer, Double) #-}
644     {-# SPECIALIZE truncate :: Double -> Integer #-}
645     {-# SPECIALIZE round    :: Double -> Integer #-}
646     {-# SPECIALIZE ceiling  :: Double -> Integer #-}
647     {-# SPECIALIZE floor    :: Double -> Integer #-}
648
649 #if defined(__UNBOXED_INSTANCES__)
650     {-# SPECIALIZE properFraction :: Double -> (Int#, Double) #-}
651     {-# SPECIALIZE truncate :: Double -> Int# #-}
652     {-# SPECIALIZE round    :: Double -> Int# #-}
653     {-# SPECIALIZE ceiling  :: Double -> Int# #-}
654     {-# SPECIALIZE floor    :: Double -> Int# #-}
655 #endif
656
657     properFraction x
658       = case (decodeFloat x)      of { (m,n) ->
659         let  b = floatRadix x     in
660         if n >= 0 then
661             (fromInteger m * fromInteger b ^ n, 0.0)
662         else
663             case (quotRem m (b^(negate n))) of { (w,r) ->
664             (fromInteger w, encodeFloat r n)
665             }
666         }
667
668     truncate x  = case properFraction x of
669                      (n,_) -> n
670
671     round x     = case properFraction x of
672                      (n,r) -> let
673                                 m         = if r < 0.0 then n - 1 else n + 1
674                                 half_down = abs r - 0.5
675                               in
676                               case (compare half_down 0.0) of
677                                 LT -> n
678                                 EQ -> if even n then n else m
679                                 GT -> m
680
681     ceiling x   = case properFraction x of
682                     (n,r) -> if r > 0.0 then n + 1 else n
683
684     floor x     = case properFraction x of
685                     (n,r) -> if r < 0.0 then n - 1 else n
686
687 instance  RealFloat Double  where
688     floatRadix _        =  FLT_RADIX        -- from float.h
689     floatDigits _       =  DBL_MANT_DIG     -- ditto
690     floatRange _        =  (DBL_MIN_EXP, DBL_MAX_EXP) -- ditto
691
692     decodeFloat (D# d#)
693       = case decodeDouble# d#   of
694           ReturnIntAndGMP exp# a# s# d# ->
695             (J# a# s# d#, I# exp#)
696
697     encodeFloat (J# a# s# d#) (I# e#)
698       = case encodeDouble# a# s# d# e#  of { dbl# -> D# dbl# }
699
700     exponent x          = case decodeFloat x of
701                             (m,n) -> if m == 0 then 0 else n + floatDigits x
702
703     significand x       = case decodeFloat x of
704                             (m,_) -> encodeFloat m (negate (floatDigits x))
705
706     scaleFloat k x      = case decodeFloat x of
707                             (m,n) -> encodeFloat m (n+k)
708     isNaN x = 
709       (0::Int) /= unsafePerformIO (_ccall_ isDoubleNaN x) {- a _pure_function! -}
710     isInfinite x =
711       (0::Int) /= unsafePerformIO (_ccall_ isDoubleInfinite x) {- ditto -}
712     isDenormalized x =
713       (0::Int) /= unsafePerformIO (_ccall_ isDoubleDenormalized x) -- ..
714     isNegativeZero x =
715       (0::Int) /= unsafePerformIO (_ccall_ isDoubleNegativeZero x) -- ...
716     isIEEE x    = True
717
718 instance  Show Double  where
719     showsPrec   x = showSigned showFloat x
720     showList = showList__ (showsPrec 0) 
721 \end{code}
722
723
724 %*********************************************************
725 %*                                                      *
726 \subsection{Common code for @Float@ and @Double@}
727 %*                                                      *
728 %*********************************************************
729
730 The @Enum@ instances for Floats and Doubles are slightly unusual.
731 The @toEnum@ function truncates numbers to Int.  The definitions
732 of @enumFrom@ and @enumFromThen@ allow floats to be used in arithmetic
733 series: [0,0.1 .. 1.0].  However, roundoff errors make these somewhat
734 dubious.  This example may have either 10 or 11 elements, depending on
735 how 0.1 is represented.
736
737 NOTE: The instances for Float and Double do not make use of the default
738 methods for @enumFromTo@ and @enumFromThenTo@, as these rely on there being
739 a `non-lossy' conversion to and from Ints. Instead we make use of the 
740 1.2 default methods (back in the days when Enum had Ord as a superclass)
741 for these (@numericEnumFromTo@ and @numericEnumFromThenTo@ below.)
742
743 \begin{code}
744 instance  Enum Float  where
745     toEnum         =  fromIntegral
746     fromEnum       =  fromInteger . truncate   -- may overflow
747     enumFrom       =  numericEnumFrom
748     enumFromThen   =  numericEnumFromThen
749     enumFromThenTo =  numericEnumFromThenTo
750
751 instance  Enum Double  where
752     toEnum         =  fromIntegral
753     fromEnum       =  fromInteger . truncate   -- may overflow
754     enumFrom       =  numericEnumFrom
755     enumFromThen   =  numericEnumFromThen
756     enumFromThenTo =  numericEnumFromThenTo
757
758 numericEnumFrom         :: (Real a) => a -> [a]
759 numericEnumFromThen     :: (Real a) => a -> a -> [a]
760 numericEnumFromThenTo   :: (Real a) => a -> a -> a -> [a]
761 numericEnumFrom         =  iterate (+1)
762 numericEnumFromThen n m =  iterate (+(m-n)) n
763 numericEnumFromThenTo n m p = takeWhile (if m >= n then (<= p) else (>= p))
764                                       (numericEnumFromThen n m)
765 \end{code}
766
767
768 %*********************************************************
769 %*                                                      *
770 \subsection{The @Ratio@ and @Rational@ types}
771 %*                                                      *
772 %*********************************************************
773
774 \begin{code}
775 data  (Eval a, Integral a)      => Ratio a = !a :% !a  deriving (Eq)
776 type  Rational          =  Ratio Integer
777 \end{code}
778
779 \begin{code}
780 {-# SPECIALISE (%) :: Integer -> Integer -> Rational #-}
781
782 (%)                     :: (Integral a) => a -> a -> Ratio a
783 numerator, denominator  :: (Integral a) => Ratio a -> a
784 approxRational          :: (RealFrac a) => a -> a -> Rational
785
786 \end{code}
787
788 \tr{reduce} is a subsidiary function used only in this module .
789 It normalises a ratio by dividing both numerator and denominator by
790 their greatest common divisor.
791
792 \begin{code}
793 reduce x 0              =  error "{Ratio.%}: zero denominator"
794 reduce x y              =  (x `quot` d) :% (y `quot` d)
795                            where d = gcd x y
796 \end{code}
797
798 \begin{code}
799 x % y                   =  reduce (x * signum y) (abs y)
800
801 numerator (x:%y)        =  x
802
803 denominator (x:%y)      =  y
804 \end{code}
805
806
807 @approxRational@, applied to two real fractional numbers x and epsilon,
808 returns the simplest rational number within epsilon of x.  A rational
809 number n%d in reduced form is said to be simpler than another n'%d' if
810 abs n <= abs n' && d <= d'.  Any real interval contains a unique
811 simplest rational; here, for simplicity, we assume a closed rational
812 interval.  If such an interval includes at least one whole number, then
813 the simplest rational is the absolutely least whole number.  Otherwise,
814 the bounds are of the form q%1 + r%d and q%1 + r'%d', where abs r < d
815 and abs r' < d', and the simplest rational is q%1 + the reciprocal of
816 the simplest rational between d'%r' and d%r.
817
818 \begin{code}
819 approxRational x eps    =  simplest (x-eps) (x+eps)
820         where simplest x y | y < x      =  simplest y x
821                            | x == y     =  xr
822                            | x > 0      =  simplest' n d n' d'
823                            | y < 0      =  - simplest' (-n') d' (-n) d
824                            | otherwise  =  0 :% 1
825                                         where xr  = toRational x
826                                               n   = numerator xr
827                                               d   = denominator xr
828                                               nd' = toRational y
829                                               n'  = numerator nd'
830                                               d'  = denominator nd'
831
832               simplest' n d n' d'       -- assumes 0 < n%d < n'%d'
833                         | r == 0     =  q :% 1
834                         | q /= q'    =  (q+1) :% 1
835                         | otherwise  =  (q*n''+d'') :% n''
836                                      where (q,r)      =  quotRem n d
837                                            (q',r')    =  quotRem n' d'
838                                            nd''       =  simplest' d' r' d r
839                                            n''        =  numerator nd''
840                                            d''        =  denominator nd''
841 \end{code}
842
843
844 \begin{code}
845 instance  (Integral a)  => Ord (Ratio a)  where
846     (x:%y) <= (x':%y')  =  x * y' <= x' * y
847     (x:%y) <  (x':%y')  =  x * y' <  x' * y
848
849 instance  (Integral a)  => Num (Ratio a)  where
850     (x:%y) + (x':%y')   =  reduce (x*y' + x'*y) (y*y')
851     (x:%y) - (x':%y')   =  reduce (x*y' - x'*y) (y*y')
852     (x:%y) * (x':%y')   =  reduce (x * x') (y * y')
853     negate (x:%y)       =  (-x) :% y
854     abs (x:%y)          =  abs x :% y
855     signum (x:%y)       =  signum x :% 1
856     fromInteger x       =  fromInteger x :% 1
857
858 instance  (Integral a)  => Real (Ratio a)  where
859     toRational (x:%y)   =  toInteger x :% toInteger y
860
861 instance  (Integral a)  => Fractional (Ratio a)  where
862     (x:%y) / (x':%y')   =  (x*y') % (y*x')
863     recip (x:%y)        =  if x < 0 then (-y) :% (-x) else y :% x
864     fromRational (x:%y) =  fromInteger x :% fromInteger y
865
866 instance  (Integral a)  => RealFrac (Ratio a)  where
867     properFraction (x:%y) = (fromIntegral q, r:%y)
868                             where (q,r) = quotRem x y
869
870 instance  (Integral a)  => Enum (Ratio a)  where
871     enumFrom            =  iterate ((+)1)
872     enumFromThen n m    =  iterate ((+)(m-n)) n
873     toEnum n            =  fromIntegral n :% 1
874     fromEnum            =  fromInteger . truncate
875
876 instance  (Integral a)  => Show (Ratio a)  where
877     showsPrec p (x:%y)  =  showParen (p > ratio_prec)
878                                (shows x . showString " % " . shows y)
879
880 -- defn. also used by the Read (Ratio a) instance PrelRead.
881 ratio_prec :: Int
882 ratio_prec = 7
883
884 \end{code}
885
886 \begin{code}
887 --Exported from std library Numeric, defined here to
888 --avoid mut. rec. between PrelNum and Numeric.
889 showSigned :: (Real a) => (a -> ShowS) -> Int -> a -> ShowS
890 showSigned showPos p x
891   | x < 0     = showParen (p > 6) (showChar '-' . showPos (-x))
892   | otherwise = showPos x
893
894 showSignedInteger :: Int -> Integer -> ShowS
895 showSignedInteger p n r
896   = -- from HBC version; support code follows
897     if n < 0 && p > 6 then '(':jtos n++(')':r) else jtos n ++ r
898
899 jtos :: Integer -> String
900 jtos n
901  | n < 0     = '-' : jtos' (-n) []
902  | otherwise = jtos' n []
903
904 jtos' :: Integer -> String -> String
905 jtos' n cs
906   | n < 10    = chr (fromInteger (n + ord_0)) : cs
907   | otherwise = jtos' q (chr (toInt r + (ord_0::Int)) : cs)
908   where
909     (q,r) = n `quotRem` 10
910
911 showFloat x  =  showString (formatRealFloat FFGeneric Nothing x)
912
913 -- These are the format types.  This type is not exported.
914
915 data FFFormat = FFExponent | FFFixed | FFGeneric --no need: deriving (Eq, Ord, Show)
916
917 formatRealFloat :: (RealFloat a) => FFFormat -> Maybe Int -> a -> String
918 formatRealFloat fmt decs x = s
919  where 
920   base = 10
921   base_i = toInteger base
922
923   s 
924    | isNaN x      = "NaN"
925    | isInfinite x = (\ str -> if x < 0 then '-':str else str) "Infinity"
926    | x < 0 || isNegativeZero x = '-' : doFmt fmt (floatToDigits base_i (-x))
927    | otherwise    = doFmt fmt (floatToDigits base_i x)
928
929   doFmt fmt (is, e) =
930     let ds = map intToDigit is in
931     case fmt of
932      FFGeneric ->
933       doFmt (if e <0 || e > 7 then FFExponent else FFFixed)
934             (is,e)
935      FFExponent ->
936       case decs of
937        Nothing ->
938         let e' = if e==0 then 0 else e-1 in
939         (case ds of
940           [d]    -> d : ".0e"
941           (d:ds) -> d : '.' : ds ++ "e") ++ show e'
942        Just dec ->
943         let dec' = max dec 1 in
944         case is of
945          [0] -> '0':'.':take dec' (repeat '0') ++ "e0"
946          _ ->
947           let
948            (ei,is') = roundTo base (dec'+1) is
949            d:ds = map intToDigit (if ei > 0 then init is' else is')
950           in
951           d:'.':ds ++ 'e':show (e-1+ei)
952      FFFixed ->
953       let
954        mk0 ls = case ls of { "" -> "0" ; _ -> ls}
955       in
956       case decs of
957        Nothing ->
958          let
959           f 0 s ds = mk0 (reverse s) ++ '.':mk0 ds
960           f n s "" = f (n-1) ('0':s) ""
961           f n s (d:ds) = f (n-1) (d:s) ds
962          in
963          f e "" ds
964        Just dec ->
965         let dec' = max dec 1 in
966         if e >= 0 then
967          let
968           (ei,is') = roundTo base (dec' + e) is
969           (ls,rs)  = splitAt (e+ei) (map intToDigit is')
970          in
971          mk0 ls ++ (if null rs then "" else '.':rs)
972         else
973          let
974           (ei,is') = roundTo base dec' (replicate (-e) 0 ++ is)
975           d:ds = map intToDigit (if ei > 0 then is' else 0:is')
976          in
977          d : '.' : ds
978          
979
980 roundTo :: Int -> Int -> [Int] -> (Int,[Int])
981 roundTo base d is =
982  let
983   v = f d is
984  in
985  case v of
986   (0,is) -> v
987   (1,is) -> (1, 1:is)
988  where
989   b2 = base `div` 2
990
991   f n [] = (0, replicate n 0)
992   f 0 (i:_) = (if i>=b2 then 1 else 0, [])
993   f d (i:is) =
994     let 
995      (c,ds) = f (d-1) is
996      i' = c + i
997     in
998     if i' == base then (1,0:ds) else (0,i':ds)
999
1000 --
1001 -- Based on "Printing Floating-Point Numbers Quickly and Accurately"
1002 -- by R.G. Burger and R.K. Dybvig in PLDI 96.
1003 -- This version uses a much slower logarithm estimator. It should be improved.
1004
1005 -- This function returns a list of digits (Ints in [0..base-1]) and an
1006 -- exponent.
1007 --floatToDigits :: (RealFloat a) => Integer -> a -> ([Int], Int)
1008 floatToDigits _ 0 = ([0], 0)
1009 floatToDigits base x =
1010  let 
1011   (f0, e0) = decodeFloat x
1012   (minExp0, _) = floatRange x
1013   p = floatDigits x
1014   b = floatRadix x
1015   minExp = minExp0 - p -- the real minimum exponent
1016   -- Haskell requires that f be adjusted so denormalized numbers
1017   -- will have an impossibly low exponent.  Adjust for this.
1018   (f, e) = 
1019    let n = minExp - e0 in
1020    if n > 0 then (f0 `div` (b^n), e0+n) else (f0, e0)
1021   (r, s, mUp, mDn) =
1022    if e >= 0 then
1023     let be = b^ e in
1024     if f == b^(p-1) then
1025       (f*be*b*2, 2*b, be*b, b)
1026     else
1027       (f*be*2, 2, be, be)
1028    else
1029     if e > minExp && f == b^(p-1) then
1030       (f*b*2, b^(-e+1)*2, b, 1)
1031     else
1032       (f*2, b^(-e)*2, 1, 1)
1033   k =
1034    let 
1035     k0 =
1036      if b == 2 && base == 10 then
1037         -- logBase 10 2 is slightly bigger than 3/10 so
1038         -- the following will err on the low side.  Ignoring
1039         -- the fraction will make it err even more.
1040         -- Haskell promises that p-1 <= logBase b f < p.
1041         (p - 1 + e0) * 3 `div` 10
1042      else
1043         ceiling ((log (fromInteger (f+1)) + fromInt e * log (fromInteger b)) /
1044                   log (fromInteger base))
1045
1046     fixup n =
1047       if n >= 0 then
1048         if r + mUp <= expt base n * s then n else fixup (n+1)
1049       else
1050         if expt base (-n) * (r + mUp) <= s then n else fixup (n+1)
1051    in
1052    fixup k0
1053
1054   gen ds rn sN mUpN mDnN =
1055    let
1056     (dn, rn') = (rn * base) `divMod` sN
1057     mUpN' = mUpN * base
1058     mDnN' = mDnN * base
1059    in
1060    case (rn' < mDnN', rn' + mUpN' > sN) of
1061     (True,  False) -> dn : ds
1062     (False, True)  -> dn+1 : ds
1063     (True,  True)  -> if rn' * 2 < sN then dn : ds else dn+1 : ds
1064     (False, False) -> gen (dn:ds) rn' sN mUpN' mDnN'
1065   
1066   rds = 
1067    if k >= 0 then
1068       gen [] r (s * expt base k) mUp mDn
1069    else
1070      let bk = expt base (-k) in
1071      gen [] (r * bk) s (mUp * bk) (mDn * bk)
1072  in
1073  (map toInt (reverse rds), k)
1074
1075 \end{code}
1076
1077 @showRational@ converts a Rational to a string that looks like a
1078 floating point number, but without converting to any floating type
1079 (because of the possible overflow).
1080
1081 From/by Lennart, 94/09/26
1082
1083 \begin{code}
1084 showRational :: Int -> Rational -> String
1085 showRational n r =
1086     if r == 0 then
1087         "0.0"
1088     else
1089         let (r', e) = normalize r
1090         in  prR n r' e
1091
1092 startExpExp = 4 :: Int
1093
1094 -- make sure 1 <= r < 10
1095 normalize :: Rational -> (Rational, Int)
1096 normalize r = if r < 1 then
1097                   case norm startExpExp (1 / r) 0 of (r', e) -> (10 / r', -e-1)
1098               else
1099                   norm startExpExp r 0
1100         where norm :: Int -> Rational -> Int -> (Rational, Int)
1101               -- Invariant: r*10^e == original r
1102               norm 0  r e = (r, e)
1103               norm ee r e =
1104                 let n = 10^ee
1105                     tn = 10^n
1106                 in  if r >= tn then norm ee (r/tn) (e+n) else norm (ee-1) r e
1107
1108 drop0 "" = ""
1109 drop0 (c:cs) = c : reverse (dropWhile (=='0') (reverse cs))
1110
1111 prR :: Int -> Rational -> Int -> String
1112 prR n r e | r <  1  = prR n (r*10) (e-1)                -- final adjustment
1113 prR n r e | r >= 10 = prR n (r/10) (e+1)
1114 prR n r e0 =
1115         let s = show ((round (r * 10^n))::Integer)
1116             e = e0+1
1117         in  if e > 0 && e < 8 then
1118                 take e s ++ "." ++ drop0 (drop e s)
1119             else if e <= 0 && e > -3 then
1120                 "0." ++ take (-e) (repeat '0') ++ drop0 s
1121             else
1122                 head s : "."++ drop0 (tail s) ++ "e" ++ show e0
1123 \end{code}
1124
1125
1126 [In response to a request for documentation of how fromRational works,
1127 Joe Fasel writes:] A quite reasonable request!  This code was added to
1128 the Prelude just before the 1.2 release, when Lennart, working with an
1129 early version of hbi, noticed that (read . show) was not the identity
1130 for floating-point numbers.  (There was a one-bit error about half the
1131 time.)  The original version of the conversion function was in fact
1132 simply a floating-point divide, as you suggest above. The new version
1133 is, I grant you, somewhat denser.
1134
1135 Unfortunately, Joe's code doesn't work!  Here's an example:
1136
1137 main = putStr (shows (1.82173691287639817263897126389712638972163e-300::Double) "\n")
1138
1139 This program prints
1140         0.0000000000000000
1141 instead of
1142         1.8217369128763981e-300
1143
1144 Lennart's code follows, and it works...
1145
1146 \begin{pseudocode}
1147 {-# GENERATE_SPECS fromRational__ a{Double#,Double} #-}
1148 fromRat :: (RealFloat a) => Rational -> a
1149 fromRat x = x'
1150         where x' = f e
1151
1152 --              If the exponent of the nearest floating-point number to x 
1153 --              is e, then the significand is the integer nearest xb^(-e),
1154 --              where b is the floating-point radix.  We start with a good
1155 --              guess for e, and if it is correct, the exponent of the
1156 --              floating-point number we construct will again be e.  If
1157 --              not, one more iteration is needed.
1158
1159               f e   = if e' == e then y else f e'
1160                       where y      = encodeFloat (round (x * (1 % b)^^e)) e
1161                             (_,e') = decodeFloat y
1162               b     = floatRadix x'
1163
1164 --              We obtain a trial exponent by doing a floating-point
1165 --              division of x's numerator by its denominator.  The
1166 --              result of this division may not itself be the ultimate
1167 --              result, because of an accumulation of three rounding
1168 --              errors.
1169
1170               (s,e) = decodeFloat (fromInteger (numerator x) `asTypeOf` x'
1171                                         / fromInteger (denominator x))
1172 \end{pseudocode}
1173
1174 Now, here's Lennart's code.
1175
1176 \begin{code}
1177 {-# SPECIALISE fromRat :: 
1178         Rational -> Double,
1179         Rational -> Float #-}
1180
1181 --fromRat :: (RealFloat a) => Rational -> a
1182 fromRat x
1183   | x == 0    = encodeFloat 0 0         -- Handle exceptional cases
1184   | x < 0     = - fromRat' (-x)         -- first.
1185   | otherwise = fromRat' x
1186
1187 -- Conversion process:
1188 -- Scale the rational number by the RealFloat base until
1189 -- it lies in the range of the mantissa (as used by decodeFloat/encodeFloat).
1190 -- Then round the rational to an Integer and encode it with the exponent
1191 -- that we got from the scaling.
1192 -- To speed up the scaling process we compute the log2 of the number to get
1193 -- a first guess of the exponent.
1194
1195 fromRat' :: (RealFloat a) => Rational -> a
1196 fromRat' x = r
1197   where b = floatRadix r
1198         p = floatDigits r
1199
1200         (minExp0, _) = floatRange r
1201
1202         minExp = minExp0 - p            -- the real minimum exponent
1203
1204         xMin = toRational (expt b (p-1))
1205         xMax = toRational (expt b p)
1206
1207         p0 = (integerLogBase b (numerator x) - integerLogBase b (denominator x) - p) `max` minExp
1208         f = if p0 < 0 then 1 % expt b (-p0) else expt b p0 % 1
1209         (x', p') = scaleRat (toRational b) minExp xMin xMax p0 (x / f)
1210
1211         r = encodeFloat (round x') p'
1212
1213 -- Scale x until xMin <= x < xMax, or p (the exponent) <= minExp.
1214 scaleRat :: Rational -> Int -> Rational -> Rational -> Int -> Rational -> (Rational, Int)
1215 scaleRat b minExp xMin xMax p x
1216     | p <= minExp = (x, p)
1217     | x >= xMax   = scaleRat b minExp xMin xMax (p+1) (x/b)
1218     | x < xMin    = scaleRat b minExp xMin xMax (p-1) (x*b)
1219     | otherwise   = (x, p)
1220
1221 -- Exponentiation with a cache for the most common numbers.
1222 minExpt = 0::Int
1223 maxExpt = 1100::Int
1224
1225 expt :: Integer -> Int -> Integer
1226 expt base n
1227  | base == 2 && n >= minExpt && n <= maxExpt = expts!n
1228  | otherwise                                 = base^n
1229
1230 expts :: Array Int Integer
1231 expts = array (minExpt,maxExpt) [(n,2^n) | n <- [minExpt .. maxExpt]]
1232
1233 -- Compute the (floor of the) log of i in base b.
1234 -- Simplest way would be just divide i by b until it's smaller then b, but that would
1235 -- be very slow!  We are just slightly more clever.
1236 integerLogBase :: Integer -> Integer -> Int
1237 integerLogBase b i
1238   | i < b     = 0
1239   | otherwise = doDiv (i `div` (b^l)) l
1240      where
1241         -- Try squaring the base first to cut down the number of divisions.
1242         l = 2 * integerLogBase (b*b) i
1243
1244         doDiv :: Integer -> Int -> Int
1245         doDiv i l = if i < b then l else doDiv (i `div` b) (l+1)
1246
1247 \end{code}
1248
1249 %*********************************************************
1250 %*                                                      *
1251 \subsection{Numeric primops}
1252 %*                                                      *
1253 %*********************************************************
1254
1255 Definitions of the boxed PrimOps; these will be
1256 used in the case of partial applications, etc.
1257
1258 \begin{code}
1259 plusFloat   (F# x) (F# y) = F# (plusFloat# x y)
1260 minusFloat  (F# x) (F# y) = F# (minusFloat# x y)
1261 timesFloat  (F# x) (F# y) = F# (timesFloat# x y)
1262 divideFloat (F# x) (F# y) = F# (divideFloat# x y)
1263 negateFloat (F# x)        = F# (negateFloat# x)
1264
1265 gtFloat     (F# x) (F# y) = gtFloat# x y
1266 geFloat     (F# x) (F# y) = geFloat# x y
1267 eqFloat     (F# x) (F# y) = eqFloat# x y
1268 neFloat     (F# x) (F# y) = neFloat# x y
1269 ltFloat     (F# x) (F# y) = ltFloat# x y
1270 leFloat     (F# x) (F# y) = leFloat# x y
1271
1272 float2Int   (F# x) = I# (float2Int# x)
1273 int2Float   (I# x) = F# (int2Float# x)
1274
1275 expFloat    (F# x) = F# (expFloat# x)
1276 logFloat    (F# x) = F# (logFloat# x)
1277 sqrtFloat   (F# x) = F# (sqrtFloat# x)
1278 sinFloat    (F# x) = F# (sinFloat# x)
1279 cosFloat    (F# x) = F# (cosFloat# x)
1280 tanFloat    (F# x) = F# (tanFloat# x)
1281 asinFloat   (F# x) = F# (asinFloat# x)
1282 acosFloat   (F# x) = F# (acosFloat# x)
1283 atanFloat   (F# x) = F# (atanFloat# x)
1284 sinhFloat   (F# x) = F# (sinhFloat# x)
1285 coshFloat   (F# x) = F# (coshFloat# x)
1286 tanhFloat   (F# x) = F# (tanhFloat# x)
1287
1288 powerFloat  (F# x) (F# y) = F# (powerFloat# x y)
1289
1290 -- definitions of the boxed PrimOps; these will be
1291 -- used in the case of partial applications, etc.
1292
1293 plusDouble   (D# x) (D# y) = D# (x +## y)
1294 minusDouble  (D# x) (D# y) = D# (x -## y)
1295 timesDouble  (D# x) (D# y) = D# (x *## y)
1296 divideDouble (D# x) (D# y) = D# (x /## y)
1297 negateDouble (D# x)        = D# (negateDouble# x)
1298
1299 gtDouble    (D# x) (D# y) = x >## y
1300 geDouble    (D# x) (D# y) = x >=## y
1301 eqDouble    (D# x) (D# y) = x ==## y
1302 neDouble    (D# x) (D# y) = x /=## y
1303 ltDouble    (D# x) (D# y) = x <## y
1304 leDouble    (D# x) (D# y) = x <=## y
1305
1306 double2Int   (D# x) = I# (double2Int#   x)
1307 int2Double   (I# x) = D# (int2Double#   x)
1308 double2Float (D# x) = F# (double2Float# x)
1309 float2Double (F# x) = D# (float2Double# x)
1310
1311 expDouble    (D# x) = D# (expDouble# x)
1312 logDouble    (D# x) = D# (logDouble# x)
1313 sqrtDouble   (D# x) = D# (sqrtDouble# x)
1314 sinDouble    (D# x) = D# (sinDouble# x)
1315 cosDouble    (D# x) = D# (cosDouble# x)
1316 tanDouble    (D# x) = D# (tanDouble# x)
1317 asinDouble   (D# x) = D# (asinDouble# x)
1318 acosDouble   (D# x) = D# (acosDouble# x)
1319 atanDouble   (D# x) = D# (atanDouble# x)
1320 sinhDouble   (D# x) = D# (sinhDouble# x)
1321 coshDouble   (D# x) = D# (coshDouble# x)
1322 tanhDouble   (D# x) = D# (tanhDouble# x)
1323
1324 powerDouble  (D# x) (D# y) = D# (x **## y)
1325 \end{code}