[commit: ghc] master: Implement phase 1 of expanded Floating (6457903)
git at git.haskell.org
git at git.haskell.org
Mon Dec 21 12:41:45 UTC 2015
Repository : ssh://git@git.haskell.org/ghc
On branch : master
Link : http://ghc.haskell.org/trac/ghc/changeset/6457903e7671b6096d2cca5d965f43daee3572a6/ghc
>---------------------------------------------------------------
commit 6457903e7671b6096d2cca5d965f43daee3572a6
Author: Dan Doel <dan.doel at gmail.com>
Date: Sun Dec 20 15:19:52 2015 +0100
Implement phase 1 of expanded Floating
- This part of the proposal is to add log1p, expm1, log1pexp and
log1mexp to the Floating class, and export the full Floating class
from Numeric
Reviewers: ekmett, #core_libraries_committee, bgamari, hvr, austin
Reviewed By: ekmett, #core_libraries_committee, bgamari
Subscribers: Phyx, RyanGlScott, ekmett, thomie
Differential Revision: https://phabricator.haskell.org/D1605
GHC Trac Issues: #11166
>---------------------------------------------------------------
6457903e7671b6096d2cca5d965f43daee3572a6
libraries/base/Data/Complex.hs | 15 ++++++++
libraries/base/GHC/Float.hs | 86 ++++++++++++++++++++++++++++++++++++++++++
libraries/base/Numeric.hs | 1 +
libraries/base/changelog.md | 3 ++
rts/RtsSymbols.c | 4 ++
5 files changed, 109 insertions(+)
diff --git a/libraries/base/Data/Complex.hs b/libraries/base/Data/Complex.hs
index 31550d5..dd831bb 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 e6180af..ddf9cf0 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 4e24bfe..51be3a1 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 33a5114..fa57556 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*
diff --git a/rts/RtsSymbols.c b/rts/RtsSymbols.c
index 8413d31..4b0a1d5 100644
--- a/rts/RtsSymbols.c
+++ b/rts/RtsSymbols.c
@@ -157,6 +157,10 @@
SymI_HasProto(erfc) \
SymI_HasProto(erff) \
SymI_HasProto(erfcf) \
+ SymI_HasProto(expm1) \
+ SymI_HasProto(expm1f) \
+ SymI_HasProto(log1p) \
+ SymI_HasProto(log1pf) \
SymI_HasProto(memcpy) \
SymI_HasProto(rts_InstallConsoleEvent) \
SymI_HasProto(rts_ConsoleHandlerDone) \
More information about the ghc-commits
mailing list