Challenge: find the first value where two functions differ
On Sat, 5 Aug 2017 06:17 am, Terry Reedy wrote:
> With 53 binary digits, all counts from 0 to 2**53 - 1 are exactly
> represented with a exponent of 0, 2**53 = 2**52 * 2, so it is exactly
> represented with an exponent of 1. Many other higher counts can be
> exactly represented with various exponents, but 2**53 + 1 requires 54
> digits. So it is the first count that does not have an exact binary64
Indeed. But although 2**53 + 1 cannot be represented as a float, there is no
difference between its integer square root and the integer square root of
2**53. Both are 94906265.
>>> We know that:
>>> - for n <= 2**53, isqrt_float(n) is exact;
> but it is not ;-). In this range, the only possibility is
> math.sqrt(m*m-1) being m instead of m - delta and that is the case for
> Chris' answer 4503599761588224 = 67108865**2 - 1.
> int(math.float(4503599761588224) is 67108865 instead of 67108864.
Indeed. I had forgotten one possibility.
I was thinking that there are only two cases when n is less than 2**53:
(1) n is a perfect square, in which case sqrt(n) is also a whole number which
can be represented exactly, and sqrt(n) is the same as integer sqrt(n);
(2) or n is not a perfect square, in which case sqrt(n) returns a number of the
form a.b, where a is the integer sqrt of n.
But in fact there's a third possibility: if the .b part of a.b is sufficiently
large that it rounds up to (a+1).0. It was this case that I had forgotten
about: the possibility of rounding up to the next whole number.
That's exactly what happens here with Chris' case. The exact mathematically
correct value of sqrt(4503599761588224) is a surd, but we can use Decimal to
get an approximate answer:
py> from decimal import Decimal
py> m = Decimal(4503599761588224)
Now 67108864.99999999254941951410... is not exactly representable as a binary
float. The two possibilities are:
67108864.99999999 # round down
67108865.00000000 # round up
Rounding down is *closer* to the correct value, which suggests that rounding up
is a bug:
py> m.sqrt() - Decimal('67108864.99999999') # round down
py> m.sqrt() - Decimal('67108865.0') # round up
Am I right that math.sqrt(4503599761588224) should round down?
?Cheer up,? they said, ?things could be worse.? So I cheered up, and sure
enough, things got worse.