[ANN] more error-free transformations

[Jeffrey Sarnoff]

ErrorFreeArith.jl offers error-free transformations not (yet?) included in the ErrorFreeTransforms package by dsiem.

These operations convey the usual arithmetic result accompanied by a residual value that is usually lost to rounding.

This gives the correct value at twice the working precision (correctly rounded for +,-,*,/; still 1/(1/x) = x or x ± ulp(x)).

[Tom Breloff]

Thanks Jeffrey. Can you expand on the specifics of the package?  What would you say are the primary use cases? How does this differ from interval arithmetic or Unums?

[Jeffrey Sarnoff]

These are distinct operators that substitute directly for (+),(-),(*),(/) in situations where one wants to obtain more of mathematically true result than is usually available:

two = 2.0; sqrt2 = sqrt(2);

residualValueRoundedAway = Float64(sqrt(big(2)) - sqrt2) # -9.667293313452913e-17

 

mostSignficantPart, leastSignificantPart = eftSqrt(two)

mostSignificantPart ==  1.4142135623730951

leastSignificantPart == -9.667293313452912e-17 # we recover the residual value, itself at Float64 precision 

so we obtain the arithmetic result at twice the 'working' precision (in two parts, the mspart == the usual result).

exp1log2 = exp(1.0)*log(2.0);                                                                       #  1.88416938536372

residualValueRoundedAway = Float64(exp(big(1))*log(big(2)) - exp1log2) #  8.146538547111741e-17

mostSignficantPart, leastSignificantPart = eftProd2( exp(1.0), log(2.0) )   # (1.88416938536372, -8.177744937186283e-17)

------

These transformations have the additional benefit that the two parts are well separated, they do not overlap in the working precision.

So, in all cases, mostSignificantPart + leastSignificantPart == mostSignificantPart.  

They are as well separated as possible, without losing information.

These functions are well-suited to assisting the implementation of extended precision Floating Point math.

Another application (that, until otherwise informed, I'll say is from me) is to accelerate inline rounding:

  (RoundFast.jl, there to see how).

Assuming one had a Float64 unum-ish capability, a double-double float would extend the precision.

(Ultimately, all these parts should meld)

[Tom Breloff]

Actually this is really smart. You can represent exact values (lsp == 0) or open intervals one ulp wide... This is a good chunk of the value of fixed sized unums. I'll be keeping a close eye on the package. Thanks Jeffrey. 

[Jeffrey Sarnoff]

Yep; and I appreciate the note.