C++ double division and return

Question

I understand that doubles are only approximations. But I was surprised that

double f(int x, int y)
{
    return double(x)/y;
}

and

double f(int x, int y)
{
    double z = double(x)/y;
    return z;
}

can return different values. Does anyone know why?

Solution 1

The main reason "why" is that the standard allows it (except that I'm not 100% sure that it does). The pragmatic reason is that on some processors (including Intel 32 bit processors), floating point registers have more precision than a double. Intermediate calculations are done in the registers, if they fit (and the standard clearly allows intermediate results to have more precision), so the results of any given expression will depend on how the compiler manages its registers, and when it spills to memory. And many (most?) compilers for such processors will return floating point values in a register, where they keep the extra precision. (I'm not sure that this is legal. When assigning a double, and I think when copy constructing a double, the value must be made to fit. Returning a value is copy construction, so I think it should force the double to fit. Regardless... most compilers don't so you have to deal with it, regardless of what the standard might or might not have to say about it.)

I think you can prevent this by explicitly casting the final results, i.e.:

return double( double(x)/y );

, but I'm not sure (neither about what the standard says about it, nor about what compilers actually do).

Note that whether the results are different or not depends somewhat on what you do with them. If you immediately assign them to a double, or pass them to another function which requires a double, then all that happens is that the rounding occurs a little later (but still on the same value). If you use the value in an expression, however, all bets are off. In your case, 2.0 * f( a, b ) might result in a different value depending on your implementation of f, even after rounding. And most noticeably, f( a, b ) == f( a, b ) may return false. (I've seen a case where this caused std::sort to crash. The comparison function used did something along the lines of return lhs.f() < rhs.f();, and was returning true when lhs and rhs were references to the same object; the compiler spilled the results of the first function it called to memory, then did the comparison with the value in the register which the second function returned.)

Solution 2

The only way I can think of that this is going to happen is with an optimizing compiler on an x87 FPU. The x87 FPU uses long double precision, so that values in registers can be held with more precision than values in double variables on the stack. It's at least conceivable that an optmising compiler given the above code might use a register for the return value in one case and the stack for the other.

But I would like to see a real example.

Comments

[Joseph Mansfield]

Do you have an example of that happening?

[LihO]

Post a complete test case with description of your problem.

[Tony The Lion]

I've helped you on your way to posting an example that shows the problem.

[fatihk]

if there is no pic, it din't happen

[JustAnotherCurious]

Hard to believe in this. What compiler, flags and OS and kernel OS and processor do you use?

[James Kanze]

@TonyTheLion Your example doesn't show the problem. The problem depends on when the compiler truncates an intermediate value with greater precision to a double. In your code, the value will be rounded to double when calling <<, even if it hadn't been rounded before.

[James Kanze]

@sftrabbit I've seen the problem happen in a comparison function for std::sort. The comparison returned something like lhs.f() < rhs.f(), and returned true when lhs and rhs were references to the same object (which ultimately caused std::sort to access out of bounds, and crashed the program).

[James Kanze]

@JustAnotherCurious I've seen it on Intel architectures, with g++, at least (but I think VC++ behaves similarly).

[Tony The Lion]

@JamesKanze, I know it doesn't show the problem, I was merely encouraging OP to create an example that creates the problem

[James Kanze]

FFIW, with VC++, you may have to specify /arch:IA32 for the problem to occur.

[john]

Or /arch:SSE2 to make it go away.

[Alexey Frunze]

@JustAnotherCurious I can reproduce the problem with x86(32-bit) MinGW gcc 4.6.2 if I compile this test with gcc -std=gnu89 fpdiff.c -o fpdiff.exe -O0. I can't seem to reproduce it with online compilers. On my PC it fails test 2.

[Pete Becker]

"doubles are only approximations" is not true. They have well-defined, precise semantics. Conversion from text to double usually results in an approximation, but that's not significantly different from converting, say, 1.2 to an int; the primary difference is that programmers have been taught from infancy how integer arithmetic works, and have internalized its quirks, but have generally not been taught how floating-point arithmetic works, so resort to "it's approximate" or "it's random" to justify seeing results they didn't expect. (That said, the problem discussed here is real)

[James Kanze]

@PeteBecker I think the reason most programmers say that "doubles are only approximations" is that they assume the abstraction for double is real numbers. (As abstractions go, double is a pretty poor abstraction of real numbers, given that most of the usual properties, like associativeness, don't hold.)

[Pete Becker]

@JamesKanze - you give programmers credit for far more understanding than I do.

[James Kanze]

It's almost certainly a case of the x87 FPU. And it's not that the compiler returns one in a floating point register, the other in an actual double; it's that the compiler enforces the rule that the returned value must be that contained in the variable (which has been rounded to double), either reading what it just wrote to return it, or explicitly forcing the rounding on the value in the FPU.