log-domain 0.1.0.1 → 0.2

raw patch · 2 files changed

+75/−3 lines, 2 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

• log-domain.cabal +1/−1
• src/Numeric/Log.hs +74/−2

@@ -1,6 +1,6 @@ name:          log-domain category:      Numeric-version:       0.1.0.1+version:       0.2 license:       BSD3 cabal-version: >= 1.8 license-file:  LICENSE

@@ -15,8 +15,10 @@ module Numeric.Log   ( Log(..)   , Precise(..)+  , sum   ) where +import Prelude hiding (maximum, sum) import Control.Applicative import Control.Comonad import Control.DeepSeq@@ -24,10 +26,13 @@ import Data.Complex import Data.Data import Data.Distributive-import Data.Foldable+import Data.Foldable as Foldable hiding (sum) import Data.Functor.Bind import Data.Functor.Extend import Data.Hashable+import Data.Int+import Data.List as List hiding (sum)+import Data.Monoid import Data.SafeCopy import Data.Semigroup.Foldable import Data.Semigroup.Traversable@@ -40,6 +45,8 @@ -- | @Log@-domain @Float@ and @Double@ values. newtype Log a = Log { runLog :: a } deriving (Eq,Ord,Data,Typeable,Generic) +deriveSafeCopy 1 'base ''Log+ instance (Floating a, Show a) => Show (Log a) where   showsPrec d (Log a) = showsPrec d (exp a) @@ -128,7 +135,23 @@   Log a >>= f = f a   {-# INLINE (>>=) #-} -deriveSafeCopy 1 'base ''Log+instance (RealFloat a, Precise a, Enum a) => Enum (Log a) where+  succ a = a + 1+  {-# INLINE succ #-}+  pred a = a - 1+  {-# INLINE pred #-}+  toEnum   = fromIntegral+  {-# INLINE toEnum #-}+  fromEnum = round . exp . runLog+  {-# INLINE fromEnum #-}+  enumFrom (Log a) = [ Log (log b) | b <- enumFrom (exp a) ]+  {-# INLINE enumFrom #-}+  enumFromThen (Log a) (Log b) = [ Log (log c) | c <- enumFromThen (exp a) (exp b) ]+  {-# INLINE enumFromThen #-}+  enumFromTo (Log a) (Log b) = [ Log (log c) | c <- enumFromTo (exp a) (exp b) ]+  {-# INLINE enumFromTo #-}+  enumFromThenTo (Log a) (Log b) (Log c) = [ Log (log d) | d <- enumFromThenTo (exp a) (exp b) (exp c) ]+  {-# INLINE enumFromThenTo #-}  -- | Negative infinity negInf :: Fractional a => a@@ -175,9 +198,58 @@   toRational (Log a) = toRational (exp a)   {-# INLINE toRational #-} +data Acc1 a = Acc1 {-# UNPACK #-} !Int64 !a++instance (Precise a, RealFloat a) => Monoid (Log a) where+  mempty  = Log negInf+  {-# INLINE mempty #-}+  mappend = (+)+  {-# INLINE mappend #-}+  mconcat [] = 0+  mconcat (Log z:zs) = Log $ case List.foldl' step1 (Acc1 0 z) zs of+    Acc1 nm1 a+      | isInfinite a -> a+      | otherwise    -> a + log1p (List.foldl' (step2 a) 0 zs + fromIntegral nm1)+    where+      step1 (Acc1 n y) (Log x) = Acc1 (n + 1) (max x y)+      step2 a r (Log x) = r + expm1 (x - a)+  {-# INLINE mconcat #-}+ logMap :: Floating a => (a -> a) -> Log a -> Log a logMap f = Log . log . f . exp . runLog {-# INLINE logMap #-}++data Acc a = Acc {-# UNPACK #-} !Int64 !a | None++-- | Efficiently and accurately compute the sum of a set of log-domain numbers+--+-- While folding with @(+)@ accomplishes the same end, it requires an+-- additional @n-2@ logarithms to sum @n@ terms. In addition,+-- here we introduce fewer opportunities for round-off error.+--+-- While for small quantities the naive sum accumulates error,+--+-- >>> let xs = replicate 40000 (Log 1e-4) :: [Log Float]+-- >>> Prelude.sum xs+-- 40001.3+--+-- This sum gives a more accurate result,+--+-- >>> Numeric.Log.sum xs+-- 40004.01+--+-- /NB:/ This does require two passes over the data.+sum :: (RealFloat a, Ord a, Precise a, Foldable f) => f (Log a) -> Log a+sum xs = Log $ case Foldable.foldl' step1 None xs of+  None -> negInf+  Acc nm1 a+    | isInfinite a -> a+    | otherwise    -> a + log1p (Foldable.foldl' (step2 a) 0 xs + fromIntegral nm1)+  where+    step1 None      (Log x) = Acc 0 x+    step1 (Acc n y) (Log x) = Acc (n + 1) (max x y)+    step2 a r (Log x) = r + expm1 (x - a)+{-# INLINE sum #-}  instance (RealFloat a, Precise a) => Floating (Log a) where   pi = Log (log pi)