[Libreoffice-bugs] [Bug 138360] calculation: basic math fail, !calc rounding wrong! dividing by inverse sometimes different from multiply, 7E10/1E-5 wrong,

Thu Nov 26 14:13:23 UTC 2020

https://bugs.documentfoundation.org/show_bug.cgi?id=138360

--- Comment #2 from b. <newbie-02 at gmx.de> ---
hello @erAck, 

generally speaking: i agree that 'fixing floating point math' is a tough
problem, but there has been progress, and i hope more is possible, 

a (workaround) solution/improvement (at least for the rounding problems) could
be to store predefined values for the three 'painpoints' (10^-5, 10^-9, 10^-15)
(respectively their reciprocal values 1/1E-5, 1/1E-9, 1/1E-15), and use them
instead of 'on the fly' wrong calculated values to avoid the evil influence
they actually have (if / as long as there are no better ways?!?),  

   --- but who wants 'better problems'   ;-) 

as far as i was able to check with calc as it calculates just now there are
only these three values for which an error becomes visible with
'=RAWSUBTRACT(1/[X];ROUND(1/[X];0))', absolute it is greatest with -0,125 for
1E-15, relative ('=RAWSUBTRACT(1/[X];ROUND(1/[X];0))/ROUND(1/[X];0)') it is
greatest for 1E-5 with (-1,45519152283669E-16), 

for all other arguments from 1E-1 to 1E-308 the results is '0', 

