Floating point trouble with x86's extended precision

[Volker Reichelt <reichelt at igpm dot rwth-aachen dot de>]

In the process of revamping the non-bugs section of the bug reporting
instructions I came across a problem with the excess precision of the
x86 FPU:

I wanted to construct a very extreme case where rounding leads to
unexpected behaviour and came up with the following example:

============================================
int main()
{
  double          x=sqrt(2.0);
  volatile double y=x;
  if (x!=y) printf("Uh oh!\n");
  return 0;
}
============================================

Compiling this with GCC 3.3.1 and using -O2 I get the output:
"Uh oh!"

This is because GCC keeps x in a FPU register all the time (thus always
using the extended 80 bit precision) whereas y is written to memory (which
involves rounding to 64 bit) and reloaded into a FPU register, so that
the last couple of bits are now missing.

I claim that this is not a bug, since GCC should be allowed to optimize
away all stores to memory locations (and therefore optimize away all
rounding to 64 bit precision).

Brad Lucier, however, argues, that this should be considered a GCC bug.
(He's not opposed to keep intermediate results in a FPU register, if I
understand him correctly, but when an assignment to a variable takes
place, according to him, the value should be rounded to 64 bits.)

See the thread "[www-patch] bugs.html rewrite, part 2: revamp floating point
section" in gcc-patches:
http://gcc.gnu.org/ml/gcc-patches/2003-08/msg01055.html

An updated patch that isolates the problem can be found below.

The question is: Is this considered a GCC bug or not?

Regards,
Volker


Index: bugs.html
===================================================================
RCS file: /cvs/gcc/wwwdocs/htdocs/bugs.html,v
retrieving revision 1.75
diff -u -p -r1.75 bugs.html
--- bugs.html	20 Aug 2003 15:33:56 -0000	1.75
+++ bugs.html	20 Aug 2003 15:45:20 -0000
@@ -472,7 +472,30 @@ to round to the nearest representable nu
 <p>This is not a bug in the compiler, but an inherent limitation of
 the floating point types. Please study
 <a href="http://www.validlab.com/goldberg/paper.ps";>this paper</a>
-for more information.</p></dd>
+for more information.</p>
+
+<p>On some architectures, calculations inside the floating point unit (FPU)
+are done with higher precision than memory operations like storing floating
+point numbers (e.g., a <code>double</code> has 80 bits in the FPU on x86,
+but only 64 bits in memory). Thus, writing data to memory involves rounding,
+and the following C program might print <code>"Uh oh!"</code>.</p>
+<blockquote><pre>
+int main()
+{
+  double          x=sqrt(2.0);
+  volatile double y=x;
+  if (x!=y) printf("Uh oh!\n");
+  return 0;
+}
+</pre></blockquote>
+
+<p>Again, this is not a bug in the compiler, but a result of rounding
+(in conjunction with the excess precision of the FPU).</p>
+
+<p>You can force the compiler to store floating point numbers after each
+computation by using the command-line option "<code>-ffloat-store</code>",
+but this is only meant as a work-around and can massively slow down
+programs.</p></dd>
 </dl>
 
 <hr />