ghc/tests/codeGen/should_run/cg034.hs

   1 -- !! fromRational woes
   2 import Ratio -- 1.3
   3
   4 main = putStr (
   5    shows tinyFloat  ( '\n'
   6  : shows t_f        ( '\n'
   7  : shows hugeFloat  ( '\n'
   8  : shows h_f        ( '\n'
   9  : shows tinyDouble ( '\n'
  10  : shows t_d        ( '\n'
  11  : shows hugeDouble ( '\n'
  12  : shows h_d        ( '\n'
  13  : shows x_f        ( '\n'
  14  : shows x_d        ( '\n'
  15  : shows y_f        ( '\n'
  16  : shows y_d        ( "\n"
  17  )))))))))))))
  18   where
  19     t_f :: Float
  20     t_d :: Double
  21     h_f :: Float
  22     h_d :: Double
  23     x_f :: Float
  24     x_d :: Double
  25     y_f :: Float
  26     y_d :: Double
  27     t_f = fromRationalX (toRational tinyFloat)
  28     t_d = fromRationalX (toRational tinyDouble)
  29     h_f = fromRationalX (toRational hugeFloat)
  30     h_d = fromRationalX (toRational hugeDouble)
  31     x_f = fromRationalX (1.82173691287639817263897126389712638972163e-300 :: Rational)
  32     x_d = fromRationalX (1.82173691287639817263897126389712638972163e-300 :: Rational)
  33     y_f = 1.82173691287639817263897126389712638972163e-300
  34     y_d = 1.82173691287639817263897126389712638972163e-300
  35
  36 fromRationalX :: (RealFloat a) => Rational -> a
  37 fromRationalX r =
  38         let
  39             h = ceiling (huge `asTypeOf` x)
  40             b = toInteger (floatRadix x)
  41             x = fromRat 0 r
  42             fromRat e0 r' =
  43                 let d = denominator r'
  44                     n = numerator r'
  45                 in  if d > h then
  46                        let e = integerLogBase b (d `div` h) + 1
  47                        in  fromRat (e0-e) (n % (d `div` (b^e)))
  48                     else if abs n > h then
  49                        let e = integerLogBase b (abs n `div` h) + 1
  50                        in  fromRat (e0+e) ((n `div` (b^e)) % d)
  51                     else
  52                        scaleFloat e0 (rationalToRealFloat {-fromRational-} r')
  53         in  x
  54
  55 {-
  56 fromRationalX r =
  57   rationalToRealFloat r
  58 {- Hmmm...
  59         let
  60             h = ceiling (huge `asTypeOf` x)
  61             b = toInteger (floatRadix x)
  62             x = fromRat 0 r
  63
  64             fromRat e0 r' =
  65 {--}            trace (shows e0 ('/' : shows r' ('/' : shows h "\n"))) (
  66                 let d = denominator r'
  67                     n = numerator r'
  68                 in  if d > h then
  69                        let e = integerLogBase b (d `div` h) + 1
  70                        in  fromRat (e0-e) (n % (d `div` (b^e)))
  71                     else if abs n > h then
  72                        let e = integerLogBase b (abs n `div` h) + 1
  73                        in  fromRat (e0+e) ((n `div` (b^e)) % d)
  74                     else
  75                        scaleFloat e0 (rationalToRealFloat r')
  76                        -- now that we know things are in-bounds,
  77                        -- we use the "old" Prelude code.
  78 {--}            )
  79         in  x
  80 -}
  81 -}
  82
  83 -- Compute the discrete log of i in base b.
  84 -- Simplest way would be just divide i by b until it's smaller then b, but that would
  85 -- be very slow!  We are just slightly more clever.
  86 integerLogBase :: Integer -> Integer -> Int
  87 integerLogBase b i =
  88      if i < b then
  89         0
  90      else
  91         -- Try squaring the base first to cut down the number of divisions.
  92         let l = 2 * integerLogBase (b*b) i
  93             doDiv :: Integer -> Int -> Int
  94             doDiv i l = if i < b then l else doDiv (i `div` b) (l+1)
  95         in  doDiv (i `div` (b^l)) l
  96
  97
  98 ------------
  99
 100 -- Compute smallest and largest floating point values.
 101 tiny :: (RealFloat a) => a
 102 tiny =
 103         let (l, _) = floatRange x
 104             x = encodeFloat 1 (l-1)
 105         in  x
 106
 107 huge :: (RealFloat a) => a
 108 huge =
 109         let (_, u) = floatRange x
 110             d = floatDigits x
 111             x = encodeFloat (floatRadix x ^ d - 1) (u - d)
 112         in  x
 113
 114 tinyDouble = tiny :: Double
 115 tinyFloat  = tiny :: Float
 116 hugeDouble = huge :: Double
 117 hugeFloat  = huge :: Float
 118
 119 {-
 120 [In response to a request by simonpj, Joe Fasel writes:]
 121
 122 A quite reasonable request!  This code was added to the Prelude just
 123 before the 1.2 release, when Lennart, working with an early version
 124 of hbi, noticed that (read . show) was not the identity for
 125 floating-point numbers.  (There was a one-bit error about half the time.)
 126 The original version of the conversion function was in fact simply
 127 a floating-point divide, as you suggest above.  The new version is,
 128 I grant you, somewhat denser.
 129
 130 How's this?
 131
 132 --Joe
 133 -}
 134
 135
 136 rationalToRealFloat :: (RealFloat a) => Rational -> a
 137
 138 rationalToRealFloat x   =  x'
 139         where x'    = f e
 140
 141 --              If the exponent of the nearest floating-point number to x
 142 --              is e, then the significand is the integer nearest xb^(-e),
 143 --              where b is the floating-point radix.  We start with a good
 144 --              guess for e, and if it is correct, the exponent of the
 145 --              floating-point number we construct will again be e.  If
 146 --              not, one more iteration is needed.
 147
 148               f e   = if e' == e then y else f e'
 149                       where y      = encodeFloat (round (x * (1%b)^^e)) e
 150                             (_,e') = decodeFloat y
 151               b     = floatRadix x'
 152
 153 --              We obtain a trial exponent by doing a floating-point
 154 --              division of x's numerator by its denominator.  The
 155 --              result of this division may not itself be the ultimate
 156 --              result, because of an accumulation of three rounding
 157 --              errors.
 158
 159               (s,e) = decodeFloat (fromInteger (numerator x) `asTypeOf` x'
 160                                         / fromInteger (denominator x))
 161