2

In Python, are either

n**0.5  # or
math.sqrt(n)

recognized when a number is a perfect square? Specifically, should I worry that when I use 

int(n**0.5)  # instead of
int(n**0.5 + 0.000000001)

I might accidentally end up with the number one less than the actual square root due to precision error?

Andy Hayden
  • 359,921
  • 101
  • 625
  • 535
Sylvester V Lowell
  • 1,243
  • 1
  • 12
  • 13

7 Answers7

5

As several answers have suggested integer arithmetic, I'll recommend the gmpy2 library. It provides functions for checking if a number is a perfect power, calculating integer square roots, and integer square root with remainder.

>>> import gmpy2
>>> gmpy2.is_power(9)
True
>>> gmpy2.is_power(10)
False
>>> gmpy2.isqrt(10)
mpz(3)
>>> gmpy2.isqrt_rem(10)
(mpz(3), mpz(1))

Disclaimer: I maintain gmpy2.

casevh
  • 11,093
  • 1
  • 24
  • 35
4

Yes, you should worry:

In [11]: int((100000000000000000000000000000000000**2) ** 0.5)
Out[11]: 99999999999999996863366107917975552L

In [12]: int(math.sqrt(100000000000000000000000000000000000**2))
Out[12]: 99999999999999996863366107917975552L

obviously adding the 0.000000001 doesn't help here either...

As @DSM points out, you can use the decimal library:

In [21]: from decimal import Decimal

In [22]: x = Decimal('100000000000000000000000000000000000')

In [23]: (x ** 2).sqrt() == x
Out[23]: True

for numbers over 10**999999999, provided you keep a check on the precision (configurable), it'll throw an error rather than an incorrect answer...

Andy Hayden
  • 359,921
  • 101
  • 625
  • 535
  • The rounding of the result into the precision of a `float` is causing the problem rather than any deficiency in `sqrt` itself. The only way to avoid it is to use an integer square root function. `>>> '%0.0f' % float(100000000000000000000000000000000000) '99999999999999996863366107917975552'` – Mark Ransom Jun 06 '13 at 18:58
  • 2
    Note that you'd still have to make sure that the precision was sufficient when using `Decimal` -- it makes getting the right answer possible, not automatic. – DSM Jun 06 '13 at 19:18
  • You're getting ridiculous numbers because the thing you're squaring then square-rooting is not actually 1036, but rather the closest `double` to it. – tmyklebu Jun 07 '13 at 01:19
  • @tmyklebu yes... and hence worth worrying about. – Andy Hayden Jun 07 '13 at 01:47
  • @AndyHayden: My point is that it's a representation problem, not an arithmetic problem. – tmyklebu Jun 07 '13 at 01:50
  • @tmyklebu yes... and hence worth worrying about. – Andy Hayden Jun 07 '13 at 01:51
  • @AndyHayden Since the asker has picked yours as the accepted answer, perhaps you could add some text to explain the reasons behind the behaviour that you have demonstrated. That way readers can gain understanding. – David Heffernan Jun 07 '13 at 15:10
  • @DavidHeffernan Good advice, I will "beef it up" :) – Andy Hayden Jun 07 '13 at 15:22
3

Both **0.5 and math.sqrt() perform the calculation using floating point arithmetic. The input is converted to float before the square root is calculated.

Do these calculations recognize when the input value is a perfect square?

No they do not. Floating arithmetic has no concept of perfect squares.

large integers may not be representable, for values where the number has more significant digits than available in the floating point mantissa. It's easy to see therefore that for non-representable input values, n**0.5 may be innaccurate. And you proposed fix by adding a small value will not in general fix the problem.

If your input is an integer then you should consider performing your calculation using integer arithmetic. That ultimately is the right way to deal with this.

