5

I would like to understand why something is happening. I needed to implement the integer square root in Python (isqrt(64) = 8 = isqrt(80)). I was convinced that the naive approach:

def isqrt(n):
    return int(math.sqrt(n))

was bound to fail occasionally when the n passed is a square number, assuming that Python converts to floats, then performs a square root calculation on floats. For instance, calling isqrt(13*13) I expected that after converting to floats and calculating sqrt, you could get something like 12.999999843, which after casting to integer would give us 12.

But I performed large loops testing values, big and small, and always get the correct result. It would seem there is no need to implement a special square root for integers, after all!

Not understanding bothers me, as much as when something that was supposed to work is failing. Why is this happening?

There is another question regarding integer square root in python: Integer square root in python

In the isqrt() defined there, a +0.5 is added to n, which I guess was included precisely to fight the issue I mentioned I was expecting, but cannot find in a specific case.

EDIT: Forgot to specify, I am using Python 2.7

Community
  • 1
  • 1
zeycus
  • 860
  • 1
  • 8
  • 20
  • _Why is this happening?_ http://stackoverflow.com/questions/588004/is-floating-point-math-broken – Łukasz Rogalski Jan 21 '16 at 08:20
  • In Python 2.7.6 I get `math.sqrt(5) == 2.23606797749979`, so `sqrt` results are always floating point even for small integer arguments. But I think that changed in Python 3. Python is dynamically typed, so can in principle change the type of the result as needed - giving an integer result (including large integers) where the argument is integer if that's what the powers that be decided. It makes a lot of sense, but I haven't checked (hence comment not answer). –  Jan 21 '16 at 08:25
  • 1
    @Rogalski My problem here is, why is this (what your link describes) *not* happening? – zeycus Jan 21 '16 at 08:35
  • In CPython 2.7, `math.sqrt()` is directly forwarded (complex number handling magic aside) to the corresponding function in `libm` (note that `Decimal.sqrt()` is different). Not surprisingly, a [pure C version](https://gist.github.com/anonymous/0c88c85b0262d025ad17) has the same behaviour – dhke Jan 21 '16 at 08:50

1 Answers1

4

Using python 2.7 on a machine with the C type long as a 64 bit integer and C type double implemented as a 64 bit IEEE floating point number yields

>>> import math
>>> x = (2<<53) + 1
>>> int(math.sqrt(x*x)) == x
False

I cheated and chose a number that 64 bit IEEE floating point can't represent exactly (but python's integer type can), (2<<53) + 1. Python 2.7 computes x*x as a python 2.7 long integer. (Note: This is distinct from a C long; python 2.7 can represent 2<<600 as an integer, but C cannot.)

Speaking of 2<<600,

>>> import math
>>> x = 2<<600
>>> int(math.sqrt(x*x)) == x
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
OverflowError: long int too large to convert to float
David Hammen
  • 32,454
  • 9
  • 60
  • 108
  • @Davin Nice example, you tested far beyond what I had tried. If I understand, what is happening is that for square numbers with an exact floating-number representation, the square root is also exactly calculated. Something I did not know! – zeycus Jan 21 '16 at 11:17
  • FWIW, assuming IEEE 754 semantics (with C's `double` matching IEEE 754's binary64), the smallest `n` for which `int(math.sqrt(n))` doesn't give the right result is `n = 2**52 + 2**27`. – Mark Dickinson Dec 30 '19 at 16:22