patternMinor
Full-precision summation in Haskell
Viewed 0 times
precisionfullsummationhaskell
Problem
I'm trying to speed up the following code for summing doubles with full precision, based on a Python recipe from Raymond Hettinger. Any thoughts on how to make this go faster?
This is already a second draft; you can see a previous version of the code in this gist.
-- Full precision summation based on http://code.activestate.com/recipes/393090/
msum :: [Double] -> Double
msum [] = 0
msum (x:xs) = sum $ foldl' (\ps x' -> reverse (partials x' ps)) [x] xs
partials :: Double -> [Double] -> [Double]
partials x = foldl' (\acc p -> case hilo x p of (hi, 0) -> hi:acc
(hi, lo) -> hi:lo:acc) []
hilo :: Double -> Double -> (Double, Double)
hilo x y | abs x < abs y = hilo y x
| otherwise = (hi, lo) where hi = x + y
lo = y - (hi - x)This is already a second draft; you can see a previous version of the code in this gist.
Solution
You could write
BTW: How could
Rounding?
hilo as hilo :: Double -> Double -> [Double], which returns [hi,lo] for lo /= 0 and [hi] for lo = 0. Then you can write partials aspartials x = foldl' (\acc p -> (hilo x p) ++ acc) -> hi:acc []BTW: How could
lo ever become non-zero consideringlo = y - (hi - x) = y - (x + y - x) = y - y = 0Rounding?
Code Snippets
partials x = foldl' (\acc p -> (hilo x p) ++ acc) -> hi:acc []lo = y - (hi - x) = y - (x + y - x) = y - y = 0Context
StackExchange Code Review Q#1834, answer score: 2
Revisions (0)
No revisions yet.