From 811bc680cd96c8b367dc14a772ff00951694c8ef Mon Sep 17 00:00:00 2001 From: Jeffrey Sarnoff Date: Mon, 13 Jul 2020 13:40:18 -0400 Subject: [PATCH] robustify `four_sum(a,b,c,d)` --- src/math/errorfree.jl | 37 ++++++++++++++++++------------------- 1 file changed, 18 insertions(+), 19 deletions(-) diff --git a/src/math/errorfree.jl b/src/math/errorfree.jl index 7130e99d..e583df74 100644 --- a/src/math/errorfree.jl +++ b/src/math/errorfree.jl @@ -28,43 +28,42 @@ end Computes `hi = fl(a+b+c)` and `md = err(a+b+c), lo = err(md)`. """ -function three_sum(a::T, b::T, c::T) where {T<:FloatWithFMA} - s, t = two_sum(b, c) - hi, u = two_sum(a, s) - md, lo = two_sum(u, t) +function three_sum(a::T,b::T,c::T) where {T<:FloatWithFMA} + t0, t1 = two_sum(a, b) + hi, t2 = two_sum(t0, c) + md, lo = two_sum(t2, t1) hi, md = two_hilo_sum(hi, md) return hi, md, lo end - """ four_sum(a, b, c, d) Computes `hi = fl(a+b+c+d)` and `hm = err(a+b+c+d), ml = err(hm), lo = err(ml)`. """ function four_sum(a::T,b::T,c::T,d::T) where {T<:FloatWithFMA} - t0, t1 = two_sum(a , b) - t0, t2 = two_sum(t0, c) - hi, t3 = two_sum(t0, d) - t0, t1 = two_sum(t1, t2) - hm, t2 = two_sum(t0, t3) # here, t0 >= t3 - ml, lo = two_sum(t1, t2) + t0, t1 = two_sum(a, b) + t2, t3 = two_sum(c, d) + hi, t4 = two_sum(t0, t2) + t5, lo = two_sum(t1, t3) + hm, ml = two_sum(t4, t5) + ml, lo = two_hilo_sum(ml, lo) + hm, ml = two_hilo_sum(hm, ml) + hi, hm = two_hilo_sum(hi,hm) return hi, hm, ml, lo end - - """ two_sum(a, b, c) Computes `hi = fl(a+b+c)` and `lo = err(a+b+c)`. """ -function two_sum(a::T, b::T, c::T) where {T<:FloatWithFMA} - s, t = two_sum(b, c) - x, u = two_sum(a, s) - y = u + t - x, y = two_sum(x, y) - return x, y +function two_sum(a::T,b::T,c::T) where {T<:FloatWithFMA} + t0, t1 = two_sum(a, b) + hi, t2 = two_sum(t0, c) + lo = t2 + t1 + hi, lo = two_hilo_sum(hi, lo) + return hi, lo end """