Re: FLOATING-POINT CONSISTENCY, -FFLOAT-STORE, AND X86

[Toon Moene <toon at moene dot indiv dot nluug dot nl>]

Geert Bosch wrote:
> 
> On Sat, 19 Dec 1998 21:37:49 +0100, Toon Moene wrote:
> 
>   In the mean time, it would be useful for the compiler to warn about
>   testing floating point variables for (in)equality.
> 
> Testing for equality is perfectly fine on systems with IEEE arithmetic
> and many algorithms would be impossible to write efficiently if one would
> regard floating-point as a fuzzy kind of real value. Your statement would
> be true in the pre-IEEE era, but fortunately fpt arithmetic is well-defined
> on the large majority of current systems.

Oh, but I wasn't suggesting that the floating point *instructions* are
not to be trusted,

> If you want to know why, I advise you to read "What Every Computer Scientist
> Should Know About Floating-Point Arithmetic", by David Goldberg, in ACM
> Computing Surveys, vol. 23 nr. 1, march 1991, available in PostScript at:
> http://swift.lanl.gov/Internal/Computing/SunOS_Compilers/common-tools/numerical_comp_guide/goldberg1.ps

[ Thanks for the reference - I read through most of it
  now; you probably refer to the section "Languages and
  compilers" ? ]

What I am thinking of is what I presented earlier, namely that in the
following root-finding function:

      SUBROUTINE FINDROOT(GUESS)
  10  FINDROOT = [some-expression-involving-GUESS]
      IF (FINDROOT .EQ. GUESS) RETURN
      GUESS = FINDROOT
      GOTO 10
      END

the .EQ. is a mistake.  Even when using exact rounding (i.e. the IEEE
Standard), it is impossible to prove that,  given some arbitrary GUESS
and arbitrary some-expression-involving-GUESS, the "solution" FINDROOT
and GUESS will *not* oscillate between two numbers one bit apart
indefinitely.

[ Note that in Fortran, the name of a function is a
  local variable which will magically turn into the
  function result upon return ]

This is the case that always fails on the i386 [we get a bug report for
something like this about once a month on fortran@gnu.org] because
FINDROOT will end up in a register (being freshly computed) and GUESS
will be retrieved from memory because some-expression will be so huge as
to drive all floating point values out of their registers.

In My Not So Humble Opinion, it is a good thing that this fails on the
i386, because it cannot be proven to always succeed on *any* floating
point implementation, so people better learn (or else ;-).

Regards,

-- 
Toon Moene (toon@moene.indiv.nluug.nl)
Saturnushof 14, 3738 XG  Maartensdijk, The Netherlands
Phone: +31 346 214290; Fax: +31 346 214286
g77 Support: fortran@gnu.org; egcs: egcs-bugs@cygnus.com