From Bugzilla Helper: User-Agent: Mozilla/5.0 (Macintosh; U; PPC Mac OS X Mach-O; en-US; rv:1.7.5) Gecko/20041201 Camino/0.8.2 Description of problem: In a Fortran program compiled with g77, the behavior differs depending on whether or not the use of things like DSIN() and DCOS() are used in-line in an expression, versus first being assigned to DOUBLE PRECISION variables and then using the variables in the expression. Version-Release number of selected component (if applicable): gcc-g77-3.2.3-49 How reproducible: Always Steps to Reproduce: 1. Save the attachment test case program as "mollweide_bug.f" (or download it from the URL supplied with this bug report) 2. Compile it with "g77 -g -o mollweide_bug mollweide_bug.f" 3. Run the resultant "mollweide_bug" binary and note that the results fo r the two routines (which should be identical) are, in fact, not at all identical Actual Results: % mollweide_bug THIS VERSION WORKS: DEL_THETA = 6.98536454E-05 DEL_THETA = -2.59862096E-05 DEL_THETA = -1.46490641E-05 DEL_THETA = -5.82396393E-06 DEL_THETA = -9.55089071E-07 DEL_THETA = -2.45142533E-08 DEL_THETA = 0. AND THIS VERSION DOES NOT: DEL_THETA = 6.95029915E-05 DEL_THETA = -2.58266687E-05 DEL_THETA = -1.45964678E-05 DEL_THETA = -5.69967789E-06 DEL_THETA = -1.07194915E-06 DEL_THETA = 9.96385755E-08 DEL_THETA = 9.91175102E-08 DEL_THETA = -5.11485558E-07 DEL_THETA = 1.01300971E-07 DEL_THETA = 1.00757952E-07 Expected Results: The output from the two routines should be the same. Additional info: Linux pauling 2.4.21-27.0.1.ELsmp #1 SMP Mon Dec 20 18:47:45 EST 2004 i686 i686 i386 GNU/Linux
Created attachment 110087 [details] Test case - Fortran source file Compile via "g77 -g -o mollweide_bug mollweide_bug.f", run resultant binary.
The reason for this is mis-designed i387 FPU unit. All FPU computation are done in one selected precision (by default 80bit extended double) and rounding to other precisions happens only when storing the values into memory. Not sure about Fortran, but in C it is certainly allowed to compute expressions in excess precision (as long as downcasts result in rounding). See e.g. http://gcc.gnu.org/ml/gcc/2005-01/msg01082.html and associated thread for details. If you have SSE2 capable CPU (Pentium4+, recent AMD chips), you can fix this by using -msse2 -mfpmath=sse, otherwise simply live with the broken hardware and count with the possibility of excess precision.