compensated-0.3: Compensated floating-point arithmetic

Portabilitynon-portable
Stabilityexperimental
MaintainerEdward Kmett <ekmett@gmail.com>
Safe HaskellNone

Numeric.Compensated

Contents

Description

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:

Synopsis

Documentation

class (RealFrac a, Precise a, Floating a) => Compensable a whereSource

Associated Types

data Compensated a Source

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

Methods

with :: Compensable a => Compensated a -> (a -> a -> r) -> rSource

This extracts both the primal and residual components of a Compensated number.

compensated :: Compensable a => a -> a -> Compensated aSource

Used internally to construct compensated values that satisfy our residual contract.

When in doubt, use add a b compensated instead of compensated a b

magic :: aSource

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

Instances

_Compensated :: Compensable a => Iso' (Compensated a) (a, a)Source

This provides the isomorphism between the compact representation we store these in internally and the naive pair of the primal and residual components.

type Overcompensated a = Compensated (Compensated a)Source

primal :: Compensable a => Lens' (Compensated a) aSource

This Lens lets us edit the primal directly, leaving the residual untouched.

residual :: Compensable a => Lens' (Compensated a) aSource

This Lens lets us edit the residual directly, leaving the primal untouched.

uncompensated :: Compensable a => Compensated a -> aSource

Extract the primal component of a compensated value, when and if compensation is no longer required.

fadd :: Num a => a -> a -> (a -> a -> r) -> rSource

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.

lifting scalars

add :: Num a => a -> a -> (a -> a -> r) -> rSource

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.

times :: Compensable a => a -> a -> (a -> a -> r) -> rSource

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.

squared :: Compensable a => a -> (a -> a -> r) -> rSource

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.

divide :: Compensable a => a -> a -> (a -> a -> r) -> rSource

split :: Compensable a => a -> (a -> a -> r) -> rSource

error-free split of a floating point number into two parts.

Note: these parts do not satisfy the compensated contract

kahan :: (Foldable f, Compensable a) => f a -> Compensated aSource

Perform Kahan summation over a list.

(+^) :: Compensable a => a -> Compensated a -> Compensated aSource

Calculate a scalar + compensated sum with Kahan summation.

(*^) :: Compensable a => a -> Compensated a -> Compensated aSource

Compute a * Compensated a

compensated operators

square :: Compensable a => Compensated a -> Compensated aSource

Calculate a fast square of a compensated number.