strange clusters of java.lang.Double with respect to multiplication

Question

Something really odd about these java.lang.Double things... Here's the [Clojure] code - same as Java, for all [execution-related] intents and purposes:

(doseq [x (range 100)
        :let [x33 (str (+ x 0.33))
              x99 (str (+ x 0.99))
              tx33 (* 100 (Double/parseDouble x33))
              tx99 (* 100 (Double/parseDouble x99))]]
  (do (prn (str x33 "->" tx33))
      (prn (str x99 "->" tx99))))

For Java folks - this goes through the range of java.lang.Long-s, creates two java.lang.String-s for each (x33 and x99), by adding 0.33 and 0.99 to it, respectively. It thin parses the Strings to Doubles, and multiplies the result by 100, and prints the mapping. Observe the results (the nil at the end is irrelevant):

"0.33->33.0"
"0.99->99.0"
"1.33->133.0"
"1.99->199.0"
"2.33->233.0"
"2.99->299.0"
"3.33->333.0"
"3.99->399.0"
"4.33->433.0"
"4.99->499.0"
"5.33->533.0"
"5.99->599.0"
"6.33->633.0"
"6.99->699.0"
"7.33->733.0"
"7.99->799.0"
"8.33->833.0"
"8.99->899.0"
"9.33->933.0"
"9.99->999.0"
"10.33->1033.0"
"10.99->1099.0"
"11.33->1133.0"
"11.99->1199.0"
"12.33->1233.0"
"12.99->1299.0"
"13.33->1333.0"
"13.99->1399.0"
"14.33->1433.0"
"14.99->1499.0"
"15.33->1533.0"
"15.99->1599.0"
"16.33->1632.9999999999998"
"16.99->1698.9999999999998"
"17.33->1732.9999999999998"
"17.99->1798.9999999999998"
"18.33->1832.9999999999998"
"18.99->1898.9999999999998"
"19.33->1932.9999999999998"
"19.99->1998.9999999999998"
"20.33->2032.9999999999998"
"20.99->2099.0"
"21.33->2133.0"
"21.99->2199.0"
"22.33->2233.0"
"22.99->2299.0"
"23.33->2333.0"
"23.99->2399.0"
"24.33->2433.0"
"24.99->2499.0"
"25.33->2533.0"
"25.99->2599.0"
"26.33->2633.0"
"26.99->2699.0"
"27.33->2733.0"
"27.99->2799.0"
"28.33->2833.0"
"28.99->2899.0"
"29.33->2933.0"
"29.99->2999.0"
"30.33->3033.0"
"30.99->3099.0"
"31.33->3133.0"
"31.99->3199.0"
"32.33->3233.0"
"32.99->3299.0"
"33.33->3333.0"
"33.99->3399.0"
"34.33->3433.0"
"34.99->3499.0"
"35.33->3533.0"
"35.99->3599.0"
"36.33->3633.0"
"36.99->3699.0"
"37.33->3733.0"
"37.99->3799.0"
"38.33->3833.0"
"38.99->3899.0"
"39.33->3933.0"
"39.99->3999.0"
"40.33->4033.0"
"40.99->4099.0"
"41.33->4133.0"
"41.99->4199.0"
"42.33->4233.0"
"42.99->4299.0"
"43.33->4333.0"
"43.99->4399.0"
"44.33->4433.0"
"44.99->4499.0"
"45.33->4533.0"
"45.99->4599.0"
"46.33->4633.0"
"46.99->4699.0"
"47.33->4733.0"
"47.99->4799.0"
"48.33->4833.0"
"48.99->4899.0"
"49.33->4933.0"
"49.99->4999.0"
"50.33->5033.0"
"50.99->5099.0"
"51.33->5133.0"
"51.99->5199.0"
"52.33->5233.0"
"52.99->5299.0"
"53.33->5333.0"
"53.99->5399.0"
"54.33->5433.0"
"54.99->5499.0"
"55.33->5533.0"
"55.99->5599.0"
"56.33->5633.0"
"56.99->5699.0"
"57.33->5733.0"
"57.99->5799.0"
"58.33->5833.0"
"58.99->5899.0"
"59.33->5933.0"
"59.99->5999.0"
"60.33->6033.0"
"60.99->6099.0"
"61.33->6133.0"
"61.99->6199.0"
"62.33->6233.0"
"62.99->6299.0"
"63.33->6333.0"
"63.99->6399.0"
"64.33->6433.0"
"64.99->6498.999999999999"
"65.33->6533.0"
"65.99->6598.999999999999"
"66.33->6633.0"
"66.99->6698.999999999999"
"67.33->6733.0"
"67.99->6798.999999999999"
"68.33->6833.0"
"68.99->6898.999999999999"
"69.33->6933.0"
"69.99->6998.999999999999"
"70.33->7033.0"
"70.99->7098.999999999999"
"71.33->7133.0"
"71.99->7198.999999999999"
"72.33->7233.0"
"72.99->7298.999999999999"
"73.33->7333.0"
"73.99->7398.999999999999"
"74.33->7433.0"
"74.99->7498.999999999999"
"75.33->7533.0"
"75.99->7598.999999999999"
"76.33->7633.0"
"76.99->7698.999999999999"
"77.33->7733.0"
"77.99->7798.999999999999"
"78.33->7833.0"
"78.99->7898.999999999999"
"79.33->7933.0"
"79.99->7998.999999999999"
"80.33->8033.0"
"80.99->8098.999999999999"
"81.33->8133.0"
"81.99->8199.0"
"82.33->8233.0"
"82.99->8299.0"
"83.33->8333.0"
"83.99->8399.0"
"84.33->8433.0"
"84.99->8499.0"
"85.33->8533.0"
"85.99->8599.0"
"86.33->8633.0"
"86.99->8699.0"
"87.33->8733.0"
"87.99->8799.0"
"88.33->8833.0"
"88.99->8899.0"
"89.33->8933.0"
"89.99->8999.0"
"90.33->9033.0"
"90.99->9099.0"
"91.33->9133.0"
"91.99->9199.0"
"92.33->9233.0"
"92.99->9299.0"
"93.33->9333.0"
"93.99->9399.0"
"94.33->9433.0"
"94.99->9499.0"
"95.33->9533.0"
"95.99->9599.0"
"96.33->9633.0"
"96.99->9699.0"
"97.33->9733.0"
"97.99->9799.0"
"98.33->9833.0"
"98.99->9899.0"
"99.33->9933.0"
"99.99->9999.0"
nil

