Here’s an interesting paragraph from Handbook of Floating Point Arithmetic by Jean-Michel Muller *et al*.

Many common transformations of a code into some

naivelyequivalent code become impossible if one takes into account special cases such as NaNs, signed zeros, and rounding modes. For instance, the expression`x + 0`

is not equivalent to the expression`x`

if`x`

is -0 and the rounding is to nearest.

Wait, what?!

First let’s see what this paragraph is *not* saying. What does the following code print?

#include <stdio.h> int main() { double t = 1e-200; double x = -1e-200 * t; // x == -0 double y = x + 0; printf( (x == y) ? "0 == 0\n" : "0 != 0\n"); }

I compiled the code above with gcc using default options, which means rounding to nearest.

It prints `0 == 0`

as you would have guessed if you weren’t spooked by the quotation above.

Let’s go back and read that paragraph more carefully. It speaks of **equivalent code**. In the code above, the variables `x`

and `y`

are equal in the sense that `x == y`

returns true, but that does not mean that `x`

and `y`

are *equivalent*. The two variables are not equal bit-for-bit and so some code could behave differently when given one rather than the other.

Let’s define

void hexdump(double x) { unsigned long* pu; pu = (unsigned long*)&x; printf("%lx\n", *pu); }

and add the lines

hexdump(x); hexdump(y);

to `main`

.

This prints the following.

8000000000000000 0

We see that `x`

has its top bit set but `y`

does not. More on why here.

You missed escaping <stdio.h>.

Thanks. I always forget that when I include C code here.