Bug in optimization under Linunx/x86

[jonathan dot desena at jhuapl dot edu]

Originator 
    Jonathan DeSena 
Organization 
    Johns Hopkins University Applied Physics Laboratory
Synopsis 
    g77 produces code which generates different results when optimized. 
Severity 
    serious 
Priority 
    medium/low
Category 
    fortran (others too??)
    optimization 
    target (Linux/x86)
Class 
    wrong-code ?
Release 
    gcc version 2.95.4 20010604
    GNU F77 version 2.95.4 20010604
    GNU Fortran Front End version 0.5.25 20010319
Environment
    Linux kernel 2.4.3
    Intel Pentium III
Description 
    Optimized code produces different results than unoptimized.  The
unoptimized code produces the expected results compared to doing the
operation by calculator, or using another computer (Solaris/SPARc with
f77 produces exact same results as unoptimized g77 linux/Intel).
    While the differences are right about at the limit for the precision
of
a 32 bit floating point number, they imply that the value in memory is
different, which I believe to be a bug.  More specifically, the lowest
bit appears to be being altered somehow.  To show this look at the
output from the test code included below.  I get:
unoptimized:
  0.381159384232887533E-05  0.274937897920608521E+00
-0.450124204158782959E+00
  0.381159384232887533E-05  0.159900173544883728E+00
-0.680199623107910156E+00
optimized:
  0.381159384232887533E-05  0.274937897920608521E+00
-0.450124233961105347E+00
  0.381159384232887533E-05  0.159900173544883728E+00
-0.680199682712554932E+00

    The first two numbers are identical as expected.  The third number,
which is the result of the arithmetic operation, varies.  

    From the unoptimized case ,the floating point,
-0.450124204158782959E+00, is represented by the hex 0xbee676ae (last
bit=0). But the optimized result, -0.450124233961105347E+00, is
0xbee676af (last bit=1).
Again -0.680199623107910156E+00 is 0xbf2e2190 (last bit=0)
and   -0.680199682712554932E+00 is 0xbf2e2191 (last bit=1)    

    This implies that optimized, since the last bit is not always
accurate, the REAL might as well be a 31 bit number rather than 32.

    Note that this problem only manifests itself under certain
conditions which are commented in the included code.

How-To-Repeat
      PROGRAM TEST
C     *This program demonstrates a difference in the results
C     *of optimized versus unoptimized code on Linux/x86.
C     *Optimized code must be compiled with -O2 or higher as
C     *written.  Code is included in comments to get it to work
C     *with only -O optimization.

	  INTEGER I1,I2
C     *the following ints will cause the problem, others may or may not
	  I1=72132
	  I2=41951
C     *uncomment the following and use B=5.0*A-3.0 below for 
C     *error with -O optimization
C      I2=154481

	  CALL TESTME(I1)
	  CALL TESTME(I2)
	  STOP
	  END

C     *A seperate subroutine is required for some reason
	  SUBROUTINE TESTME(I)

C     *The following parameter does not have to be the same as used
C     *here, but not all numbers will result in the problem given a
C     *certain integer, I, either.
	  PARAMETER (S=3.8115939E-6)
	  INTEGER I
	  REAL A,B

C     *The next line seems to be required for the problem to occur.
C     *Note that the result, A, is the same optimized or not.
	  A=FLOAT(I)*S

C     *Now do some simple arithmetic. this should not let the compiler
C     *do any reordering of operations: it must multiply and then 
C     *subtract. This step is required.  Multiplication or 
C     *addition/subtraction alone are not enough. The choice of 
C     *constants (2.0,1.0) is somewhat arbitrary, but changing them will
C     *cause certain integers to "work" or not.
	  B=2.0*A-1.0
C	  *Comment previous line and uncomment following line as well as 
C	  *I2=154481 above to show error with -O optimization
C      B=5.0*A-1.0

C     *Print out the results
	  PRINT 1, S,A,B

C     *Format with more digits than implied by DOUBLE PRECISION
 1    FORMAT (3E26.18)
	  RETURN
	  END
 
Fix 
    Only fix is to also compile using -ffloat-store, or to use DOUBLE
PRECISION instead of REAL.  This, however, is a somewhat unsatisfactory
solution to me, since I am unsure if this fixes all similar errors which
I have yet to isolate.

--
Jonathan DeSena
Johns Hopkins University Applied Physics Laboratory
Power Projection Systems Department
Aviation Systems Group