logfloat 0.11.1 → 0.12
raw patch · 5 files changed
+194/−121 lines, 5 files
Files
- logfloat.cabal +20/−5
- src/Data/Number/LogFloat.hs +115/−112
- src/Data/Number/RealToFrac.hs +2/−0
- src/Data/Number/Transfinite.hs +45/−3
- src/Hugs/RealFloat.hs +12/−1
logfloat.cabal view
@@ -1,9 +1,9 @@ ------------------------------------------------------------------- wren ng thornton <wren@community.haskell.org> ~ 2009.01.29+-- wren ng thornton <wren@community.haskell.org> ~ 2009.03.17 ---------------------------------------------------------------- Name: logfloat-Version: 0.11.1+Version: 0.12 Cabal-Version: >= 1.2 Build-Type: Simple Stability: experimental@@ -26,6 +26,11 @@ Description: base-3.0 broke out array and other packages Default: False +Flag useFFI+ Description: Use FFI to link C's log1p. Improves accuracy and+ same performance. Disable only if your compiler+ doesn't support FFI.+ Default: True Library Hs-Source-Dirs: src@@ -39,9 +44,19 @@ else Build-depends: base < 3.0 - Hugs-Options: -98 +o -F'cpp -P -traditional -D__HUGS__=200609'- if impl(ghc < 6.10)- GHC-Options: -fno-warn-orphans+ if flag(useFFI)+ -- BUG: (Cabal 1.2 + Haddock) See the INSTALL file.+ --GHC-Options: -D__USE_FFI__+ CPP-Options: -D__USE_FFI__+ includes: math.h+ extra-libraries: m+ + -- BUG: (Cabal <= 1.6 + Hugs) See the INSTALL file.+ Hugs-Options: -98 +o+ if impl(hugs)+ -- BUG: (Cabal 1.2 + Haddock) See the INSTALL file.+ --GHC-Options: -D__HUGS__=200609+ CPP-Options: -D__HUGS__=200609 ---------------------------------------------------------------- ----------------------------------------------------------- fin.
src/Data/Number/LogFloat.hs view
@@ -1,11 +1,14 @@ -- FlexibleContexts needed by our RealToFrac contexts -- CPP needed for IArray UArray instance+-- FFI is for log1p {-# LANGUAGE FlexibleContexts- , CPP #-}+ , CPP+ , ForeignFunctionInterface+ #-} -- Removed -Wall because -fno-warn-orphans was removed in GHC 6.10-{-# OPTIONS_GHC -fwarn-tabs #-}+{-# OPTIONS_GHC -Wall -fwarn-tabs #-} -- Unfortunately we need -fglasgow-exts in order to actually pick -- up on the rules (see -ddump-rules). The -frewrite-rules flag@@ -14,37 +17,15 @@ -- cf <http://www.mail-archive.com/glasgow-haskell-users@haskell.org/msg14313.html> {-# OPTIONS_GHC -O2 -fvia-C -optc-O3 -fexcess-precision -fglasgow-exts #-} --- Version History--- (v0.11.1) Added IArray UArray instance--- (v0.11) Broke Data.Number.RealToFrac out--- (v0.10) Fixed bugs in Hugs for PartialOrd and Transfinite.--- Also added maxPO, minPO, comparingPO--- (v0.9.1) Fixed some PartialOrd stuff and sanitized documentation--- (v0.9.0) s/toFractional/realToFrac/g.--- Also moved realToFrac and log to Transfinite--- (v0.8.6) Removed buggy RULES--- (v0.8.5) Gave up and converted from lhs to hs so Hackage docs work--- (v0.8.4) Broke out Transfinite--- (v0.8.3) Documentation updates--- (v0.8.2) Announced release--- (v0.8) Did a bunch of tweaking. Things should be decent now--- (v0.7) Haddockified--- (v0.6) Fixed monomorphism.--- (v0.5) Added optimization rules.--- (v0.4) Translated to Haskell at revision 2007.12.20.--- (v0.3) Converted extensive comments to POD format.--- (v0.2) Did a bunch of profiling, optimizing, and debugging.--- (v0.1) Initial version created for hw5 for NLP with Jason Eisner.--- ------------------------------------------------------------------- ~ 2009.03.07+-- ~ 2009.03.10 -- | -- Module : Data.Number.LogFloat -- Copyright : Copyright (c) 2007--2009 wren ng thornton -- License : BSD3 -- Maintainer : wren@community.haskell.org -- Stability : stable--- Portability : portable (with CPP)+-- Portability : portable (with CPP, FFI) -- -- This module presents a type for storing numbers in the log-domain. -- The main reason for doing this is to prevent underflow when@@ -61,7 +42,7 @@ -- The 'LogFloat' of this module is restricted to non-negative -- numbers for efficiency's sake, see the forthcoming -- "Data.Number.LogFloat.Signed" for doing signed log-domain--- calculations.+-- calculations. (Or harass the maintainer to write it already.) ---------------------------------------------------------------- module Data.Number.LogFloat@@ -70,10 +51,17 @@ module Data.Number.Transfinite , module Data.Number.RealToFrac - -- * @LogFloat@ data type and conversion functions+ -- * @LogFloat@ data type , LogFloat- , logFloat, logToLogFloat- , fromLogFloat, logFromLogFloat+ -- ** Isomorphism to normal-domain+ , logFloat+ , fromLogFloat+ -- ** Isomorphism to log-domain+ , logToLogFloat+ , logFromLogFloat+ + -- * Accurate versions of logarithm\/exponentiation+ , log1p, expm1 ) where import Prelude hiding (log, realToFrac, isInfinite, isNaN)@@ -96,73 +84,9 @@ import NonStdUnsafeCoerce (unsafeCoerce) #endif ----------------------------------------------------------------------- Try to add in some optimizations. Why these need to be--- down here and localized to the module, I don't know. We don't--- do anything foolish like this, but our clients might, or they--- might be generated by other code transformations. Note that due--- to the fuzz, these equations are not strictly true, even though--- they are mathematically correct.--{-# RULES-"log/exp" forall x. log (exp x) = x-"log.exp" log . exp = id--"exp/log" forall x. exp (log x) = x-"exp.log" exp . log = id- #-}---- We'd like to be able to take advantage of general rule versions--- of our operators for 'LogFloat', with rules like @log x + log y--- = log (x * y)@ and @log x - log y = log (x / y)@. However the--- problem is that those equations could be driven in either direction--- depending on whether we think time performance or non-underflow--- performance is more important, and the answers may be different--- at every call site.------ Since we implore users to do normal-domain computations whenever--- it would not degenerate accuracy, we should not rewrite their--- decisions in any way. The log\/exp fusion strictly improves both--- time and accuracy, so those are safe. But the buck stops with--- them.----- These should only fire when it's type-safe--- This should already happen, but...--- TODO: Check the logs to see if it ever fires--- N.B. these are orphaned-{-# RULES-"toRational/fromRational" forall x. toRational (fromRational x) = x-"toRational.fromRational" toRational . fromRational = id- #-}----------------------------------------------------------------------- | Reduce the number of constant string literals we need to store.-errorOutOfRange :: String -> a-errorOutOfRange fun = error $! "Data.Number.LogFloat."++fun- ++ ": argument out of range"----- | We need these guards in order to ensure some invariants.-guardNonNegative :: String -> Double -> Double-guardNonNegative fun x | x >= 0 = x- | otherwise = errorOutOfRange fun----- TODO: since we're using Hugs.RealFloat instead of Prelude now,--- is it still non-portable?------ | It's unfortunate that 'notANumber' is not equal to itself, but--- we can hack around that. GHC gives NaN for the log of negatives--- and so we could ideally take advantage of @log . guardNonNegative--- fun = guardIsANumber fun . log@ to simplify things, but Hugs--- raises an error so that's non-portable.-guardIsANumber :: String -> Double -> Double-guardIsANumber fun x | isNaN x = errorOutOfRange fun- | otherwise = x+#ifdef __GLASGOW_HASKELL__+import Foreign.Storable (Storable)+#endif ---------------------------------------------------------------- --@@ -196,15 +120,17 @@ ( Eq , Ord -- Should we really perpetuate the Ord lie? #ifdef __GLASGOW_HASKELL__- , IArray UArray- -- At least GHC 6.8.2 can derive IArray UArray (without+ -- At least GHC 6.8.2 can derive these (without -- GeneralizedNewtypeDeriving). The H98 Report doesn't include- -- that among the options for automatic derivation though.+ -- them among the options for automatic derivation though.+ , IArray UArray+ , Storable #endif ) #if __HUGS__ || __NHC__+-- TODO: Storable instance. Though Foreign.Storable isn't in Hugs(Sept06) -- These two operators make it much easier to read the instance. -- Hopefully inlining everything will get rid of the eta overhead.@@ -235,7 +161,7 @@ logToLFFunc = ($. unsafeLogToLogFloat ~> id ~> logFromLogFloat) -- | Remove the extranious 'isNaN' test of 'logToLogFloat', when--- we know we can.+-- we know it's safe. {-# INLINE unsafeLogToLogFloat #-} unsafeLogToLogFloat :: Double -> LogFloat unsafeLogToLogFloat = LogFloat@@ -281,8 +207,31 @@ ------------------------------------------------------------------- | A constructor which does semantic conversion from normal-domain--- to log-domain.+-- | Reduce the number of constant string literals we need to store.+errorOutOfRange :: String -> a+{-# NOINLINE errorOutOfRange #-}+errorOutOfRange fun = error $! "Data.Number.LogFloat."++fun+ ++ ": argument out of range"+++-- | We need these guards in order to ensure some invariants.+guardNonNegative :: String -> Double -> Double+guardNonNegative fun x | x >= 0 = x+ | otherwise = errorOutOfRange fun+++-- | In general @log . guardNonNegative fun == guardIsANumber fun . log@+-- This even holds on Hugs now. However, the latter issues an error+-- from 'Data.Number.Transfinite.log' whereas the former issues an+-- error from the calling function and is therefore more helpful+-- for debugging.+guardIsANumber :: String -> Double -> Double+guardIsANumber fun x | isNaN x = errorOutOfRange fun+ | otherwise = x++----------------------------------------------------------------+-- | Constructor which does semantic conversion from normal-domain+-- to log-domain. Throws errors on negative input. logFloat :: (Real a, RealToFrac a Double) => a -> LogFloat {-# SPECIALIZE logFloat :: Double -> LogFloat #-} logFloat = LogFloat . log . guardNonNegative "logFloat" . realToFrac@@ -296,7 +245,7 @@ -- constructors\/destructors. -- -- | Constructor which assumes the argument is already in the--- log-domain.+-- log-domain. Throws errors on @notANumber@ input. logToLogFloat :: (Real a, RealToFrac a Double) => a -> LogFloat {-# SPECIALIZE logToLogFloat :: Double -> LogFloat #-} logToLogFloat = LogFloat . guardIsANumber "logToLogFloat" . realToFrac@@ -310,7 +259,7 @@ fromLogFloat (LogFloat x) = realToFrac (exp x) --- | Return the log-domain value itself without costly conversion+-- | Return the log-domain value itself without conversion. logFromLogFloat :: (Fractional a, Transfinite a, RealToFrac Double a) => LogFloat -> a {-# SPECIALIZE logFromLogFloat :: LogFloat -> Double #-}@@ -358,6 +307,56 @@ ----------------------------------------------------------------+#ifdef __USE_FFI__+#define LOG1P_WHICH_VERSION specialized version.+#else+#define LOG1P_WHICH_VERSION naive version! \+ Contact the maintainer with any FFI difficulties.+#endif+-- | Definition: @log1p == log . (1+)@. The C language provides a+-- special definition for 'log1p' which is more accurate than doing+-- the naive thing, especially for very small arguments. For example,+-- the naive version underflows around @2 ** -53@, whereas the+-- specialized version underflows around @2 ** -1074@. This function+-- is used by ('+') and ('-') on @LogFloat@.+--+-- /This installation was compiled to use the LOG1P_WHICH_VERSION/++#ifdef __USE_FFI__+foreign import ccall unsafe "math.h log1p"+ log1p :: Double -> Double++-- Technically we should use 'Foreign.C.CDouble' however there's+-- no isomorphism provided to normal 'Double'. The former is+-- documented as being a newtype of the later, and so this should+-- be safe.++#else++log1p :: Double -> Double+{-# INLINE log1p #-}+log1p x = log (1 + x)+#endif+++-- | Definition: @expm1 == (subtract 1) . exp@. The C language+-- provides a special definition for 'expm1' which is more accurate+-- than doing the naive thing, especially for very small arguments.+-- This function isn't needed internally, but is provided for+-- symmetry with 'log1p'.+--+-- /This installation was compiled to use the LOG1P_WHICH_VERSION/++#ifdef __USE_FFI__+foreign import ccall unsafe "math.h expm1"+ expm1 :: Double -> Double+#else+expm1 :: Double -> Double+{-# INLINE expm1 #-}+expm1 x = exp x - 1+#endif++---------------------------------------------------------------- -- These all work without causing underflow. However, do note that -- they tend to induce more of the floating-point fuzz than using -- regular floating numbers because @exp . log@ doesn't really equal@@ -376,16 +375,16 @@ (*) (LogFloat x) (LogFloat y) = LogFloat (x+y)-+ (+) (LogFloat x) (LogFloat y)- | x >= y = LogFloat (x + log (1 + exp (y - x)))- | otherwise = LogFloat (y + log (1 + exp (x - y)))-+ | x >= y = LogFloat (x + log1p (exp (y - x)))+ | otherwise = LogFloat (y + log1p (exp (x - y)))+ -- Without the guard this would return NaN instead of error (-) (LogFloat x) (LogFloat y)- | x >= y = LogFloat (x + log (1 - exp (y - x)))+ | x >= y = LogFloat (x + log1p (negate (exp (y - x)))) | otherwise = errorOutOfRange "(-)"-+ signum (LogFloat x) | x == negativeInfinity = 0 | x > negativeInfinity = 1@@ -417,7 +416,11 @@ -- Just for fun. The more coersion functions the better. Though -- Rationals are very buggy when it comes to transfinite values instance Real LogFloat where- toRational (LogFloat x) = toRational (exp x)+ toRational (LogFloat x)+ | isInfinite ex || isNaN ex = errorOutOfRange "toRational"+ | otherwise = toRational ex+ where+ ex = exp x {- -- Commented out because I'm not sure about requiring MPTCs. Of course, those are already required by "Data.Number.Transfinite" so it's pretty moot...
src/Data/Number/RealToFrac.hs view
@@ -97,10 +97,12 @@ instance RealToFrac Integer Float where -- TODO: is there a more primitive way?+ {-# INLINE realToFrac #-} realToFrac j = Prelude.realToFrac j instance RealToFrac Integer Double where -- TODO: is there a more primitive way?+ {-# INLINE realToFrac #-} realToFrac j = Prelude.realToFrac j
src/Data/Number/Transfinite.hs view
@@ -1,7 +1,14 @@ {-# OPTIONS_GHC -Wall -fwarn-tabs #-} +-- Unfortunately we need -fglasgow-exts in order to actually pick+-- up on the rules (see -ddump-rules). The -frewrite-rules flag+-- doesn't do what you want.+-- <http://hackage.haskell.org/trac/ghc/ticket/2213>+-- <http://www.mail-archive.com/glasgow-haskell-users@haskell.org/msg14313.html>+{-# OPTIONS_GHC -fglasgow-exts #-}+ ------------------------------------------------------------------- ~ 2009.01.29+-- ~ 2009.03.09 -- | -- Module : Data.Number.Transfinite -- Copyright : Copyright (c) 2007--2009 wren ng thornton@@ -91,10 +98,14 @@ -- operations must return @notANumber@, where @inf@ is any value -- which @isInfinite@: --- -- * @inf + inf@+ -- * @infinity + negativeInfinity@ --- -- * @inf - inf@+ -- * @negativeInfinity + infinity@ --+ -- * @infinity - infinity@+ --+ -- * @negativeInfinity - negativeInfinity@+ -- -- * @inf * 0@ -- -- * @0 * inf@@@ -163,6 +174,14 @@ -- for your @0@ and @negativeInfinity@. If it doesn't, then you -- should avoid importing our @log@ and will probably want converters -- to handle the discrepancy.+--+-- For GHC, this version of @log@ has rules for fusion with @exp@.+-- These can give different behavior by preventing overflow to+-- @infinity@ and preventing errors for taking the logarithm of+-- negative values. For 'Double' and 'Float' they can also give+-- different answers due to eliminating floating point fuzz. The+-- rules strictly improve mathematical accuracy, however they should+-- be noted in case your code depends on the implementation details. log :: (Floating a, Transfinite a) => a -> a {-# SPECIALIZE log :: Double -> Double #-}@@ -180,5 +199,28 @@ -- hack our way around without using PartialOrd by using isNaN, (== -- 0), ((>0).signum) but that would be less efficient. +----------------------------------------------------------------+-- These rules moved here from "LogFloat" in v0.11.2+{-# RULES+"log/exp" forall x. log (exp x) = x+"log.exp" log . exp = id++"exp/log" forall x. exp (log x) = x+"exp.log" exp . log = id+ #-}++-- We'd like to be able to take advantage of general rule versions+-- of our operators for 'LogFloat', with rules like @log x + log y+-- = log (x * y)@ and @log x - log y = log (x / y)@. However the+-- problem is that those equations could be driven in either direction+-- depending on whether we think time performance or non-underflow+-- performance is more important, and the answers may be different+-- at every call site.+--+-- Since we implore users to do normal-domain computations whenever+-- it would not degenerate accuracy, we should not rewrite their+-- decisions in any way. The log\/exp fusion strictly improves both+-- time and accuracy, so those are safe. But the buck stops with+-- them. ---------------------------------------------------------------- ----------------------------------------------------------- fin.
src/Hugs/RealFloat.hs view
@@ -3,8 +3,15 @@ {-# OPTIONS_GHC -Wall -fwarn-tabs #-} +#if defined(__HUGS__) && (__HUGS__ <= 200609)+#define REALFLOAT_VERSION corrected Hugs version.+#elif defined(__GLASGOW_HASKELL__) || defined(__NHC__)+#define REALFLOAT_VERSION normal Prelude version. This should be correct.+#else+#define REALFLOAT_VERSION normal Prelude version. This could be buggy.+#endif ------------------------------------------------------------------- ~ 2009.01.29+-- ~ 2009.03.10 -- | -- Module : Hugs.RealFloat -- Copyright : Copyright (c) 2007--2009 wren ng thornton@@ -25,6 +32,8 @@ -- N.B. The corrected definitions have only been tested to work for -- 'Float' and 'Double'. These definitions should probably not be -- used for other 'RealFloat' types.+--+-- /This installation was compiled with the REALFLOAT_VERSION/ ---------------------------------------------------------------- module Hugs.RealFloat ( isInfinite@@ -36,6 +45,7 @@ ---------------------------------------------------------------- isInfinite :: (RealFloat a) => a -> Bool+{-# SPECIALIZE isInfinite :: Double -> Bool #-} {-# INLINE isInfinite #-} #if defined(__HUGS__) && (__HUGS__ <= 200609) isInfinite x = (1/0) == abs x@@ -45,6 +55,7 @@ isNaN :: (RealFloat a) => a -> Bool+{-# SPECIALIZE isNaN :: Double -> Bool #-} {-# INLINE isNaN #-} #if defined(__HUGS__) && (__HUGS__ <= 200609) isNaN x = compareEQ x 0 && compareEQ x 1