[Rd] bug? (PR#8074)

[Thomas Lumley]

On Tue, 16 Aug 2005, Paul Mosquin wrote:
> I guess that I expect R to act pretty much as C or C++ would do if I were to 
> program the same code.  It's a bit of a surprise that assignment of 
> rationals, well within precision, followed by multiplication leading to a 
> result well within precision picks up those extra bits along the way. 
> Something to watch out for, to be sure.

But those rationals are *not* well within precision. 0.2 is a infinite 
repeating binary fraction (in base 16 it is 0.333....) so it is not stored 
precisely. 0.04 is also not stored precisely, and it so happens that the 
error in representing 0.04 is not the same as the error in representing 
0.02*0.02.

Of course this will still happen in C: R is written in C.
For example, on my computer the following C program
-----------------------
#include <stdio.h>

int main(){
   double d=0.2;
   double dd;

   dd=d*d;
   if(dd==d)
      printf("Equal\n");
   else
      printf("Difference=%20.18f\n",0.04-dd);
}
----------------------
prints
[al:~] thomas% ./a.out
Difference=-0.000000000000000007

which happens to agree with the result R gives, though this isn't 
guaranteed.  You simply cannot rely on floating point equality unless you 
know how the last bit rounding errors are handled.

 	-thomas