Re: a-b+b != a

[James Sherman Jr.]

On Mon, Sep 4, 2017 at 10:53 AM, stn021 <address@hidden> wrote:
> Hi,
>
> this is weird:
>
> r1 = round( 1e3*randn(5) ) / 1e3 ;
> r2 = round( 1e3*randn(5) ) / 1e3 ;
> r3 = r1 - r2 + r2 ;
> not_zero = r1 - r3
>
> not_zero =
>
>    0.0000e+00   0.0000e+00   5.5511e-17   0.0000e+00   1.1102e-16
>    5.5511e-17   0.0000e+00   0.0000e+00   0.0000e+00   0.0000e+00
>    0.0000e+00   0.0000e+00   0.0000e+00   0.0000e+00  -2.0817e-16
>   -5.5511e-17   0.0000e+00   0.0000e+00  -1.1102e-16   0.0000e+00
>    2.7756e-17   0.0000e+00   0.0000e+00   0.0000e+00   0.0000e+00
>
> So a-b+b != a
>
> The difference to zero is small, around 1e-17. But iterations can
> cause this error to increase.
>
> I use leasqr() and in each iteration the last line is
>   retval = retval - someval + someval
>
> With that additional line I get quite different results compared to
> the same program without this line even though they should be
> identical.
>
> Is there a way to avoid this phenomenon ?
>
> THX
> Stefan
>
> _______________________________________________
> Help-octave mailing list
> address@hidden
> https://lists.gnu.org/mailman/listinfo/help-octave

This is fundamental issue caused by the use of floating point numbers.
I'd suggest reading the article to get a better background in floating
point, but the key section is
https://en.wikipedia.org/wiki/Floating-point_arithmetic#Accuracy_problems

One way to avoid this issue to never test equality (or not equality)
when dealing with floating point numbers.  So instead of
>if ( a == b )
>  do stuff
>end
instead use something like:
>if ( abs(a-b) < eps(a) )
> do stuff
>end

Hope this helps,
James Sherman