David Heffernan
  • 601,492
  • 42
  • 1,072
  • 1,490
  • Do you have a counterexample that doesn't involve subnormals, underflow, or overflow? – tmyklebu Jun 07 '13 at 01:17
  • @tmyklebu Counter example to what. I'm addressing whether to code detects perfect squares. It does not because it uses floating point arithmetic which has no concept of perfect squares. – David Heffernan Jun 07 '13 at 06:33
  • The first sentence of the question asks about recognition. The second asks about getting a wrong answer. So the question about a counter example is asking whether you can show a value that is a perfect square for which `sqrt` does not return the exact square root. – Eric Postpischil Jun 07 '13 at 06:39
  • @eric that's trivially easy to show. All you need is a number that is not representable. – David Heffernan Jun 07 '13 at 06:44
  • @tmyklebu Forget about subnormals. And overflow/underflow don't ever come into this. I mean, for which integer n does n**0.5 either overflow or underflow. All you need is a non-representable number. – David Heffernan Jun 07 '13 at 06:53
  • @DavidHeffernan: You should address your answer to the question to tmyklebu, not me. I just explained the question to you. I did not pose it. – Eric Postpischil Jun 07 '13 at 13:44
  • @EricPostpischil Yes, I did that too. – David Heffernan Jun 07 '13 at 13:49
  • `x*x` can overflow or underflow. I'm asking for an example of a number `x`, either integer or floating-point, such that `sqrt(x*x) != x`. – tmyklebu Jun 07 '13 at 22:50
  • @tmyklebu Then feel free to ask your own question on that topic. The issue here is that large integers are not representable. – David Heffernan Jun 08 '13 at 07:03
1

You can use the round(number, significant_figures) before converting to an int, I cannot recall if python truncs or rounds when doing a float-to-integer conversion.

In any case, since python uses floating point arithmetic, all the pitfalls apply. See:
http://docs.python.org/2/tutorial/floatingpoint.html

NeonMan
  • 623
  • 10
  • 24
1

Perfect-square values will have no fractional components, so your main worry would be very large values, and for such values a difference of 1 or 2 being significant means you're going to want a specific numerical library that supports such high precision (as DSM mentions, the Decimal library, standard since Python 2.4, should be able to do what you want as it supports arbitrary precision.

http://docs.python.org/library/decimal.html

JAB
  • 20,783
  • 6
  • 71
  • 80
  • +1. `Decimal('10000000000000000000000000000000000000000000000000000000000000000000000').sqrt() Decimal('1.000000000000000000000000000E+35')` – Mark Ransom Jun 06 '13 at 19:03
1

sqrt is one of the easier math library functions to implement, and any math library of reasonable quality will implement it with faithful rounding (sub-ULP accuracy). If the input is a perfect square, its square root is representable (in a reasonable floating-point format). In this case, faithful rounding guarantees the result is exact.

This addresses only the value actually passed to sqrt. Whether a number can be converted without error from another format to the floating-point input for sqrt is a separate issue.

Eric Postpischil
  • 195,579
  • 13
  • 168
  • 312
  • I think you can say something stronger, namely that, for reasonable exponents, `sqrt(x*x) == x` because `(x+ulp/2)*(x+ulp/2)`, which is where `sqrt` starts giving the wrong value, is always more than an ulp higher and `(x-ulp/2)*(x-ulp/2)` is more than an ulp smaller except when `x` is within a few ulps of a power of two. – tmyklebu Jun 07 '13 at 01:29
  • @tmyklebu: That seems like an argument using correct rounding (error is at most 1/2 ULP). Faithful rounding allows up to one ULP (exclusive). Anytime the result is not exactly representable, there are two values that may be returned. So you cannot guarantee `sqrt(x*x) == x` unless `x*x` is exact. – Eric Postpischil Jun 07 '13 at 06:34
  • Doesn't IEEE guarantee a correctly rounded `sqrt`? (I'm aware that Python probably has no reason to care about what IEEE says. But I've never seen a math library get `sqrt` wrong, so, for practical purposes, I feel OK with relying on correctly rounded `sqrt` and making recommendations based on correctly rounded `sqrt`. I think my argument might also work for faithfully-rounded `sqrt`, but there are more edge cases, notably around `sqrt(2)`, that would require special attention.) – tmyklebu Jun 07 '13 at 22:51
  • @tmyklebu: No. IEEE 754-2008 recommends a correctly round `sqrt`. It does not mandate it. And, as you note, Python does not guarantee IEEE 754. I can state definitely that your argument does not work for faithful rounding. Given any `x`, if `sqrt(x*x)` returns `x`, I can replace it with a `sqrt` that is identical except that `sqrt(x*x)` returns the other value (different from x) that faithful rounding allows. – Eric Postpischil Jun 07 '13 at 23:01
0

This question is very old and a lot has changed over the years. Now there is a math.isqrt function in python that does almost the same as int(math.sqrt(...)), but does not loose in precision:

>>> int(math.sqrt(100000000000000000000000000000000000**2))
99999999999999996863366107917975552
>>> math.isqrt(100000000000000000000000000000000000**2)
100000000000000000000000000000000000
Uladzislau_Khamkou
  • 150
  • 7