{-# LANGUAGE BangPatterns, DeriveDataTypeable, FlexibleContexts,
    MultiParamTypeClasses, TypeFamilies, CPP #-}
{-# OPTIONS_GHC -fno-warn-name-shadowing #-}
-- |
-- Module    : Numeric.Sum
-- Copyright : (c) 2014 Bryan O'Sullivan
-- License   : BSD3
--
-- Maintainer  : bos@serpentine.com
-- Stability   : experimental
-- Portability : portable
--
-- Functions for summing floating point numbers more accurately than
-- the naive 'Prelude.sum' function and its counterparts in the
-- @vector@ package and elsewhere.
--
-- When used with floating point numbers, in the worst case, the
-- 'Prelude.sum' function accumulates numeric error at a rate
-- proportional to the number of values being summed. The algorithms
-- in this module implement different methods of /compensated
-- summation/, which reduce the accumulation of numeric error so that
-- it either grows much more slowly than the number of inputs
-- (e.g. logarithmically), or remains constant.
module Numeric.Sum (
    -- * Summation type class
      Summation(..)
    , sumVector
    -- ** Usage
    -- $usage

    -- * Kahan-Babuška-Neumaier summation
    , KBNSum(..)
    , kbn

    -- * Order-2 Kahan-Babuška summation
    , KB2Sum(..)
    , kb2

    -- * Less desirable approaches

    -- ** Kahan summation
    , KahanSum(..)
    , kahan

    -- ** Pairwise summation
    , pairwiseSum

    -- * References
    -- $references
    ) where

import Control.Arrow ((***))
import Control.DeepSeq (NFData(..))
import Data.Bits (shiftR)
import Data.Data (Typeable, Data)
import Data.Monoid                  (Monoid(..))
#if MIN_VERSION_base(4,9,0)
import Data.Semigroup               (Semigroup(..))
#endif
import Data.Vector.Generic          (Vector(..), foldl')
-- Needed for GHC 7.2 & 7.4 to derive Unbox instances
import Control.Monad (liftM)
import Data.Vector.Generic.Mutable (MVector(..))

import qualified Data.Foldable as F
import qualified Data.Vector as V
import qualified Data.Vector.Generic as G
import qualified Data.Vector.Generic.Mutable as GM
import qualified Data.Vector.Unboxed as U

-- | A class for summation of floating point numbers.
class Summation s where
    -- | The identity for summation.
    zero :: s

    -- | Add a value to a sum.
    add  :: s -> Double -> s

    -- | Sum a collection of values.
    --
    -- Example:
    -- @foo = 'Numeric.Sum.sum' 'kbn' [1,2,3]@
    sum  :: (F.Foldable f) => (s -> Double) -> f Double -> Double
    sum  s -> Double
f = s -> Double
f (s -> Double) -> (f Double -> s) -> f Double -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (s -> Double -> s) -> s -> f Double -> s
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
F.foldl' s -> Double -> s
forall s. Summation s => s -> Double -> s
add s
forall s. Summation s => s
zero
    {-# INLINE sum #-}

instance Summation Double where
    zero :: Double
zero = Double
0
    add :: Double -> Double -> Double
add = Double -> Double -> Double
forall a. Num a => a -> a -> a
(+)

-- | Kahan summation. This is the least accurate of the compensated
-- summation methods.  In practice, it only beats naive summation for
-- inputs with large magnitude.  Kahan summation can be /less/
-- accurate than naive summation for small-magnitude inputs.
--
-- This summation method is included for completeness. Its use is not
-- recommended.  In practice, 'KBNSum' is both 30% faster and more
-- accurate.
data KahanSum = KahanSum {-# UNPACK #-} !Double {-# UNPACK #-} !Double
              deriving (KahanSum -> KahanSum -> Bool
(KahanSum -> KahanSum -> Bool)
-> (KahanSum -> KahanSum -> Bool) -> Eq KahanSum
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: KahanSum -> KahanSum -> Bool
$c/= :: KahanSum -> KahanSum -> Bool
== :: KahanSum -> KahanSum -> Bool
$c== :: KahanSum -> KahanSum -> Bool
Eq, Int -> KahanSum -> ShowS
[KahanSum] -> ShowS
KahanSum -> String
(Int -> KahanSum -> ShowS)
-> (KahanSum -> String) -> ([KahanSum] -> ShowS) -> Show KahanSum
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [KahanSum] -> ShowS
$cshowList :: [KahanSum] -> ShowS
show :: KahanSum -> String
$cshow :: KahanSum -> String
showsPrec :: Int -> KahanSum -> ShowS
$cshowsPrec :: Int -> KahanSum -> ShowS
Show, Typeable, Typeable KahanSum
DataType
Constr
Typeable KahanSum
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> KahanSum -> c KahanSum)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c KahanSum)
-> (KahanSum -> Constr)
-> (KahanSum -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c KahanSum))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KahanSum))
-> ((forall b. Data b => b -> b) -> KahanSum -> KahanSum)
-> (forall r r'.
    (r -> r' -> r)
    -> r -> (forall d. Data d => d -> r') -> KahanSum -> r)
-> (forall r r'.
    (r' -> r -> r)
    -> r -> (forall d. Data d => d -> r') -> KahanSum -> r)
-> (forall u. (forall d. Data d => d -> u) -> KahanSum -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> KahanSum -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> KahanSum -> m KahanSum)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> KahanSum -> m KahanSum)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> KahanSum -> m KahanSum)
-> Data KahanSum
KahanSum -> DataType
KahanSum -> Constr
(forall b. Data b => b -> b) -> KahanSum -> KahanSum
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KahanSum -> c KahanSum
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KahanSum
forall a.
Typeable a
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> a -> c a)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c a)
-> (a -> Constr)
-> (a -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c a))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c a))
-> ((forall b. Data b => b -> b) -> a -> a)
-> (forall r r'.
    (r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall r r'.
    (r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall u. (forall d. Data d => d -> u) -> a -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> a -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> Data a
forall u. Int -> (forall d. Data d => d -> u) -> KahanSum -> u
forall u. (forall d. Data d => d -> u) -> KahanSum -> [u]
forall r r'.
(r -> r' -> r)
-> r -> (forall d. Data d => d -> r') -> KahanSum -> r
forall r r'.
(r' -> r -> r)
-> r -> (forall d. Data d => d -> r') -> KahanSum -> r
forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KahanSum
forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KahanSum -> c KahanSum
forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c KahanSum)
forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KahanSum)
$cKahanSum :: Constr
$tKahanSum :: DataType
gmapMo :: (forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
$cgmapMo :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
gmapMp :: (forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
$cgmapMp :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
gmapM :: (forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
$cgmapM :: forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> KahanSum -> m KahanSum
gmapQi :: Int -> (forall d. Data d => d -> u) -> KahanSum -> u
$cgmapQi :: forall u. Int -> (forall d. Data d => d -> u) -> KahanSum -> u
gmapQ :: (forall d. Data d => d -> u) -> KahanSum -> [u]
$cgmapQ :: forall u. (forall d. Data d => d -> u) -> KahanSum -> [u]
gmapQr :: (r' -> r -> r)
-> r -> (forall d. Data d => d -> r') -> KahanSum -> r
$cgmapQr :: forall r r'.
(r' -> r -> r)
-> r -> (forall d. Data d => d -> r') -> KahanSum -> r
gmapQl :: (r -> r' -> r)
-> r -> (forall d. Data d => d -> r') -> KahanSum -> r
$cgmapQl :: forall r r'.
(r -> r' -> r)
-> r -> (forall d. Data d => d -> r') -> KahanSum -> r
gmapT :: (forall b. Data b => b -> b) -> KahanSum -> KahanSum
$cgmapT :: (forall b. Data b => b -> b) -> KahanSum -> KahanSum
dataCast2 :: (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KahanSum)
$cdataCast2 :: forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KahanSum)
dataCast1 :: (forall d. Data d => c (t d)) -> Maybe (c KahanSum)
$cdataCast1 :: forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c KahanSum)
dataTypeOf :: KahanSum -> DataType
$cdataTypeOf :: KahanSum -> DataType
toConstr :: KahanSum -> Constr
$ctoConstr :: KahanSum -> Constr
gunfold :: (forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KahanSum
$cgunfold :: forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KahanSum
gfoldl :: (forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KahanSum -> c KahanSum
$cgfoldl :: forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KahanSum -> c KahanSum
$cp1Data :: Typeable KahanSum
Data)

instance U.Unbox KahanSum
newtype instance U.MVector s KahanSum = MV_KahanSum (U.MVector s (Double, Double))
instance MVector U.MVector KahanSum where
  {-# INLINE GM.basicLength #-}
  {-# INLINE GM.basicUnsafeSlice #-}
  {-# INLINE basicOverlaps #-}
  {-# INLINE basicUnsafeNew #-}
  {-# INLINE basicInitialize #-}
  {-# INLINE basicUnsafeReplicate #-}
  {-# INLINE basicUnsafeRead #-}
  {-# INLINE basicUnsafeWrite #-}
  {-# INLINE basicClear #-}
  {-# INLINE basicSet #-}
  {-# INLINE GM.basicUnsafeCopy #-}
  {-# INLINE basicUnsafeMove #-}
  {-# INLINE basicUnsafeGrow #-}
  basicLength :: MVector s KahanSum -> Int
basicLength (MV_KahanSum mvec) = MVector s (Double, Double) -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
GM.basicLength MVector s (Double, Double)
mvec
  basicUnsafeSlice :: Int -> Int -> MVector s KahanSum -> MVector s KahanSum
basicUnsafeSlice Int
idx Int
len (MV_KahanSum mvec) = MVector s (Double, Double) -> MVector s KahanSum
forall s. MVector s (Double, Double) -> MVector s KahanSum
MV_KahanSum (Int
-> Int -> MVector s (Double, Double) -> MVector s (Double, Double)
forall (v :: * -> * -> *) a s.
MVector v a =>
Int -> Int -> v s a -> v s a
GM.basicUnsafeSlice Int
idx Int
len MVector s (Double, Double)
mvec)
  basicOverlaps :: MVector s KahanSum -> MVector s KahanSum -> Bool
basicOverlaps (MV_KahanSum mvec) (MV_KahanSum mvec') = MVector s (Double, Double) -> MVector s (Double, Double) -> Bool
forall (v :: * -> * -> *) a s.
MVector v a =>
v s a -> v s a -> Bool
basicOverlaps MVector s (Double, Double)
mvec MVector s (Double, Double)
mvec'
  basicUnsafeNew :: Int -> m (MVector (PrimState m) KahanSum)
basicUnsafeNew Int
len = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KahanSum
forall s. MVector s (Double, Double) -> MVector s KahanSum
MV_KahanSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KahanSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KahanSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Int -> m (MVector (PrimState m) (Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> m (v (PrimState m) a)
basicUnsafeNew Int
len
  basicInitialize :: MVector (PrimState m) KahanSum -> m ()
basicInitialize (MV_KahanSum mvec) = MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
basicInitialize MVector (PrimState m) (Double, Double)
mvec
  basicUnsafeReplicate :: Int -> KahanSum -> m (MVector (PrimState m) KahanSum)
basicUnsafeReplicate Int
len KahanSum
val = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KahanSum
forall s. MVector s (Double, Double) -> MVector s KahanSum
MV_KahanSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KahanSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KahanSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Int
-> (Double, Double) -> m (MVector (PrimState m) (Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> a -> m (v (PrimState m) a)
basicUnsafeReplicate Int
len ((\ (KahanSum Double
a Double
b) -> (Double
a, Double
b)) KahanSum
val)
  basicUnsafeRead :: MVector (PrimState m) KahanSum -> Int -> m KahanSum
basicUnsafeRead (MV_KahanSum mvec) Int
idx = (\ (Double
a, Double
b) -> Double -> Double -> KahanSum
KahanSum Double
a Double
b) ((Double, Double) -> KahanSum) -> m (Double, Double) -> m KahanSum
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` MVector (PrimState m) (Double, Double) -> Int -> m (Double, Double)
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m a
basicUnsafeRead MVector (PrimState m) (Double, Double)
mvec Int
idx
  basicUnsafeWrite :: MVector (PrimState m) KahanSum -> Int -> KahanSum -> m ()
basicUnsafeWrite (MV_KahanSum mvec) Int
idx KahanSum
val = MVector (PrimState m) (Double, Double)
-> Int -> (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> a -> m ()
basicUnsafeWrite MVector (PrimState m) (Double, Double)
mvec Int
idx ((\ (KahanSum Double
a Double
b) -> (Double
a, Double
b)) KahanSum
val)
  basicClear :: MVector (PrimState m) KahanSum -> m ()
basicClear (MV_KahanSum mvec) = MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
basicClear MVector (PrimState m) (Double, Double)
mvec
  basicSet :: MVector (PrimState m) KahanSum -> KahanSum -> m ()
basicSet (MV_KahanSum mvec) KahanSum
val = MVector (PrimState m) (Double, Double) -> (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> a -> m ()
basicSet MVector (PrimState m) (Double, Double)
mvec ((\ (KahanSum Double
a Double
b) -> (Double
a, Double
b)) KahanSum
val)
  basicUnsafeCopy :: MVector (PrimState m) KahanSum
-> MVector (PrimState m) KahanSum -> m ()
basicUnsafeCopy (MV_KahanSum mvec) (MV_KahanSum mvec') = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
GM.basicUnsafeCopy MVector (PrimState m) (Double, Double)
mvec MVector (PrimState m) (Double, Double)
mvec'
  basicUnsafeMove :: MVector (PrimState m) KahanSum
-> MVector (PrimState m) KahanSum -> m ()
basicUnsafeMove (MV_KahanSum mvec) (MV_KahanSum mvec') = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
basicUnsafeMove MVector (PrimState m) (Double, Double)
mvec MVector (PrimState m) (Double, Double)
mvec'
  basicUnsafeGrow :: MVector (PrimState m) KahanSum
-> Int -> m (MVector (PrimState m) KahanSum)
basicUnsafeGrow (MV_KahanSum mvec) Int
len = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KahanSum
forall s. MVector s (Double, Double) -> MVector s KahanSum
MV_KahanSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KahanSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KahanSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` MVector (PrimState m) (Double, Double)
-> Int -> m (MVector (PrimState m) (Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m (v (PrimState m) a)
basicUnsafeGrow MVector (PrimState m) (Double, Double)
mvec Int
len

newtype instance U.Vector KahanSum = V_KahanSum (U.Vector (Double, Double))
instance Vector U.Vector KahanSum where
  {-# INLINE basicUnsafeFreeze #-}
  {-# INLINE basicUnsafeThaw #-}
  {-# INLINE G.basicLength #-}
  {-# INLINE G.basicUnsafeSlice #-}
  {-# INLINE basicUnsafeIndexM #-}
  {-# INLINE G.basicUnsafeCopy #-}
  {-# INLINE elemseq #-}
  basicUnsafeFreeze :: Mutable Vector (PrimState m) KahanSum -> m (Vector KahanSum)
basicUnsafeFreeze (MV_KahanSum mvec) = Vector (Double, Double) -> Vector KahanSum
V_KahanSum (Vector (Double, Double) -> Vector KahanSum)
-> m (Vector (Double, Double)) -> m (Vector KahanSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Mutable Vector (PrimState m) (Double, Double)
-> m (Vector (Double, Double))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> m (v a)
basicUnsafeFreeze MVector (PrimState m) (Double, Double)
Mutable Vector (PrimState m) (Double, Double)
mvec
  basicUnsafeThaw :: Vector KahanSum -> m (Mutable Vector (PrimState m) KahanSum)
basicUnsafeThaw (V_KahanSum vec) = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KahanSum
forall s. MVector s (Double, Double) -> MVector s KahanSum
MV_KahanSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KahanSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KahanSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Vector (Double, Double)
-> m (Mutable Vector (PrimState m) (Double, Double))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
v a -> m (Mutable v (PrimState m) a)
basicUnsafeThaw Vector (Double, Double)
vec
  basicLength :: Vector KahanSum -> Int
basicLength (V_KahanSum vec) = Vector (Double, Double) -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.basicLength Vector (Double, Double)
vec
  basicUnsafeSlice :: Int -> Int -> Vector KahanSum -> Vector KahanSum
basicUnsafeSlice Int
idx Int
len (V_KahanSum vec) = Vector (Double, Double) -> Vector KahanSum
V_KahanSum (Int -> Int -> Vector (Double, Double) -> Vector (Double, Double)
forall (v :: * -> *) a. Vector v a => Int -> Int -> v a -> v a
G.basicUnsafeSlice Int
idx Int
len Vector (Double, Double)
vec)
  basicUnsafeIndexM :: Vector KahanSum -> Int -> m KahanSum
basicUnsafeIndexM (V_KahanSum vec) Int
idx = (\ (Double
a, Double
b) -> Double -> Double -> KahanSum
KahanSum Double
a Double
b) ((Double, Double) -> KahanSum) -> m (Double, Double) -> m KahanSum
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Vector (Double, Double) -> Int -> m (Double, Double)
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, Monad m) =>
v a -> Int -> m a
basicUnsafeIndexM Vector (Double, Double)
vec Int
idx
  basicUnsafeCopy :: Mutable Vector (PrimState m) KahanSum -> Vector KahanSum -> m ()
basicUnsafeCopy (MV_KahanSum mvec) (V_KahanSum vec) = Mutable Vector (PrimState m) (Double, Double)
-> Vector (Double, Double) -> m ()
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> v a -> m ()
G.basicUnsafeCopy MVector (PrimState m) (Double, Double)
Mutable Vector (PrimState m) (Double, Double)
mvec Vector (Double, Double)
vec
  elemseq :: Vector KahanSum -> KahanSum -> b -> b
elemseq (V_KahanSum vec) KahanSum
val = Vector (Double, Double) -> (Double, Double) -> b -> b
forall (v :: * -> *) a b. Vector v a => v a -> a -> b -> b
elemseq Vector (Double, Double)
vec ((\ (KahanSum Double
a Double
b) -> (Double
a, Double
b)) KahanSum
val)


instance Summation KahanSum where
    zero :: KahanSum
zero = Double -> Double -> KahanSum
KahanSum Double
0 Double
0
    add :: KahanSum -> Double -> KahanSum
add  = KahanSum -> Double -> KahanSum
kahanAdd

instance NFData KahanSum where
    rnf :: KahanSum -> ()
rnf !KahanSum
_ = ()

-- | @since 0.3.0.0
instance Monoid KahanSum where
  mempty :: KahanSum
mempty = KahanSum
forall s. Summation s => s
zero
  KahanSum
s mappend :: KahanSum -> KahanSum -> KahanSum
`mappend` KahanSum Double
s' Double
_ = KahanSum -> Double -> KahanSum
forall s. Summation s => s -> Double -> s
add KahanSum
s Double
s'

#if MIN_VERSION_base(4,9,0)
-- | @since 0.3.0.0
instance Semigroup KahanSum where
  <> :: KahanSum -> KahanSum -> KahanSum
(<>) = KahanSum -> KahanSum -> KahanSum
forall a. Monoid a => a -> a -> a
mappend
#endif

kahanAdd :: KahanSum -> Double -> KahanSum
kahanAdd :: KahanSum -> Double -> KahanSum
kahanAdd (KahanSum Double
sum Double
c) Double
x = Double -> Double -> KahanSum
KahanSum Double
sum' Double
c'
  where sum' :: Double
sum' = Double
sum Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
y
        c' :: Double
c'   = (Double
sum' Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
sum) Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
y
        y :: Double
y    = Double
x Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
c

-- | Return the result of a Kahan sum.
kahan :: KahanSum -> Double
kahan :: KahanSum -> Double
kahan (KahanSum Double
sum Double
_) = Double
sum

-- | Kahan-Babuška-Neumaier summation. This is a little more
-- computationally costly than plain Kahan summation, but is /always/
-- at least as accurate.
data KBNSum = KBNSum {-# UNPACK #-} !Double {-# UNPACK #-} !Double
            deriving (KBNSum -> KBNSum -> Bool
(KBNSum -> KBNSum -> Bool)
-> (KBNSum -> KBNSum -> Bool) -> Eq KBNSum
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: KBNSum -> KBNSum -> Bool
$c/= :: KBNSum -> KBNSum -> Bool
== :: KBNSum -> KBNSum -> Bool
$c== :: KBNSum -> KBNSum -> Bool
Eq, Int -> KBNSum -> ShowS
[KBNSum] -> ShowS
KBNSum -> String
(Int -> KBNSum -> ShowS)
-> (KBNSum -> String) -> ([KBNSum] -> ShowS) -> Show KBNSum
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [KBNSum] -> ShowS
$cshowList :: [KBNSum] -> ShowS
show :: KBNSum -> String
$cshow :: KBNSum -> String
showsPrec :: Int -> KBNSum -> ShowS
$cshowsPrec :: Int -> KBNSum -> ShowS
Show, Typeable, Typeable KBNSum
DataType
Constr
Typeable KBNSum
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> KBNSum -> c KBNSum)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c KBNSum)
-> (KBNSum -> Constr)
-> (KBNSum -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c KBNSum))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KBNSum))
-> ((forall b. Data b => b -> b) -> KBNSum -> KBNSum)
-> (forall r r'.
    (r -> r' -> r)
    -> r -> (forall d. Data d => d -> r') -> KBNSum -> r)
-> (forall r r'.
    (r' -> r -> r)
    -> r -> (forall d. Data d => d -> r') -> KBNSum -> r)
-> (forall u. (forall d. Data d => d -> u) -> KBNSum -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> KBNSum -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> KBNSum -> m KBNSum)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> KBNSum -> m KBNSum)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> KBNSum -> m KBNSum)
-> Data KBNSum
KBNSum -> DataType
KBNSum -> Constr
(forall b. Data b => b -> b) -> KBNSum -> KBNSum
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KBNSum -> c KBNSum
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KBNSum
forall a.
Typeable a
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> a -> c a)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c a)
-> (a -> Constr)
-> (a -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c a))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c a))
-> ((forall b. Data b => b -> b) -> a -> a)
-> (forall r r'.
    (r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall r r'.
    (r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall u. (forall d. Data d => d -> u) -> a -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> a -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> Data a
forall u. Int -> (forall d. Data d => d -> u) -> KBNSum -> u
forall u. (forall d. Data d => d -> u) -> KBNSum -> [u]
forall r r'.
(r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> KBNSum -> r
forall r r'.
(r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> KBNSum -> r
forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KBNSum
forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KBNSum -> c KBNSum
forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c KBNSum)
forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KBNSum)
$cKBNSum :: Constr
$tKBNSum :: DataType
gmapMo :: (forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
$cgmapMo :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
gmapMp :: (forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
$cgmapMp :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
gmapM :: (forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
$cgmapM :: forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> KBNSum -> m KBNSum
gmapQi :: Int -> (forall d. Data d => d -> u) -> KBNSum -> u
$cgmapQi :: forall u. Int -> (forall d. Data d => d -> u) -> KBNSum -> u
gmapQ :: (forall d. Data d => d -> u) -> KBNSum -> [u]
$cgmapQ :: forall u. (forall d. Data d => d -> u) -> KBNSum -> [u]
gmapQr :: (r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> KBNSum -> r
$cgmapQr :: forall r r'.
(r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> KBNSum -> r
gmapQl :: (r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> KBNSum -> r
$cgmapQl :: forall r r'.
(r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> KBNSum -> r
gmapT :: (forall b. Data b => b -> b) -> KBNSum -> KBNSum
$cgmapT :: (forall b. Data b => b -> b) -> KBNSum -> KBNSum
dataCast2 :: (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KBNSum)
$cdataCast2 :: forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KBNSum)
dataCast1 :: (forall d. Data d => c (t d)) -> Maybe (c KBNSum)
$cdataCast1 :: forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c KBNSum)
dataTypeOf :: KBNSum -> DataType
$cdataTypeOf :: KBNSum -> DataType
toConstr :: KBNSum -> Constr
$ctoConstr :: KBNSum -> Constr
gunfold :: (forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KBNSum
$cgunfold :: forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KBNSum
gfoldl :: (forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KBNSum -> c KBNSum
$cgfoldl :: forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KBNSum -> c KBNSum
$cp1Data :: Typeable KBNSum
Data)

instance U.Unbox KBNSum
newtype instance U.MVector s KBNSum = MV_KBNSum (U.MVector s (Double, Double))
instance MVector U.MVector KBNSum where
  {-# INLINE GM.basicLength #-}
  {-# INLINE GM.basicUnsafeSlice #-}
  {-# INLINE basicOverlaps #-}
  {-# INLINE basicUnsafeNew #-}
  {-# INLINE basicInitialize #-}
  {-# INLINE basicUnsafeReplicate #-}
  {-# INLINE basicUnsafeRead #-}
  {-# INLINE basicUnsafeWrite #-}
  {-# INLINE basicClear #-}
  {-# INLINE basicSet #-}
  {-# INLINE GM.basicUnsafeCopy #-}
  {-# INLINE basicUnsafeMove #-}
  {-# INLINE basicUnsafeGrow #-}
  basicLength :: MVector s KBNSum -> Int
basicLength (MV_KBNSum mvec) = MVector s (Double, Double) -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
GM.basicLength MVector s (Double, Double)
mvec
  basicUnsafeSlice :: Int -> Int -> MVector s KBNSum -> MVector s KBNSum
basicUnsafeSlice Int
idx Int
len (MV_KBNSum mvec) = MVector s (Double, Double) -> MVector s KBNSum
forall s. MVector s (Double, Double) -> MVector s KBNSum
MV_KBNSum (Int
-> Int -> MVector s (Double, Double) -> MVector s (Double, Double)
forall (v :: * -> * -> *) a s.
MVector v a =>
Int -> Int -> v s a -> v s a
GM.basicUnsafeSlice Int
idx Int
len MVector s (Double, Double)
mvec)
  basicOverlaps :: MVector s KBNSum -> MVector s KBNSum -> Bool
basicOverlaps (MV_KBNSum mvec) (MV_KBNSum mvec') = MVector s (Double, Double) -> MVector s (Double, Double) -> Bool
forall (v :: * -> * -> *) a s.
MVector v a =>
v s a -> v s a -> Bool
basicOverlaps MVector s (Double, Double)
mvec MVector s (Double, Double)
mvec'
  basicUnsafeNew :: Int -> m (MVector (PrimState m) KBNSum)
basicUnsafeNew Int
len = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KBNSum
forall s. MVector s (Double, Double) -> MVector s KBNSum
MV_KBNSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KBNSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KBNSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Int -> m (MVector (PrimState m) (Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> m (v (PrimState m) a)
basicUnsafeNew Int
len
  basicInitialize :: MVector (PrimState m) KBNSum -> m ()
basicInitialize (MV_KBNSum mvec) = MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
basicInitialize MVector (PrimState m) (Double, Double)
mvec
  basicUnsafeReplicate :: Int -> KBNSum -> m (MVector (PrimState m) KBNSum)
basicUnsafeReplicate Int
len KBNSum
val = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KBNSum
forall s. MVector s (Double, Double) -> MVector s KBNSum
MV_KBNSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KBNSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KBNSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Int
-> (Double, Double) -> m (MVector (PrimState m) (Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> a -> m (v (PrimState m) a)
basicUnsafeReplicate Int
len ((\ (KBNSum Double
a Double
b) -> (Double
a, Double
b)) KBNSum
val)
  basicUnsafeRead :: MVector (PrimState m) KBNSum -> Int -> m KBNSum
basicUnsafeRead (MV_KBNSum mvec) Int
idx = (\ (Double
a, Double
b) -> Double -> Double -> KBNSum
KBNSum Double
a Double
b) ((Double, Double) -> KBNSum) -> m (Double, Double) -> m KBNSum
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` MVector (PrimState m) (Double, Double) -> Int -> m (Double, Double)
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m a
basicUnsafeRead MVector (PrimState m) (Double, Double)
mvec Int
idx
  basicUnsafeWrite :: MVector (PrimState m) KBNSum -> Int -> KBNSum -> m ()
basicUnsafeWrite (MV_KBNSum mvec) Int
idx KBNSum
val = MVector (PrimState m) (Double, Double)
-> Int -> (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> a -> m ()
basicUnsafeWrite MVector (PrimState m) (Double, Double)
mvec Int
idx ((\ (KBNSum Double
a Double
b) -> (Double
a, Double
b)) KBNSum
val)
  basicClear :: MVector (PrimState m) KBNSum -> m ()
basicClear (MV_KBNSum mvec) = MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
basicClear MVector (PrimState m) (Double, Double)
mvec
  basicSet :: MVector (PrimState m) KBNSum -> KBNSum -> m ()
basicSet (MV_KBNSum mvec) KBNSum
val = MVector (PrimState m) (Double, Double) -> (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> a -> m ()
basicSet MVector (PrimState m) (Double, Double)
mvec ((\ (KBNSum Double
a Double
b) -> (Double
a, Double
b)) KBNSum
val)
  basicUnsafeCopy :: MVector (PrimState m) KBNSum
-> MVector (PrimState m) KBNSum -> m ()
basicUnsafeCopy (MV_KBNSum mvec) (MV_KBNSum mvec') = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
GM.basicUnsafeCopy MVector (PrimState m) (Double, Double)
mvec MVector (PrimState m) (Double, Double)
mvec'
  basicUnsafeMove :: MVector (PrimState m) KBNSum
-> MVector (PrimState m) KBNSum -> m ()
basicUnsafeMove (MV_KBNSum mvec) (MV_KBNSum mvec') = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) (Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
basicUnsafeMove MVector (PrimState m) (Double, Double)
mvec MVector (PrimState m) (Double, Double)
mvec'
  basicUnsafeGrow :: MVector (PrimState m) KBNSum
-> Int -> m (MVector (PrimState m) KBNSum)
basicUnsafeGrow (MV_KBNSum mvec) Int
len = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KBNSum
forall s. MVector s (Double, Double) -> MVector s KBNSum
MV_KBNSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KBNSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KBNSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` MVector (PrimState m) (Double, Double)
-> Int -> m (MVector (PrimState m) (Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m (v (PrimState m) a)
basicUnsafeGrow MVector (PrimState m) (Double, Double)
mvec Int
len

newtype instance U.Vector KBNSum = V_KBNSum (U.Vector (Double, Double))
instance Vector U.Vector KBNSum where
  {-# INLINE basicUnsafeFreeze #-}
  {-# INLINE basicUnsafeThaw #-}
  {-# INLINE G.basicLength #-}
  {-# INLINE G.basicUnsafeSlice #-}
  {-# INLINE basicUnsafeIndexM #-}
  {-# INLINE G.basicUnsafeCopy #-}
  {-# INLINE elemseq #-}
  basicUnsafeFreeze :: Mutable Vector (PrimState m) KBNSum -> m (Vector KBNSum)
basicUnsafeFreeze (MV_KBNSum mvec) = Vector (Double, Double) -> Vector KBNSum
V_KBNSum (Vector (Double, Double) -> Vector KBNSum)
-> m (Vector (Double, Double)) -> m (Vector KBNSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Mutable Vector (PrimState m) (Double, Double)
-> m (Vector (Double, Double))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> m (v a)
basicUnsafeFreeze MVector (PrimState m) (Double, Double)
Mutable Vector (PrimState m) (Double, Double)
mvec
  basicUnsafeThaw :: Vector KBNSum -> m (Mutable Vector (PrimState m) KBNSum)
basicUnsafeThaw (V_KBNSum vec) = MVector (PrimState m) (Double, Double)
-> MVector (PrimState m) KBNSum
forall s. MVector s (Double, Double) -> MVector s KBNSum
MV_KBNSum (MVector (PrimState m) (Double, Double)
 -> MVector (PrimState m) KBNSum)
-> m (MVector (PrimState m) (Double, Double))
-> m (MVector (PrimState m) KBNSum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Vector (Double, Double)
-> m (Mutable Vector (PrimState m) (Double, Double))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
v a -> m (Mutable v (PrimState m) a)
basicUnsafeThaw Vector (Double, Double)
vec
  basicLength :: Vector KBNSum -> Int
basicLength (V_KBNSum vec) = Vector (Double, Double) -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.basicLength Vector (Double, Double)
vec
  basicUnsafeSlice :: Int -> Int -> Vector KBNSum -> Vector KBNSum
basicUnsafeSlice Int
idx Int
len (V_KBNSum vec) = Vector (Double, Double) -> Vector KBNSum
V_KBNSum (Int -> Int -> Vector (Double, Double) -> Vector (Double, Double)
forall (v :: * -> *) a. Vector v a => Int -> Int -> v a -> v a
G.basicUnsafeSlice Int
idx Int
len Vector (Double, Double)
vec)
  basicUnsafeIndexM :: Vector KBNSum -> Int -> m KBNSum
basicUnsafeIndexM (V_KBNSum vec) Int
idx = (\ (Double
a, Double
b) -> Double -> Double -> KBNSum
KBNSum Double
a Double
b) ((Double, Double) -> KBNSum) -> m (Double, Double) -> m KBNSum
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Vector (Double, Double) -> Int -> m (Double, Double)
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, Monad m) =>
v a -> Int -> m a
basicUnsafeIndexM Vector (Double, Double)
vec Int
idx
  basicUnsafeCopy :: Mutable Vector (PrimState m) KBNSum -> Vector KBNSum -> m ()
basicUnsafeCopy (MV_KBNSum mvec) (V_KBNSum vec) = Mutable Vector (PrimState m) (Double, Double)
-> Vector (Double, Double) -> m ()
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> v a -> m ()
G.basicUnsafeCopy MVector (PrimState m) (Double, Double)
Mutable Vector (PrimState m) (Double, Double)
mvec Vector (Double, Double)
vec
  elemseq :: Vector KBNSum -> KBNSum -> b -> b
elemseq (V_KBNSum vec) KBNSum
val = Vector (Double, Double) -> (Double, Double) -> b -> b
forall (v :: * -> *) a b. Vector v a => v a -> a -> b -> b
elemseq Vector (Double, Double)
vec ((\ (KBNSum Double
a Double
b) -> (Double
a, Double
b)) KBNSum
val)


instance Summation KBNSum where
    zero :: KBNSum
zero = Double -> Double -> KBNSum
KBNSum Double
0 Double
0
    add :: KBNSum -> Double -> KBNSum
add  = KBNSum -> Double -> KBNSum
kbnAdd

instance NFData KBNSum where
    rnf :: KBNSum -> ()
rnf !KBNSum
_ = ()

-- | @since 0.3.0.0
instance Monoid KBNSum where
  mempty :: KBNSum
mempty = KBNSum
forall s. Summation s => s
zero
  KBNSum
s mappend :: KBNSum -> KBNSum -> KBNSum
`mappend` KBNSum Double
s' Double
c' = KBNSum -> Double -> KBNSum
forall s. Summation s => s -> Double -> s
add (KBNSum -> Double -> KBNSum
forall s. Summation s => s -> Double -> s
add KBNSum
s Double
s') Double
c'

#if MIN_VERSION_base(4,9,0)
-- | @since 0.3.0.0
instance Semigroup KBNSum where
  <> :: KBNSum -> KBNSum -> KBNSum
(<>) = KBNSum -> KBNSum -> KBNSum
forall a. Monoid a => a -> a -> a
mappend
#endif

kbnAdd :: KBNSum -> Double -> KBNSum
kbnAdd :: KBNSum -> Double -> KBNSum
kbnAdd (KBNSum Double
sum Double
c) Double
x = Double -> Double -> KBNSum
KBNSum Double
sum' Double
c'
  where c' :: Double
c' | Double -> Double
forall a. Num a => a -> a
abs Double
sum Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
>= Double -> Double
forall a. Num a => a -> a
abs Double
x = Double
c Double -> Double -> Double
forall a. Num a => a -> a -> a
+ ((Double
sum Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
sum') Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
x)
           | Bool
otherwise        = Double
c Double -> Double -> Double
forall a. Num a => a -> a -> a
+ ((Double
x Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
sum') Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
sum)
        sum' :: Double
sum'                  = Double
sum Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
x

-- | Return the result of a Kahan-Babuška-Neumaier sum.
kbn :: KBNSum -> Double
kbn :: KBNSum -> Double
kbn (KBNSum Double
sum Double
c) = Double
sum Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
c

-- | Second-order Kahan-Babuška summation.  This is more
-- computationally costly than Kahan-Babuška-Neumaier summation,
-- running at about a third the speed.  Its advantage is that it can
-- lose less precision (in admittedly obscure cases).
--
-- This method compensates for error in both the sum and the
-- first-order compensation term, hence the use of \"second order\" in
-- the name.
data KB2Sum = KB2Sum {-# UNPACK #-} !Double
                     {-# UNPACK #-} !Double
                     {-# UNPACK #-} !Double
            deriving (KB2Sum -> KB2Sum -> Bool
(KB2Sum -> KB2Sum -> Bool)
-> (KB2Sum -> KB2Sum -> Bool) -> Eq KB2Sum
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: KB2Sum -> KB2Sum -> Bool
$c/= :: KB2Sum -> KB2Sum -> Bool
== :: KB2Sum -> KB2Sum -> Bool
$c== :: KB2Sum -> KB2Sum -> Bool
Eq, Int -> KB2Sum -> ShowS
[KB2Sum] -> ShowS
KB2Sum -> String
(Int -> KB2Sum -> ShowS)
-> (KB2Sum -> String) -> ([KB2Sum] -> ShowS) -> Show KB2Sum
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [KB2Sum] -> ShowS
$cshowList :: [KB2Sum] -> ShowS
show :: KB2Sum -> String
$cshow :: KB2Sum -> String
showsPrec :: Int -> KB2Sum -> ShowS
$cshowsPrec :: Int -> KB2Sum -> ShowS
Show, Typeable, Typeable KB2Sum
DataType
Constr
Typeable KB2Sum
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> KB2Sum -> c KB2Sum)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c KB2Sum)
-> (KB2Sum -> Constr)
-> (KB2Sum -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c KB2Sum))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KB2Sum))
-> ((forall b. Data b => b -> b) -> KB2Sum -> KB2Sum)
-> (forall r r'.
    (r -> r' -> r)
    -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r)
-> (forall r r'.
    (r' -> r -> r)
    -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r)
-> (forall u. (forall d. Data d => d -> u) -> KB2Sum -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> KB2Sum -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum)
-> Data KB2Sum
KB2Sum -> DataType
KB2Sum -> Constr
(forall b. Data b => b -> b) -> KB2Sum -> KB2Sum
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KB2Sum -> c KB2Sum
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KB2Sum
forall a.
Typeable a
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> a -> c a)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c a)
-> (a -> Constr)
-> (a -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c a))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c a))
-> ((forall b. Data b => b -> b) -> a -> a)
-> (forall r r'.
    (r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall r r'.
    (r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall u. (forall d. Data d => d -> u) -> a -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> a -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> Data a
forall u. Int -> (forall d. Data d => d -> u) -> KB2Sum -> u
forall u. (forall d. Data d => d -> u) -> KB2Sum -> [u]
forall r r'.
(r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r
forall r r'.
(r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r
forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KB2Sum
forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KB2Sum -> c KB2Sum
forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c KB2Sum)
forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KB2Sum)
$cKB2Sum :: Constr
$tKB2Sum :: DataType
gmapMo :: (forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
$cgmapMo :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
gmapMp :: (forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
$cgmapMp :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
gmapM :: (forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
$cgmapM :: forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> KB2Sum -> m KB2Sum
gmapQi :: Int -> (forall d. Data d => d -> u) -> KB2Sum -> u
$cgmapQi :: forall u. Int -> (forall d. Data d => d -> u) -> KB2Sum -> u
gmapQ :: (forall d. Data d => d -> u) -> KB2Sum -> [u]
$cgmapQ :: forall u. (forall d. Data d => d -> u) -> KB2Sum -> [u]
gmapQr :: (r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r
$cgmapQr :: forall r r'.
(r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r
gmapQl :: (r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r
$cgmapQl :: forall r r'.
(r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> KB2Sum -> r
gmapT :: (forall b. Data b => b -> b) -> KB2Sum -> KB2Sum
$cgmapT :: (forall b. Data b => b -> b) -> KB2Sum -> KB2Sum
dataCast2 :: (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KB2Sum)
$cdataCast2 :: forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c KB2Sum)
dataCast1 :: (forall d. Data d => c (t d)) -> Maybe (c KB2Sum)
$cdataCast1 :: forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c KB2Sum)
dataTypeOf :: KB2Sum -> DataType
$cdataTypeOf :: KB2Sum -> DataType
toConstr :: KB2Sum -> Constr
$ctoConstr :: KB2Sum -> Constr
gunfold :: (forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KB2Sum
$cgunfold :: forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c KB2Sum
gfoldl :: (forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KB2Sum -> c KB2Sum
$cgfoldl :: forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> KB2Sum -> c KB2Sum
$cp1Data :: Typeable KB2Sum
Data)

instance U.Unbox KB2Sum
newtype instance U.MVector s KB2Sum = MV_KB2Sum (U.MVector s (Double, Double, Double))
instance MVector U.MVector KB2Sum where
  {-# INLINE GM.basicLength #-}
  {-# INLINE GM.basicUnsafeSlice #-}
  {-# INLINE basicOverlaps #-}
  {-# INLINE basicUnsafeNew #-}
  {-# INLINE basicInitialize #-}
  {-# INLINE basicUnsafeReplicate #-}
  {-# INLINE basicUnsafeRead #-}
  {-# INLINE basicUnsafeWrite #-}
  {-# INLINE basicClear #-}
  {-# INLINE basicSet #-}
  {-# INLINE GM.basicUnsafeCopy #-}
  {-# INLINE basicUnsafeMove #-}
  {-# INLINE basicUnsafeGrow #-}
  basicLength :: MVector s KB2Sum -> Int
basicLength (MV_KB2Sum mvec) = MVector s (Double, Double, Double) -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
GM.basicLength MVector s (Double, Double, Double)
mvec
  basicUnsafeSlice :: Int -> Int -> MVector s KB2Sum -> MVector s KB2Sum
basicUnsafeSlice Int
idx Int
len (MV_KB2Sum mvec) = MVector s (Double, Double, Double) -> MVector s KB2Sum
forall s. MVector s (Double, Double, Double) -> MVector s KB2Sum
MV_KB2Sum (Int
-> Int
-> MVector s (Double, Double, Double)
-> MVector s (Double, Double, Double)
forall (v :: * -> * -> *) a s.
MVector v a =>
Int -> Int -> v s a -> v s a
GM.basicUnsafeSlice Int
idx Int
len MVector s (Double, Double, Double)
mvec)
  basicOverlaps :: MVector s KB2Sum -> MVector s KB2Sum -> Bool
basicOverlaps (MV_KB2Sum mvec) (MV_KB2Sum mvec') = MVector s (Double, Double, Double)
-> MVector s (Double, Double, Double) -> Bool
forall (v :: * -> * -> *) a s.
MVector v a =>
v s a -> v s a -> Bool
basicOverlaps MVector s (Double, Double, Double)
mvec MVector s (Double, Double, Double)
mvec'
  basicUnsafeNew :: Int -> m (MVector (PrimState m) KB2Sum)
basicUnsafeNew Int
len = MVector (PrimState m) (Double, Double, Double)
-> MVector (PrimState m) KB2Sum
forall s. MVector s (Double, Double, Double) -> MVector s KB2Sum
MV_KB2Sum (MVector (PrimState m) (Double, Double, Double)
 -> MVector (PrimState m) KB2Sum)
-> m (MVector (PrimState m) (Double, Double, Double))
-> m (MVector (PrimState m) KB2Sum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Int -> m (MVector (PrimState m) (Double, Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> m (v (PrimState m) a)
basicUnsafeNew Int
len
  basicInitialize :: MVector (PrimState m) KB2Sum -> m ()
basicInitialize (MV_KB2Sum mvec) = MVector (PrimState m) (Double, Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
basicInitialize MVector (PrimState m) (Double, Double, Double)
mvec
  basicUnsafeReplicate :: Int -> KB2Sum -> m (MVector (PrimState m) KB2Sum)
basicUnsafeReplicate Int
len KB2Sum
val = MVector (PrimState m) (Double, Double, Double)
-> MVector (PrimState m) KB2Sum
forall s. MVector s (Double, Double, Double) -> MVector s KB2Sum
MV_KB2Sum (MVector (PrimState m) (Double, Double, Double)
 -> MVector (PrimState m) KB2Sum)
-> m (MVector (PrimState m) (Double, Double, Double))
-> m (MVector (PrimState m) KB2Sum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Int
-> (Double, Double, Double)
-> m (MVector (PrimState m) (Double, Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> a -> m (v (PrimState m) a)
basicUnsafeReplicate Int
len ((\ (KB2Sum Double
a Double
b Double
c) -> (Double
a, Double
b, Double
c)) KB2Sum
val)
  basicUnsafeRead :: MVector (PrimState m) KB2Sum -> Int -> m KB2Sum
basicUnsafeRead (MV_KB2Sum mvec) Int
idx = (\ (Double
a, Double
b, Double
c) -> Double -> Double -> Double -> KB2Sum
KB2Sum Double
a Double
b Double
c) ((Double, Double, Double) -> KB2Sum)
-> m (Double, Double, Double) -> m KB2Sum
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` MVector (PrimState m) (Double, Double, Double)
-> Int -> m (Double, Double, Double)
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m a
basicUnsafeRead MVector (PrimState m) (Double, Double, Double)
mvec Int
idx
  basicUnsafeWrite :: MVector (PrimState m) KB2Sum -> Int -> KB2Sum -> m ()
basicUnsafeWrite (MV_KB2Sum mvec) Int
idx KB2Sum
val = MVector (PrimState m) (Double, Double, Double)
-> Int -> (Double, Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> a -> m ()
basicUnsafeWrite MVector (PrimState m) (Double, Double, Double)
mvec Int
idx ((\ (KB2Sum Double
a Double
b Double
c) -> (Double
a, Double
b, Double
c)) KB2Sum
val)
  basicClear :: MVector (PrimState m) KB2Sum -> m ()
basicClear (MV_KB2Sum mvec) = MVector (PrimState m) (Double, Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
basicClear MVector (PrimState m) (Double, Double, Double)
mvec
  basicSet :: MVector (PrimState m) KB2Sum -> KB2Sum -> m ()
basicSet (MV_KB2Sum mvec) KB2Sum
val = MVector (PrimState m) (Double, Double, Double)
-> (Double, Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> a -> m ()
basicSet MVector (PrimState m) (Double, Double, Double)
mvec ((\ (KB2Sum Double
a Double
b Double
c) -> (Double
a, Double
b, Double
c)) KB2Sum
val)
  basicUnsafeCopy :: MVector (PrimState m) KB2Sum
-> MVector (PrimState m) KB2Sum -> m ()
basicUnsafeCopy (MV_KB2Sum mvec) (MV_KB2Sum mvec') = MVector (PrimState m) (Double, Double, Double)
-> MVector (PrimState m) (Double, Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
GM.basicUnsafeCopy MVector (PrimState m) (Double, Double, Double)
mvec MVector (PrimState m) (Double, Double, Double)
mvec'
  basicUnsafeMove :: MVector (PrimState m) KB2Sum
-> MVector (PrimState m) KB2Sum -> m ()
basicUnsafeMove (MV_KB2Sum mvec) (MV_KB2Sum mvec') = MVector (PrimState m) (Double, Double, Double)
-> MVector (PrimState m) (Double, Double, Double) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
basicUnsafeMove MVector (PrimState m) (Double, Double, Double)
mvec MVector (PrimState m) (Double, Double, Double)
mvec'
  basicUnsafeGrow :: MVector (PrimState m) KB2Sum
-> Int -> m (MVector (PrimState m) KB2Sum)
basicUnsafeGrow (MV_KB2Sum mvec) Int
len = MVector (PrimState m) (Double, Double, Double)
-> MVector (PrimState m) KB2Sum
forall s. MVector s (Double, Double, Double) -> MVector s KB2Sum
MV_KB2Sum (MVector (PrimState m) (Double, Double, Double)
 -> MVector (PrimState m) KB2Sum)
-> m (MVector (PrimState m) (Double, Double, Double))
-> m (MVector (PrimState m) KB2Sum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` MVector (PrimState m) (Double, Double, Double)
-> Int -> m (MVector (PrimState m) (Double, Double, Double))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m (v (PrimState m) a)
basicUnsafeGrow MVector (PrimState m) (Double, Double, Double)
mvec Int
len

newtype instance U.Vector KB2Sum = V_KB2Sum (U.Vector (Double, Double, Double))
instance Vector U.Vector KB2Sum where
  {-# INLINE basicUnsafeFreeze #-}
  {-# INLINE basicUnsafeThaw #-}
  {-# INLINE G.basicLength #-}
  {-# INLINE G.basicUnsafeSlice #-}
  {-# INLINE basicUnsafeIndexM #-}
  {-# INLINE G.basicUnsafeCopy #-}
  {-# INLINE elemseq #-}
  basicUnsafeFreeze :: Mutable Vector (PrimState m) KB2Sum -> m (Vector KB2Sum)
basicUnsafeFreeze (MV_KB2Sum mvec) = Vector (Double, Double, Double) -> Vector KB2Sum
V_KB2Sum (Vector (Double, Double, Double) -> Vector KB2Sum)
-> m (Vector (Double, Double, Double)) -> m (Vector KB2Sum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Mutable Vector (PrimState m) (Double, Double, Double)
-> m (Vector (Double, Double, Double))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> m (v a)
basicUnsafeFreeze MVector (PrimState m) (Double, Double, Double)
Mutable Vector (PrimState m) (Double, Double, Double)
mvec
  basicUnsafeThaw :: Vector KB2Sum -> m (Mutable Vector (PrimState m) KB2Sum)
basicUnsafeThaw (V_KB2Sum vec) = MVector (PrimState m) (Double, Double, Double)
-> MVector (PrimState m) KB2Sum
forall s. MVector s (Double, Double, Double) -> MVector s KB2Sum
MV_KB2Sum (MVector (PrimState m) (Double, Double, Double)
 -> MVector (PrimState m) KB2Sum)
-> m (MVector (PrimState m) (Double, Double, Double))
-> m (MVector (PrimState m) KB2Sum)
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Vector (Double, Double, Double)
-> m (Mutable Vector (PrimState m) (Double, Double, Double))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
v a -> m (Mutable v (PrimState m) a)
basicUnsafeThaw Vector (Double, Double, Double)
vec
  basicLength :: Vector KB2Sum -> Int
basicLength (V_KB2Sum vec) = Vector (Double, Double, Double) -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.basicLength Vector (Double, Double, Double)
vec
  basicUnsafeSlice :: Int -> Int -> Vector KB2Sum -> Vector KB2Sum
basicUnsafeSlice Int
idx Int
len (V_KB2Sum vec) = Vector (Double, Double, Double) -> Vector KB2Sum
V_KB2Sum (Int
-> Int
-> Vector (Double, Double, Double)
-> Vector (Double, Double, Double)
forall (v :: * -> *) a. Vector v a => Int -> Int -> v a -> v a
G.basicUnsafeSlice Int
idx Int
len Vector (Double, Double, Double)
vec)
  basicUnsafeIndexM :: Vector KB2Sum -> Int -> m KB2Sum
basicUnsafeIndexM (V_KB2Sum vec) Int
idx = (\ (Double
a, Double
b, Double
c) -> Double -> Double -> Double -> KB2Sum
KB2Sum Double
a Double
b Double
c) ((Double, Double, Double) -> KB2Sum)
-> m (Double, Double, Double) -> m KB2Sum
forall (m :: * -> *) a1 r. Monad m => (a1 -> r) -> m a1 -> m r
`liftM` Vector (Double, Double, Double)
-> Int -> m (Double, Double, Double)
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, Monad m) =>
v a -> Int -> m a
basicUnsafeIndexM Vector (Double, Double, Double)
vec Int
idx
  basicUnsafeCopy :: Mutable Vector (PrimState m) KB2Sum -> Vector KB2Sum -> m ()
basicUnsafeCopy (MV_KB2Sum mvec) (V_KB2Sum vec) = Mutable Vector (PrimState m) (Double, Double, Double)
-> Vector (Double, Double, Double) -> m ()
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> v a -> m ()
G.basicUnsafeCopy MVector (PrimState m) (Double, Double, Double)
Mutable Vector (PrimState m) (Double, Double, Double)
mvec Vector (Double, Double, Double)
vec
  elemseq :: Vector KB2Sum -> KB2Sum -> b -> b
elemseq (V_KB2Sum vec) KB2Sum
val = Vector (Double, Double, Double)
-> (Double, Double, Double) -> b -> b
forall (v :: * -> *) a b. Vector v a => v a -> a -> b -> b
elemseq Vector (Double, Double, Double)
vec ((\ (KB2Sum Double
a Double
b Double
c) -> (Double
a, Double
b, Double
c)) KB2Sum
val)

instance Summation KB2Sum where
    zero :: KB2Sum
zero = Double -> Double -> Double -> KB2Sum
KB2Sum Double
0 Double
0 Double
0
    add :: KB2Sum -> Double -> KB2Sum
add  = KB2Sum -> Double -> KB2Sum
kb2Add

instance NFData KB2Sum where
    rnf :: KB2Sum -> ()
rnf !KB2Sum
_ = ()

-- | @since 0.3.0.0
instance Monoid KB2Sum where
  mempty :: KB2Sum
mempty = KB2Sum
forall s. Summation s => s
zero
  KB2Sum
s mappend :: KB2Sum -> KB2Sum -> KB2Sum
`mappend` KB2Sum Double
s' Double
c' Double
cc' = KB2Sum -> Double -> KB2Sum
forall s. Summation s => s -> Double -> s
add (KB2Sum -> Double -> KB2Sum
forall s. Summation s => s -> Double -> s
add (KB2Sum -> Double -> KB2Sum
forall s. Summation s => s -> Double -> s
add KB2Sum
s Double
s') Double
c') Double
cc'

#if MIN_VERSION_base(4,9,0)
-- | @since 0.3.0.0
instance Semigroup KB2Sum where
  <> :: KB2Sum -> KB2Sum -> KB2Sum
(<>) = KB2Sum -> KB2Sum -> KB2Sum
forall a. Monoid a => a -> a -> a
mappend
#endif


kb2Add :: KB2Sum -> Double -> KB2Sum
kb2Add :: KB2Sum -> Double -> KB2Sum
kb2Add (KB2Sum Double
sum Double
c Double
cc) Double
x = Double -> Double -> Double -> KB2Sum
KB2Sum Double
sum' Double
c' Double
cc'
  where sum' :: Double
sum'                 = Double
sum Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
x
        c' :: Double
c'                   = Double
c Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
k
        cc' :: Double
cc' | Double -> Double
forall a. Num a => a -> a
abs Double
c Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
>= Double -> Double
forall a. Num a => a -> a
abs Double
k = Double
cc Double -> Double -> Double
forall a. Num a => a -> a -> a
+ ((Double
c Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
c') Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
k)
            | Bool
otherwise      = Double
cc Double -> Double -> Double
forall a. Num a => a -> a -> a
+ ((Double
k Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
c') Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
c)
        k :: Double
k | Double -> Double
forall a. Num a => a -> a
abs Double
sum Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
>= Double -> Double
forall a. Num a => a -> a
abs Double
x = (Double
sum Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
sum') Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
x
          | Bool
otherwise        = (Double
x Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
sum') Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
sum

-- | Return the result of an order-2 Kahan-Babuška sum.
kb2 :: KB2Sum -> Double
kb2 :: KB2Sum -> Double
kb2 (KB2Sum Double
sum Double
c Double
cc) = Double
sum Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
c Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
cc

-- | /O(n)/ Sum a vector of values.
sumVector :: (Vector v Double, Summation s) =>
             (s -> Double) -> v Double -> Double
sumVector :: (s -> Double) -> v Double -> Double
sumVector s -> Double
f = s -> Double
f (s -> Double) -> (v Double -> s) -> v Double -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (s -> Double -> s) -> s -> v Double -> s
forall (v :: * -> *) b a.
Vector v b =>
(a -> b -> a) -> a -> v b -> a
foldl' s -> Double -> s
forall s. Summation s => s -> Double -> s
add s
forall s. Summation s => s
zero
{-# INLINE sumVector #-}

-- | /O(n)/ Sum a vector of values using pairwise summation.
--
-- This approach is perhaps 10% faster than 'KBNSum', but has poorer
-- bounds on its error growth.  Instead of having roughly constant
-- error regardless of the size of the input vector, in the worst case
-- its accumulated error grows with /O(log n)/.
pairwiseSum :: (Vector v Double) => v Double -> Double
pairwiseSum :: v Double -> Double
pairwiseSum v Double
v
  | Int
len Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
256 = v Double -> Double
forall (v :: * -> *) a. (Vector v a, Num a) => v a -> a
G.sum v Double
v
  | Bool
otherwise  = (Double -> Double -> Double) -> (Double, Double) -> Double
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry Double -> Double -> Double
forall a. Num a => a -> a -> a
(+) ((Double, Double) -> Double)
-> (v Double -> (Double, Double)) -> v Double -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (v Double -> Double
forall (v :: * -> *). Vector v Double => v Double -> Double
pairwiseSum (v Double -> Double)
-> (v Double -> Double) -> (v Double, v Double) -> (Double, Double)
forall (a :: * -> * -> *) b c b' c'.
Arrow a =>
a b c -> a b' c' -> a (b, b') (c, c')
*** v Double -> Double
forall (v :: * -> *). Vector v Double => v Double -> Double
pairwiseSum) ((v Double, v Double) -> (Double, Double))
-> (v Double -> (v Double, v Double))
-> v Double
-> (Double, Double)
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
                 Int -> v Double -> (v Double, v Double)
forall (v :: * -> *) a. Vector v a => Int -> v a -> (v a, v a)
G.splitAt (Int
len Int -> Int -> Int
forall a. Bits a => a -> Int -> a
`shiftR` Int
1) (v Double -> Double) -> v Double -> Double
forall a b. (a -> b) -> a -> b
$ v Double
v
  where len :: Int
len = v Double -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v Double
v
{-# SPECIALIZE pairwiseSum :: V.Vector Double -> Double #-}
{-# SPECIALIZE pairwiseSum :: U.Vector Double -> Double #-}

-- $usage
--
-- Most of these summation algorithms are intended to be used via the
-- 'Summation' typeclass interface. Explicit type annotations should
-- not be necessary, as the use of a function such as 'kbn' or 'kb2'
-- to extract the final sum out of a 'Summation' instance gives the
-- compiler enough information to determine the precise type of
-- summation algorithm to use.
--
-- As an example, here is a (somewhat silly) function that manually
-- computes the sum of elements in a list.
--
-- @
-- sillySumList :: [Double] -> Double
-- sillySumList = loop 'zero'
--   where loop s []     = 'kbn' s
--         loop s (x:xs) = 'seq' s' loop s' xs
--           where s'    = 'add' s x
-- @
--
-- In most instances, you can simply use the much more general 'Numeric.Sum.sum'
-- function instead of writing a summation function by hand.
--
-- @
-- -- Avoid ambiguity around which sum function we are using.
-- import Prelude hiding (sum)
-- --
-- betterSumList :: [Double] -> Double
-- betterSumList xs = 'Numeric.Sum.sum' 'kbn' xs
-- @

-- Note well the use of 'seq' in the example above to force the
-- evaluation of intermediate values.  If you must write a summation
-- function by hand, and you forget to evaluate the intermediate
-- values, you are likely to incur a space leak.
--
-- Here is an example of how to compute a prefix sum in which the
-- intermediate values are as accurate as possible.
--
-- @
-- prefixSum :: [Double] -> [Double]
-- prefixSum xs = map 'kbn' . 'scanl' 'add' 'zero' $ xs
-- @

-- $references
--
-- * Kahan, W. (1965), Further remarks on reducing truncation
--   errors. /Communications of the ACM/ 8(1):40.
--
-- * Neumaier, A. (1974), Rundungsfehleranalyse einiger Verfahren zur
--   Summation endlicher Summen.
--   /Zeitschrift für Angewandte Mathematik und Mechanik/ 54:39–51.
--
-- * Klein, A. (2006), A Generalized
--   Kahan-Babuška-Summation-Algorithm. /Computing/ 76(3):279-293.
--
-- * Higham, N.J. (1993), The accuracy of floating point
--   summation. /SIAM Journal on Scientific Computing/ 14(4):783–799.