Accuracy of floating point calculations

[H. S. Teoh]

On Tue, Oct 29, 2019 at 07:10:08PM +0000, Twilight via Digitalmars-d-learn wrote:
> On Tuesday, 29 October 2019 at 16:11:45 UTC, Daniel Kozak wrote:
> > On Tue, Oct 29, 2019 at 5:09 PM Daniel Kozak <kozzi11 at gmail.com> wrote:
> > > If you use gdc or ldc you will get same results as c++, or you can
> > > use C log directly:
> > > 
> > > import std.stdio;
> > > import std.math : pow;
> > > import core.stdc.math;
> > > 
> > > void main()
> > > {
> > >      writefln("%12.3F",log(1-0.9999)/log(1-(1-0.6)^^20));
> > > }
> > 
> > AFAIK dmd use real  for floating point operations instead of double
> 
> Thanks for the clarification. It appears then that because of dmd's
> real calculations, it produces more accurate results, but maybe
> slower.

Yes, it will be somewhat more accurate, depending on the exact
calculation you're performing. But it depends on the x87 coprocessor,
which hasn't been improved for many years now, and not much attention
has been paid to it, so it would appear that 64-bit double arithmetic
using SIMD or MMX instructions would probably run faster.  (I'd profile
it just to be sure, though. Sometimes performance predictions can be
very wrong.)

So roughly speaking, if you want accuracy, use real, if you want speed,
use float or double.


> (Calculating the result with the high precision calculator at
> https://keisan.casio.com/calculator agrees with dmd.)

To verify accuracy, it's usually safer to use an arbitrary-precision
calculator instead of assuming that the most common answer is the right
one (it may be far off, depending on what exactly is being computed and
how, e.g., due to catastrophic cancellation and things like that). Like
`bc -l` if you're running *nix, e.g. the input:

	scale=32
	l(1 - 0.9999) / l(1 - (1 - 0.6)^20)

gives:

	837675572.37859373067028812966306043501772

So it appears that the C++ answer is less accurate. Note that not all of
the digits above are trustworthy; running the same calculation with
scale=100 gives:

	837675572.3785937306702880546932327627909527172023597021486261165664\
	994508853029795054669322261827298817174322

which shows a divergence in digits after the 15th decimal, meaning that
the subsequent digits are probably garbage values. This is probably
because the logarithm function near 0 is poorly-conditioned, so you
could potentially be getting complete garbage from your floating-point
operations if you're not careful.

Floating-point is a bear. Every programmer should learn to tame it lest
they get mauled:

	https://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html

:-D


T

-- 
May you live all the days of your life. -- Jonathan Swift