{-# LANGUAGE CPP #-}
{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE StandaloneDeriving #-}
{-# LANGUAGE DeriveDataTypeable #-}
{-# LANGUAGE PatternGuards #-}
{-# LANGUAGE Trustworthy #-}

#ifndef MIN_VERSION_vector
#define MIN_VERSION_vector(x,y,z) 1
#endif

--------------------------------------------------------------------
-- |
-- Copyright :  (c) Edward Kmett 2013-2015
-- License   :  BSD3
-- Maintainer:  Edward Kmett <ekmett@gmail.com>
-- Stability :  experimental
-- Portability: non-portable
--
-- This module provides a fairly extensive API for compensated
-- floating point arithmetic based on Knuth's error free
-- transformation, various algorithms by Ogita, Rump and Oishi,
-- Hida, Li and Bailey, Kahan summation, etc. with custom compensated
-- arithmetic circuits to do multiplication, division, etc. of compensated
-- numbers.
--
-- In general if @a@ has x bits of significand, @Compensated a@ gives
-- you twice that. You can iterate this construction for arbitrary
-- precision.
--
-- References:
--
-- * <http://web.mit.edu/tabbott/Public/quaddouble-debian/qd-2.3.4-old/docs/qd.pdf>
--
-- * <http://www.ti3.tuhh.de/paper/rump/OgRuOi05.pdf>
--
-- * Donald Knuth's \"The Art of Computer Programming, Volume 2: Seminumerical Algorithms\"
--
-- * <http://en.wikipedia.org/wiki/Kahan_summation_algorithm>
--------------------------------------------------------------------
module Numeric.Compensated
  ( Compensable(..)
  , _Compensated
  , Overcompensated
  , primal
  , residual
  , uncompensated
  , fadd
  -- * lifting scalars
  , add, times, squared, divide, split
  , kahan, (+^), (*^)
  -- * compensated operators
  , square
  ) where

#if __GLASGOW_HASKELL__ < 710
import Control.Applicative
#endif
import Control.Lens as L
import Control.DeepSeq
import Data.Binary as Binary
import Data.Data
import Data.Foldable as Foldable
import Data.Function (on)
import Data.Hashable
import Data.Ratio
import Data.SafeCopy
import Data.Serialize as Serialize
import Data.Bytes.Serial as Bytes
#if !(MIN_VERSION_base(4,11,0))
import Data.Semigroup
#endif
import Data.Vector.Unboxed as U
import Data.Vector.Generic as G
import Data.Vector.Generic.Mutable as M
import Foreign.Ptr
import Foreign.Storable
import Numeric (Floating(..))
import Text.Read as T
import Text.Show as T

-- $setup
-- >>> import Numeric.Compensated

-- | @'add' a b k@ computes @k x y@ such that
--
-- > x + y = a + b
-- > x = fl(a + b)
--
-- Which is to say that @x@ is the floating point image of @(a + b)@ and
-- @y@ stores the residual error term.
add :: Num a => a -> a -> (a -> a -> r) -> r
add :: a -> a -> (a -> a -> r) -> r
add a
a a
b a -> a -> r
k = a -> a -> r
k a
x a
y where
  x :: a
x = a
a a -> a -> a
forall a. Num a => a -> a -> a
+ a
b
  z :: a
z = a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
a
  y :: a
y = (a
a a -> a -> a
forall a. Num a => a -> a -> a
- (a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
z)) a -> a -> a
forall a. Num a => a -> a -> a
+ (a
b a -> a -> a
forall a. Num a => a -> a -> a
- a
z)
{-# INLINE add #-}

-- | @'fadd' a b k@ computes @k x y@ such that
--
-- > x + y = a + b
-- > x = fl(a + b)
--
-- but only under the assumption that @'abs' a '>=' 'abs' b@. If you
-- aren't sure, use 'add'.
--
-- Which is to say that @x@ is the floating point image of @(a + b)@ and
-- @y@ stores the residual error term.
fadd :: Num a => a -> a -> (a -> a -> r) -> r
fadd :: a -> a -> (a -> a -> r) -> r
fadd a
a a
b a -> a -> r
k = a -> a -> r
k a
x (a
b a -> a -> a
forall a. Num a => a -> a -> a
- (a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
a)) where
  x :: a
x = a
a a -> a -> a
forall a. Num a => a -> a -> a
+ a
b
{-# INLINE fadd #-}

-- | @'times' a b k@ computes @k x y@ such that
--
-- > x + y = a * b
-- > x = fl(a * b)
--
-- Which is to say that @x@ is the floating point image of @(a * b)@ and
-- @y@ stores the residual error term.
--
-- This could be nicer if we had access to a hardware fused multiply-add.
times :: Compensable a => a -> a -> (a -> a -> r) -> r
times :: a -> a -> (a -> a -> r) -> r
times a
a a
b a -> a -> r
k =
  a -> (a -> a -> r) -> r
forall a r. Compensable a => a -> (a -> a -> r) -> r
split a
a ((a -> a -> r) -> r) -> (a -> a -> r) -> r
forall a b. (a -> b) -> a -> b
$ \a
a1 a
a2 ->
  a -> (a -> a -> r) -> r
forall a r. Compensable a => a -> (a -> a -> r) -> r
split a
b ((a -> a -> r) -> r) -> (a -> a -> r) -> r
forall a b. (a -> b) -> a -> b
$ \a
b1 a
b2 ->
  let x :: a
x = a
a a -> a -> a
forall a. Num a => a -> a -> a
* a
b in a -> a -> r
k a
x (a
a2a -> a -> a
forall a. Num a => a -> a -> a
*a
b2 a -> a -> a
forall a. Num a => a -> a -> a
- (((a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
a1a -> a -> a
forall a. Num a => a -> a -> a
*a
b1) a -> a -> a
forall a. Num a => a -> a -> a
- a
a2a -> a -> a
forall a. Num a => a -> a -> a
*a
b1) a -> a -> a
forall a. Num a => a -> a -> a
- a
a1a -> a -> a
forall a. Num a => a -> a -> a
*a
b2))
  -- let x = a * b in k x (((a1*b1-x)+a1*b2+b2*b1)+a2*b2)
{-# INLINEABLE times #-}

-- this is a variant on a division algorithm by Liddicoat and Flynn
divide :: Compensable a => a -> a -> (a -> a -> r) -> r
divide :: a -> a -> (a -> a -> r) -> r
divide a
a a
b = Compensated a -> (a -> a -> r) -> r
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with (Compensated a
aX Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
* Compensated a
ms) where
  x0 :: a
x0   = a -> a
forall a. Fractional a => a -> a
recip a
b
  aX :: Compensated a
aX   = a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
a a
x0 a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated -- calculate aX
  m :: Compensated a
m    = a
1 a -> Compensated a -> Compensated a
forall a. Compensable a => a -> Compensated a -> Compensated a
+^ Compensated a -> Compensated a
forall a. Num a => a -> a
negate (a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
b a
x0 a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated)
  mm :: Compensated a
mm   = Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
*Compensated a
m
  ms :: Compensated a
ms   = Compensated a
1Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+((Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+Compensated a
mm)Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
*Compensated a
mm)
{-# INLINEABLE divide #-}

-- | Priest's renormalization algorithm
--
-- @renorm a b c@ generates a 'Compensated' number assuming @a '>=' b '>=' c@.
renorm :: Compensable a => a -> a -> a -> Compensated a
renorm :: a -> a -> a -> Compensated a
renorm a
a a
b a
c =
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
fadd a
b a
c ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
fadd a
a a
x1 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
fadd a
x2 (a
y1 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y2) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated
{-# INLINE renorm #-}

-- | @'squared' a k@ computes @k x y@ such that
--
-- > x + y = a * a
-- > x = fl(a * a)
--
-- Which is to say that @x@ is the floating point image of @(a * a)@ and
-- @y@ stores the residual error term.
squared :: Compensable a => a -> (a -> a -> r) -> r
squared :: a -> (a -> a -> r) -> r
squared a
a a -> a -> r
k =
  a -> (a -> a -> r) -> r
forall a r. Compensable a => a -> (a -> a -> r) -> r
split a
a ((a -> a -> r) -> r) -> (a -> a -> r) -> r
forall a b. (a -> b) -> a -> b
$ \a
a1 a
a2 ->
  let x :: a
x = a
a a -> a -> a
forall a. Num a => a -> a -> a
* a
a in a -> a -> r
k a
x (a
a2a -> a -> a
forall a. Num a => a -> a -> a
*a
a2 a -> a -> a
forall a. Num a => a -> a -> a
- ((a
x a -> a -> a
forall a. Num a => a -> a -> a
- a
a1a -> a -> a
forall a. Num a => a -> a -> a
*a
a1) a -> a -> a
forall a. Num a => a -> a -> a
- a
2a -> a -> a
forall a. Num a => a -> a -> a
*(a
a2a -> a -> a
forall a. Num a => a -> a -> a
*a
a1)))
{-# INLINE squared #-}

-- | Calculate a fast square of a compensated number.
square :: Compensable a => Compensated a -> Compensated a
square :: Compensated a -> Compensated a
square Compensated a
m =
  Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
  a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> (a -> a -> r) -> r
squared a
a ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
a a
b ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1 (a
x2a -> a -> a
forall a. Num a => a -> a -> a
*a
2) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
  a -> a -> a -> Compensated a
forall a. Compensable a => a -> a -> a -> Compensated a
renorm a
x1 a
x3 (a
ba -> a -> a
forall a. Num a => a -> a -> a
*a
b a -> a -> a
forall a. Num a => a -> a -> a
+ a
2a -> a -> a
forall a. Num a => a -> a -> a
*a
y2 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y3)
{-# INLINE square #-}

-- | error-free split of a floating point number into two parts.
--
-- Note: these parts do not satisfy the `compensated` contract
split :: Compensable a => a -> (a -> a -> r) -> r
split :: a -> (a -> a -> r) -> r
split a
a a -> a -> r
k = a -> a -> r
k a
x a
y where
  c :: a
c = a
forall a. Compensable a => a
magica -> a -> a
forall a. Num a => a -> a -> a
*a
a
  x :: a
x = a
c a -> a -> a
forall a. Num a => a -> a -> a
- (a
c a -> a -> a
forall a. Num a => a -> a -> a
- a
a)
  y :: a
y = a
a a -> a -> a
forall a. Num a => a -> a -> a
- a
x
{-# INLINEABLE split #-}

-- | Calculate a scalar + compensated sum with Kahan summation.
(+^) :: Compensable a => a -> Compensated a -> Compensated a
a
a +^ :: a -> Compensated a -> Compensated a
+^ Compensated a
m = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
b a
c -> let y :: a
y = a
a a -> a -> a
forall a. Num a => a -> a -> a
- a
c; t :: a
t = a
b a -> a -> a
forall a. Num a => a -> a -> a
+ a
y in a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
t ((a
t a -> a -> a
forall a. Num a => a -> a -> a
- a
b) a -> a -> a
forall a. Num a => a -> a -> a
- a
y)
{-# INLINE (+^) #-}

-- | Compute @a * 'Compensated' a@
(*^) :: Compensable a => a -> Compensated a -> Compensated a
a
c *^ :: a -> Compensated a -> Compensated a
*^ Compensated a
m =
  Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \ a
a a
b ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
c a
a ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
c a
b ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
fadd a
x1 a
x2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
  a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1 a
y3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x4 a
y4 ->
  a -> a -> a -> Compensated a
forall a. Compensable a => a -> a -> a -> Compensated a
renorm a
x3 a
x4 (a
y4 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y2)
{-# INLINE (*^) #-}

class (RealFrac a, Floating a) => Compensable a where
  -- | This provides a numeric data type with effectively doubled precision by
  -- using Knuth's error free transform and a number of custom compensated
  -- arithmetic circuits.
  --
  -- This construction can be iterated, doubling precision each time.
  --
  -- >>> round (Prelude.product [2..100] :: Compensated (Compensated (Compensated Double)))
  -- 93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000
  --
  -- >>> Prelude.product [2..100]
  -- 93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000
  data Compensated a

  -- | This extracts both the 'primal' and 'residual' components of a 'Compensated'
  -- number.
  with :: Compensable a => Compensated a -> (a -> a -> r) -> r

  -- | Used internally to construct 'compensated' values that satisfy our residual contract.
  --
  -- When in doubt, use @'add' a b 'compensated'@ instead of @'compensated' a b@
  compensated :: Compensable a => a -> a -> Compensated a

  -- | This 'magic' number is used to 'split' the significand in half, so we can multiply
  -- them separately without losing precision in 'times'.
  magic :: a

instance Compensable Double where
  data Compensated Double = CD {-# UNPACK #-} !Double {-# UNPACK #-} !Double
  with :: Compensated Double -> (Double -> Double -> r) -> r
with (CD a b) Double -> Double -> r
k = Double -> Double -> r
k Double
a Double
b
  {-# INLINE with #-}
  compensated :: Double -> Double -> Compensated Double
compensated = Double -> Double -> Compensated Double
CD
  {-# INLINE compensated #-}
  magic :: Double
magic = Double
134217729
  {-# INLINE magic #-}

instance Compensable Float where
  data Compensated Float = CF {-# UNPACK #-} !Float {-# UNPACK #-} !Float
  with :: Compensated Float -> (Float -> Float -> r) -> r
with (CF a b) Float -> Float -> r
k = Float -> Float -> r
k Float
a Float
b
  {-# INLINE with #-}
  compensated :: Float -> Float -> Compensated Float
compensated = Float -> Float -> Compensated Float
CF
  {-# INLINE compensated #-}
  magic :: Float
magic = Float
4097
  {-# INLINE magic #-}

instance Compensable a => Compensable (Compensated a) where
  data Compensated (Compensated a) = CC !(Compensated a) !(Compensated a)
  with :: Compensated (Compensated a)
-> (Compensated a -> Compensated a -> r) -> r
with (CC a b) Compensated a -> Compensated a -> r
k = Compensated a -> Compensated a -> r
k Compensated a
a Compensated a
b
  {-# INLINE with #-}
  compensated :: Compensated a -> Compensated a -> Compensated (Compensated a)
compensated = Compensated a -> Compensated a -> Compensated (Compensated a)
forall a.
Compensated a -> Compensated a -> Compensated (Compensated a)
CC
  {-# INLINE compensated #-}
  magic :: Compensated a
magic = a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a
forall a. Compensable a => a
magic a -> a -> a
forall a. Num a => a -> a -> a
- a
1) (a
forall a. Compensable a => a
magic a -> a -> a
forall a. Num a => a -> a -> a
- a
1) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \ a
x a
y -> a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
x (a
y a -> a -> a
forall a. Num a => a -> a -> a
+ a
1)
  {-# INLINE magic #-}

#if __GLASGOW_HASKELL__ < 707
instance Typeable1 Compensated where
  typeOf1 _ = mkTyConApp (mkTyCon3 "analytics" "Data.Analytics.Numeric.Compensated" "Compensated") []
#else
deriving instance Typeable Compensated
#endif

instance (Compensable a, Hashable a) => Hashable (Compensated a) where
  hashWithSalt :: Int -> Compensated a -> Int
hashWithSalt Int
n Compensated a
m = Compensated a -> (a -> a -> Int) -> Int
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Int) -> Int) -> (a -> a -> Int) -> Int
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Int -> (a, a) -> Int
forall a. Hashable a => Int -> a -> Int
hashWithSalt Int
n (a
a,a
b)

instance (Compensable a, Data a) => Data (Compensated a) where
  gfoldl :: (forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> Compensated a -> c (Compensated a)
gfoldl forall d b. Data d => c (d -> b) -> d -> c b
f forall g. g -> c g
z Compensated a
m = Compensated a -> (a -> a -> c (Compensated a)) -> c (Compensated a)
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> c (Compensated a)) -> c (Compensated a))
-> (a -> a -> c (Compensated a)) -> c (Compensated a)
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> (a -> a -> Compensated a) -> c (a -> a -> Compensated a)
forall g. g -> c g
z a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated c (a -> a -> Compensated a) -> a -> c (a -> Compensated a)
forall d b. Data d => c (d -> b) -> d -> c b
`f` a
a c (a -> Compensated a) -> a -> c (Compensated a)
forall d b. Data d => c (d -> b) -> d -> c b
`f` a
b
  toConstr :: Compensated a -> Constr
toConstr Compensated a
_ = Constr
compensatedConstr
  gunfold :: (forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c (Compensated a)
gunfold forall b r. Data b => c (b -> r) -> c r
k forall r. r -> c r
z Constr
c = case Constr -> Int
constrIndex Constr
c of
    Int
1 -> c (a -> Compensated a) -> c (Compensated a)
forall b r. Data b => c (b -> r) -> c r
k (c (a -> a -> Compensated a) -> c (a -> Compensated a)
forall b r. Data b => c (b -> r) -> c r
k ((a -> a -> Compensated a) -> c (a -> a -> Compensated a)
forall r. r -> c r
z a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated))
    Int
_ -> [Char] -> c (Compensated a)
forall a. HasCallStack => [Char] -> a
error [Char]
"gunfold"
  dataTypeOf :: Compensated a -> DataType
dataTypeOf Compensated a
_ = DataType
compensatedDataType
  dataCast1 :: (forall d. Data d => c (t d)) -> Maybe (c (Compensated a))
dataCast1 forall d. Data d => c (t d)
f = c (t a) -> Maybe (c (Compensated a))
forall k1 k2 (c :: k1 -> *) (t :: k2 -> k1) (t' :: k2 -> k1)
       (a :: k2).
(Typeable t, Typeable t') =>
c (t a) -> Maybe (c (t' a))
gcast1 c (t a)
forall d. Data d => c (t d)
f

compensatedConstr :: Constr
compensatedConstr :: Constr
compensatedConstr = DataType -> [Char] -> [[Char]] -> Fixity -> Constr
mkConstr DataType
compensatedDataType [Char]
"compensated" [] Fixity
Prefix
{-# NOINLINE compensatedConstr #-}

compensatedDataType :: DataType
compensatedDataType :: DataType
compensatedDataType = [Char] -> [Constr] -> DataType
mkDataType [Char]
"Data.Analytics.Numeric.Compensated" [Constr
compensatedConstr]
{-# NOINLINE compensatedDataType #-}

instance (Compensable a, NFData a) => NFData (Compensated a) where
  rnf :: Compensated a -> ()
rnf Compensated a
m = Compensated a -> (a -> a -> ()) -> ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> ()) -> ()) -> (a -> a -> ()) -> ()
forall a b. (a -> b) -> a -> b
$ \a
x a
y -> a -> ()
forall a. NFData a => a -> ()
rnf a
x () -> () -> ()
`seq` a -> ()
forall a. NFData a => a -> ()
rnf a
y
  {-# INLINE rnf #-}

instance (Compensable a, Show a) => Show (Compensated a) where
  showsPrec :: Int -> Compensated a -> ShowS
showsPrec Int
d Compensated a
m = Compensated a -> (a -> a -> ShowS) -> ShowS
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> ShowS) -> ShowS) -> (a -> a -> ShowS) -> ShowS
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Bool -> ShowS -> ShowS
showParen (Int
d Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
> Int
10) (ShowS -> ShowS) -> ShowS -> ShowS
forall a b. (a -> b) -> a -> b
$
    [Char] -> ShowS
showString [Char]
"compensated " ShowS -> ShowS -> ShowS
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> a -> ShowS
forall a. Show a => Int -> a -> ShowS
T.showsPrec Int
11 a
a ShowS -> ShowS -> ShowS
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Char -> ShowS
showChar Char
' ' ShowS -> ShowS -> ShowS
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> a -> ShowS
forall a. Show a => Int -> a -> ShowS
T.showsPrec Int
11 a
b

instance (Compensable a, Read a) => Read (Compensated a) where
  readPrec :: ReadPrec (Compensated a)
readPrec = ReadPrec (Compensated a) -> ReadPrec (Compensated a)
forall a. ReadPrec a -> ReadPrec a
parens (ReadPrec (Compensated a) -> ReadPrec (Compensated a))
-> ReadPrec (Compensated a) -> ReadPrec (Compensated a)
forall a b. (a -> b) -> a -> b
$ Int -> ReadPrec (Compensated a) -> ReadPrec (Compensated a)
forall a. Int -> ReadPrec a -> ReadPrec a
prec Int
10 (ReadPrec (Compensated a) -> ReadPrec (Compensated a))
-> ReadPrec (Compensated a) -> ReadPrec (Compensated a)
forall a b. (a -> b) -> a -> b
$ do
    Ident [Char]
"compensated" <- ReadPrec Lexeme
lexP
    a
a <- ReadPrec a -> ReadPrec a
forall a. ReadPrec a -> ReadPrec a
step ReadPrec a
forall a. Read a => ReadPrec a
T.readPrec
    a
b <- ReadPrec a -> ReadPrec a
forall a. ReadPrec a -> ReadPrec a
step ReadPrec a
forall a. Read a => ReadPrec a
T.readPrec
    Compensated a -> ReadPrec (Compensated a)
forall (m :: * -> *) a. Monad m => a -> m a
return (Compensated a -> ReadPrec (Compensated a))
-> Compensated a -> ReadPrec (Compensated a)
forall a b. (a -> b) -> a -> b
$ a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
a a
b

type Overcompensated a = Compensated (Compensated a)

-- | This provides the isomorphism between the compact representation we store these in internally
-- and the naive pair of the 'primal' and 'residual' components.
_Compensated :: Compensable a => Iso' (Compensated a) (a, a)
_Compensated :: Iso' (Compensated a) (a, a)
_Compensated = (Compensated a -> (a, a))
-> ((a, a) -> Compensated a) -> Iso' (Compensated a) (a, a)
forall s a b t. (s -> a) -> (b -> t) -> Iso s t a b
iso (Compensated a -> (a -> a -> (a, a)) -> (a, a)
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
`with` (,)) ((a -> a -> Compensated a) -> (a, a) -> Compensated a
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated)
{-# INLINE _Compensated #-}

-- | This 'Lens' lets us edit the 'primal' directly, leaving the 'residual' untouched.
primal :: Compensable a => Lens' (Compensated a) a
primal :: Lens' (Compensated a) a
primal a -> f a
f Compensated a
c = Compensated a -> (a -> a -> f (Compensated a)) -> f (Compensated a)
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
c ((a -> a -> f (Compensated a)) -> f (Compensated a))
-> (a -> a -> f (Compensated a)) -> f (Compensated a)
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> a -> f a
f a
a f a -> (a -> Compensated a) -> f (Compensated a)
forall (f :: * -> *) a b. Functor f => f a -> (a -> b) -> f b
<&> \a
a' -> a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
a' a
b
{-# INLINE primal #-}

-- | This 'Lens' lets us edit the 'residual' directly, leaving the 'primal' untouched.
residual :: Compensable a => Lens' (Compensated a) a
residual :: Lens' (Compensated a) a
residual a -> f a
f Compensated a
c = Compensated a -> (a -> a -> f (Compensated a)) -> f (Compensated a)
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
c ((a -> a -> f (Compensated a)) -> f (Compensated a))
-> (a -> a -> f (Compensated a)) -> f (Compensated a)
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
a (a -> Compensated a) -> f a -> f (Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> a -> f a
f a
b
{-# INLINE residual #-}

-- | Extract the 'primal' component of a 'compensated' value, when and if compensation
-- is no longer required.
uncompensated :: Compensable a => Compensated a -> a
uncompensated :: Compensated a -> a
uncompensated Compensated a
c = Compensated a -> (a -> a -> a) -> a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
c a -> a -> a
forall a b. a -> b -> a
const
{-# INLINE uncompensated #-}

type instance Index (Compensated a) = Int
instance (Compensable a, Compensable b) => Each (Compensated a) (Compensated b) a b where
  each :: (a -> f b) -> Compensated a -> f (Compensated b)
each a -> f b
f Compensated a
m = Compensated a -> (a -> a -> f (Compensated b)) -> f (Compensated b)
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> f (Compensated b)) -> f (Compensated b))
-> (a -> a -> f (Compensated b)) -> f (Compensated b)
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> b -> b -> Compensated b
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (b -> b -> Compensated b) -> f b -> f (b -> Compensated b)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> a -> f b
f a
a f (b -> Compensated b) -> f b -> f (Compensated b)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> a -> f b
f a
b
  {-# INLINE each #-}

instance Compensable a => Eq (Compensated a) where
  Compensated a
m == :: Compensated a -> Compensated a -> Bool
== Compensated a
n = Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
c a
d -> a
a a -> a -> Bool
forall a. Eq a => a -> a -> Bool
== a
c Bool -> Bool -> Bool
&& a
b a -> a -> Bool
forall a. Eq a => a -> a -> Bool
== a
d
  Compensated a
m /= :: Compensated a -> Compensated a -> Bool
/= Compensated a
n = Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
c a
d -> a
a a -> a -> Bool
forall a. Eq a => a -> a -> Bool
/= a
c Bool -> Bool -> Bool
|| a
b a -> a -> Bool
forall a. Eq a => a -> a -> Bool
/= a
d
  {-# INLINE (==) #-}

instance Compensable a => Ord (Compensated a) where
  compare :: Compensated a -> Compensated a -> Ordering
compare Compensated a
m Compensated a
n = Compensated a -> (a -> a -> Ordering) -> Ordering
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Ordering) -> Ordering)
-> (a -> a -> Ordering) -> Ordering
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Compensated a -> (a -> a -> Ordering) -> Ordering
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Ordering) -> Ordering)
-> (a -> a -> Ordering) -> Ordering
forall a b. (a -> b) -> a -> b
$ \a
c a
d -> a -> a -> Ordering
forall a. Ord a => a -> a -> Ordering
compare a
a a
c Ordering -> Ordering -> Ordering
forall a. Semigroup a => a -> a -> a
<> a -> a -> Ordering
forall a. Ord a => a -> a -> Ordering
compare a
b a
d
  {-# INLINE compare #-}
  Compensated a
m <= :: Compensated a -> Compensated a -> Bool
<= Compensated a
n = Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
c a
d -> case a -> a -> Ordering
forall a. Ord a => a -> a -> Ordering
compare a
a a
c of
    Ordering
LT -> Bool
True
    Ordering
EQ -> a
b a -> a -> Bool
forall a. Ord a => a -> a -> Bool
<= a
d
    Ordering
GT -> Bool
False
  {-# INLINE (<=) #-}
  Compensated a
m >= :: Compensated a -> Compensated a -> Bool
>= Compensated a
n = Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
c a
d -> case a -> a -> Ordering
forall a. Ord a => a -> a -> Ordering
compare a
a a
c of
    Ordering
LT -> Bool
False
    Ordering
EQ -> a
b a -> a -> Bool
forall a. Ord a => a -> a -> Bool
>= a
d
    Ordering
GT -> a
a a -> a -> Bool
forall a. Ord a => a -> a -> Bool
>= a
c -- @compare x NaN@ and @compare naN x@ always return 'GT', but @m >= n@ should be 'False'
  {-# INLINE (>=) #-}
  Compensated a
m > :: Compensated a -> Compensated a -> Bool
> Compensated a
n = Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
c a
d -> case a -> a -> Ordering
forall a. Ord a => a -> a -> Ordering
compare a
a a
c of
    Ordering
LT -> Bool
False
    Ordering
EQ -> a
b a -> a -> Bool
forall a. Ord a => a -> a -> Bool
> a
d
    Ordering
GT -> a
a a -> a -> Bool
forall a. Ord a => a -> a -> Bool
> a
c -- @compare x NaN@ and @compare naN x@ always return 'GT', but @m >= n@ should be 'False'
  {-# INLINE (>) #-}
  Compensated a
m < :: Compensated a -> Compensated a -> Bool
< Compensated a
n = Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> Compensated a -> (a -> a -> Bool) -> Bool
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Bool) -> Bool) -> (a -> a -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \a
c a
d -> case a -> a -> Ordering
forall a. Ord a => a -> a -> Ordering
compare a
a a
c of
    Ordering
LT -> Bool
True
    Ordering
EQ -> a
b a -> a -> Bool
forall a. Ord a => a -> a -> Bool
< a
d
    Ordering
GT -> Bool
False
  {-# INLINE (<) #-}

instance Compensable a => Semigroup (Compensated a) where
  <> :: Compensated a -> Compensated a -> Compensated a
(<>) = Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
(+)
  {-# INLINE (<>) #-}

instance Compensable a => Monoid (Compensated a) where
  mempty :: Compensated a
mempty = a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
0 a
0
  {-# INLINE mempty #-}
  mappend :: Compensated a -> Compensated a -> Compensated a
mappend = Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
(+)
  {-# INLINE mappend #-}

-- | Perform Kahan summation over a list.
kahan :: (Foldable f, Compensable a) => f a -> Compensated a
kahan :: f a -> Compensated a
kahan = (a -> Compensated a -> Compensated a)
-> Compensated a -> f a -> Compensated a
forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
Foldable.foldr a -> Compensated a -> Compensated a
forall a. Compensable a => a -> Compensated a -> Compensated a
(+^) Compensated a
forall a. Monoid a => a
mempty
{-# INLINE kahan #-}

instance Compensable a => Num (Compensated a) where
  Compensated a
m + :: Compensated a -> Compensated a -> Compensated a
+ Compensated a
n =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a  a
b  ->
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
c a
d ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add  a
a  a
c ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1  a
d ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add  a
b a
x2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x1 a
x3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x4 a
y4 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x4 (a
y2 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y3 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y4) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated
  {-# INLINE (+) #-}

{-
  m + n =
    with m $ \a b ->
    with n $ \c d ->
    add a c $ \x1 y1 ->
    add b d $ \x2 y2 ->
    renorm x1 x2 (y1 + y2)
  {-# INLINE (+) #-}
-}


  Compensated a
m * :: Compensated a -> Compensated a -> Compensated a
* Compensated a
n =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
n ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
c a
d ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
a a
c ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
b a
c ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times a
a a
d ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x1 a
x2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x4 a
y4 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x3 a
x4 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x5 a
y5 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1 a
y4 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x6 a
y6 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y5 a
x6 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x7 a
y7 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x5 a
x7 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x8 a
y8 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x8 (a
ba -> a -> a
forall a. Num a => a -> a -> a
*a
d a -> a -> a
forall a. Num a => a -> a -> a
+ a
y2 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y3 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y6 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y7 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y8) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated
  {-# INLINE (*) #-}

{-
  m * n =
    with m $ \a b ->
    with n $ \c d ->
    times a c $ \x1 y1 ->
    times b c $ \x2 y2 ->
    times a d $ \x3 y3 ->
    add y1 x2 $ \x4 y4 ->
    add x3 x4 $ \x5 y5 ->
    renorm x1 x5 (b * d + y2 + y4 + y3 + y5)
  {-# INLINE (*) #-}
-}

  negate :: Compensated a -> Compensated a
negate Compensated a
m = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> (a -> a) -> a -> a -> Compensated a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
on a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a -> a
forall a. Num a => a -> a
negate)
  -- {-# INLINE negate #-}

  Compensated a
x - :: Compensated a -> Compensated a -> Compensated a
- Compensated a
y = Compensated a
x Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a -> Compensated a
forall a. Num a => a -> a
negate Compensated a
y
  {-# INLINE (-) #-}

  signum :: Compensated a -> Compensated a
signum Compensated a
m = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
_ -> a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a
forall a. Num a => a -> a
signum a
a) a
0
  {-# INLINE signum #-}

  abs :: Compensated a -> Compensated a
abs Compensated a
m = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    if a
a a -> a -> Bool
forall a. Ord a => a -> a -> Bool
< a
0
    then a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a
forall a. Num a => a -> a
negate a
a) (a -> a
forall a. Num a => a -> a
negate a
b)
    else a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
a a
b
  {-# INLINE abs #-}

  fromInteger :: Integer -> Compensated a
fromInteger Integer
i = a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x (Integer -> a
forall a. Num a => Integer -> a
fromInteger (Integer
i Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
- a -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round a
x)) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated where
    x :: a
x = Integer -> a
forall a. Num a => Integer -> a
fromInteger Integer
i
  {-# INLINE fromInteger #-}

instance Compensable a => Enum (Compensated a) where
  succ :: Compensated a -> Compensated a
succ Compensated a
a = Compensated a
a Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
1
  {-# INLINE succ #-}
  pred :: Compensated a -> Compensated a
pred Compensated a
a = Compensated a
a Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
- Compensated a
1
  {-# INLINE pred #-}
  toEnum :: Int -> Compensated a
toEnum Int
i = a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x (Int -> a
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
- a -> Int
forall a b. (RealFrac a, Integral b) => a -> b
round a
x)) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated where
    x :: a
x = Int -> a
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
i
  {-# INLINE toEnum #-}
  fromEnum :: Compensated a -> Int
fromEnum = Compensated a -> Int
forall a b. (RealFrac a, Integral b) => a -> b
round
  {-# INLINE fromEnum #-}
  enumFrom :: Compensated a -> [Compensated a]
enumFrom Compensated a
a = Compensated a
a Compensated a -> [Compensated a] -> [Compensated a]
forall a. a -> [a] -> [a]
: Compensated a -> [Compensated a]
forall a. Enum a => a -> [a]
Prelude.enumFrom (Compensated a
a Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
1)
  {-# INLINE enumFrom #-}
  enumFromThen :: Compensated a -> Compensated a -> [Compensated a]
enumFromThen Compensated a
a Compensated a
b = Compensated a
a Compensated a -> [Compensated a] -> [Compensated a]
forall a. a -> [a] -> [a]
: Compensated a -> Compensated a -> [Compensated a]
forall a. Enum a => a -> a -> [a]
Prelude.enumFromThen Compensated a
b (Compensated a
b Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
- Compensated a
a Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
b)
  {-# INLINE enumFromThen #-}
  enumFromTo :: Compensated a -> Compensated a -> [Compensated a]
enumFromTo Compensated a
a Compensated a
b
    | Compensated a
a Compensated a -> Compensated a -> Bool
forall a. Ord a => a -> a -> Bool
<= Compensated a
b = Compensated a
a Compensated a -> [Compensated a] -> [Compensated a]
forall a. a -> [a] -> [a]
: Compensated a -> Compensated a -> [Compensated a]
forall a. Enum a => a -> a -> [a]
Prelude.enumFromTo (Compensated a
a Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
1) Compensated a
b
    | Bool
otherwise = []
  {-# INLINE enumFromTo #-}
  enumFromThenTo :: Compensated a -> Compensated a -> Compensated a -> [Compensated a]
enumFromThenTo Compensated a
a Compensated a
b Compensated a
c
    | Compensated a
a Compensated a -> Compensated a -> Bool
forall a. Ord a => a -> a -> Bool
<= Compensated a
b    = Compensated a -> [Compensated a]
up Compensated a
a
    | Bool
otherwise = Compensated a -> [Compensated a]
down Compensated a
a
    where
     delta :: Compensated a
delta = Compensated a
b Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
- Compensated a
a
     up :: Compensated a -> [Compensated a]
up Compensated a
x | Compensated a
x Compensated a -> Compensated a -> Bool
forall a. Ord a => a -> a -> Bool
<= Compensated a
c    = Compensated a
x Compensated a -> [Compensated a] -> [Compensated a]
forall a. a -> [a] -> [a]
: Compensated a -> [Compensated a]
up (Compensated a
x Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
delta)
          | Bool
otherwise = []
     down :: Compensated a -> [Compensated a]
down Compensated a
x | Compensated a
c Compensated a -> Compensated a -> Bool
forall a. Ord a => a -> a -> Bool
<= Compensated a
x    = Compensated a
x Compensated a -> [Compensated a] -> [Compensated a]
forall a. a -> [a] -> [a]
: Compensated a -> [Compensated a]
down (Compensated a
x Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
delta)
            | Bool
otherwise = []
  {-# INLINE enumFromThenTo #-}

instance Compensable a => Fractional (Compensated a) where
  recip :: Compensated a -> Compensated a
recip Compensated a
m = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add (a -> a
forall a. Fractional a => a -> a
recip a
a) (-a
b a -> a -> a
forall a. Fractional a => a -> a -> a
/ (a
a a -> a -> a
forall a. Num a => a -> a -> a
* a
a)) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated
  {-# INLINE recip #-}

  -- | A variant on a hardware division algorithm by Liddicoat and Flynn
  Compensated a
a / :: Compensated a -> Compensated a -> Compensated a
/ Compensated a
b = (Compensated a
aCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
*Compensated a
x0) Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
* (Compensated a
1Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+((Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+Compensated a
mm)Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
*Compensated a
mm)) where
    x0 :: Compensated a
x0  = Compensated a -> Compensated a
forall a. Fractional a => a -> a
recip Compensated a
b
    m :: Compensated a
m   = Compensated a
1 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
- Compensated a
bCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
*Compensated a
x0
    mm :: Compensated a
mm  = Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
*Compensated a
m
  -- {-# INLINE (/) #-}

  fromRational :: Rational -> Compensated a
fromRational Rational
r = Integer -> Compensated a
forall a. Num a => Integer -> a
fromInteger (Rational -> Integer
forall a. Ratio a -> a
numerator Rational
r) Compensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/ Integer -> Compensated a
forall a. Num a => Integer -> a
fromInteger (Rational -> Integer
forall a. Ratio a -> a
denominator Rational
r)
  -- {-# INLINE fromRational #-}

instance Compensable a => Real (Compensated a) where
  toRational :: Compensated a -> Rational
toRational Compensated a
m = Compensated a -> (a -> a -> Rational) -> Rational
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((Rational -> Rational -> Rational)
-> (a -> Rational) -> a -> a -> Rational
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
on Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
(+) a -> Rational
forall a. Real a => a -> Rational
toRational)
  -- {-# INLINE toRational #-}

instance Compensable a => RealFrac (Compensated a) where
  properFraction :: Compensated a -> (b, Compensated a)
properFraction Compensated a
m = Compensated a
-> (a -> a -> (b, Compensated a)) -> (b, Compensated a)
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> (b, Compensated a)) -> (b, Compensated a))
-> (a -> a -> (b, Compensated a)) -> (b, Compensated a)
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> case a -> (b, a)
forall a b. (RealFrac a, Integral b) => a -> (b, a)
properFraction a
a of
    (b
w, a
p) -> a -> a -> (a -> a -> (b, Compensated a)) -> (b, Compensated a)
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
p a
b ((a -> a -> (b, Compensated a)) -> (b, Compensated a))
-> (a -> a -> (b, Compensated a)) -> (b, Compensated a)
forall a b. (a -> b) -> a -> b
$ \ a
x a
y -> case a -> (b, a)
forall a b. (RealFrac a, Integral b) => a -> (b, a)
properFraction a
x of
      (b
w',a
q) -> (b
w b -> b -> b
forall a. Num a => a -> a -> a
+ b
w', a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
q a
y a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated)
  -- {-# INLINE properFraction #-}

instance (Compensable a, Binary a) => Binary (Compensated a) where
  get :: Get (Compensated a)
get = a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> Compensated a) -> Get a -> Get (a -> Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Get a
forall t. Binary t => Get t
Binary.get Get (a -> Compensated a) -> Get a -> Get (Compensated a)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Get a
forall t. Binary t => Get t
Binary.get
  put :: Compensated a -> Put
put Compensated a
m = Compensated a -> (a -> a -> Put) -> Put
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Put) -> Put) -> (a -> a -> Put) -> Put
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> do
    a -> Put
forall t. Binary t => t -> Put
Binary.put a
a
    a -> Put
forall t. Binary t => t -> Put
Binary.put a
b

instance (Compensable a, Serialize a) => Serialize (Compensated a) where
  get :: Get (Compensated a)
get = a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> Compensated a) -> Get a -> Get (a -> Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Get a
forall t. Serialize t => Get t
Serialize.get Get (a -> Compensated a) -> Get a -> Get (Compensated a)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Get a
forall t. Serialize t => Get t
Serialize.get
  put :: Putter (Compensated a)
put Compensated a
m = Compensated a -> (a -> a -> PutM ()) -> PutM ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> PutM ()) -> PutM ()) -> (a -> a -> PutM ()) -> PutM ()
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> do
    a -> PutM ()
forall t. Serialize t => Putter t
Serialize.put a
a
    a -> PutM ()
forall t. Serialize t => Putter t
Serialize.put a
b

instance (Compensable a, Serial a) => Serial (Compensated a) where
  deserialize :: m (Compensated a)
deserialize = a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> Compensated a) -> m a -> m (a -> Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> m a
forall a (m :: * -> *). (Serial a, MonadGet m) => m a
Bytes.deserialize m (a -> Compensated a) -> m a -> m (Compensated a)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> m a
forall a (m :: * -> *). (Serial a, MonadGet m) => m a
Bytes.deserialize
  serialize :: Compensated a -> m ()
serialize Compensated a
m = Compensated a -> (a -> a -> m ()) -> m ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> m ()) -> m ()) -> (a -> a -> m ()) -> m ()
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> do
    a -> m ()
forall a (m :: * -> *). (Serial a, MonadPut m) => a -> m ()
Bytes.serialize a
a
    a -> m ()
forall a (m :: * -> *). (Serial a, MonadPut m) => a -> m ()
Bytes.serialize a
b

-- ಠ_ಠ this unnecessarily expects that the format won't change, because I can't derive a better instance.
instance (Compensable a, Serialize a, Typeable a) => SafeCopy (Compensated a) where
  -- safecopy-0.10.0 changed its default implementations for these methods.
  -- The implementations below are copied from the pre-0.10.0 defaults.
  errorTypeName :: Proxy (Compensated a) -> [Char]
errorTypeName Proxy (Compensated a)
_ = [Char]
"<unknown type>"
  getCopy :: Contained (Get (Compensated a))
getCopy = Get (Compensated a) -> Contained (Get (Compensated a))
forall a. a -> Contained a
contain Get (Compensated a)
forall t. Serialize t => Get t
Serialize.get
  putCopy :: Compensated a -> Contained (PutM ())
putCopy = PutM () -> Contained (PutM ())
forall a. a -> Contained a
contain (PutM () -> Contained (PutM ()))
-> (Compensated a -> PutM ())
-> Compensated a
-> Contained (PutM ())
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Compensated a -> PutM ()
forall t. Serialize t => Putter t
Serialize.put

instance (Compensable a, Storable a) => Storable (Compensated a) where
  sizeOf :: Compensated a -> Int
sizeOf Compensated a
_ = a -> Int
forall a. Storable a => a -> Int
sizeOf (a
forall a. HasCallStack => a
undefined :: a) Int -> Int -> Int
forall a. Num a => a -> a -> a
* Int
2
  -- {-# INLINE sizeOf #-}
  alignment :: Compensated a -> Int
alignment Compensated a
_ = a -> Int
forall a. Storable a => a -> Int
alignment (a
forall a. HasCallStack => a
undefined :: a)
  -- {-# INLINE alignment #-}
  peekElemOff :: Ptr (Compensated a) -> Int -> IO (Compensated a)
peekElemOff Ptr (Compensated a)
p Int
o | Ptr a
q <- Ptr (Compensated a) -> Ptr a
forall a b. Ptr a -> Ptr b
castPtr Ptr (Compensated a)
p, Int
o2 <- Int
o Int -> Int -> Int
forall a. Num a => a -> a -> a
* Int
2 =
    a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> Compensated a) -> IO a -> IO (a -> Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Ptr a -> Int -> IO a
forall a. Storable a => Ptr a -> Int -> IO a
peekElemOff Ptr a
q Int
o2 IO (a -> Compensated a) -> IO a -> IO (Compensated a)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Ptr a -> Int -> IO a
forall a. Storable a => Ptr a -> Int -> IO a
peekElemOff Ptr a
q (Int
o2Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
  -- {-# INLINE peekElemOff #-}
  pokeElemOff :: Ptr (Compensated a) -> Int -> Compensated a -> IO ()
pokeElemOff Ptr (Compensated a)
p Int
o Compensated a
m | Ptr a
q <- Ptr (Compensated a) -> Ptr a
forall a b. Ptr a -> Ptr b
castPtr Ptr (Compensated a)
p, Int
o2 <- Int
o Int -> Int -> Int
forall a. Num a => a -> a -> a
* Int
2 = Compensated a -> (a -> a -> IO ()) -> IO ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> IO ()) -> IO ()) -> (a -> a -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> do
    Ptr a -> Int -> a -> IO ()
forall a. Storable a => Ptr a -> Int -> a -> IO ()
pokeElemOff Ptr a
q Int
o2 a
a
    Ptr a -> Int -> a -> IO ()
forall a. Storable a => Ptr a -> Int -> a -> IO ()
pokeElemOff Ptr a
q (Int
o2Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1) a
b
  -- {-# INLINE pokeElemOff #-}
  peekByteOff :: Ptr b -> Int -> IO (Compensated a)
peekByteOff Ptr b
p Int
o | Ptr Any
q <- Ptr b -> Ptr Any
forall a b. Ptr a -> Ptr b
castPtr Ptr b
p =
    a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> Compensated a) -> IO a -> IO (a -> Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Ptr Any -> Int -> IO a
forall a b. Storable a => Ptr b -> Int -> IO a
peekByteOff Ptr Any
q Int
o IO (a -> Compensated a) -> IO a -> IO (Compensated a)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Ptr Any -> Int -> IO a
forall a b. Storable a => Ptr b -> Int -> IO a
peekByteOff Ptr Any
q (Int
o Int -> Int -> Int
forall a. Num a => a -> a -> a
+ a -> Int
forall a. Storable a => a -> Int
sizeOf (a
forall a. HasCallStack => a
undefined :: a))
  -- {-# INLINE peekByteOff #-}
  pokeByteOff :: Ptr b -> Int -> Compensated a -> IO ()
pokeByteOff Ptr b
p Int
o Compensated a
m | Ptr Any
q <- Ptr b -> Ptr Any
forall a b. Ptr a -> Ptr b
castPtr Ptr b
p = Compensated a -> (a -> a -> IO ()) -> IO ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> IO ()) -> IO ()) -> (a -> a -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> do
    Ptr Any -> Int -> a -> IO ()
forall a b. Storable a => Ptr b -> Int -> a -> IO ()
pokeByteOff Ptr Any
q Int
o a
a
    Ptr Any -> Int -> a -> IO ()
forall a b. Storable a => Ptr b -> Int -> a -> IO ()
pokeByteOff Ptr Any
q (Int
oInt -> Int -> Int
forall a. Num a => a -> a -> a
+a -> Int
forall a. Storable a => a -> Int
sizeOf (a
forall a. HasCallStack => a
undefined :: a)) a
b
  -- {-# INLINE pokeByteOff #-}
  peek :: Ptr (Compensated a) -> IO (Compensated a)
peek Ptr (Compensated a)
p | Ptr a
q <- Ptr (Compensated a) -> Ptr a
forall a b. Ptr a -> Ptr b
castPtr Ptr (Compensated a)
p = a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> Compensated a) -> IO a -> IO (a -> Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Ptr a -> IO a
forall a. Storable a => Ptr a -> IO a
peek Ptr a
q IO (a -> Compensated a) -> IO a -> IO (Compensated a)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Ptr a -> Int -> IO a
forall a. Storable a => Ptr a -> Int -> IO a
peekElemOff Ptr a
q Int
1
  -- {-# INLINE peek #-}
  poke :: Ptr (Compensated a) -> Compensated a -> IO ()
poke Ptr (Compensated a)
p Compensated a
m | Ptr a
q <- Ptr (Compensated a) -> Ptr a
forall a b. Ptr a -> Ptr b
castPtr Ptr (Compensated a)
p = Compensated a -> (a -> a -> IO ()) -> IO ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> IO ()) -> IO ()) -> (a -> a -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \a
a a
b -> do
    Ptr a -> a -> IO ()
forall a. Storable a => Ptr a -> a -> IO ()
poke Ptr a
q a
a
    Ptr a -> Int -> a -> IO ()
forall a. Storable a => Ptr a -> Int -> a -> IO ()
pokeElemOff Ptr a
q Int
1 a
b
  -- {-# INLINE poke #-}

newtype instance U.MVector s (Compensated a) = MV_Compensated (U.MVector s (a,a))
newtype instance U.Vector (Compensated a) = V_Compensated (U.Vector (a, a))

instance (Compensable a, Unbox a) => M.MVector U.MVector (Compensated a) where
  basicLength :: MVector s (Compensated a) -> Int
basicLength (MV_Compensated v) = MVector s (a, a) -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
M.basicLength MVector s (a, a)
v
  {-# INLINE basicLength #-}
  basicUnsafeSlice :: Int
-> Int -> MVector s (Compensated a) -> MVector s (Compensated a)
basicUnsafeSlice Int
i Int
n (MV_Compensated v) = MVector s (a, a) -> MVector s (Compensated a)
forall s a. MVector s (a, a) -> MVector s (Compensated a)
MV_Compensated (MVector s (a, a) -> MVector s (Compensated a))
-> MVector s (a, a) -> MVector s (Compensated a)
forall a b. (a -> b) -> a -> b
$ Int -> Int -> MVector s (a, a) -> MVector s (a, a)
forall (v :: * -> * -> *) a s.
MVector v a =>
Int -> Int -> v s a -> v s a
M.basicUnsafeSlice Int
i Int
n MVector s (a, a)
v
  {-# INLINE basicUnsafeSlice #-}
  basicOverlaps :: MVector s (Compensated a) -> MVector s (Compensated a) -> Bool
basicOverlaps (MV_Compensated v1) (MV_Compensated v2) = MVector s (a, a) -> MVector s (a, a) -> Bool
forall (v :: * -> * -> *) a s.
MVector v a =>
v s a -> v s a -> Bool
M.basicOverlaps MVector s (a, a)
v1 MVector s (a, a)
v2
  {-# INLINE basicOverlaps #-}
  basicUnsafeNew :: Int -> m (MVector (PrimState m) (Compensated a))
basicUnsafeNew Int
n = MVector (PrimState m) (a, a)
-> MVector (PrimState m) (Compensated a)
forall s a. MVector s (a, a) -> MVector s (Compensated a)
MV_Compensated (MVector (PrimState m) (a, a)
 -> MVector (PrimState m) (Compensated a))
-> m (MVector (PrimState m) (a, a))
-> m (MVector (PrimState m) (Compensated a))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Int -> m (MVector (PrimState m) (a, a))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> m (v (PrimState m) a)
M.basicUnsafeNew Int
n
  {-# INLINE basicUnsafeNew #-}
  basicUnsafeReplicate :: Int -> Compensated a -> m (MVector (PrimState m) (Compensated a))
basicUnsafeReplicate Int
n Compensated a
m = Compensated a
-> (a -> a -> m (MVector (PrimState m) (Compensated a)))
-> m (MVector (PrimState m) (Compensated a))
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> m (MVector (PrimState m) (Compensated a)))
 -> m (MVector (PrimState m) (Compensated a)))
-> (a -> a -> m (MVector (PrimState m) (Compensated a)))
-> m (MVector (PrimState m) (Compensated a))
forall a b. (a -> b) -> a -> b
$ \a
x a
y -> MVector (PrimState m) (a, a)
-> MVector (PrimState m) (Compensated a)
forall s a. MVector s (a, a) -> MVector s (Compensated a)
MV_Compensated (MVector (PrimState m) (a, a)
 -> MVector (PrimState m) (Compensated a))
-> m (MVector (PrimState m) (a, a))
-> m (MVector (PrimState m) (Compensated a))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Int -> (a, a) -> m (MVector (PrimState m) (a, a))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
Int -> a -> m (v (PrimState m) a)
M.basicUnsafeReplicate Int
n (a
x,a
y)
  {-# INLINE basicUnsafeReplicate #-}
  basicUnsafeRead :: MVector (PrimState m) (Compensated a) -> Int -> m (Compensated a)
basicUnsafeRead (MV_Compensated v) Int
i = (a -> a -> Compensated a) -> (a, a) -> Compensated a
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated ((a, a) -> Compensated a) -> m (a, a) -> m (Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> MVector (PrimState m) (a, a) -> Int -> m (a, a)
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m a
M.basicUnsafeRead MVector (PrimState m) (a, a)
v Int
i
  {-# INLINE basicUnsafeRead #-}
  basicUnsafeWrite :: MVector (PrimState m) (Compensated a)
-> Int -> Compensated a -> m ()
basicUnsafeWrite (MV_Compensated v) Int
i Compensated a
m = Compensated a -> (a -> a -> m ()) -> m ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> m ()) -> m ()) -> (a -> a -> m ()) -> m ()
forall a b. (a -> b) -> a -> b
$ \ a
x a
y -> MVector (PrimState m) (a, a) -> Int -> (a, a) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> a -> m ()
M.basicUnsafeWrite MVector (PrimState m) (a, a)
v Int
i (a
x,a
y)
  {-# INLINE basicUnsafeWrite #-}
  basicClear :: MVector (PrimState m) (Compensated a) -> m ()
basicClear (MV_Compensated v) = MVector (PrimState m) (a, a) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
M.basicClear MVector (PrimState m) (a, a)
v
  {-# INLINE basicClear #-}
  basicSet :: MVector (PrimState m) (Compensated a) -> Compensated a -> m ()
basicSet (MV_Compensated v) Compensated a
m = Compensated a -> (a -> a -> m ()) -> m ()
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> m ()) -> m ()) -> (a -> a -> m ()) -> m ()
forall a b. (a -> b) -> a -> b
$ \ a
x a
y -> MVector (PrimState m) (a, a) -> (a, a) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> a -> m ()
M.basicSet MVector (PrimState m) (a, a)
v (a
x,a
y)
  {-# INLINE basicSet #-}
  basicUnsafeCopy :: MVector (PrimState m) (Compensated a)
-> MVector (PrimState m) (Compensated a) -> m ()
basicUnsafeCopy (MV_Compensated v1) (MV_Compensated v2) = MVector (PrimState m) (a, a)
-> MVector (PrimState m) (a, a) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
M.basicUnsafeCopy MVector (PrimState m) (a, a)
v1 MVector (PrimState m) (a, a)
v2
  {-# INLINE basicUnsafeCopy #-}
  basicUnsafeMove :: MVector (PrimState m) (Compensated a)
-> MVector (PrimState m) (Compensated a) -> m ()
basicUnsafeMove (MV_Compensated v1) (MV_Compensated v2) = MVector (PrimState m) (a, a)
-> MVector (PrimState m) (a, a) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> v (PrimState m) a -> m ()
M.basicUnsafeMove MVector (PrimState m) (a, a)
v1 MVector (PrimState m) (a, a)
v2
  {-# INLINE basicUnsafeMove #-}
  basicUnsafeGrow :: MVector (PrimState m) (Compensated a)
-> Int -> m (MVector (PrimState m) (Compensated a))
basicUnsafeGrow (MV_Compensated v) Int
n = MVector (PrimState m) (a, a)
-> MVector (PrimState m) (Compensated a)
forall s a. MVector s (a, a) -> MVector s (Compensated a)
MV_Compensated (MVector (PrimState m) (a, a)
 -> MVector (PrimState m) (Compensated a))
-> m (MVector (PrimState m) (a, a))
-> m (MVector (PrimState m) (Compensated a))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> MVector (PrimState m) (a, a)
-> Int -> m (MVector (PrimState m) (a, a))
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> Int -> m (v (PrimState m) a)
M.basicUnsafeGrow MVector (PrimState m) (a, a)
v Int
n
  {-# INLINE basicUnsafeGrow #-}
#if MIN_VERSION_vector(0,11,0)
  basicInitialize :: MVector (PrimState m) (Compensated a) -> m ()
basicInitialize (MV_Compensated v) = MVector (PrimState m) (a, a) -> m ()
forall (v :: * -> * -> *) a (m :: * -> *).
(MVector v a, PrimMonad m) =>
v (PrimState m) a -> m ()
M.basicInitialize MVector (PrimState m) (a, a)
v
  {-# INLINE basicInitialize #-}
#endif

instance (Compensable a, Unbox a) => G.Vector U.Vector (Compensated a) where
  basicUnsafeFreeze :: Mutable Vector (PrimState m) (Compensated a)
-> m (Vector (Compensated a))
basicUnsafeFreeze (MV_Compensated v) = Vector (a, a) -> Vector (Compensated a)
forall a. Vector (a, a) -> Vector (Compensated a)
V_Compensated (Vector (a, a) -> Vector (Compensated a))
-> m (Vector (a, a)) -> m (Vector (Compensated a))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Mutable Vector (PrimState m) (a, a) -> m (Vector (a, a))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> m (v a)
G.basicUnsafeFreeze MVector (PrimState m) (a, a)
Mutable Vector (PrimState m) (a, a)
v
  {-# INLINE basicUnsafeFreeze #-}
  basicUnsafeThaw :: Vector (Compensated a)
-> m (Mutable Vector (PrimState m) (Compensated a))
basicUnsafeThaw (V_Compensated v) = MVector (PrimState m) (a, a)
-> MVector (PrimState m) (Compensated a)
forall s a. MVector s (a, a) -> MVector s (Compensated a)
MV_Compensated (MVector (PrimState m) (a, a)
 -> MVector (PrimState m) (Compensated a))
-> m (MVector (PrimState m) (a, a))
-> m (MVector (PrimState m) (Compensated a))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Vector (a, a) -> m (Mutable Vector (PrimState m) (a, a))
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
v a -> m (Mutable v (PrimState m) a)
G.basicUnsafeThaw Vector (a, a)
v
  {-# INLINE basicUnsafeThaw #-}
  basicLength :: Vector (Compensated a) -> Int
basicLength (V_Compensated v) = Vector (a, a) -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.basicLength Vector (a, a)
v
  {-# INLINE basicLength #-}
  basicUnsafeSlice :: Int -> Int -> Vector (Compensated a) -> Vector (Compensated a)
basicUnsafeSlice Int
i Int
n (V_Compensated v) = Vector (a, a) -> Vector (Compensated a)
forall a. Vector (a, a) -> Vector (Compensated a)
V_Compensated (Vector (a, a) -> Vector (Compensated a))
-> Vector (a, a) -> Vector (Compensated a)
forall a b. (a -> b) -> a -> b
$ Int -> Int -> Vector (a, a) -> Vector (a, a)
forall (v :: * -> *) a. Vector v a => Int -> Int -> v a -> v a
G.basicUnsafeSlice Int
i Int
n Vector (a, a)
v
  {-# INLINE basicUnsafeSlice #-}
  basicUnsafeIndexM :: Vector (Compensated a) -> Int -> m (Compensated a)
basicUnsafeIndexM (V_Compensated v) Int
i
                = (a -> a -> Compensated a) -> (a, a) -> Compensated a
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated ((a, a) -> Compensated a) -> m (a, a) -> m (Compensated a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Vector (a, a) -> Int -> m (a, a)
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, Monad m) =>
v a -> Int -> m a
G.basicUnsafeIndexM Vector (a, a)
v Int
i
  {-# INLINE basicUnsafeIndexM #-}
  basicUnsafeCopy :: Mutable Vector (PrimState m) (Compensated a)
-> Vector (Compensated a) -> m ()
basicUnsafeCopy (MV_Compensated mv) (V_Compensated v)
                = Mutable Vector (PrimState m) (a, a) -> Vector (a, a) -> m ()
forall (v :: * -> *) a (m :: * -> *).
(Vector v a, PrimMonad m) =>
Mutable v (PrimState m) a -> v a -> m ()
G.basicUnsafeCopy MVector (PrimState m) (a, a)
Mutable Vector (PrimState m) (a, a)
mv Vector (a, a)
v
  {-# INLINE basicUnsafeCopy #-}
  elemseq :: Vector (Compensated a) -> Compensated a -> b -> b
elemseq Vector (Compensated a)
_ Compensated a
m b
z = Compensated a -> (a -> a -> b) -> b
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> b) -> b) -> (a -> a -> b) -> b
forall a b. (a -> b) -> a -> b
$ \a
x a
y -> Vector a -> a -> b -> b
forall (v :: * -> *) a b. Vector v a => v a -> a -> b -> b
G.elemseq (Vector a
forall a. HasCallStack => a
undefined :: U.Vector a) a
x
                                 (b -> b) -> b -> b
forall a b. (a -> b) -> a -> b
$ Vector a -> a -> b -> b
forall (v :: * -> *) a b. Vector v a => v a -> a -> b -> b
G.elemseq (Vector a
forall a. HasCallStack => a
undefined :: U.Vector a) a
y b
z
  {-# INLINE elemseq #-}

-- | /NB:/ Experimental and partially implemented.
--
-- Other than sqrt, the accuracy of these is basically uncalculated! In fact many of these are known to be wrong! Patches and improvements are welcome.
instance Compensable a => Floating (Compensated a) where
#ifdef SPECIALIZE_INSTANCES
  {-# SPECIALIZE instance Floating (Compensated Double) #-}
  {-# SPECIALIZE instance Floating (Compensated Float) #-}
  {-# SPECIALIZE instance Compensable a => Floating (Compensated (Compensated a)) #-}
#endif

  exp :: Compensated a -> Compensated a
exp Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
exp a
a) (a -> a
forall a. Floating a => a -> a
exp a
b) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated

  sin :: Compensated a -> Compensated a
sin Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
sin a
a) (a -> a
forall a. Floating a => a -> a
cos a
b) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
sin a
b) (a -> a
forall a. Floating a => a -> a
cos a
a) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x1 a
x2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1 a
y2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x4 a
y4 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x4 a
y3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x5 a
y5 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x5 a
x3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x6 a
y6 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add (a
y4 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y5 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y6) a
x6 a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated

  cos :: Compensated a -> Compensated a
cos Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
cos a
a) (a -> a
forall a. Floating a => a -> a
cos a
b) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (-a -> a
forall a. Floating a => a -> a
sin a
b) (a -> a
forall a. Floating a => a -> a
sin a
a) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x1 a
x2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1 a
y2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x4 a
y4 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x4 a
y3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x5 a
y5 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x5 a
x3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x6 a
y6 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add (a
y4 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y5 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y6) a
x6 a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated

  tan :: Compensated a -> Compensated a
tan Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add (a -> a
forall a. Floating a => a -> a
tan a
a) (a -> a
forall a. Floating a => a -> a
tan a
b) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated Compensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/
    (a
1 a -> Compensated a -> Compensated a
forall a. Compensable a => a -> Compensated a -> Compensated a
+^ a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
tan a
a) (a -> a
forall a. Floating a => a -> a
tan a
b) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated)

  sinh :: Compensated a -> Compensated a
sinh Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
sinh a
a) (a -> a
forall a. Floating a => a -> a
cosh a
b) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
cosh a
a) (a -> a
forall a. Floating a => a -> a
sinh a
b) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x1 a
x2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1 a
y2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x4 a
y4 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x4 a
y3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x5 a
y5 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x5 a
x3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x6 a
y6 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add (a
y4 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y5 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y6) a
x6 a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated

  cosh :: Compensated a -> Compensated a
cosh Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
cosh a
a) (a -> a
forall a. Floating a => a -> a
cosh a
b) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x1 a
y1 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
sinh a
b) (a -> a
forall a. Floating a => a -> a
sinh a
a) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x2 a
y2 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x1 a
x2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x3 a
y3 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
y1 a
y2 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x4 a
y4 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x4 a
y3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x5 a
y5 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add a
x5 a
x3 ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
x6 a
y6 ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add (a
y4 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y5 a -> a -> a
forall a. Num a => a -> a -> a
+ a
y6) a
x6 a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated

  tanh :: Compensated a -> Compensated a
tanh Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \a
a a
b ->
    a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
fadd (a -> a
forall a. Floating a => a -> a
tanh a
a) (a -> a
forall a. Floating a => a -> a
tanh a
b) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated Compensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/
    (a
1 a -> Compensated a -> Compensated a
forall a. Compensable a => a -> Compensated a -> Compensated a
+^ a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Compensable a => a -> a -> (a -> a -> r) -> r
times (a -> a
forall a. Floating a => a -> a
tanh a
a) (a -> a
forall a. Floating a => a -> a
tanh a
b) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated)

  -- This requires an accurate 'exp', which we currently lack.
  log :: Compensated a -> Compensated a
log Compensated a
m =
    Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with Compensated a
m ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ \ a
a a
b -> let
      xy1 :: Compensated a
xy1 = a -> a -> (a -> a -> Compensated a) -> Compensated a
forall a r. Num a => a -> a -> (a -> a -> r) -> r
add (a -> a
forall a. Floating a => a -> a
log a
a) (a
ba -> a -> a
forall a. Fractional a => a -> a -> a
/a
a) a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated
      xy2 :: Compensated a
xy2 = Compensated a
xy1 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
m Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
* Compensated a -> Compensated a
forall a. Floating a => a -> a
exp (-Compensated a
xy1) Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
- Compensated a
1 -- Newton Raphson step 1
    in Compensated a
xy2 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
m Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
* Compensated a -> Compensated a
forall a. Floating a => a -> a
exp (-Compensated a
xy2) Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
- Compensated a
1      -- Newton Raphson step 2

  -- | Hardware sqrt improved by the Babylonian algorithm (Newton Raphson)
  sqrt :: Compensated a -> Compensated a
sqrt Compensated a
m = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with (Compensated a
z4 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/Compensated a
z4) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ (a -> a -> Compensated a) -> (a -> a) -> a -> a -> Compensated a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
on a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> a
forall a. Fractional a => a -> a -> a
/a
2) where
    z0 :: a
z0 = a -> a
forall a. Floating a => a -> a
sqrt (Compensated a
mCompensated a -> Getting a (Compensated a) a -> a
forall s a. s -> Getting a s a -> a
^.Getting a (Compensated a) a
forall a. Compensable a => Lens' (Compensated a) a
primal)
    z1 :: Compensated a
z1 = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with (a
z0 a -> Compensated a -> Compensated a
forall a. Compensable a => a -> Compensated a -> Compensated a
+^ (Compensated a
m Compensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/ a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated a
z0 a
0)) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ (a -> a -> Compensated a) -> (a -> a) -> a -> a -> Compensated a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
on a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> a
forall a. Fractional a => a -> a -> a
/a
2)
    z2 :: Compensated a
z2 = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with (Compensated a
z1 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/Compensated a
z1) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ (a -> a -> Compensated a) -> (a -> a) -> a -> a -> Compensated a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
on a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> a
forall a. Fractional a => a -> a -> a
/a
2)
    z3 :: Compensated a
z3 = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with (Compensated a
z2 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/Compensated a
z2) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ (a -> a -> Compensated a) -> (a -> a) -> a -> a -> Compensated a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
on a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> a
forall a. Fractional a => a -> a -> a
/a
2)
    z4 :: Compensated a
z4 = Compensated a -> (a -> a -> Compensated a) -> Compensated a
forall a r.
(Compensable a, Compensable a) =>
Compensated a -> (a -> a -> r) -> r
with (Compensated a
z3 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
mCompensated a -> Compensated a -> Compensated a
forall a. Fractional a => a -> a -> a
/Compensated a
z3) ((a -> a -> Compensated a) -> Compensated a)
-> (a -> a -> Compensated a) -> Compensated a
forall a b. (a -> b) -> a -> b
$ (a -> a -> Compensated a) -> (a -> a) -> a -> a -> Compensated a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
on a -> a -> Compensated a
forall a. (Compensable a, Compensable a) => a -> a -> Compensated a
compensated (a -> a -> a
forall a. Fractional a => a -> a -> a
/a
2)

  -- TODO: do log1p, expm1, log1mexp, log1pexp right!
  log1p :: Compensated a -> Compensated a
log1p Compensated a
a = Compensated a -> Compensated a
forall a. Floating a => a -> a
log (Compensated a
1 Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a
a)
  {-# INLINE log1p #-}
  expm1 :: Compensated a -> Compensated a
expm1 Compensated a
a = Compensated a -> Compensated a
forall a. Floating a => a -> a
exp Compensated a
a Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
- Compensated a
1
  {-# INLINE expm1 #-}
  log1mexp :: Compensated a -> Compensated a
log1mexp Compensated a
a | Compensated a
a Compensated a -> Compensated a -> Bool
forall a. Ord a => a -> a -> Bool
<= Compensated a -> Compensated a
forall a. Floating a => a -> a
log Compensated a
2 = Compensated a -> Compensated a
forall a. Floating a => a -> a
log (Compensated a -> Compensated a
forall a. Num a => a -> a
negate (Compensated a -> Compensated a
forall a. Floating a => a -> a
expm1 (Compensated a -> Compensated a
forall a. Num a => a -> a
negate Compensated a
a)))
             | Bool
otherwise  = Compensated a -> Compensated a
forall a. Floating a => a -> a
log1p (Compensated a -> Compensated a
forall a. Num a => a -> a
negate (Compensated a -> Compensated a
forall a. Floating a => a -> a
exp (Compensated a -> Compensated a
forall a. Num a => a -> a
negate Compensated a
a)))
  {-# INLINE log1mexp #-}
  log1pexp :: Compensated a -> Compensated a
log1pexp Compensated a
a
    | Compensated a
a Compensated a -> Compensated a -> Bool
forall a. Ord a => a -> a -> Bool
<= Compensated a
18   = Compensated a -> Compensated a
forall a. Floating a => a -> a
log1p (Compensated a -> Compensated a
forall a. Floating a => a -> a
exp Compensated a
a)
    | Compensated a
a Compensated a -> Compensated a -> Bool
forall a. Ord a => a -> a -> Bool
<= Compensated a
100  = Compensated a
a Compensated a -> Compensated a -> Compensated a
forall a. Num a => a -> a -> a
+ Compensated a -> Compensated a
forall a. Floating a => a -> a
exp (Compensated a -> Compensated a
forall a. Num a => a -> a
negate Compensated a
a)
    | Bool
otherwise = Compensated a
a
  {-# INLINE log1pexp #-}

  -- (**)    = error "TODO"
  pi :: Compensated a
pi      = [Char] -> Compensated a
forall a. HasCallStack => [Char] -> a
error [Char]
"TODO"
  asin :: Compensated a -> Compensated a
asin    = [Char] -> Compensated a -> Compensated a
forall a. HasCallStack => [Char] -> a
error [Char]
"TODO"
  atan :: Compensated a -> Compensated a
atan    = [Char] -> Compensated a -> Compensated a
forall a. HasCallStack => [Char] -> a
error [Char]
"TODO"
  acos :: Compensated a -> Compensated a
acos    = [Char] -> Compensated a -> Compensated a
forall a. HasCallStack => [Char] -> a
error [Char]
"TODO"
  asinh :: Compensated a -> Compensated a
asinh   = [Char] -> Compensated a -> Compensated a
forall a. HasCallStack => [Char] -> a
error [Char]
"TODO"
  atanh :: Compensated a -> Compensated a
atanh   = [Char] -> Compensated a -> Compensated a
forall a. HasCallStack => [Char] -> a
error [Char]
"TODO"
  acosh :: Compensated a -> Compensated a
acosh   = [Char] -> Compensated a -> Compensated a
forall a. HasCallStack => [Char] -> a
error [Char]
"TODO"