11

Recently, I wrote a small program and compiled it using mingw32(on Windows8) of 2 different versions. Surprisingly, I got two different reusults. I tried disassmbling it but found nothing special. Could anyone help me? Thank you.

the exe files: https://www.dropbox.com/s/69sq1ttjgwv1qm3/asm.7z

results: 720720（gcc version 4.5.2）, 720719（gcc version 4.7.0）

compiler flags: -lstdc++ -static

Code snipped as following:

#include <iostream>
#include <cmath>

using namespace std;

int main()
{
    int a = 55440, b = 13;
    a *= pow(b, 1);
    cout << a << endl;
    return 0;
}

Assembly output(4.5.2):

http://pastebin.com/EJAkVAaH

Assembly output(4.7.0):

http://pastebin.com/kzbbFGs6

Alexey Frunze
  • 61,140
  • 12
  • 83
  • 180
Zhe Chen
  • 2,918
  • 4
  • 23
  • 39
  • What are the results? And what are the compiler flags? – John Smith Mar 10 '13 at 06:17
  • Also, what are the gcc versions? It's most likely a compiler bug. I'm getting the correct answer (720720) with gcc 4.7.2 – John Smith Mar 10 '13 at 06:20
  • 1
    Why did you write such strange code? – Pubby Mar 10 '13 at 06:21
  • 1
    It's a problem found in a project, and for test purpose I write this code. @Pubby – Zhe Chen Mar 10 '13 at 06:24
  • results: 720720(gcc version 4.5.2), 720719(gcc version 4.7.0). I really want to know the reason behind the problem. @JohnSmith – Zhe Chen Mar 10 '13 at 06:26
  • Can you also give us the assembly outputs for main from 4.5.2 and 4.7.0? You can use the -S flag for this. – John Smith Mar 10 '13 at 06:38
  • 2
    In general, you should avoid using floating point arithmetics for integer calculations. Though, if you really need to mix them, consider rounding when converting to int: `a = int (0.5 + a*pow(b, 1));`. I suppose that `pow(b,1)` is returning 12.9999... instead of 13. – comocomocomocomo Mar 10 '13 at 07:03
  • @comocomocomocomo But I traced them, finding the results in FPU register are both correct. After excuting `fistp` command, the values stored into memory become different. Thus, I feel confused. – Zhe Chen Mar 10 '13 at 07:09
  • Sorry. All I know is that floating point arithmetics use to provoke small errors here and there. You have to try to minimize their effect instead of amplifying it. Regarding the exact reason in this case, I have no idea :-( – comocomocomocomo Mar 10 '13 at 07:28
  • 2
    @chenzhekl, how did you verify that the contents of the register were `13`? Were you able to see the extended precision bits? Perhaps your debugger/output was rounded to a standard `double` type showing `13` instead of showing the extended precision bits. – edA-qa mort-ora-y Mar 10 '13 at 10:23

1 Answers1

9

I've been able to reproduce the problem with a single version of the compiler.

Mine is MinGW g++ 4.6.2.

When I compile the program as g++ -g -O2 bugflt.cpp -o bugflt.exe, I get 720720.

This is the disassembly of main():

_main:
        pushl   %ebp
        movl    %esp, %ebp
        andl    $-16, %esp
        subl    $16, %esp
        call    ___main
        movl    $720720, 4(%esp)
        movl    $__ZSt4cout, (%esp)
        call    __ZNSolsEi
        movl    %eax, (%esp)
        call    __ZSt4endlIcSt11char_traitsIcEERSt13basic_ostreamIT_T0_ES6_
        xorl    %eax, %eax
        leave
        ret

As you can see, the value is calculated at compile time.

When I compile it as g++ -g -O2 -fno-inline bugflt.cpp -o bugflt.exe, I get 720719.

This is the disassembly of main():

_main:
        pushl   %ebp
        movl    %esp, %ebp
        andl    $-16, %esp
        subl    $32, %esp
        call    ___main
        movl    $1, 4(%esp)
        movl    $13, (%esp)
        call    __ZSt3powIiiEN9__gnu_cxx11__promote_2INS0_11__enable_ifIXaasrSt15__is_arithmeticIT_E7__valuesrS3_IT0_E7__valueES4_E6__typeES6_E6__typeES4_S6_
        fmuls   LC1
        fnstcw  30(%esp)
        movw    30(%esp), %ax
        movb    $12, %ah
        movw    %ax, 28(%esp)
        fldcw   28(%esp)
        fistpl  4(%esp)
        fldcw   30(%esp)
        movl    $__ZSt4cout, (%esp)
        call    __ZNSolsEi
        movl    $__ZSt4endlIcSt11char_traitsIcEERSt13basic_ostreamIT_T0_ES6_, 4(%esp)
        movl    %eax, (%esp)
        call    __ZNSolsEPFRSoS_E
        xorl    %eax, %eax
        leave
        ret
...
LC1:
        .long   1196986368 // 55440.0 exactly

If I replace the call to exp() with loading 13.0 like this:

_main:
        pushl   %ebp
        movl    %esp, %ebp
        andl    $-16, %esp
        subl    $32, %esp
        call    ___main
        movl    $1, 4(%esp)
        movl    $13, (%esp)

//        call    __ZSt3powIiiEN9__gnu_cxx11__promote_2INS0_11__enable_ifIXaasrSt15__is_arithmeticIT_E7__valuesrS3_IT0_E7__valueES4_E6__typeES6_E6__typeES4_S6_
        fildl    (%esp)

        fmuls   LC1
        fnstcw  30(%esp)
        movw    30(%esp), %ax
        movb    $12, %ah
        movw    %ax, 28(%esp)
        fldcw   28(%esp)
        fistpl  4(%esp)
        fldcw   30(%esp)
        movl    $__ZSt4cout, (%esp)
        call    __ZNSolsEi
        movl    $__ZSt4endlIcSt11char_traitsIcEERSt13basic_ostreamIT_T0_ES6_, 4(%esp)
        movl    %eax, (%esp)
        call    __ZNSolsEPFRSoS_E
        xorl    %eax, %eax
        leave
        ret

I get 720720.

If I set the same rounding and precision control fields of the x87 FPU control word for the duration of exp() as for the fistpl 4(%esp) instruction like this:

_main:
        pushl   %ebp
        movl    %esp, %ebp
        andl    $-16, %esp
        subl    $32, %esp
        call    ___main
        movl    $1, 4(%esp)
        movl    $13, (%esp)

        fnstcw  30(%esp)
        movw    30(%esp), %ax
        movb    $12, %ah
        movw    %ax, 28(%esp)
        fldcw   28(%esp)

        call    __ZSt3powIiiEN9__gnu_cxx11__promote_2INS0_11__enable_ifIXaasrSt15__is_arithmeticIT_E7__valuesrS3_IT0_E7__valueES4_E6__typeES6_E6__typeES4_S6_

        fldcw   30(%esp)

        fmuls   LC1
        fnstcw  30(%esp)
        movw    30(%esp), %ax
        movb    $12, %ah
        movw    %ax, 28(%esp)
        fldcw   28(%esp)
        fistpl  4(%esp)
        fldcw   30(%esp)
        movl    $__ZSt4cout, (%esp)
        call    __ZNSolsEi
        movl    $__ZSt4endlIcSt11char_traitsIcEERSt13basic_ostreamIT_T0_ES6_, 4(%esp)
        movl    %eax, (%esp)
        call    __ZNSolsEPFRSoS_E
        xorl    %eax, %eax
        leave
        ret

I get 720720 as well.

From this I can only conclude that exp() isn't calculating 131 precisely as 13.0.

It may be worth looking at the source code of that __gnu_cxx::__promote_2<__gnu_cxx::__enable_if<(std::__is_arithmetic<int>::__value)&&(std::__is_arithmetic<int>::__value), int>::__type, int>::__type std::pow<int, int>(int, int) to see exactly how it manages to screw up exponentiation with integers (see, unlike C's exp() it takes two ints instead of two doubles).

But I wouldn't blame exp() for that. C++11 defines float pow(float, float) and long double pow(long double, long double) in addition to C's double pow(double, double). But there's no double pow(int, int) in the standard.

The fact that the compiler provides a version for integer arguments does not make any additional guarantee about the precision of the result. If exp() calculates ab as

    ab = 2b * log2(a)

or as

    ab = eb * ln(a)

for floating-point values, there definitely can be rounding errors in the process.

If the "integer" version of exp() does something similar and incurs a similar loss of precision due to rounding errors, it still does its job right. And it does it even if the loss of precision is due to some silly bug and not because of the normal rounding errors.

However surprising this behavior may seem, it's correct. Or so I believe until proven wrong.

Alexey Frunze
  • 61,140
  • 12
  • 83
  • 180