Re: fuzzy rounding, roundoff, roundup(), et cetera.

[John Denker <jsd av8n com>]

On 7/8/21 10:52 AM, b. via gnumeric-list wrote:

> i wanted to ask for some concrete help, not to start again long discussions 
> about the inaccuracy of fp-math, IEEE, gnumeric, Excel and Calc.

Sometimes it's better to get the answer that is correct rather
than the answer you wanted or expected.

============

There are at least three things going on here:

1) Gnumeric uses fuzzy rounding and fuzzy comparisons. This is
 deeply embedded in the code. This is (mostly) not a bug in the
 code; the code is (mostly) doing what it was designed to do.

 It's a bad design, but that's not something we need to discuss
 right now. For the next level of detail on this, see here:
   https://www.av8n.com/physics/spreadsheet-tips.htm#sec-fuzzy

 The bottom line is that decent algorithms will not be affected
 by the depredations of fuzzy rounding.

2) On top of that, there are bugs in the code. Fuzzy rounding is
 not implemented entirely correctly. I submitted a bug report with
 a patch here:
   https://gitlab.gnome.org/GNOME/gnumeric/-/issues/589

 The fact that these bugs have persisted for years is strong
 evidence that nobody depends on fuzzy rounding. Nobody really
 benefits from it.

3a) Fuzzy rounding definitely causes roundup(x) to give wrong
 answers. An example of this is x='0.24999999999999997

3b) Even if you use correct (unfuzzy) rounding, it seems likely
 that roundup()-- as currently implemented -- will be vulnerable
 to rare, small errors due to floating-point roundoff.

 I don't have an acid test of this, because it gets jumbled up
 with the fuzzy business, but I have preliminary indications
 that x='0.1250000001000003 is an example of this.

 It's super-easy to defend against such roundoff errors.
 A spreadsheet that does most of the work is here:
   https://www.av8n.com/hack/roundup.gnumeric

===================
Here are some preliminary thoughts on what roundup() ought to be doing:

        **Specification for roundup(x, places)**
        
0       For simplicity, assume x>=0 and places>=0
        The generalization is straightforward.
        
1       The input (x) is a floating point number.
        This is a nontrivial statement, for the following reason:
        you can't claim that '0.750000000000000111 is bigger than '0.7500000000000001
        and therefore it should allegedly round up to 0.7500000000000002
        That's because in floating point, '0.750000000000000111 and '0.7500000000000001
        are the same number. It will round up to itself.
        
2       The output will be of the form nume/denom,
        where both nume and denom are whole numbers.
        I trust the FPU to do the division as accurately as possible,
        so all we have to do is choose the proper nume and denom.
        
3       The denominator is a power of 10. It will be whole number.
        Floating point can represent 10^16 or even 10^22 exactly.
        
4       The numerator is the smallest whole number such that
        nume/denom >= x.
        
5       The numerator can be /estimated/ by taking the ceil() of
        x*10^places. This will be very close, but due to inevitable
        floating point imperfections, it might be off by plus or minus 1.
        
6       The notions of fuzzy rounding and fuzzy comparison must not be used here.
        For details, see here:
        https://www.av8n.com/physics/spreadsheet-tips.htm#sec-fuzzy