Arithmetic involving a float and a "lesser" type should be done in float in C89, at least when optimization is turned off. In some cases (shown in the attached program), the arithmetic is done in double, carrying extra precision, which is *not* wanted, and makes the program results different from those on other machines.
Created attachment 5842 [details] Source code showing the bug Compile with -O0 with or without -DVF to show the bug, producing the (incorrect) result 12345.123456 The correct result seems to be obtained when compiling with -O3 and without -DVF, producing 12345.123047
Created attachment 5843 [details] Plain text description of the bug Full description of the bug and the possible test cases
Created attachment 5853 [details] Source code showing that float arithmetic is done in double The result of this code shows that gcc does float arithmetic in double, which it is not allowed to do when compiling with -O0 optimization. Both SGI and DEC C compilers produce the correct result, i.e. that addition of two floats and of two doubles should produce different results ... (whatever the optimization levels for them, by the way)
Created attachment 5854 [details] Source code showing that extra precision is kept also from functions This source code shows that gcc also keeps extra double precision out of functions declared as float, not even honoring the explicit cast (which is not needed here by the way). The code prints twice 12345.123456 whereas 12345.123047 should be printed for a machine running with 32 bits IEEE floats. Once again, even if it can be discussed whether this is a bug or not when optimizing the code, it is *not* allowed at -O0 or the compiler doesn't conform (and SGI and DEC compilers produce the right result whatever the optimization level)
I do not think this is a dup though. Where in the C standard says that this must work this way? Also what target are you are testing on?
Subject: Re: Invalid use of extra precision floating-point with -O0 optimization anz at obs-nice dot fr wrote: > This source code shows that gcc also keeps extra double precision out of > functions declared as float, not even honoring the explicit cast (which is > not needed here by the way). The gccbug.txt attachment states that this is an x86 target. So this is the classic x86 excess precision problem which has been known for over a decade. This is partly the fault of the x86 FPU, and partly the fault of the gcc x86 backend. The classic x86 FPU has no float or double arithmetic. It has only long double arithmetic. This means the only way to get the result you want is to emit extra instructions to convert results to float or double, but this makes the code run slower, and introduces double-rounding errors. So there is no way to win here. Gcc by design emits fast excess-precision results. This is compounded by problems in the gcc backend where it sometimes accidently loses the excess precision, resulting in unpredictable rounding errors. There are some things that can be done about this. 1) Use -ffloat-store. This forces user variables to be rounded but not intermediates, so it solves some but not all problems. This is the simplest solution. 2) Set the processor rounding more to float or double instead of the default long double. However, once you set the processor rounding mode, you can no longer use larger precisions, so this is useful only if your program uses only one precision. Changing this may break glibc math routines. 3) Do arithmetic in the SSE registers instead of the classic FP register stack. The SSE registers have float and double arithmetic instructions, and hence do not have the excess precision problem. This requires recent gcc versions and recent x86 processors. This is probably the best solution. 4) Don't use x86 processors for FP work. This is probably impractical for most people though. 5) Fix the gcc backend so that slow double-rounded code is an option. This is hard.
Subject: Re: New: Invalid use of extra precision floating-point with -O0 optimization anz at obs-nice dot fr wrote: > Arithmetic involving a float and a "lesser" type should be done in float in C89, > at least when optimization is turned off. This is false. K&R C required float operations to be performed in double. C89 allows float operations to be performed as float, but does not require it. In my ANSI C 89 manual, this is section 3.2.1.5 Usual Arithmetic Conversions, last sentence "The values of floating operands and of the results of floating expressions may be represented in greater precision and range than that required by the type; the types are not changed thereby." This gives the compiler the option of using float/double/long double for operations on float types. This is unfortunately necessary, because some hardware, the classic x86 FPU in particular, does not have a full set of operations. It only has one, which for linux, is long double by default.
This is a dup of bug 323 then, thanks Jim for the explanation. *** This bug has been marked as a duplicate of 323 ***
Created attachment 5867 [details] ANSWERS: this is not a duplicate of bug 323 and it is not fixed Answer to the answers I got: Even if arithmetic involving floats *can* be done with higher accuracy, casts and assignments have to be obeyed and cannot be optimized away.
(In reply to comment #8) > This is a dup of bug 323 then, thanks Jim for the explanation. > > *** This bug has been marked as a duplicate of 323 ***
I do not see any problems with this: tin:~/src/gnu/gcctest>gcc pr14384.c tin:~/src/gnu/gcctest>./a.out 12345.123456 tin:~/src/gnu/gcctest>gcc pr14384.c -O tin:~/src/gnu/gcctest>./a.out 12345.123047 tin:~/src/gnu/gcctest>gcc pr14384.c -O3 tin:~/src/gnu/gcctest>./a.out 12345.123047 tin:~/src/gnu/gcctest>gcc pr14384.c -mfpmath=sse pr14384.c:1: warning: SSE instruction set disabled, using 387 arithmetics ^[[Atin:~/src/gnu/gcctest>gcc pr14384.c -mfpmath=sse -msse tin:~/src/gnu/gcctest>./a.out 12345.123047
(In reply to comment #11) The problem doesn't lie with the source gccbug.c, which I suppose is the one you renamed to pr14384.c It lies with the source code gccbug1c.c, in the fourth attachment. This program, when run on my machine (see the 2nd attachment for the version of gcc and system info), prints 12345.123456 12345.123456 showing that the result of the function has *not* been narrowed to float, even though this *is* mandated by the C standard. If the narrowing casts are obeyed, as they should, at least the second result must be 12345.123047 on a machine running with 32-bit IEEE float, and this is what is obtained if compiled with -mfpmath=sse -msse As Jim Wilson said, the C standard *allows* for float expressions to be computed with higher accuracy, but if a (float) cast (or an assignment) is present, it *must* be obeyed (see my 5th attachment), and this is not the case with gcc (whatever the optimization level in that case). *IF* you don't obey this rule, then programs that determine themselve the machine accuracy in floating-point (like the famous paranoia program, see http://gcc.gnu.org/ml/gcc/2002-10/msg00950.html http://gcc.gnu.org/ml/gcc/2001-05/msg01583.html could fail in obscure ways. To summarize once more, even if f1+f2 *may* be computed with higher accuracy than float if f1 and f2 are float (float)(f1+f2) cannot, and the extra precision of the x87 processor *must* be discarded.
no you are missing the point that f1+f2 is not done in double precission but only done in single. d1 = ((double)f1)+offset; Will give you the addition done in double precission.
(In reply to comment #13) > no you are missing the point that f1+f2 is not done in double precission but only done in single. > > d1 = ((double)f1)+offset; > > Will give you the addition done in double precission. NO, if the addition was done in single, 12345.123047 would be printed (twice) (check with -mfpmath=sse -msse). The fact that 12345.123456 is printed twice proves that BOTH additions are done in double precision, because otherwise 12345.123456 would have been rounded to 12345.123047 because float has not enough accuracy to hold 12345.123456 exactly. And if you use the line you wrote > d1 = ((double)f1)+offset; you will see that it changes nothing because gcc computes with (double) extra precision anyway. And if it is OK to compute with extra precision in the first expression d1 = f1+offset; it is not in the second d2 = addf(f1); because there addf() has the return type float (and I even added an explicit cast to to float), so the result *HAS* to be rounded to single precision and produce 12345.123047 as result.
Well ICC 6.0 gives exactly the same result (even at -O0).
Created attachment 5893 [details] An even simpler source file demonstrating the problem This even simpler code shows that gcc doesn't round an expression to (float) and keeps the extra precision even when an explicit cast is present, where it shouldn't be allowed to do so. The same expression is computed twice, in (double) and in forced (float) with an explicit cast. They must differ, unless float and double are the same.
Subject: Re: Invalid use of extra precision floating-point with -O0 optimization anz at obs-nice dot fr wrote: > This even simpler code shows that gcc doesn't round an expression to (float) > and > keeps the extra precision even when an explicit cast is present, where it > shouldn't be allowed to do so. Yes. As I explained a week ago, this is the excess-precision problem that has been known for over a decade. There is no easy solution in the compiler because the design of the x86 FP register stack is flawed. You either have to use a compiler workaround, or else you have to complain to Intel, or else you have to stop doing FP work on the x86 FP register stack. Or you could try writing a fix for gcc yourself, but this will be hard, and will hurt performance, and will introduce double-rounding errors. So you may still have problems even with a fix. That is why we haven't already done it. Continuing to complain to us will only annoy us. We can't fix it.
Read up on the mailing lists for more information on why this is considered a non-bug and a dup of bug 323. The main problem is x87 (the FP part of the x86) is not suited for IEEE fp work unless you change some stuff, that is one reason why sse improved things in this area and one reason why in x86_64 x87 is not recomened any more. *** This bug has been marked as a duplicate of 323 ***