3

In adherence with IEEE-754, the following

#include <iostream>
#include <cmath>

int main() {
  std::cout << std::exp(-1/0.) << "\n";
  return 0;
}

yields 0 as the result. The reason being, -1/0.0 gives -INF, and exp of that gives 0.

Now, I need my code to run with -ffast-math. I understand that this option makes, in particular, infinities non-IEEE-conformant. So, I don't think I can rely on the intermediate -INF result.

However, the exponential function is (at least, mathematically) very stable for negative arguments: any finite value < -1e+3 already yields 0 as the result. So, I would be completely ok if the compiler just replaced the infinity with a large finite value.

At least in g++-5.4, this seems to work reliably, but is it prudent to rely on this?


I'm not yet convinced by the comments and answer saying that the original code is unreliable, but here's what I've decided to go with to be sure:

#include <iostream>
#include <cmath>
#include <limits>

int main() {
  double a = 0.;   // In practice, this is actually a small, nonnegative,
                   // usually-but-not-always >0 number
  double tiny = 1e-50;
  double b = std::exp(-1/(a+tiny));
  std::cout << b << "\n";
  return 0;
}

This is ugly, but it ensures that the divisor is always greater than zero whilst avoiding any conditional.

I considered using epsilon instead of tiny, but it doesn't really express the intent any better. (In fact it would skew the results; what I actually want to express is a + tiny == a for any a such that exp(-1/a) is defined.)

leftaroundabout
  • 117,950
  • 5
  • 174
  • 319
  • As for your problem, what *is* your problem? Why would you ever want to do division with zero? Don't do bad things, and you should not have to worry about these kind of things. – Some programmer dude May 17 '18 at 11:01
  • @Someprogrammerdude I don't _want_ to do division with zero, but I need the neg-inverse-exponential of a value that can be very small. Mathematically, it's always `>0`, but I can't guarantee that it doesn't come out exactly zero numerically. So I wonder if I need to explicitly catch that case, if IEEE already guarantees that the zero case is handled _as if_ it were finite. – leftaroundabout May 17 '18 at 11:05
  • *"The result of the / operator is the quotient from the division of the first operand by the second; the result of the % operator is the remainder. In both operations, **if the value of the second operand is zero, the behavior is undefined.**"* -- [C11/6.5.5p5](https://port70.net/~nsz/c/c11/n1570.html#6.5.5p5)... You should probably write an example that doesn't invoke undefined behavior. – autistic May 17 '18 at 11:06
  • 1
    @Sebivor that's the integer division operator you're talking about. – leftaroundabout May 17 '18 at 11:07
  • @leftaroundabout Citation please... – autistic May 17 '18 at 11:07
  • Values very close to zero are still not zero and therefore valid in a division, it's the division with *exactly* zero that's invalid. So why not simply have a check for *exactly* zero? – Some programmer dude May 17 '18 at 11:08
  • @Someprogrammerdude well, an explicit check would clunk up the code and conceivably introduce a branch at runtime. Unneeded code is always bad... thus my question if it _is_ needed. – leftaroundabout May 17 '18 at 11:18
  • I just want to make it clear that there is no such thing as "the integer division operator" in C. There's "the arithmetic division operator", which is what I linked to, and that covers both integer and floating point. You can see this in [C11/6.5.5p2](https://port70.net/~nsz/c/c11/n1570.html#6.5.5p2), just three points above what I linked to: *"Each of the operands shall have arithmetic type. The operands of the % operator shall have integer type."* Your floating point division by zero *is* undefined behavior. – autistic May 17 '18 at 11:18
  • @Sebivor ok, interesting. I didn't ask about C though. – leftaroundabout May 17 '18 at 11:19
  • @leftaroundabout The C++ standard has similar clauses in [expr.mul](http://cpp14.centaur.ath.cx/expr.mul.html). – autistic May 17 '18 at 11:24
  • If you're targeting GCC exclusively, then [macros like `likely()` or `unlikely()`](https://kernelnewbies.org/FAQ/LikelyUnlikely) can help you avoid branch problems, if it's unlikely that you will have exactly zero. But whatever you do, don't rely on the behavior of UB. – Some programmer dude May 17 '18 at 11:25
  • `-ffast-math optimisations` is particular to a given compiler. I recommend to add that compiler to the tags list. You may have to sacrifice a tag though, perhaps floating-point or ieee-754. – chux - Reinstate Monica May 17 '18 at 15:33
  • The days of `-ffast-math` being a 'win' for performant code are pretty much over. They're most likely to be a source of bugs that are subtle and hard to isolate. Floating point resources and performance on modern CPUs make this flag a bad decision - if for no other reason in that it generates incorrect code. – Brett Hale May 20 '18 at 07:44
  • @BrettHale that may well be the case for some domains, but we're doing computational plasma physics. The last time we did a full-program comparison for one of our simulators (that was an old Fortran one, but compiled with modern GCC), `-ffast-math` did make a significant difference, about 20% if I recall correctly. You're right that it would probably be a good idea to profile this more consequently, but I definitely can't introduce code _knowing_ it may break with `-ffast-math`. – leftaroundabout May 20 '18 at 10:10

1 Answers1

4

No, it is not prudent to rely on this. As mentioned in the stackoverflow question you linked to, if you rely in any way on infinities being handled correctly, and you still want to use -ffast-math, you should also use the -fno-finite-math-only flag.

G. Sliepen
  • 7,637
  • 1
  • 15
  • 31