Re: checking for tiny inequalities (was: Re: puzzling Octave behaviour)

[Przemek Klosowski]

   Which I understand as "any number smaller than eps is undistinguishable
   from 0". However, AFAIK, 64-bit double precision floating point permits
   representing something as small as 10^-323 or so, and although I've done
   comparisons to some multiple of eps as a threshold to check for 0s,
   I'm unsure on how would that work if I really needed to work with values
   around 10^-20.

'eps' is the smaller number undistinguishable from zero IN
CALCULATIONS.  Consider the expression (1+a)-1, which should result in
the value of 'a', but look here: (1+eps)-1 gives 2.2204e-16 (which is
equal to 'eps'), but (1+eps/2)-1 results in a 0.

True, you can represent numbers as small as 1e-323, which would be
0.0000.....0001 with 322 zeros after decimal point --- but the
floating point representation DOES NOT specify all 323 significant
digits in this number.  The 'missing' digits are implicitly specified
by the exponent: the 15 or so significant digits are shifted left or
right by the exponent, and zeros are implicitly introduced in between
the decimal point and the position of the significant digits.

In fact, log10(1/eps) pretty much tells you how many significant
digits you have in your floating point representation.

When you try to add or compare numbers whose relative magnitudes
(exponents) differ more than the number of significant digits, the
floating point representation must equalize their exponents first,
resulting in completely shifting out the significant digits of one of
the numbers, with the resulting loss of precision. 

Some further trickiness: the numerical value of eps/2 is
1.11022302462516e-16, but if we actually use it in calculation, we get
this:

(1+eps)-1                    ans =  2.22044604925031e-16   (OK, 'eps')
(1+1.11022302462516e-16)-1   ans =  2.22044604925031e-16   ( huh?)
(1+1.11022302462515e-16)-1   ans = 0                       (OK again)

So, we not only have issues of rounding in calculations, but also
rounding in I/O: the numerical value of eps/2 (1.11022302462516e-16)
is read in and rounded in such a way that the calculation returns
eps. Any tweaks to that value even in the least digit (last line),
result in complete loss of precision.