ghc/lib/ghc/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 -H20m #-}
  13 #include "../includes/ieee-flpt.h"
  14 \end{code}
  15
  16 \begin{code}
  17 {-# OPTIONS -fno-implicit-prelude #-}
  18
  19 module PrelNum where
  20
  21 import {-# SOURCE #-}   IOBase  ( error )
  22 import PrelList
  23 import PrelBase
  24 import GHC
  25
  26 infixr 8  ^, ^^, **
  27 infixl 7  %, `quot`, `rem`, `div`, `mod`
  28 \end{code}
  29
  30
  31 %*********************************************************
  32 %*                                                      *
  33 \subsection{Standard numeric classes}
  34 %*                                                      *
  35 %*********************************************************
  36
  37 \begin{code}
  38 class  (Num a, Ord a) => Real a  where
  39     toRational          ::  a -> Rational
  40
  41 class  (Real a, Enum a) => Integral a  where
  42     quot, rem, div, mod :: a -> a -> a
  43     quotRem, divMod     :: a -> a -> (a,a)
  44     toInteger           :: a -> Integer
  45     toInt               :: a -> Int -- partain: Glasgow extension
  46
  47     n `quot` d          =  q  where (q,r) = quotRem n d
  48     n `rem` d           =  r  where (q,r) = quotRem n d
  49     n `div` d           =  q  where (q,r) = divMod n d
  50     n `mod` d           =  r  where (q,r) = divMod n d
  51     divMod n d          =  if signum r == negate (signum d) then (q-1, r+d) else qr
  52                            where qr@(q,r) = quotRem n d
  53
  54 class  (Num a) => Fractional a  where
  55     (/)                 :: a -> a -> a
  56     recip               :: a -> a
  57     fromRational        :: Rational -> a
  58
  59     recip x             =  1 / x
  60
  61 class  (Fractional a) => Floating a  where
  62     pi                  :: a
  63     exp, log, sqrt      :: a -> a
  64     (**), logBase       :: a -> a -> a
  65     sin, cos, tan       :: a -> a
  66     asin, acos, atan    :: a -> a
  67     sinh, cosh, tanh    :: a -> a
  68     asinh, acosh, atanh :: a -> a
  69
  70     x ** y              =  exp (log x * y)
  71     logBase x y         =  log y / log x
  72     sqrt x              =  x ** 0.5
  73     tan  x              =  sin  x / cos  x
  74     tanh x              =  sinh x / cosh x
  75
  76 class  (Real a, Fractional a) => RealFrac a  where
  77     properFraction      :: (Integral b) => a -> (b,a)
  78     truncate, round     :: (Integral b) => a -> b
  79     ceiling, floor      :: (Integral b) => a -> b
  80
  81     truncate x          =  m  where (m,_) = properFraction x
  82
  83     round x             =  let (n,r) = properFraction x
  84                                m     = if r < 0 then n - 1 else n + 1
  85                            in case signum (abs r - 0.5) of
  86                                 -1 -> n
  87                                 0  -> if even n then n else m
  88                                 1  -> m
  89
  90     ceiling x           =  if r > 0 then n + 1 else n
  91                            where (n,r) = properFraction x
  92
  93     floor x             =  if r < 0 then n - 1 else n
  94                            where (n,r) = properFraction x
  95
  96 class  (RealFrac a, Floating a) => RealFloat a  where
  97     floatRadix          :: a -> Integer
  98     floatDigits         :: a -> Int
  99     floatRange          :: a -> (Int,Int)
 100     decodeFloat         :: a -> (Integer,Int)
 101     encodeFloat         :: Integer -> Int -> a
 102     exponent            :: a -> Int
 103     significand         :: a -> a
 104     scaleFloat          :: Int -> a -> a
 105     isNaN, isInfinite, isDenormalized, isNegativeZero, isIEEE
 106                         :: a -> Bool
 107
 108     exponent x          =  if m == 0 then 0 else n + floatDigits x
 109                            where (m,n) = decodeFloat x
 110
 111     significand x       =  encodeFloat m (negate (floatDigits x))
 112                            where (m,_) = decodeFloat x
 113
 114     scaleFloat k x      =  encodeFloat m (n+k)
 115                            where (m,n) = decodeFloat x
 116 \end{code}
 117
 118 %*********************************************************
 119 %*                                                      *
 120 \subsection{Overloaded numeric functions}
 121 %*                                                      *
 122 %*********************************************************
 123
 124 \begin{code}
 125 even, odd       :: (Integral a) => a -> Bool
 126 even n          =  n `rem` 2 == 0
 127 odd             =  not . even
 128
 129 {-# GENERATE_SPECS gcd a{Int#,Int,Integer} #-}
 130 gcd             :: (Integral a) => a -> a -> a
 131 gcd 0 0         =  error "Prelude.gcd: gcd 0 0 is undefined"
 132 gcd x y         =  gcd' (abs x) (abs y)
 133                    where gcd' x 0  =  x
 134                          gcd' x y  =  gcd' y (x `rem` y)
 135
 136 {-# GENERATE_SPECS lcm a{Int#,Int,Integer} #-}
 137 lcm             :: (Integral a) => a -> a -> a
 138 lcm _ 0         =  0
 139 lcm 0 _         =  0
 140 lcm x y         =  abs ((x `quot` (gcd x y)) * y)
 141
 142 (^)             :: (Num a, Integral b) => a -> b -> a
 143 x ^ 0           =  1
 144 x ^ n | n > 0   =  f x (n-1) x
 145                    where f _ 0 y = y
 146                          f x n y = g x n  where
 147                                    g x n | even n  = g (x*x) (n `quot` 2)
 148                                          | otherwise = f x (n-1) (x*y)
 149 _ ^ _           = error "Prelude.^: negative exponent"
 150
 151 (^^)            :: (Fractional a, Integral b) => a -> b -> a
 152 x ^^ n          =  if n >= 0 then x^n else recip (x^(negate n))
 153
 154 fromIntegral    :: (Integral a, Num b) => a -> b
 155 fromIntegral    =  fromInteger . toInteger
 156
 157 fromRealFrac    :: (RealFrac a, Fractional b) => a -> b
 158 fromRealFrac    =  fromRational . toRational
 159
 160 atan2           :: (RealFloat a) => a -> a -> a
 161 atan2 y x       =  case (signum y, signum x) of
 162                         ( 0, 1) ->  0
 163                         ( 1, 0) ->  pi/2
 164                         ( 0,-1) ->  pi
 165                         (-1, 0) ->  (negate pi)/2
 166                         ( _, 1) ->  atan (y/x)
 167                         ( _,-1) ->  atan (y/x) + pi
 168                         ( 0, 0) ->  error "Prelude.atan2: atan2 of origin"
 169 \end{code}
 170
 171
 172 %*********************************************************
 173 %*                                                      *
 174 \subsection{Instances for @Int@}
 175 %*                                                      *
 176 %*********************************************************
 177
 178 \begin{code}
 179 instance  Real Int  where
 180     toRational x        =  toInteger x % 1
 181
 182 instance  Integral Int  where
 183     a@(I# _) `quotRem` b@(I# _) = (a `quotInt` b, a `remInt` b)
 184     -- OK, so I made it a little stricter.  Shoot me.  (WDP 94/10)
 185
 186     -- following chks for zero divisor are non-standard (WDP)
 187     a `quot` b          =  if b /= 0
 188                            then a `quotInt` b
 189                            else error "Integral.Int.quot{PreludeCore}: divide by 0\n"
 190     a `rem` b           =  if b /= 0
 191                            then a `remInt` b
 192                            else error "Integral.Int.rem{PreludeCore}: divide by 0\n"
 193
 194     x `div` y = if x > 0 && y < 0       then quotInt (x-y-1) y
 195                 else if x < 0 && y > 0  then quotInt (x-y+1) y
 196                 else quotInt x y
 197     x `mod` y = if x > 0 && y < 0 || x < 0 && y > 0 then
 198                     if r/=0 then r+y else 0
 199                 else
 200                     r
 201               where r = remInt x y
 202
 203     divMod x@(I# _) y@(I# _) = (x `div` y, x `mod` y)
 204     -- Stricter.  Sorry if you don't like it.  (WDP 94/10)
 205
 206 --OLD:   even x = eqInt (x `mod` 2) 0
 207 --OLD:   odd x  = neInt (x `mod` 2) 0
 208
 209     toInteger (I# n#) = int2Integer# n#  -- give back a full-blown Integer
 210     toInt x           = x
 211
 212 \end{code}
 213
 214
 215 %*********************************************************
 216 %*                                                      *
 217 \subsection{Type @Integer@}
 218 %*                                                      *
 219 %*********************************************************
 220
 221 These types are used to return from integer primops
 222
 223 \begin{code}
 224 data Return2GMPs     = Return2GMPs     Int# Int# ByteArray# Int# Int# ByteArray#
 225 data ReturnIntAndGMP = ReturnIntAndGMP Int# Int# Int# ByteArray#
 226 \end{code}
 227
 228 Instances
 229
 230 \begin{code}
 231 instance  Eq Integer  where
 232     (J# a1 s1 d1) == (J# a2 s2 d2)
 233       = (cmpInteger# a1 s1 d1 a2 s2 d2) ==# 0#
 234
 235     (J# a1 s1 d1) /= (J# a2 s2 d2)
 236       = (cmpInteger# a1 s1 d1 a2 s2 d2) /=# 0#
 237
 238 instance  Ord Integer  where
 239     (J# a1 s1 d1) <= (J# a2 s2 d2)
 240       = (cmpInteger# a1 s1 d1 a2 s2 d2) <=# 0#
 241
 242     (J# a1 s1 d1) <  (J# a2 s2 d2)
 243       = (cmpInteger# a1 s1 d1 a2 s2 d2) <# 0#
 244
 245     (J# a1 s1 d1) >= (J# a2 s2 d2)
 246       = (cmpInteger# a1 s1 d1 a2 s2 d2) >=# 0#
 247
 248     (J# a1 s1 d1) >  (J# a2 s2 d2)
 249       = (cmpInteger# a1 s1 d1 a2 s2 d2) ># 0#
 250
 251     x@(J# a1 s1 d1) `max` y@(J# a2 s2 d2)
 252       = if ((cmpInteger# a1 s1 d1 a2 s2 d2) ># 0#) then x else y
 253
 254     x@(J# a1 s1 d1) `min` y@(J# a2 s2 d2)
 255       = if ((cmpInteger# a1 s1 d1 a2 s2 d2) <# 0#) then x else y
 256
 257     compare (J# a1 s1 d1) (J# a2 s2 d2)
 258        = case cmpInteger# a1 s1 d1 a2 s2 d2 of { res# ->
 259          if res# <# 0# then LT else
 260          if res# ># 0# then GT else EQ
 261          }
 262
 263 instance  Num Integer  where
 264     (+) (J# a1 s1 d1) (J# a2 s2 d2)
 265       = plusInteger# a1 s1 d1 a2 s2 d2
 266
 267     (-) (J# a1 s1 d1) (J# a2 s2 d2)
 268       = minusInteger# a1 s1 d1 a2 s2 d2
 269
 270     negate (J# a s d) = negateInteger# a s d
 271
 272     (*) (J# a1 s1 d1) (J# a2 s2 d2)
 273       = timesInteger# a1 s1 d1 a2 s2 d2
 274
 275     -- ORIG: abs n = if n >= 0 then n else -n
 276
 277     abs n@(J# a1 s1 d1)
 278       = case 0 of { J# a2 s2 d2 ->
 279         if (cmpInteger# a1 s1 d1 a2 s2 d2) >=# 0#
 280         then n
 281         else negateInteger# a1 s1 d1
 282         }
 283
 284     signum n@(J# a1 s1 d1)
 285       = case 0 of { J# a2 s2 d2 ->
 286         let
 287             cmp = cmpInteger# a1 s1 d1 a2 s2 d2
 288         in
 289         if      cmp >#  0# then 1
 290         else if cmp ==# 0# then 0
 291         else                    (negate 1)
 292         }
 293
 294     fromInteger x       =  x
 295
 296     fromInt (I# n#)     =  int2Integer# n# -- gives back a full-blown Integer
 297
 298 instance  Real Integer  where
 299     toRational x        =  x % 1
 300
 301 instance  Integral Integer where
 302     quotRem (J# a1 s1 d1) (J# a2 s2 d2)
 303       = case (quotRemInteger# a1 s1 d1 a2 s2 d2) of
 304           Return2GMPs a3 s3 d3 a4 s4 d4
 305             -> (J# a3 s3 d3, J# a4 s4 d4)
 306
 307 {- USING THE UNDERLYING "GMP" CODE IS DUBIOUS FOR NOW:
 308
 309     divMod (J# a1 s1 d1) (J# a2 s2 d2)
 310       = case (divModInteger# a1 s1 d1 a2 s2 d2) of
 311           Return2GMPs a3 s3 d3 a4 s4 d4
 312             -> (J# a3 s3 d3, J# a4 s4 d4)
 313 -}
 314     toInteger n      = n
 315     toInt (J# a s d) = case (integer2Int# a s d) of { n# -> I# n# }
 316
 317     -- the rest are identical to the report default methods;
 318     -- you get slightly better code if you let the compiler
 319     -- see them right here:
 320     n `quot` d  =  q  where (q,r) = quotRem n d
 321     n `rem` d   =  r  where (q,r) = quotRem n d
 322     n `div` d   =  q  where (q,r) = divMod n d
 323     n `mod` d   =  r  where (q,r) = divMod n d
 324
 325     divMod n d  =  case (quotRem n d) of { qr@(q,r) ->
 326                    if signum r == negate (signum d) then (q - 1, r+d) else qr }
 327                    -- Case-ified by WDP 94/10
 328
 329 instance  Enum Integer  where
 330     enumFrom n           =  n : enumFrom (n + 1)
 331     enumFromThen m n     =  en' m (n - m)
 332                             where en' m n = m : en' (m + n) n
 333     enumFromTo n m       =  takeWhile (<= m) (enumFrom n)
 334     enumFromThenTo n m p =  takeWhile (if m >= n then (<= p) else (>= p))
 335                                       (enumFromThen n m)
 336
 337 instance  Show Integer  where
 338     showsPrec   x = showSignedInteger x
 339     showList = showList__ (showsPrec 0)
 340
 341 integer_0, integer_1, integer_2, integer_m1 :: Integer
 342 integer_0 = 0; integer_1 = 1; integer_2 = 2; integer_m1 = -1
 343 \end{code}
 344
 345
 346 %*********************************************************
 347 %*                                                      *
 348 \subsection{Type @Float@}
 349 %*                                                      *
 350 %*********************************************************
 351
 352 \begin{code}
 353 instance  Num Float  where
 354     (+)         x y     =  plusFloat x y
 355     (-)         x y     =  minusFloat x y
 356     negate      x       =  negateFloat x
 357     (*)         x y     =  timesFloat x y
 358     abs x | x >= 0.0    =  x
 359           | otherwise   =  negateFloat x
 360     signum x | x == 0.0  = 0
 361              | x > 0.0   = 1
 362              | otherwise = negate 1
 363     fromInteger n       =  encodeFloat n 0
 364     fromInt i           =  int2Float i
 365
 366 instance  Real Float  where
 367     toRational x        =  (m%1)*(b%1)^^n
 368                            where (m,n) = decodeFloat x
 369                                  b     = floatRadix  x
 370
 371 instance  Fractional Float  where
 372     (/) x y             =  divideFloat x y
 373     fromRational x      =  fromRational__ x
 374     recip x             =  1.0 / x
 375
 376 instance  Floating Float  where
 377     pi                  =  3.141592653589793238
 378     exp x               =  expFloat x
 379     log x               =  logFloat x
 380     sqrt x              =  sqrtFloat x
 381     sin x               =  sinFloat x
 382     cos x               =  cosFloat x
 383     tan x               =  tanFloat x
 384     asin x              =  asinFloat x
 385     acos x              =  acosFloat x
 386     atan x              =  atanFloat x
 387     sinh x              =  sinhFloat x
 388     cosh x              =  coshFloat x
 389     tanh x              =  tanhFloat x
 390     (**) x y            =  powerFloat x y
 391     logBase x y         =  log y / log x
 392
 393     asinh x = log (x + sqrt (1.0+x*x))
 394     acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))
 395     atanh x = log ((x+1.0) / sqrt (1.0-x*x))
 396
 397 instance  RealFrac Float  where
 398
 399     {-# SPECIALIZE properFraction :: Float -> (Int, Float) #-}
 400     {-# SPECIALIZE truncate :: Float -> Int #-}
 401     {-# SPECIALIZE round    :: Float -> Int #-}
 402     {-# SPECIALIZE ceiling  :: Float -> Int #-}
 403     {-# SPECIALIZE floor    :: Float -> Int #-}
 404
 405     {-# SPECIALIZE properFraction :: Float -> (Integer, Float) #-}
 406     {-# SPECIALIZE truncate :: Float -> Integer #-}
 407     {-# SPECIALIZE round    :: Float -> Integer #-}
 408     {-# SPECIALIZE ceiling  :: Float -> Integer #-}
 409     {-# SPECIALIZE floor    :: Float -> Integer #-}
 410
 411     properFraction x
 412       = case (decodeFloat x)      of { (m,n) ->
 413         let  b = floatRadix x     in
 414         if n >= 0 then
 415             (fromInteger m * fromInteger b ^ n, 0.0)
 416         else
 417             case (quotRem m (b^(negate n))) of { (w,r) ->
 418             (fromInteger w, encodeFloat r n)
 419             }
 420         }
 421
 422     truncate x  = case properFraction x of
 423                      (n,_) -> n
 424
 425     round x     = case properFraction x of
 426                      (n,r) -> let
 427                                 m         = if r < 0.0 then n - 1 else n + 1
 428                                 half_down = abs r - 0.5
 429                               in
 430                               case (compare half_down 0.0) of
 431                                 LT -> n
 432                                 EQ -> if even n then n else m
 433                                 GT -> m
 434
 435     ceiling x   = case properFraction x of
 436                     (n,r) -> if r > 0.0 then n + 1 else n
 437
 438     floor x     = case properFraction x of
 439                     (n,r) -> if r < 0.0 then n - 1 else n
 440
 441 instance  RealFloat Float  where
 442     floatRadix _        =  FLT_RADIX        -- from float.h
 443     floatDigits _       =  FLT_MANT_DIG     -- ditto
 444     floatRange _        =  (FLT_MIN_EXP, FLT_MAX_EXP) -- ditto
 445
 446     decodeFloat (F# f#)
 447       = case decodeFloat# f#    of
 448           ReturnIntAndGMP exp# a# s# d# ->
 449             (J# a# s# d#, I# exp#)
 450
 451     encodeFloat (J# a# s# d#) (I# e#)
 452       = case encodeFloat# a# s# d# e# of { flt# -> F# flt# }
 453
 454     exponent x          = case decodeFloat x of
 455                             (m,n) -> if m == 0 then 0 else n + floatDigits x
 456
 457     significand x       = case decodeFloat x of
 458                             (m,_) -> encodeFloat m (negate (floatDigits x))
 459
 460     scaleFloat k x      = case decodeFloat x of
 461                             (m,n) -> encodeFloat m (n+k)
 462
 463 instance  Show Float  where
 464     showsPrec   x = showSigned showFloat x
 465     showList = showList__ (showsPrec 0)
 466 \end{code}
 467
 468 %*********************************************************
 469 %*                                                      *
 470 \subsection{Type @Double@}
 471 %*                                                      *
 472 %*********************************************************
 473
 474 \begin{code}
 475 instance  Num Double  where
 476     (+)         x y     =  plusDouble x y
 477     (-)         x y     =  minusDouble x y
 478     negate      x       =  negateDouble x
 479     (*)         x y     =  timesDouble x y
 480     abs x | x >= 0.0    =  x
 481           | otherwise   =  negateDouble x
 482     signum x | x == 0.0  = 0
 483              | x > 0.0   = 1
 484              | otherwise = negate 1
 485     fromInteger n       =  encodeFloat n 0
 486     fromInt (I# n#)     =  case (int2Double# n#) of { d# -> D# d# }
 487
 488 instance  Real Double  where
 489     toRational x        =  (m%1)*(b%1)^^n
 490                            where (m,n) = decodeFloat x
 491                                  b     = floatRadix  x
 492
 493 instance  Fractional Double  where
 494     (/) x y             =  divideDouble x y
 495     fromRational x      =  fromRational__ x
 496     recip x             =  1.0 / x
 497
 498 instance  Floating Double  where
 499     pi                  =  3.141592653589793238
 500     exp x               =  expDouble x
 501     log x               =  logDouble x
 502     sqrt x              =  sqrtDouble x
 503     sin  x              =  sinDouble x
 504     cos  x              =  cosDouble x
 505     tan  x              =  tanDouble x
 506     asin x              =  asinDouble x
 507     acos x              =  acosDouble x
 508     atan x              =  atanDouble x
 509     sinh x              =  sinhDouble x
 510     cosh x              =  coshDouble x
 511     tanh x              =  tanhDouble x
 512     (**) x y            =  powerDouble x y
 513     logBase x y         =  log y / log x
 514
 515     asinh x = log (x + sqrt (1.0+x*x))
 516     acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))
 517     atanh x = log ((x+1.0) / sqrt (1.0-x*x))
 518
 519 instance  RealFrac Double  where
 520
 521     {-# SPECIALIZE properFraction :: Double -> (Int, Double) #-}
 522     {-# SPECIALIZE truncate :: Double -> Int #-}
 523     {-# SPECIALIZE round    :: Double -> Int #-}
 524     {-# SPECIALIZE ceiling  :: Double -> Int #-}
 525     {-# SPECIALIZE floor    :: Double -> Int #-}
 526
 527     {-# SPECIALIZE properFraction :: Double -> (Integer, Double) #-}
 528     {-# SPECIALIZE truncate :: Double -> Integer #-}
 529     {-# SPECIALIZE round    :: Double -> Integer #-}
 530     {-# SPECIALIZE ceiling  :: Double -> Integer #-}
 531     {-# SPECIALIZE floor    :: Double -> Integer #-}
 532
 533 #if defined(__UNBOXED_INSTANCES__)
 534     {-# SPECIALIZE properFraction :: Double -> (Int#, Double) #-}
 535     {-# SPECIALIZE truncate :: Double -> Int# #-}
 536     {-# SPECIALIZE round    :: Double -> Int# #-}
 537     {-# SPECIALIZE ceiling  :: Double -> Int# #-}
 538     {-# SPECIALIZE floor    :: Double -> Int# #-}
 539 #endif
 540
 541     properFraction x
 542       = case (decodeFloat x)      of { (m,n) ->
 543         let  b = floatRadix x     in
 544         if n >= 0 then
 545             (fromInteger m * fromInteger b ^ n, 0.0)
 546         else
 547             case (quotRem m (b^(negate n))) of { (w,r) ->
 548             (fromInteger w, encodeFloat r n)
 549             }
 550         }
 551
 552     truncate x  = case properFraction x of
 553                      (n,_) -> n
 554
 555     round x     = case properFraction x of
 556                      (n,r) -> let
 557                                 m         = if r < 0.0 then n - 1 else n + 1
 558                                 half_down = abs r - 0.5
 559                               in
 560                               case (compare half_down 0.0) of
 561                                 LT -> n
 562                                 EQ -> if even n then n else m
 563                                 GT -> m
 564
 565     ceiling x   = case properFraction x of
 566                     (n,r) -> if r > 0.0 then n + 1 else n
 567
 568     floor x     = case properFraction x of
 569                     (n,r) -> if r < 0.0 then n - 1 else n
 570
 571 instance  RealFloat Double  where
 572     floatRadix _        =  FLT_RADIX        -- from float.h
 573     floatDigits _       =  DBL_MANT_DIG     -- ditto
 574     floatRange _        =  (DBL_MIN_EXP, DBL_MAX_EXP) -- ditto
 575
 576     decodeFloat (D# d#)
 577       = case decodeDouble# d#   of
 578           ReturnIntAndGMP exp# a# s# d# ->
 579             (J# a# s# d#, I# exp#)
 580
 581     encodeFloat (J# a# s# d#) (I# e#)
 582       = case encodeDouble# a# s# d# e#  of { dbl# -> D# dbl# }
 583
 584     exponent x          = case decodeFloat x of
 585                             (m,n) -> if m == 0 then 0 else n + floatDigits x
 586
 587     significand x       = case decodeFloat x of
 588                             (m,_) -> encodeFloat m (negate (floatDigits x))
 589
 590     scaleFloat k x      = case decodeFloat x of
 591                             (m,n) -> encodeFloat m (n+k)
 592
 593 instance  Show Double  where
 594     showsPrec   x = showSigned showFloat x
 595     showList = showList__ (showsPrec 0)
 596 \end{code}
 597
 598
 599 %*********************************************************
 600 %*                                                      *
 601 \subsection{Common code for @Float@ and @Double@}
 602 %*                                                      *
 603 %*********************************************************
 604
 605 The Enum instances for Floats and Doubles are slightly unusual.
 606 The `toEnum' function truncates numbers to Int.  The definitions
 607 of enumFrom and enumFromThen allow floats to be used in arithmetic
 608 series: [0,0.1 .. 1.0].  However, roundoff errors make these somewhat
 609 dubious.  This example may have either 10 or 11 elements, depending on
 610 how 0.1 is represented.
 611
 612 \begin{code}
 613 instance  Enum Float  where
 614     toEnum              =  fromIntegral
 615     fromEnum            =  fromInteger . truncate   -- may overflow
 616     enumFrom            =  numericEnumFrom
 617     enumFromThen        =  numericEnumFromThen
 618
 619 instance  Enum Double  where
 620     toEnum              =  fromIntegral
 621     fromEnum            =  fromInteger . truncate   -- may overflow
 622     enumFrom            =  numericEnumFrom
 623     enumFromThen        =  numericEnumFromThen
 624
 625 numericEnumFrom         :: (Real a) => a -> [a]
 626 numericEnumFromThen     :: (Real a) => a -> a -> [a]
 627 numericEnumFrom         =  iterate (+1)
 628 numericEnumFromThen n m =  iterate (+(m-n)) n
 629 \end{code}
 630
 631
 632 %*********************************************************
 633 %*                                                      *
 634 \subsection{The @Ratio@ and @Rational@ types}
 635 %*                                                      *
 636 %*********************************************************
 637
 638 \begin{code}
 639 data  (Integral a)      => Ratio a = a :% a  deriving (Eq)
 640 type  Rational          =  Ratio Integer
 641 \end{code}
 642
 643 \begin{code}
 644 (%)                     :: (Integral a) => a -> a -> Ratio a
 645 numerator, denominator  :: (Integral a) => Ratio a -> a
 646 approxRational          :: (RealFrac a) => a -> a -> Rational
 647
 648
 649 reduce _ 0              =  error "{Ratio.%}: zero denominator"
 650 reduce x y              =  (x `quot` d) :% (y `quot` d)
 651                            where d = gcd x y
 652
 653 x % y                   =  reduce (x * signum y) (abs y)
 654
 655 numerator (x:%y)        =  x
 656
 657 denominator (x:%y)      =  y
 658 \end{code}
 659
 660
 661 @approxRational@, applied to two real fractional numbers x and epsilon,
 662 returns the simplest rational number within epsilon of x.  A rational
 663 number n%d in reduced form is said to be simpler than another n'%d' if
 664 abs n <= abs n' && d <= d'.  Any real interval contains a unique
 665 simplest rational; here, for simplicity, we assume a closed rational
 666 interval.  If such an interval includes at least one whole number, then
 667 the simplest rational is the absolutely least whole number.  Otherwise,
 668 the bounds are of the form q%1 + r%d and q%1 + r'%d', where abs r < d
 669 and abs r' < d', and the simplest rational is q%1 + the reciprocal of
 670 the simplest rational between d'%r' and d%r.
 671
 672 \begin{code}
 673 approxRational x eps    =  simplest (x-eps) (x+eps)
 674         where simplest x y | y < x      =  simplest y x
 675                            | x == y     =  xr
 676                            | x > 0      =  simplest' n d n' d'
 677                            | y < 0      =  - simplest' (-n') d' (-n) d
 678                            | otherwise  =  0 :% 1
 679                                         where xr@(n:%d) = toRational x
 680                                               (n':%d')  = toRational y
 681
 682               simplest' n d n' d'       -- assumes 0 < n%d < n'%d'
 683                         | r == 0     =  q :% 1
 684                         | q /= q'    =  (q+1) :% 1
 685                         | otherwise  =  (q*n''+d'') :% n''
 686                                      where (q,r)      =  quotRem n d
 687                                            (q',r')    =  quotRem n' d'
 688                                            (n'':%d'') =  simplest' d' r' d r
 689 \end{code}
 690
 691
 692 \begin{code}
 693 instance  (Integral a)  => Ord (Ratio a)  where
 694     (x:%y) <= (x':%y')  =  x * y' <= x' * y
 695     (x:%y) <  (x':%y')  =  x * y' <  x' * y
 696
 697 instance  (Integral a)  => Num (Ratio a)  where
 698     (x:%y) + (x':%y')   =  reduce (x*y' + x'*y) (y*y')
 699     (x:%y) * (x':%y')   =  reduce (x * x') (y * y')
 700     negate (x:%y)       =  (-x) :% y
 701     abs (x:%y)          =  abs x :% y
 702     signum (x:%y)       =  signum x :% 1
 703     fromInteger x       =  fromInteger x :% 1
 704
 705 instance  (Integral a)  => Real (Ratio a)  where
 706     toRational (x:%y)   =  toInteger x :% toInteger y
 707
 708 instance  (Integral a)  => Fractional (Ratio a)  where
 709     (x:%y) / (x':%y')   =  (x*y') % (y*x')
 710     recip (x:%y)        =  if x < 0 then (-y) :% (-x) else y :% x
 711     fromRational (x:%y) =  fromInteger x :% fromInteger y
 712
 713 instance  (Integral a)  => RealFrac (Ratio a)  where
 714     properFraction (x:%y) = (fromIntegral q, r:%y)
 715                             where (q,r) = quotRem x y
 716
 717 instance  (Integral a)  => Enum (Ratio a)  where
 718     enumFrom            =  iterate ((+)1)
 719     enumFromThen n m    =  iterate ((+)(m-n)) n
 720     toEnum n            =  fromIntegral n :% 1
 721     fromEnum            =  fromInteger . truncate
 722
 723 ratio_prec :: Int
 724 ratio_prec = 7
 725
 726 instance  (Integral a)  => Show (Ratio a)  where
 727     showsPrec p (x:%y)  =  showParen (p > ratio_prec)
 728                                (shows x . showString " % " . shows y)
 729 \end{code}
 730
 731 {-
 732 [In response to a request by simonpj, Joe Fasel writes:]
 733
 734 A quite reasonable request!  This code was added to the Prelude just
 735 before the 1.2 release, when Lennart, working with an early version
 736 of hbi, noticed that (read . show) was not the identity for
 737 floating-point numbers.  (There was a one-bit error about half the time.)
 738 The original version of the conversion function was in fact simply
 739 a floating-point divide, as you suggest above.  The new version is,
 740 I grant you, somewhat denser.
 741
 742 How's this?
 743
 744 Joe
 745 -}
 746
 747 \begin{code}
 748 {-# GENERATE_SPECS fromRational__ a{Double#,Double} #-}
 749 fromRational__ :: (RealFloat a) => Rational -> a
 750 fromRational__ x = x'
 751         where x' = f e
 752
 753 --              If the exponent of the nearest floating-point number to x
 754 --              is e, then the significand is the integer nearest xb^(-e),
 755 --              where b is the floating-point radix.  We start with a good
 756 --              guess for e, and if it is correct, the exponent of the
 757 --              floating-point number we construct will again be e.  If
 758 --              not, one more iteration is needed.
 759
 760               f e   = if e' == e then y else f e'
 761                       where y      = encodeFloat (round (x * (1 % b)^^e)) e
 762                             (_,e') = decodeFloat y
 763               b     = floatRadix x'
 764
 765 --              We obtain a trial exponent by doing a floating-point
 766 --              division of x's numerator by its denominator.  The
 767 --              result of this division may not itself be the ultimate
 768 --              result, because of an accumulation of three rounding
 769 --              errors.
 770
 771               (s,e) = decodeFloat (fromInteger (numerator x) `asTypeOf` x'
 772                                         / fromInteger (denominator x))
 773 \end{code}
 774
 775
 776 %*********************************************************
 777 %*                                                      *
 778 \subsection{Showing numbers}
 779 %*                                                      *
 780 %*********************************************************
 781
 782 \begin{code}
 783 showInteger n r
 784   = case quotRem n 10 of                     { (n', d) ->
 785     case (chr (ord_0 + fromIntegral d)) of { C# c# -> -- stricter than necessary
 786     let
 787         r' = C# c# : r
 788     in
 789     if n' == 0 then r' else showInteger n' r'
 790     }}
 791
 792 showSigned :: (Real a) => (a -> ShowS) -> Int -> a -> ShowS
 793 showSigned showPos p x = if x < 0 then showParen (p > 6)
 794                                                  (showChar '-' . showPos (-x))
 795                                   else showPos x
 796
 797 showSignedInteger :: Int -> Integer -> ShowS
 798 showSignedInteger p n r
 799   = -- from HBC version; support code follows
 800     if n < 0 && p > 6 then '(':jtos n++(')':r) else jtos n ++ r
 801
 802 jtos :: Integer -> String
 803 jtos n
 804   = if n < 0 then
 805         '-' : jtos' (-n) []
 806     else
 807         jtos' n []
 808
 809 jtos' :: Integer -> String -> String
 810 jtos' n cs
 811   = if n < 10 then
 812         chr (fromInteger (n + ord_0)) : cs
 813     else
 814         jtos' (n `quot` 10) (chr (fromInteger (n `rem` 10 + ord_0)) : cs)
 815 \end{code}
 816
 817 The functions showFloat below uses rational arithmetic
 818 to insure correct conversion between the floating-point radix and
 819 decimal.  It is often possible to use a higher-precision floating-
 820 point type to obtain the same results.
 821
 822 \begin{code}
 823 {-# GENERATE_SPECS showFloat a{Double#,Double} #-}
 824 showFloat:: (RealFloat a) => a -> ShowS
 825 showFloat x =
 826     if x == 0 then showString ("0." ++ take (m-1) zeros)
 827               else if e >= m-1 || e < 0 then showSci else showFix
 828     where
 829     showFix     = showString whole . showChar '.' . showString frac
 830                   where (whole,frac) = splitAt (e+1) (show sig)
 831     showSci     = showChar d . showChar '.' . showString frac
 832                       . showChar 'e' . shows e
 833                   where (d:frac) = show sig
 834     (m, sig, e) = if b == 10 then (w,   s,   n+w-1)
 835                              else (m', sig', e'   )
 836     m'          = ceiling
 837                       ((fromInt w * log (fromInteger b)) / log 10 :: Double)
 838                   + 1
 839     (sig', e')  = if      sig1 >= 10^m'     then (round (t/10), e1+1)
 840                   else if sig1 <  10^(m'-1) then (round (t*10), e1-1)
 841                                             else (sig1,          e1  )
 842     sig1        = round t
 843     t           = s%1 * (b%1)^^n * 10^^(m'-e1-1)
 844     e1          = floor (logBase 10 x)
 845     (s, n)      = decodeFloat x
 846     b           = floatRadix x
 847     w           = floatDigits x
 848
 849 zeros = repeat '0'
 850 \end{code}
 851
 852 @showRational@ converts a Rational to a string that looks like a
 853 floating point number, but without converting to any floating type
 854 (because of the possible overflow).
 855
 856 From/by Lennart, 94/09/26
 857
 858 \begin{code}
 859 showRational :: Int -> Rational -> String
 860 showRational n r =
 861     if r == 0 then
 862         "0.0"
 863     else
 864         let (r', e) = normalize r
 865         in  prR n r' e
 866
 867 startExpExp = 4 :: Int
 868
 869 -- make sure 1 <= r < 10
 870 normalize :: Rational -> (Rational, Int)
 871 normalize r = if r < 1 then
 872                   case norm startExpExp (1 / r) 0 of (r', e) -> (10 / r', -e-1)
 873               else
 874                   norm startExpExp r 0
 875         where norm :: Int -> Rational -> Int -> (Rational, Int)
 876               -- Invariant: r*10^e == original r
 877               norm 0  r e = (r, e)
 878               norm ee r e =
 879                 let n = 10^ee
 880                     tn = 10^n
 881                 in  if r >= tn then norm ee (r/tn) (e+n) else norm (ee-1) r e
 882
 883 drop0 "" = ""
 884 drop0 (c:cs) = c : reverse (dropWhile (=='0') (reverse cs))
 885
 886 prR :: Int -> Rational -> Int -> String
 887 prR n r e | r <  1  = prR n (r*10) (e-1)                -- final adjustment
 888 prR n r e | r >= 10 = prR n (r/10) (e+1)
 889 prR n r e0 =
 890         let s = show ((round (r * 10^n))::Integer)
 891             e = e0+1
 892         in  if e > 0 && e < 8 then
 893                 take e s ++ "." ++ drop0 (drop e s)
 894             else if e <= 0 && e > -3 then
 895                 "0." ++ take (-e) (repeat '0') ++ drop0 s
 896             else
 897                 head s : "."++ drop0 (tail s) ++ "e" ++ show e0
 898 \end{code}
 899
 900 %*********************************************************
 901 %*                                                      *
 902 \subsection{Numeric primops}
 903 %*                                                      *
 904 %*********************************************************
 905
 906 Definitions of the boxed PrimOps; these will be
 907 used in the case of partial applications, etc.
 908
 909 \begin{code}
 910 plusFloat   (F# x) (F# y) = F# (plusFloat# x y)
 911 minusFloat  (F# x) (F# y) = F# (minusFloat# x y)
 912 timesFloat  (F# x) (F# y) = F# (timesFloat# x y)
 913 divideFloat (F# x) (F# y) = F# (divideFloat# x y)
 914 negateFloat (F# x)        = F# (negateFloat# x)
 915
 916 gtFloat     (F# x) (F# y) = gtFloat# x y
 917 geFloat     (F# x) (F# y) = geFloat# x y
 918 eqFloat     (F# x) (F# y) = eqFloat# x y
 919 neFloat     (F# x) (F# y) = neFloat# x y
 920 ltFloat     (F# x) (F# y) = ltFloat# x y
 921 leFloat     (F# x) (F# y) = leFloat# x y
 922
 923 float2Int   (F# x) = I# (float2Int# x)
 924 int2Float   (I# x) = F# (int2Float# x)
 925
 926 expFloat    (F# x) = F# (expFloat# x)
 927 logFloat    (F# x) = F# (logFloat# x)
 928 sqrtFloat   (F# x) = F# (sqrtFloat# x)
 929 sinFloat    (F# x) = F# (sinFloat# x)
 930 cosFloat    (F# x) = F# (cosFloat# x)
 931 tanFloat    (F# x) = F# (tanFloat# x)
 932 asinFloat   (F# x) = F# (asinFloat# x)
 933 acosFloat   (F# x) = F# (acosFloat# x)
 934 atanFloat   (F# x) = F# (atanFloat# x)
 935 sinhFloat   (F# x) = F# (sinhFloat# x)
 936 coshFloat   (F# x) = F# (coshFloat# x)
 937 tanhFloat   (F# x) = F# (tanhFloat# x)
 938
 939 powerFloat  (F# x) (F# y) = F# (powerFloat# x y)
 940
 941 -- definitions of the boxed PrimOps; these will be
 942 -- used in the case of partial applications, etc.
 943
 944 plusDouble   (D# x) (D# y) = D# (x +## y)
 945 minusDouble  (D# x) (D# y) = D# (x -## y)
 946 timesDouble  (D# x) (D# y) = D# (x *## y)
 947 divideDouble (D# x) (D# y) = D# (x /## y)
 948 negateDouble (D# x)        = D# (negateDouble# x)
 949
 950 gtDouble    (D# x) (D# y) = x >## y
 951 geDouble    (D# x) (D# y) = x >=## y
 952 eqDouble    (D# x) (D# y) = x ==## y
 953 neDouble    (D# x) (D# y) = x /=## y
 954 ltDouble    (D# x) (D# y) = x <## y
 955 leDouble    (D# x) (D# y) = x <=## y
 956
 957 double2Int   (D# x) = I# (double2Int#   x)
 958 int2Double   (I# x) = D# (int2Double#   x)
 959 double2Float (D# x) = F# (double2Float# x)
 960 float2Double (F# x) = D# (float2Double# x)
 961
 962 expDouble    (D# x) = D# (expDouble# x)
 963 logDouble    (D# x) = D# (logDouble# x)
 964 sqrtDouble   (D# x) = D# (sqrtDouble# x)
 965 sinDouble    (D# x) = D# (sinDouble# x)
 966 cosDouble    (D# x) = D# (cosDouble# x)
 967 tanDouble    (D# x) = D# (tanDouble# x)
 968 asinDouble   (D# x) = D# (asinDouble# x)
 969 acosDouble   (D# x) = D# (acosDouble# x)
 970 atanDouble   (D# x) = D# (atanDouble# x)
 971 sinhDouble   (D# x) = D# (sinhDouble# x)
 972 coshDouble   (D# x) = D# (coshDouble# x)
 973 tanhDouble   (D# x) = D# (tanhDouble# x)
 974
 975 powerDouble  (D# x) (D# y) = D# (x **## y)
 976 \end{code}