As you can see, there are two clusters of such Doubles which behave differently from the rest. A mystery to me... God knows how many more there are down the path to infinity. Can anyone shed some light on this? Thanks! I am aware of the fact that not all rational numbers can be represented in a floating-point notation, but I think this is something else..

@RedDeckWins - so this is just it, "not all rationals can be represented as a floating points"? I get that - but, visually, it would appear that some number are, um... "especially" impossible to represent? I guess, my question is about this fact - why some are "more impossible" than others. — alexakarpov, Sep 02 '14 at 19:07
So, you're asking for the precise mathematical properties of the input numbers that account for the observed behavior here? Extending the bounds of your range to 1000, I observe that the clusters you describe occur at x = 16, 64, 128, and 512. — Alex, Sep 02 '14 at 19:13
ran it for range of 0..2050, and this is what I see: 1) 16 to 20, but not 20.99, are 'special'. So it's not just 16 by itself. 2) 64 to 80, but only their .99-s, not .33-s, are special. 3) 128 to 163, this time .33-s only. 4) 512 to 655, only .33-s. 5) 2048 to 2098, .99-s only. So, certain powers of 2 don't have this problem, while others do. Time to read this, I guess: http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html — alexakarpov, Sep 02 '14 at 21:00

Pascal Cuoq · Accepted Answer · 2014-09-04T21:45:30.507

IEEE 754 binary floating-point numbers should be thought of as organized in binades, intervals ranging from a power of two to the next, in which all numbers are separated by the same minimum interval (ULP).

"15.99->1599.0"
"16.33->1632.9999999999998"
...
"20.33->2032.9999999999998"
"20.99->2099.0"

The above cluster corresponds to a more-inaccurate-than-usual computation when the origin value is in the 16-32 binade and the result is in the 1024-2048 binade.