if you look at other enumerators than '1' (i did for the powers 10^0 to 10^100)
1E9/1E-5 has the max relative deviation of -1,5625E-16, absolute -0,015625, and
an 'uneven' value near it is 955037695 with a relative deviation of
-1,63606107714942E-16, absolute -0,015625, while an ULP in that range is
0,010603048383739, thus a claim "result of a fp calculation shouldn't be more
than 1 ULP off from 'standard math'" (0,5 ULP from each argument) doesn't hold
... (this claim isn't valid for all operations, i know), 

i haven't checked it extensively yet, but i think as long as theese three
painpoint values are allowed to inject deviations larger than one ULP into the
calculations (which are visible with rawsubtract), wrong results will appear
from time to time and bring up complaints, 

imho this proposal is a 'crutch', but not completely stupid, in other areas
such things are used as best practice, 

and now forgive my cumbersome thinking and persistent arguing, progress would
be slower if not again and again someone would try, and i just have to work my
way into the problems first ... 

and forgive my long writing and the repetitions below, i'm just limited in time
to make it short and concise :-( 

> You can't apply school math on IEEE-754 binary floating-point values over all ranges. Period.

accepted, but imho all efforts to get as much as possible are 'good', and there
have been improvements over the years, thus progress isn't impossible, 

> 0.00001 has no exact binary representation, hence

yes, if i'm right its the one value of 1E-1 ... 1E-308 neccessary for rounding
to decimals which is most 'off' from what it should be, 

> 70000000000 / 0.00001 => 6999999999999999

yes, see 'better' below, 

> They share IEEE-754 binary floating-point.

but there are variations how you can use it? options? or might even be flaws in
it? 

> There's no extra rounding involved in multiplication and division. Also, the formula expression is calculated as entered from left to right, there's no "use the better".

my idea was that the programmer / developer has at hand how a rounding is
calculated behind the scenes, and if he - not the user in the sheet who isn't
aware of such problems - would steer '=round(7000000000000000;-5)' to
calculate: 
'=7000000000000000/100000 (+/- [0..0,9] to steer the rounding) 'take int'
*100000', 
instead of: 
'=7000000000000000*0,00001 (+/- [0..0,9] to steer the rounding) 'take int'
/0,00001', 
calc would gain better results, 

couldn't test yet if it holds for all cases as tie to even and '4-bit cut'
stepped in, 

> There's no Calc error in this case. And this case of inaccuracy can't be circumvented unless the final value is rounded but to what criteria?

i'm thinking about such problems and testing in this direction, at the moment
stuck by - see above -, 

- at least for some cases fixing the arguments could be better than rounding
the results, e.g. 
> 70000000000 / 0.00001 => 6999999999999999
'=70000000000 * (1/0,00001)' => 6999999999999999
'=70000000000 * ROUND((1/0,00001);0) => 7000000000000000

- for many operations one can calculate arguments for a meaningful rounding
from the operands, unfortunately not for divisions?, 

if you calculate a rounded value by '=7000000000000000*0,00001 ~ ~ *100000'
then you can legitimately round to -15 decimal places (-15 from 7000.., +5 from
0,00001, -5 from *100000), 

- if you calculate 70000000000 * 100000 instead of 70000000000 / 0.00001 no
error will occur, but if you can legitimately round to '-15' decimal places, 

> Which version yielded the result 7000000000000000 for this very example? I tried 5.3.7 and it delivers the same 6999999999999999 result.

Version: 4.1.6.2
Build ID: 40ff705089295be5be0aae9b15123f687c05b0a, 

was astonished myself as it was 'better than Weitz' ... did a recheck, even on
neighboured values, no error for '/0,00001' in the range '70000000000 - 35' to
'70000000000 + 36', not 'error outside that range' but 'outside that range not
yet tested', 

> Also doing that math in a raw C++ compiled program yields 6999999999999999.
Also Gnumeric agrees.

i suspect a mistake or 'something special' at my side, but what? 

and of course you can't 'embarrass' LO or calc if / as long as basic
functionalities like IEEE or C-compiler have bugs or 'limitations', but we can
/ should insist that these basic things are improved (it can, should and
mustn't be that 50 years after a computer controlled flight to the moon
personal computers do worse calculations than simple pocket calculators for
1,50 EUR), and calc or LO can - and should - try to work around the weak points
as far as they are known, 

if i spotted a weakness in IEEE, c-compilers, gnumeric or even intel FPU it
schould be corrected there, what i dislike is 'fp-math has to have fails' and
no further care ... if there are unavoidable errors we need to learn to live
with it, 

>> didn't do extensive testing, but on a first glance even ex$el looks better, 
> Maybe they apply some magic. Ask them to reveal their tricks.

my connection to Micro$oft is a bit weak and burdened with a lot of trouble,
maybe better someone else could ask ... 

is that magic? instead of '=80000000000 / 0,00001' calculate '=80000000000
* round(1/0,00001;x)' with values for 'x' from 0 to -4, 

!!! don't try 'round(1/0,00001;-5)' ... or do!, it is funny, did it in discrete
steps in separate cells, had a minute to think why things obviously looking
right!!! failed, then remembered 'rounding to -5 fails in calc', 'even
100000,00000000000000000000 shown as result for '=round(1/1E-5;-5)' may deviate
from 100000', it doe's!, i'm still smiling about my silly attempt ... ;-)   !!!

> Yes, IEEE-754 binary floating-point is not perfect. Get over it.

i proposed IEEE-754 2008, @Mike Kaganski: 'NO!', i'm trying user space macros,
not so bad but slow, still have fun to work on it ... i'm in hope that i can
produce ideas which - if coded - will enhance calc, 

> I'll keep this bug for the rounding thing, we maybe could do better to not round integer values beyond the decimal. 

that would avoid some problems, but solving the /1E-5 issue would be better as
it also heals other calculations where it's not so easy to spot deviations, 

> Apart from that Gnumeric also delivers the same 5000000000000002 result.

school math and user expectations won't change because Gnumeric and / or ex$el
fail, 
'=((5000000000000000+1)*1000000000+0,5)/1000000000' -> 500000000000000**2**, 
'=((5000000000000000+1)/0,000000001+0,5)*0,000000001' -> 500000000000000**1**, 
not a solution for all cases, but a proof that theese cases can be calculated
correct, 

>> 8. try the same steps in calc ver. 4.1.6.2, it is not! unavoidable, 
> Oh, that version does it different? Maybe it did other things wrong. 

yes, of course, maybe someone can try to combine the strong points of both
versions?  

> It delivered better results for both

5000000000000002
6999999999999999

'm not fully satisfied with 4.1's '5000000000000000' for '=ROUND(5e15+1;9)',
but again, wrong results are not carved in stone, 

>> 9. call it a regression, 
> If so, maybe a bibisect could shed some light. Requesting. But it could also be that compilers' IEEE-754 handling or libraries changed somewhat. And that even both on Windows and Linux.

i remember that @Mike Kaganski changed something recently, tdf#130725, acc. my
limited understanding it shouldn't play into this issue, but maybe this or
similar changes imported weaker behaviour for some cases ... in combination
with former patches to circumvent former weaknesses? such is 'devs stage', 

one actual idea noted here for 'not to get lost': 0,5 ULP precision at
granularity borders ... produces overlapping 'definition ranges'? as e.g. 0,5
ULP of 2^53 reaches 'down' to including '2^53 - 1' while 2^53-1 itself covers
the range from 2^53-1-0,5 to 2^53-1+0,4999~, that leads to some ambiguity as
e.g. 2^53-1+0,4 is falling into both 'ULP-precision-ranges', is that handeled
in a meaningful way?

-- 
You are receiving this mail because:
You are the assignee for the bug.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.freedesktop.org/archives/libreoffice-bugs/attachments/20201126/295d2a2d/attachment-0001.htm>