2

Consider

int main() {
   double d = 1.0e+308;
   std::cout << (d*d)/1.0e+308;
}

compiled with g++ version 4.8.5 and the two compiler options -std=c++11 and -mfpmath=387 on a Linux CentOS system (more precisely, Linux 3.10.0-693.21.1.el7.x86_64 with Intel Xeon X5690 CPUs). As expected, this outputs the value 1e+308 because the multiplication d*d is evaluated in 80-bit extended precision and thus no overflow occurs (FLT_EVAL_METHOD returns 2 with these compiler settings). Now consider:

int main() {
   double d = 1.0e+308;
   double e;
   e = d*d;
   e = e/1.0e+308;
   std::cout << e;
}

This outputs inf - again as expected - since d*d is assigned to the double parameter e and in my understanding of the ISO/IEC 14882:2011 standard Sec. 5, 11, footnote 60, the value of e must be (or at least must behave as) a true double value, i.e. must be inf when used subsequently. However, and this is the point, when I add -O1 to the compiler options, then the output of the program is 1e+308. This violates - in my mind - the requirement that assignments (and casts) "perform their specific conversions" - see the passage in the c++11-standard document mentioned above. Did I misinterpret something here or is gcc not standard compliant (in this regard) when optimization level O1 (or higher) is used?

Deduplicator
  • 44,692
  • 7
  • 66
  • 118
Teilhart
  • 119
  • 6
  • Related: https://stackoverflow.com/questions/20869904/c-handling-of-excess-precision –  Jul 17 '18 at 20:21
  • I was under the impression that the standards do not specify any limit to the effective size of floating point values. As long as `double` meets minimum standard values it can have an infinite bit representation. Am I wrong about that? If so how did Intel get away with 80 bit floats for years? – Zan Lynx Jul 17 '18 at 20:24
  • @ZanLynx In C, C90 didn't specify much so it was implicitly allowed. C99 tightened the rules a bit so that excess precision is forbidden in some cases. In C++, I do not know enough to say whether it implicitly uses C99's rules. –  Jul 17 '18 at 20:29
  • Of course, I've always argued that if your algorithms rely on _less precision_ from floating point operations it is and always has been horribly broken. – Zan Lynx Jul 17 '18 at 20:32
  • GCC has a `-ffast-math` switch that licenses it to perform floating-point transformations not normally permitted by C/C++, and there may be similar switches such as -fsafe-math. I do not know what their default settings are in the system and version you are using, so you might check those. – Eric Postpischil Jul 18 '18 at 01:51
  • What did you expect? Consistency? That part of the std is gibberish anyway. – curiousguy Jul 22 '18 at 01:45

1 Answers1

1

This is the mentioned paragraph:

The values of the floating operands and 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.

And the footnote:

The cast and assignment operators must still perform their specific conversions as described in 5.4, 5.2.9 and 5.17.

I don't think that GCC violates this. e is represented in greater precision, and this is allowed by the standard (there is no conversion here, so the footnote doesn't apply).

If you don't like this behavior, use the -ffloat-store option, which removes this excess precision, and your program prints inf (but this option makes your program slower).

Note, this specification is a little bit weird.

Consider this. This is the same example as yours, but using float:

#include <iostream>

int main() {
   float d = 2.0e+38f;
   float e;
   e = d*d;
   e = e/2e+38f;
   std::cout << e;
}

This prints 2e+38, instead of inf. Now, if you change d's type to double: there is a conversion at e = d*d, footnote applies, and your program should print inf. And GCC behaves standard conformant, and it prints inf indeed (tested with gcc 5.4.1 and gcc 8.1.0).

geza
  • 28,403
  • 6
  • 61
  • 135
  • When I got this right, your argument is that the footnote from the standard addresses (or might address) only non-trivial type conversions (performed via explicit casts or assignments). However, when you replace the line `e = d*d;` in the 2nd program from my question by `e = d*d*1.0L`(and leave everything else as it is), then this requires the compiler to perform an implicit `long double` to `double` conversion on assigment, so the value on the right hand side of the assignment should be mapped to a double. But still, the output of the program I get on my system is finite. – Teilhart Jul 18 '18 at 17:30
  • @Teilhart: it seems a bug in GCC then. Try to change `d`'s type to long double, or use a non 1.0L constant, like 1.001L (it seems that 1.0L is handled specially, and it doesn't change `d*d*1.0L`'s type to `long double`). – geza Jul 18 '18 at 17:37
  • @Teilhart: the strange thing is, if I use `1.0L*d*d`, it works as intended, it only have problems with `d*d*1.0L`. – geza Jul 18 '18 at 17:39
  • On my system, 1.0L in the front gives `inf`, too; interestingly, if I replace `e = d*d;` by `e = static_cast(d*d);` then I again obtain a finite output value, so the explicit cast apparently remains without effect (`e = static_cast(d)*d;` gives `ìnf`, though); apart from this I observed that when I start with the 2nd example in my question and add just `<< std::endl ` right after `std::cout`, then this causes the output to turn to `inf` – Teilhart Jul 18 '18 at 18:38
  • @Teilhart: hmm, that doesn't make a difference on my machine. Are you sure that causes it? Can you reproduce it, if you toggle the presence of `std::endl`? – geza Jul 18 '18 at 19:11
  • Yes, it is reproducible - toggling the presence of `std::endl` on and off gives `inf`and the finite value alternatingly; my gcc-version, however, is older than the 2 versions you have mentioned in your answer above; but of course, also some other (subtle) difference in our system's configurations could explain our different observations; looking at the assembler code might reveal what's goinig on here but I'm not an expert in this field – Teilhart Jul 18 '18 at 19:44
  • @Teilhart: that's strange. Btw, why are you concerned about this issue? You could just use SSE, couldn't you? – geza Jul 18 '18 at 19:55
  • I work on a code simulating ice accretion on aircraft wings and tested (parts of) it with various compiler settings; rounding issues play an important role in this area since the PDEs are very stiff (up to 4th order derivatives!); that's why I came across these issues; but of course, SSE is the means of choice to avoid annomalies caused by spilling parameter values from registers to memory back and forth in an unpredicable way – Teilhart Jul 20 '18 at 05:36
  • "_this is allowed by the standard_" This contradicts everything in the definition of C++ and is "allowed" by a getaway clause that doesn't mean anything. It's such a mess that nobody wants to attack it frankly. Shame! – curiousguy Jul 22 '18 at 01:47