# a-b+b != a

13 messages
Open this post in threaded view
|

## a-b+b != a

 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 [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

 On Mon, Sep 4, 2017 at 10:53 AM, stn021 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 [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave try thisr3 = r1+ (- r2 + r2) ;use brackets.-- DAS _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

 In reply to this post by stn021 On Mon, 4 Sep 2017 16:53:29 +0200 stn021 <[hidden email]> wrote: > Is there a way to avoid this phenomenon ? Learn about errors in floating point arithmetic? This might be a start:  https://en.wikipedia.org/wiki/Kahan_summation_algorithmBut in general you should not expect (a-b)+b to equal a. Gord _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

 In reply to this post by stn021 On Mon, Sep 4, 2017 at 10:53 AM, stn021 <[hidden email]> 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 > [hidden email] > https://lists.gnu.org/mailman/listinfo/help-octaveThis 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_problemsOne 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 _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

 Well, then it is easy to solve. Simply round by some power of 2: mul = 2^20 ; r1 = round( mul*randn(5) ) / mul ; r2 = round( mul*randn(5) ) / mul ; r3 = r1 - r2 + r2 ; not_zero = r1 - r3 not_zero =    0   0   0   0   0    0   0   0   0   0    0   0   0   0   0    0   0   0   0   0    0   0   0   0   0 2017-09-04 17:58 GMT+02:00 James Sherman Jr. <[hidden email]>: > On Mon, Sep 4, 2017 at 10:53 AM, stn021 <[hidden email]> 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 >> [hidden email] >> 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 _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

 In reply to this post by Doug Stewart-4 Hi Doug,the idea here is to add weights to the result.res = y + weights .* ( res-y )If weights==1 then res should remain unchanged. Only it doesn't.So (-r2+r2) does not really work here.I did not include the weights in my original post because it did not seem relevant to the issue, sorry about that.Thx, Stefan2017-09-04 17:46 GMT+02:00 Doug Stewart <[hidden email]>:> ...> try this>> r3 = r1+ (- r2 + r2) ;>> use brackets.2017-09-04 17:46 GMT+02:00 Doug Stewart :On Mon, Sep 4, 2017 at 10:53 AM, stn021 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 [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave try thisr3 = r1+ (- r2 + r2) ;use brackets.-- DAS _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

 In reply to this post by stn021 No way! It's correctin pure mathematics b - a + a = bin computing it is not true that b -a + a = bfloating points and doubles and rounding mean that it is definitely NOT the case that mathematical identities holdOn Mon, Sep 4, 2017 at 3:53 PM, stn021 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 [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

 In reply to this post by stn021 On Mon, Sep 04, 2017 at 18:50:29 +0200, stn021 wrote: > the idea here is to add weights to the result. > > res = y + weights .* ( res-y ) > > If weights==1 then res should remain unchanged. Only it doesn't. If this property is important, then you should probably just add a special case     if (weights != 1)       res = y + weights .* (res - y);     endif -- mike _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave signature.asc (849 bytes) Download Attachment
Open this post in threaded view
|

## Re: a-b+b != a

 2017-09-04 19:00 GMT+02:00 Mike Miller <[hidden email]>: > On Mon, Sep 04, 2017 at 18:50:29 +0200, stn021 wrote: >> the idea here is to add weights to the result. >> >> res = y + weights .* ( res-y ) >> >> If weights==1 then res should remain unchanged. Only it doesn't. > > If this property is important, then you should probably just add a > special case > >     if (weights != 1) >       res = y + weights .* (res - y); >     endif > > -- > mike Hi Mike, not sure about this. The important property is that the weights are applied correctly. So weights==1/2 should lead to half the value just as weights==1 should lead to an unchanged value. Just to be sure: weights is a vector, not a scalar. So some elements can be 1 and others can be something else at the same time. These are supposed to be statistical weights assigning more or less importance to observed data. Stefan _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave
Open this post in threaded view
|

## Re: a-b+b != a

Open this post in threaded view
|

## Re: a-b+b != a

 In reply to this post by Tim Pierce TL DR; addition, subtraction, mult and div do not hold with mathematical identitiesAddition, subtraction, mult and div all run the risk of 'losing accuracy' So every one of those operations can 'break' an identity contrasted with pure mathsPut another way, Try this to illustrate to yourself        > 1e-20 + 1e20or         > 2e-20 + 2e20 - 2e20or if you really want to 'blow your mind'        > -0        > 0This is worth knowing, so if you need to perform          ( x /10000 )* 8 for example it will 9 times out of 10 be more 'accurate' to compute        ( x * 8 ) / 10000because div by a large amount is more likely to 'truncate' detail than mult ...On Mon, Sep 4, 2017 at 5:56 PM, Tim Pierce wrote:No way! It's correctin pure mathematics b - a + a = bin computing it is not true that b -a + a = bfloating points and doubles and rounding mean that it is definitely NOT the case that mathematical identities holdOn Mon, Sep 4, 2017 at 3:53 PM, stn021 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 [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave _______________________________________________ Help-octave mailing list [hidden email] https://lists.gnu.org/mailman/listinfo/help-octave