"63.99->6399.0"
"64.33->6433.0"
"64.99->6498.999999999999"
"65.33->6533.0"
"65.99->6598.999999999999"
"66.33->6633.0"
"66.99->6698.999999999999"
...
"80.33->8033.0"
"80.99->8098.999999999999"
"81.33->8133.0"
"81.99->8199.0"

For some reason, there is another cluster where one computation out of two is more inaccurate than usual, when the origin value is in the 64-128 binade and the result in the 4096-8192 binade.

You can observe what happens more precisely in each of these two clusters by looking at the computations at the binary level. Computations for values where the input value and the output value are respectively in the same binades will look very similar (and that's how they are either all accurate or all inaccurate at the same time).

The C99 program below makes the computation appear as regular as it is, within the first cluster (this program should be trivial to translate to Clojure, as long as you can make use of Java's printing functions and of the %a format for hexadecimal):

#include <stdio.h>

int main(int argc, char** argv) 
{  
  printf("%a -> %.13a (exact result %a)\n", 15.33, 15.33 * 100.0, 1533.0);
  printf("%a -> %.13a (exact result %a)\n", 16.33, 16.33 * 100.0, 1633.0);
  printf("%a -> %.13a (exact result %a)\n", 17.33, 17.33 * 100.0, 1733.0);
  printf("%a -> %.13a (exact result %a)\n", 18.33, 18.33 * 100.0, 1833.0);
  printf("...\n");
  printf("%a -> %.13a (exact result %a)\n", 20.33, 20.33 * 100.0, 2033.0);
  printf("%a -> %.13a (exact result %a)\n", 21.33, 21.33 * 100.0, 2133.0);
}

The output is:

0x1.ea8f5c28f5c29p+3 -> 0x1.7f40000000000p+10 (exact result 0x1.7f4p+10)
0x1.0547ae147ae14p+4 -> 0x1.983ffffffffffp+10 (exact result 0x1.984p+10)
0x1.1547ae147ae14p+4 -> 0x1.b13ffffffffffp+10 (exact result 0x1.b14p+10)
0x1.2547ae147ae14p+4 -> 0x1.ca3ffffffffffp+10 (exact result 0x1.ca4p+10)
...
0x1.4547ae147ae14p+4 -> 0x1.fc3ffffffffffp+10 (exact result 0x1.fc4p+10)
0x1.5547ae147ae14p+4 -> 0x1.0aa0000000000p+11 (exact result 0x1.0aap+11)

Above, a p+4 exponent in the hexadecimal representation of a number means the number is in the 16-32 binade. As you can see, all X.33 numbers in this binade have similar representations. The computations on them lead to similar exact results, but when the result is below 2048, it has to be rounded differently than when the exact result is above 2048.

Multiplying a number that ends in …ae147ae14 by 100 produces an exact result that ends in …fffffffd, or …11111111101 in binary.

When the number to multiply by 100 is 16.33, the last two binary digits need to be cut off for rounding. The exact result of the computation, ending in binary …11111111101, is closer to the previous representable number than to the next representable number, so the result of the multiplication is rounded to something ending in …111111111 in binary, or …fff in hexadecimal (or …99998 in decimal).
When the number to multiply by 100 is 21.33, because this time the result of the multiplication falls in the next binade, the last three binary digits need to be cut off for rounding. The last three binary digits being cut off are 101. They represent more than half the ULP, so the nearest representable value is the representable value above the exact value. This value, which is chosen by the FPU as result of the multiplication, ends in …00000 and is exactly the floating-point representation of 2133.0.

thank you for a very detailed an informative answer! – alexakarpov Sep 04 '14 at 21:33 — alexakarpov, Sep 04 '14 at 21:33

score -1 · Answer 2 · answered Sep 02 '14 at 18:56

-1

You're hitting a common problem with floating point numbers. If you need exact results you should use BigDecimal.

answered Sep 02 '14 at 18:56

Simone Gianni

11,426
40
49

3

Or Clojure's ratios in some cases. – Arthur Ulfeldt Sep 02 '14 at 19:04
1

And how is BigDecimal better at giving exact results? – tmyklebu Sep 03 '14 at 09:29
I'm not against downvotes, but leaving a comment with the motivation helps improve. – Simone Gianni Sep 03 '14 at 15:42

strange clusters of java.lang.Double with respect to multiplication

2 Answers2