nly on Why doesn't GCC optimize a*a*a*a*a*a to (a*a*a)*(a*a*a)? (2013)

Comment on Why doesn't GCC optimize a*a*a*a*a*a to (a*a*a)*(a*a*a)? (2013)parent

Just to make your significant digits clearer:

    0.1 + 0.2 + 0.3   = 0.6000000000000000900
    0.1 + (0.2 + 0.3) = 0.5999999999999999800

    0.1 * 0.2 * 0.3   = 0.0060000000000000010
    0.1 * (0.2 * 0.3) = 0.0060000000000000001

Using -ffast-math doesn't effect these results.

haberman 12y

I find it illuminating to look at the exact values of the floating point numbers. This takes away an air of mystery: floating point numbers do in fact have very specific values, which we can inspect exactly.

Here's 0.1 + 0.2 + 0.3:

        0.1 (0.1000000000000000055511151231257827021181583404541015625)
      + 0.2 (0.200000000000000011102230246251565404236316680908203125)
           = 0.3000000000000000166533453693773481063544750213623046875
    rounded: 0.3000000000000000444089209850062616169452667236328125

      + 0.3 (0.299999999999999988897769753748434595763683319091796875)
           = 0.600000000000000033306690738754696212708950042724609375
    rounded: 0.600000000000000088817841970012523233890533447265625

The operands can't be represented exactly as floating point so they are rounded to the nearest floating-point number before the addition even happens. And once the addition does happen, it is rounded to the nearest representable floating point value.

Now we can look at the same for 0.1 + (0.2 + 0.3):

        0.2 (0.200000000000000011102230246251565404236316680908203125)
      + 0.3 (0.299999999999999988897769753748434595763683319091796875)
           = 0.500000000000000000000000000000000000000000000000000000
    rounded: 0.500000000000000000000000000000000000000000000000000000

      + 0.1 (0.1000000000000000055511151231257827021181583404541015625)
           = 0.6000000000000000055511151231257827021181583404541015625
    rounded: 0.59999999999999997779553950749686919152736663818359375

It's interesting to note that the first addition ends up being exact! This is a bit lucky, because 0.5 is exactly representable in floating point and the initial operands got rounded in opposite directions. It's also interesting that 0.6... only turns into 0.599999... at the point that we round the answer.

You can play around with this stuff conveniently using the Python "Decimal" module.

    $ python
    Python 2.7.5 (default, Aug 25 2013, 00:04:04)-
    [GCC 4.2.1 Compatible Apple LLVM 5.0 (clang-500.0.68)] on darwin
    Type "help", "copyright", "credits" or "license" for more information.
    >>> from decimal import Decimal, getcontext
    >>> getcontext().prec = 1000  # To prevent truncation/rounding
    >>> Decimal(0.1)
    Decimal('0.1000000000000000055511151231257827021181583404541015625')
    >>> a = Decimal(0.1) + Decimal(0.2)
    >>> a
    Decimal('0.3000000000000000166533453693773481063544750213623046875')
    >>> Decimal(float(a))  # get the value rounded to the nearest fp number
    Decimal('0.3000000000000000444089209850062616169452667236328125')

You can also use my tool "dumpfp" (https://github.com/haberman/dumpfp) to inspect the value of a floating-point number in detail:

    $ ./dumpfp 0.1
    Single Precision (IEEE 32-bit):
               raw = 0x3dcccccd
              sign = 0x0
          exponent = 0x7b (-4)
       significand = 0x4ccccd

       VALUE CALCULATION =
           significand   (1 + 5033165/2^23  (1.60000002384185791015625))
         * 2^exponent    (2^-4)
         = VALUE         (13421773/2^27  (0.100000001490116119384765625))

    Double Precision (IEEE 64-bit):
               raw = 0x3fb999999999999a
              sign = 0x0
          exponent = 0x3fb (-4)
       significand = 0x999999999999a

       VALUE CALCULATION =
           significand   (1 + 1351079888211149/2^51  (1.600000000000000088817841970012523233890533447265625))
         * 2^exponent    (2^-4)
         = VALUE         (3602879701896397/2^55  (0.1000000000000000055511151231257827021181583404541015625))

Scaevolus 12y

Note that recent Python versions use Gay's correctly-rounded conversions, so .1+(.2+.3) displays as .6, since that's equivalent when converted back to binary floating point. http://bugs.python.org/issue1580

whoopdedo 12y

If you want an exact representation, C99 has the %a format specifier for base-16 floating point. It seems that Python 3 doesn't support it but Lua 5.2 does.

    Lua 5.2.2  Copyright (C) 1994-2013 Lua.org, PUC-Rio
    > print(string.format("  0.1 (%a)\n+ 0.2 (%a)\n= %a\n+ 0.3 (%a)\n= %a", 0.1,0.2,0.1+0.2, 0.3, 0.1+0.2+0.3))
      0.1 (0x1.999999999999ap-4)
    + 0.2 (0x1.999999999999ap-3)
    = 0x1.3333333333334p-2
    + 0.3 (0x1.3333333333333p-2)
    = 0x1.3333333333334p-1
    > print(string.format("  0.2 (%a)\n+ 0.3 (%a)\n= %a\n+ 0.1 (%a)\n= %a", 0.2,0.3,0.2+0.3, 0.1, 0.2+0.3+0.1))
      0.2 (0x1.999999999999ap-3)
    + 0.3 (0x1.3333333333333p-2)
    = 0x1p-1
    + 0.1 (0x1.999999999999ap-4)
    = 0x1.3333333333333p-1

enupten 12y

I don't know if you realise that by converting a double "0.1" into the python decimal of a much larger size, you're bound to get garbage.

Edit: I see what you're doing now. I didn't appreciate enough how converting to decimals had such unpleasant effects :) Thanks.

haberman 12y

I'm not sure what you're trying to say. 0.1000000000000000055511151231257827021181583404541015625 is the precise value of float(0.1). It is not "garbage" or an approximation in any way. It is exact.

Comments