diff options
Diffstat (limited to 'libraries/base')
| -rw-r--r-- | libraries/base/Data/Complex.hs | 15 | ||||
| -rw-r--r-- | libraries/base/GHC/Float.hs | 86 | ||||
| -rw-r--r-- | libraries/base/Numeric.hs | 1 | ||||
| -rw-r--r-- | libraries/base/changelog.md | 3 |
4 files changed, 105 insertions, 0 deletions
diff --git a/libraries/base/Data/Complex.hs b/libraries/base/Data/Complex.hs index 31550d5ac7..dd831bbb91 100644 --- a/libraries/base/Data/Complex.hs +++ b/libraries/base/Data/Complex.hs @@ -37,6 +37,7 @@ module Data.Complex ) where import GHC.Generics (Generic, Generic1) +import GHC.Float (Floating(..)) import Data.Data (Data) import Foreign (Storable, castPtr, peek, poke, pokeElemOff, peekElemOff, sizeOf, alignment) @@ -195,6 +196,20 @@ instance (RealFloat a) => Floating (Complex a) where acosh z = log (z + (z+1) * sqrt ((z-1)/(z+1))) atanh z = 0.5 * log ((1.0+z) / (1.0-z)) + log1p x@(a :+ b) + | abs a < 0.5 && abs b < 0.5 + , u <- 2*a + a*a + b*b = log1p (u/(1 + sqrt(u+1))) :+ atan2 (1 + a) b + | otherwise = log (1 + x) + {-# INLINE log1p #-} + + expm1 x@(a :+ b) + | a*a + b*b < 1 + , u <- expm1 a + , v <- sin (b/2) + , w <- -2*v*v = (u*w + u + w) :+ (u+1)*sin b + | otherwise = exp x - 1 + {-# INLINE expm1 #-} + instance Storable a => Storable (Complex a) where sizeOf a = 2 * sizeOf (realPart a) alignment a = alignment (realPart a) diff --git a/libraries/base/GHC/Float.hs b/libraries/base/GHC/Float.hs index e6180afcf6..ddf9cf01ca 100644 --- a/libraries/base/GHC/Float.hs +++ b/libraries/base/GHC/Float.hs @@ -4,6 +4,7 @@ , MagicHash , UnboxedTuples #-} +{-# LANGUAGE CApiFFI #-} -- We believe we could deorphan this module, by moving lots of things -- around, but we haven't got there yet: {-# OPTIONS_GHC -Wno-orphans #-} @@ -61,6 +62,46 @@ class (Fractional a) => Floating a where sinh, cosh, tanh :: a -> a asinh, acosh, atanh :: a -> a + -- | @'log1p' x@ computes @'log' (1 + x)@, but provides more precise + -- results for small (absolute) values of @x@ if possible. + -- + -- @since 4.9.0.0 + log1p :: a -> a + + -- | @'expm1' x@ computes @'exp' x - 1@, but provides more precise + -- results for small (absolute) values of @x@ if possible. + -- + -- @since 4.9.0.0 + expm1 :: a -> a + + -- | @'log1pexp' x@ computes @'log' (1 + 'exp' x)@, but provides more + -- precise results if possible. + -- + -- Examples: + -- + -- * if @x@ is a large negative number, @'log' (1 + 'exp' x)@ will be + -- imprecise for the reasons given in 'log1p'. + -- + -- * if @'exp' x@ is close to @-1@, @'log' (1 + 'exp' x)@ will be + -- imprecise for the reasons given in 'expm1'. + -- + -- @since 4.9.0.0 + log1pexp :: a -> a + + -- | @'log1mexp' x@ computes @'log' (1 - 'exp' x)@, but provides more + -- precise results if possible. + -- + -- Examples: + -- + -- * if @x@ is a large negative number, @'log' (1 - 'exp' x)@ will be + -- imprecise for the reasons given in 'log1p'. + -- + -- * if @'exp' x@ is close to @1@, @'log' (1 - 'exp' x)@ will be + -- imprecise for the reasons given in 'expm1'. + -- + -- @since 4.9.0.0 + log1mexp :: a -> a + {-# INLINE (**) #-} {-# INLINE logBase #-} {-# INLINE sqrt #-} @@ -72,6 +113,15 @@ class (Fractional a) => Floating a where tan x = sin x / cos x tanh x = sinh x / cosh x + {-# INLINE log1p #-} + {-# INLINE expm1 #-} + {-# INLINE log1pexp #-} + {-# INLINE log1mexp #-} + log1p x = log (1 + x) + expm1 x = exp x - 1 + log1pexp x = log1p (exp x) + log1mexp x = log1p (negate (exp x)) + -- | Efficient, machine-independent access to the components of a -- floating-point number. class (RealFrac a, Floating a) => RealFloat a where @@ -307,6 +357,19 @@ instance Floating Float where acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0))) atanh x = 0.5 * log ((1.0+x) / (1.0-x)) + log1p = log1pFloat + expm1 = expm1Float + + log1mexp a + | a <= log 2 = log (negate (expm1Float a)) + | otherwise = log1pFloat (negate (exp a)) + {-# INLINE log1mexp #-} + log1pexp a + | a <= 18 = log1pFloat (exp a) + | a <= 100 = a + exp (negate a) + | otherwise = a + {-# INLINE log1pexp #-} + instance RealFloat Float where floatRadix _ = FLT_RADIX -- from float.h floatDigits _ = FLT_MANT_DIG -- ditto @@ -415,6 +478,19 @@ instance Floating Double where acosh x = log (x + (x+1.0) * sqrt ((x-1.0)/(x+1.0))) atanh x = 0.5 * log ((1.0+x) / (1.0-x)) + log1p = log1pDouble + expm1 = expm1Double + + log1mexp a + | a <= log 2 = log (negate (expm1Double a)) + | otherwise = log1pDouble (negate (exp a)) + {-# INLINE log1mexp #-} + log1pexp a + | a <= 18 = log1pDouble (exp a) + | a <= 100 = a + exp (negate a) + | otherwise = a + {-# INLINE log1pexp #-} + -- RULES for Integer and Int {-# RULES "properFraction/Double->Integer" properFraction = properFractionDoubleInteger @@ -1069,6 +1145,16 @@ foreign import ccall unsafe "isDoubleDenormalized" isDoubleDenormalized :: Doubl foreign import ccall unsafe "isDoubleNegativeZero" isDoubleNegativeZero :: Double -> Int foreign import ccall unsafe "isDoubleFinite" isDoubleFinite :: Double -> Int + +------------------------------------------------------------------------ +-- libm imports for extended floating +------------------------------------------------------------------------ +foreign import capi unsafe "math.h log1p" log1pDouble :: Double -> Double +foreign import capi unsafe "math.h expm1" expm1Double :: Double -> Double +foreign import capi unsafe "math.h log1pf" log1pFloat :: Float -> Float +foreign import capi unsafe "math.h expm1f" expm1Float :: Float -> Float + + ------------------------------------------------------------------------ -- Coercion rules ------------------------------------------------------------------------ diff --git a/libraries/base/Numeric.hs b/libraries/base/Numeric.hs index 4e24bfe7a8..51be3a1698 100644 --- a/libraries/base/Numeric.hs +++ b/libraries/base/Numeric.hs @@ -56,6 +56,7 @@ module Numeric ( -- * Miscellaneous fromRat, + Floating(..) ) where diff --git a/libraries/base/changelog.md b/libraries/base/changelog.md index 33a51143a8..fa57556de2 100644 --- a/libraries/base/changelog.md +++ b/libraries/base/changelog.md @@ -112,6 +112,9 @@ * Re-export `Const` from `Control.Applicative` for backwards compatibility. + * Expand `Floating` class to include operations that allow for better + precision: `log1p`, `expm1`, `log1pexp` and `log1mexp`. These are not + available from `Prelude`, but the full class is exported from `Numeric`. ## 4.8.2.0 *Oct 2015* |