float f = 1.0;
while (f != 0.0) f = f / 2.0;

This loop runs 150 times using 32-bit precision. Why is that so? Is it getting rounded to zero?

CodePudding user response：

In common C implementations, the IEEE-754 binary32 format is used for float. It is also called “single precision.” It is a binary based format where finite numbers are represented as ±f•2^e, where f is a 24-bit binary numeral in [1, 2) and e is an integer in [−126, 127].

In this format, 1 is represented as 1.00000000000000000000000₂•2⁰. Dividing that by 2 yields ½, which is represented as 1.00000000000000000000000₂•2⁻¹. Dividing that by 2 yields 1.00000000000000000000000₂•2⁻², then 1.00000000000000000000000₂•2⁻³, and so on until we reach 1.00000000000000000000000₂•2⁻¹²⁶.

When that is divided by two, the mathematical result is 1.00000000000000000000000₂•2⁻¹²⁷, but −127 is below the normal exponent range, [−126, 127]. Instead, the significand becomes denormalized; 2⁻¹²⁷ is represented by 0.10000000000000000000000₂•2⁻¹²⁶. Dividing that by 2 yields 0.01000000000000000000000₂•2⁻¹²⁶, then 0.00100000000000000000000₂•2⁻¹²⁶, 0.00010000000000000000000₂•2⁻¹²⁶, and so on until we get to 0.00000000000000000000001₂•2⁻¹²⁶.

At this point, we have done 149 divisions by 2; 0.00000000000000000000001₂•2⁻¹²⁶ is 2⁻¹⁴⁹.

When the next division is performed, the result would be 2⁻¹⁵⁰, but that is not representable in this format. Even with the lowest non-zero significand, 0.00000000000000000000001₂, and the lowest exponent, −126, we cannot get to 2⁻¹⁵⁰. The next lower representable number is 0.00000000000000000000000₂•2⁻¹²⁶, which equals 0.

So, the real-number-arithmetic result of the division would be 2⁻¹⁵⁰, but we cannot represent that in this format. The two nearest representable numbers are 0.00000000000000000000001₂•2⁻¹²⁶ just above it and 0.00000000000000000000000₂•2⁻¹²⁶ just below it. They are equally near 2⁻¹⁵⁰. The default rounding method is to take the nearest representable number and, in case of ties, to take the number with the even low digit. So 0.00000000000000000000000₂•2⁻¹²⁶ wins the tie, and that is produced as the result for the 150^th division.

CodePudding user response：

What happens is simply that your system has only a limited number of bits available for a variable, and hence limited precision; even though, mathematically, you can halve a number (!= 0) indefinitely without ever reaching zero, in a computer implementation that has a limited precision for a float variable, that variable will inevitably, at some stage, become indistinguishable from zero. The more bits your system uses, the more precision it has and the later this will happen, but at some stage it will.

Since I suppose this is meant to be C, I just implemented it in C (with a counter counting each iteration), and indeed it ran for 150 rounds until the loop ended. I also implemented it with a double, where it ran for 1075 iterations. Keep in mind, however, that the C standard does not define the exact precision of a float variable. In most implementations it's 32 bits for a float and 64 for a double. With a long double, I get 16,446 iterations.