12

If you take a number, take its square root, drop the decimal, and then raise it to the second power, the result should always be less than or equal to the original number.

This seems to hold true in python until you try it on 99999999999999975425 for some reason.

import math

def check(n):
    assert math.pow(math.floor(math.sqrt(n)), 2) <= n

check(99999999999999975424)  # No exception.
check(99999999999999975425)  # Throws AssertionError.

It looks like math.pow(math.floor(math.sqrt(99999999999999975425)), 2) returns 1e+20.

I assume this has something to do with the way we store values in python... something related to floating point arithmetic, but I can't reason about specifically how that affects this case.

3
  • Similar question: Why is math.sqrt() incorrect for large numbers?
    – wjandrea
    Jan 20, 2022 at 3:14
  • Note that you'll also get an exception for n = 4503599761588224 on most machines, despite the fact that both 4503599761588224 and the floor of its square root (which is 67108864) are small enough to be exactly representable in IEEE 754 binary64 floating-point. So it's a bit more subtle than just using numbers larger than floating-point can represent. Feb 6, 2022 at 19:19
  • 1
    (BTW, I don't think the duplicate is appropriate here. There's a lot more that can usefully be said about this specific problem that isn't contained in the answers to that duplicate.) Feb 6, 2022 at 19:24

2 Answers 2

21

The problem is not really about sqrt or pow, the problem is you're using numbers larger than floating point can represent precisely. Standard IEEE 64 bit floating point arithmetic can't represent every integer value beyond 52 bits (plus one sign bit).

Try just converting your inputs to float and back again:

>>> int(float(99999999999999975424))
99999999999999967232
>>> int(float(99999999999999975425))
99999999999999983616

As you can see, the representable value skipped by 16384. The first step in math.sqrt is converting to float (C double), and at that moment, your value increased by enough to ruin the end result.

Short version: float can't represent large integers precisely. Use decimal if you need greater precision. Or if you don't care about the fractional component, as of 3.8, you can use math.isqrt, which works entirely in integer space (so you never experience precision loss, only the round down loss you expect), giving you the guarantee you're looking for, that the result is "the greatest integer a such that a² ≤ n".

3
  • I'd suggest you to add '64-bit float' to specify that larger floats would handle the larger value.
    – mike3996
    Jan 22, 2018 at 9:13
  • @progo: I made the note. Not relevant to base Python, where there is only float (64 bit C double), but yeah, if you had access to other floating point types (e.g. through gmpy2.mpfr with specified precision) you could perform more precise computations. Or you use the decimal module to get fixed decimal digit precision. Jan 23, 2018 at 23:38
  • Understood. People in general have so many misconceptions about floats vs decimals I think every little bit of clarity helps
    – mike3996
    Jan 24, 2018 at 12:08
8

Unlike Evan Rose's (now-deleted) answer claims, this is not due to an epsilon value in the sqrt algorithm.

Most math module functions cast their inputs to float, and math.sqrt is one of them.

99999999999999975425 cannot be represented as a float. For this input, the cast produces a float with exact numeric value 99999999999999983616, which repr shows as 9.999999999999998e+19:

>>> float(99999999999999975425)
9.999999999999998e+19
>>> int(_)
99999999999999983616L

The closest float to the square root of this number is 10000000000.0, and that's what math.sqrt returns.

Not the answer you're looking for? Browse other questions tagged or ask your own question.