FIX #4381

[ghc-base.git] / GHC / Float.lhs

diff --git a/GHC/Float.lhs b/GHC/Float.lhs
index 61d009f..fc52dd9 100644 (file)
--- a/GHC/Float.lhs
+++ b/GHC/Float.lhs
@@ -131,9 +131,20 @@ class  (RealFrac a, Floating a) => RealFloat a  where
     significand x       =  encodeFloat m (negate (floatDigits x))
                            where (m,_) = decodeFloat x
 
-    scaleFloat k x      =  encodeFloat m (n+k)
+    scaleFloat k x      =  encodeFloat m (n + clamp b k)
                            where (m,n) = decodeFloat x
-                           
+                                 (l,h) = floatRange x
+                                 d     = floatDigits x
+                                 b     = h - l + 4*d
+                                 -- n+k may overflow, which would lead
+                                 -- to wrong results, hence we clamp the
+                                 -- scaling parameter.
+                                 -- If n + k would be larger than h,
+                                 -- n + clamp b k must be too, simliar
+                                 -- for smaller than l - d.
+                                 -- Add a little extra to keep clear
+                                 -- from the boundary cases.
+
     atan2 y x
       | x > 0            =  atan (y/x)
       | x == 0 && y > 0  =  pi/2
@@ -171,10 +182,7 @@ instance  Num Float  where
     fromInteger i = F# (floatFromInteger i)
 
 instance  Real Float  where
-    toRational x | isInfinite x     = if x < 0 then -infinity else infinity
-                 | isNaN x          = notANumber
-                 | isNegativeZero x = negativeZero
-                 | otherwise        = (m%1)*(b%1)^^n
+    toRational x        =  (m%1)*(b%1)^^n
                            where (m,n) = decodeFloat x
                                  b     = floatRadix  x
 
@@ -252,7 +260,7 @@ instance  Floating Float  where
 
     asinh x = log (x + sqrt (1.0+x*x))
     acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))
-    atanh x = log ((x+1.0) / sqrt (1.0-x*x))
+    atanh x = 0.5 * log ((1.0+x) / (1.0-x))
 
 instance  RealFloat Float  where
     floatRadix _        =  FLT_RADIX        -- from float.h
@@ -271,7 +279,9 @@ instance  RealFloat Float  where
                             (m,_) -> encodeFloat m (negate (floatDigits x))
 
     scaleFloat k x      = case decodeFloat x of
-                            (m,n) -> encodeFloat m (n+k)
+                            (m,n) -> encodeFloat m (n + clamp bf k)
+                        where bf = FLT_MAX_EXP - (FLT_MIN_EXP) + 4*FLT_MANT_DIG
+
     isNaN x          = 0 /= isFloatNaN x
     isInfinite x     = 0 /= isFloatInfinite x
     isDenormalized x = 0 /= isFloatDenormalized x
@@ -306,10 +316,7 @@ instance  Num Double  where
 
 
 instance  Real Double  where
-    toRational x | isInfinite x     = if x < 0 then -infinity else infinity
-                 | isNaN x          = notANumber
-                 | isNegativeZero x = negativeZero
-                 | otherwise        = (m%1)*(b%1)^^n
+    toRational x        =  (m%1)*(b%1)^^n
                            where (m,n) = decodeFloat x
                                  b     = floatRadix  x
 
@@ -337,7 +344,7 @@ instance  Floating Double  where
 
     asinh x = log (x + sqrt (1.0+x*x))
     acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0)))
-    atanh x = log ((x+1.0) / sqrt (1.0-x*x))
+    atanh x = 0.5 * log ((1.0+x) / (1.0-x))
 
 {-# RULES "truncate/Double->Int" truncate = double2Int #-}
 instance  RealFrac Double  where
@@ -401,7 +408,8 @@ instance  RealFloat Double  where
                             (m,_) -> encodeFloat m (negate (floatDigits x))
 
     scaleFloat k x      = case decodeFloat x of
-                            (m,n) -> encodeFloat m (n+k)
+                            (m,n) -> encodeFloat m (n + clamp bd k)
+                        where bd = DBL_MAX_EXP - (DBL_MIN_EXP) + 4*DBL_MANT_DIG
 
     isNaN x             = 0 /= isDoubleNaN x
     isInfinite x        = 0 /= isDoubleInfinite x
@@ -719,7 +727,6 @@ fromRat (n :% 0) | n > 0     =  1/0        -- +Infinity
 
 fromRat (n :% d) | n > 0     = fromRat' (n :% d)
                  | n < 0     = - fromRat' ((-n) :% d)
-                 | d < 0     = 0/(-1)      -- -0.0
                  | otherwise = encodeFloat 0 0             -- Zero
 
 -- Conversion process:
@@ -978,3 +985,12 @@ showSignedFloat showPos p x
        = showParen (p > 6) (showChar '-' . showPos (-x))
    | otherwise = showPos x
 \end{code}
+
+We need to prevent over/underflow of the exponent in encodeFloat when
+called from scaleFloat, hence we clamp the scaling parameter.
+We must have a large enough range to cover the maximum difference of
+exponents returned by decodeFloat.
+\begin{code}
+clamp :: Int -> Int -> Int
+clamp bd k = max (-bd) (min bd k)
+\end{code}