VERSION KNOWN TO BE AFFECTED: 3.4.3, 3.2.2
My gcc gives the following version/compilation info:
~/Projects/BUGS $ g++ -v
Reading specs from /usr/local/gcc-3.4.3/lib/gcc/i686-pc-linux-gnu/3.4.3/specs
Configured with: ../gcc-3.4.3/configure --enable-threads=posix
--prefix=/usr/local/gcc-3.4.3 --enable-languages=c,c++,f77 --with-system-zlib
Thread model: posix
gcc version 3.4.3
OK. Here is the origin of this bug: I have a formula:
alpha = 1.0 + (x * x) - x * sqrt(x * x + 2.0);
where alpha and x are double-precision numbers. Now, given a very large
value of x (for example, x >= 5000.0), then g++ (and also gcc) would
evaluate wrongly the formula above. Somehow, the outcome of the value
depends on the optimization level.
Let's take x = 20000.0 as a severe testcase. If you do Taylor expansion on
alpha for large x values (thus, expanding in terms of 1/x^2), the leading
term of alpha is
0.5 / (x * x)
So, the leading term is 1.25e-09, which is essentially a an accurate
answer. Precise numerical evaluation using Waterloo Maple (version 9, in
this case) yields
1.24999999687500000976562496582032e-09 ("exact result", 33-digits precision)
What happens if I evaluate the function above using gcc?
gcc command result (alpha)
g++ testcase.cpp -O0 2.526212483644485e-08
g++ testcase.cpp -O1 1.251464709639549e-09
g++ testcase.cpp -O2 1.251464709639549e-09
g++ testcase.cpp -O3 0
Use the snippet below as the testcase.
What is this? The -O0 compiled code is exceedingly wrong by predicting a
number 20 times bigger.
The -O1 and -O2 results are essentially correct. Why they disagree with
the exact answer is understandable, since the subtraction is between two
VERY LARGE numbers:
alpha = 400000001 - 20000*sqrt(400000002)
If the sqrt thing can be accurate, then we should be left with a small
number of residual precision, such as predicted by gcc-compiled code above
with -O1 and -O2.
NOTE that the same formula is used in another program (my own code, not
included here). There, the -O3 executable code gives the same number as
-O1 and -O2 above. So, there is no consistency here, although the gist of
the problem remains.
I'm sorry if this bug is probably a result of pushing numerics to its
edge. I should evaluate the expression using a series expansion instead of
naively using the formula like above. But I report this problem anyway, at
least to document the problem to the community.
Although this bug is (originally) rather specific to a formula, I think
that other like this one, which computes the difference between two large
numbers, could suffer from a similar problem.
REPRODUCING THE ERROR
Use the stock g++ 3.4.3. Run on an x86 machine. I'm not sure if AMD chips
would yield the same error. I bet it would, because of the executable
binary compatibility. But you have to test it.
This may have been a specific problem to the Intel x86 machines only (I
have tested this problem using Pentium M and Pentium 4 chips, see below).
I have NOT tested using AMD chips. This may not apply to other machines.
I have tried the `-mieee-fp' switch, it does not fix the problem.
This error is also reproducible if you compile the testcase using
`gcc -x c' (the C compiler) instead of g++ (the C++ compiler).
g++ version 3.2.2 running on an Intel Pentium 4 also reproduces the same
results above, except that the -O3 binary gives the same answer as -O2
above, instead of zero.
The g++'s identification string is:
[wirawan@snowflake wirawan]$ g++ -v
Reading specs from /usr/lib/gcc-lib/i386-redhat-linux/3.2.2/specs
Configured with: ../configure --prefix=/usr --mandir=/usr/share/man
--infodir=/usr/share/info --enable-shared --enable-threads=posix
--disable-checking --with-system-zlib --enable-__cxa_atexit
Thread model: posix
gcc version 3.2.2 20030222 (Red Hat Linux 3.2.2-5)
g++ version 3.1 on Alpha ev67 gives zero regardless of the optimization
level (-O0, -O1, -O2, -O3).
g++ version 3.3 on an UltraSPARC III box also gives zero regardless of the
optimization level (-O0, -O1, -O2, -O3).
SUN C++ (version "5.4 Patch 111715-08 2003/06/08") also gives zero
regardless of the optimization level (no switch, or -xO3).
Intel C++ version 7.1 (Build 20030617Z) gives 1.251464709639549e-09
regardless of the optimization level.
Created attachment 8147 [details]
A testcase to reproduce the error reported.
x86 is werid in that it almost always uses 80 bit float point and not IEEE 64 bit floating point.
This is a duplicate of bug 323.
*** This bug has been marked as a duplicate of 323 ***