在python整数平方根

在python或标准库中是否有一个整数平方根？ 我希望它是确切的（即返回一个整数），并且如果没有解决scheme就叫。

此刻，我推出了自己的天真的一个：

def isqrt(n): i = int(math.sqrt(n) + 0.5) if i**2 == n: return i raise ValueError('input was not a perfect square') 

但是这很丑，我不太相信大整数。 如果我超过了这个值，我可以遍历这些方块并放弃，但是我认为这样做会有点慢。 另外我想我可能会重新发明轮子，像这样的东西一定存在于Python已经…

Answer 1

牛顿的方法在整数上工作得很好：

 def isqrt(n): x = n y = (x + 1) // 2 while y < x: x = y y = (x + n // x) // 2 return x 

这将返回x * x不超过n的最大整数x 。 如果你想检查结果是否只是平方根，只需执行乘法来检查n是否是一个完美的正方形。

我在我的博客上讨论这个algorithm，以及另外三个计算平方根的algorithm。

Answer 2

对不起，很晚的回应; 我偶然发现了这个页面。 如果有人在将来访问这个页面，python模块gmpy2被devise用来处理非常大的input，并且包括一个整数平方根函数。

例：

 >>> import gmpy2 >>> gmpy2.isqrt((10**100+1)**2) mpz(10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001L) >>> gmpy2.isqrt((10**100+1)**2 - 1) mpz(10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000L) 

当然，一切都会有“mpz”标签，但是mpz与int的兼容：

 >>> gmpy2.mpz(3)*4 mpz(12) >>> int(gmpy2.mpz(12)) 12 

请参阅我的其他答案 ，讨论此方法相对于此问题的其他答案的性能。

下载： https ： //code.google.com/p/gmpy/

Answer 3

长手平方根algorithm

事实certificate，有一个计算平方根的algorithm，你可以用手来计算，就像长分治疗一样。 该algorithm的每次迭代都会生成所产生的平方根的一个数字，同时会消耗您寻找的平方根的两位数字。 虽然该algorithm的“长手”版本是用十进制指定的，但它可以在任何基础上工作，二进制最简单，也许执行得最快（取决于底层的bignum表示）。

由于此algorithm对数字逐位进行操作，因此它可以为任意大的完美平方生成精确的结果，对于非完美平方，可以根据需要生成精确的小数位数（精确到小数点右侧）。

在“Dr. Math”网站上有两个很好的文字说明algorithm：

• 二进制中的平方根
• Longhand广场根

以下是Python中的一个实现：

 def exact_sqrt(x): """Calculate the square root of an arbitrarily large integer. The result of exact_sqrt(x) is a tuple (a, r) such that a**2 + r = x, where a is the largest integer such that a**2 <= x, and r is the "remainder". If x is a perfect square, then r will be zero. The algorithm used is the "long-hand square root" algorithm, as described at http://mathforum.org/library/drmath/view/52656.html Tobin Fricke 2014-04-23 Max Planck Institute for Gravitational Physics Hannover, Germany """ N = 0 # Problem so far a = 0 # Solution so far # We'll process the number two bits at a time, starting at the MSB L = x.bit_length() L += (L % 2) # Round up to the next even number for i in xrange(L, -1, -1): # Get the next group of two bits n = (x >> (2*i)) & 0b11 # Check whether we can reduce the remainder if ((N - a*a) << 2) + n >= (a<<2) + 1: b = 1 else: b = 0 a = (a << 1) | b # Concatenate the next bit of the solution N = (N << 2) | n # Concatenate the next bit of the problem return (a, Na*a) 

你可以很容易地修改这个函数来进行额外的迭代来计算平方根的小数部分。 我最感兴趣的是计算大型完美广场的根源。

我不确定这与“整数牛顿法”algorithm相比如何。 我怀疑牛顿的方法是更快的，因为它原则上可以在一次迭代中产生多位解，而“长手”algorithm每次迭代只产生一位解。

源回购： https ： //gist.github.com/tobin/11233492

Answer 4

一种select是使用decimal模块，并用足够精确的浮点数来完成：

 import decimal def isqrt(n): nd = decimal.Decimal(n) with decimal.localcontext() as ctx: ctx.prec = n.bit_length() i = int(nd.sqrt()) if i**2 != n: raise ValueError('input was not a perfect square') return i 

我认为应该工作：

 >>> isqrt(1) 1 >>> isqrt(7**14) == 7**7 True >>> isqrt(11**1000) == 11**500 True >>> isqrt(11**1000+1) Traceback (most recent call last): File "<ipython-input-121-e80953fb4d8e>", line 1, in <module> isqrt(11**1000+1) File "<ipython-input-100-dd91f704e2bd>", line 10, in isqrt raise ValueError('input was not a perfect square') ValueError: input was not a perfect square 

Answer 5

看起来像你可以这样检查：

 if int(math.sqrt(n))**2 == n: print n, 'is a perfect square' 

更新：

正如你所指出的，上面的n值很大。 对于那些看起来很有希望的，这是由Martin Guy @ UKC于1985年6月改编的示例C代码，用于维基百科文章中提到的相对简单的二进制数字逐位计算方法。计算平方根 ：

 from math import ceil, log def isqrt(n): res = 0 bit = 4**int(ceil(log(n, 4))) if n else 0 # smallest power of 4 >= the argument while bit: if n >= res + bit: n -= res + bit res = (res >> 1) + bit else: res >>= 1 bit >>= 2 return res if __name__ == '__main__': from math import sqrt # for comparison purposes for i in range(17)+[2**53, (10**100+1)**2]: is_perfect_sq = isqrt(i)**2 == i print '{:21,d}: math.sqrt={:12,.7G}, isqrt={:10,d} {}'.format( i, sqrt(i), isqrt(i), '(perfect square)' if is_perfect_sq else '') 

输出：

  0: math.sqrt= 0, isqrt= 0 (perfect square) 1: math.sqrt= 1, isqrt= 1 (perfect square) 2: math.sqrt= 1.414214, isqrt= 1 3: math.sqrt= 1.732051, isqrt= 1 4: math.sqrt= 2, isqrt= 2 (perfect square) 5: math.sqrt= 2.236068, isqrt= 2 6: math.sqrt= 2.44949, isqrt= 2 7: math.sqrt= 2.645751, isqrt= 2 8: math.sqrt= 2.828427, isqrt= 2 9: math.sqrt= 3, isqrt= 3 (perfect square) 10: math.sqrt= 3.162278, isqrt= 3 11: math.sqrt= 3.316625, isqrt= 3 12: math.sqrt= 3.464102, isqrt= 3 13: math.sqrt= 3.605551, isqrt= 3 14: math.sqrt= 3.741657, isqrt= 3 15: math.sqrt= 3.872983, isqrt= 3 16: math.sqrt= 4, isqrt= 4 (perfect square) 9,007,199,254,740,992: math.sqrt=9.490627E+07, isqrt=94,906,265 100,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,020,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,001: math.sqrt= 1E+100, isqrt=10,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,000,001 (perfect square) 

Answer 6

这是一个非常简单的实现：

 def i_sqrt(n): i = n.bit_length() >> 1 # i = floor( (1 + floor(log_2(n))) / 2 ) m = 1 << i # m = 2^i # # Fact: (2^(i + 1))^2 > n, so m has at least as many bits # as the floor of the square root of n. # # Proof: (2^(i+1))^2 = 2^(2i + 2) >= 2^(floor(log_2(n)) + 2) # >= 2^(ceil(log_2(n) + 1) >= 2^(log_2(n) + 1) > 2^(log_2(n)) = n. QED. # while m*m > n: m >>= 1 i -= 1 for k in xrange(i-1, -1, -1): x = m | (1 << k) if x*x <= n: m = x return m 

这只是一个二进制search。 将值m初始化为不超过平方根的2的最大次幂，然后检查是否可以设置每个较小的位，同时保持不大于平方根的结果。 （一次检查一个位，按降序排列。）

对于相当大的n值（比如大约10**6000或大约20000比特），这似乎是：

• 比user448810描述的牛顿方法实现更快 。
• 很多，比我的其他答案 gmpy2内置方法慢得多。
• 与nibot 描述的Longhand Square Root相比，但稍慢一些。

所有这些方法在这个尺寸的input上都成功了，但是在我的机器上，这个函数需要大约1.5秒，而@Nibot需要大约0.9秒，@ user448810需要大约19秒，而gmpy2内置方法需要不到1毫秒（！）。 例：

 >>> import random >>> import timeit >>> import gmpy2 >>> r = random.getrandbits >>> t = timeit.timeit >>> t('i_sqrt(r(20000))', 'from __main__ import *', number = 5)/5. # This function 1.5102493192883117 >>> t('exact_sqrt(r(20000))', 'from __main__ import *', number = 5)/5. # Nibot 0.8952787937686366 >>> t('isqrt(r(20000))', 'from __main__ import *', number = 5)/5. # user448810 19.326695976676184 >>> t('gmpy2.isqrt(r(20000))', 'from __main__ import *', number = 5)/5. # gmpy2 0.0003599147067689046 >>> all(i_sqrt(n)==isqrt(n)==exact_sqrt(n)[0]==int(gmpy2.isqrt(n)) for n in (r(1500) for i in xrange(1500))) True 

这个函数可以很容易地推广，虽然它不是很好，因为我没有m的初始猜测那么精确：

 def i_root(num, root, report_exactness = True): i = num.bit_length() / root m = 1 << i while m ** root < num: m <<= 1 i += 1 while m ** root > num: m >>= 1 i -= 1 for k in xrange(i-1, -1, -1): x = m | (1 << k) if x ** root <= num: m = x if report_exactness: return m, m ** root == num return m 

但是请注意， gmpy2也有一个i_root方法。

事实上，这种方法可以适用于任何（非负，递增）函数f来确定“ f整数倒数”。 但是，要selectm的有效初始值，您仍然想知道关于f一些信息。

编辑：感谢@Greggo指出，可以重写i_sqrt函数，以避免使用任何乘法。 这产生令人印象深刻的性能提升

 def improved_i_sqrt(n): assert n >= 0 if n == 0: return 0 i = n.bit_length() >> 1 # i = floor( (1 + floor(log_2(n))) / 2 ) m = 1 << i # m = 2^i # # Fact: (2^(i + 1))^2 > n, so m has at least as many bits # as the floor of the square root of n. # # Proof: (2^(i+1))^2 = 2^(2i + 2) >= 2^(floor(log_2(n)) + 2) # >= 2^(ceil(log_2(n) + 1) >= 2^(log_2(n) + 1) > 2^(log_2(n)) = n. QED. # while (m << i) > n: # (m<<i) = m*(2^i) = m*m m >>= 1 i -= 1 d = n - (m << i) # d = nm^2 for k in xrange(i-1, -1, -1): j = 1 << k new_diff = d - (((m<<1) | j) << k) # n-(m+2^k)^2 = nm^2-2*m*2^k-2^(2k) if new_diff >= 0: d = new_diff m |= j return m 

注意，通过构造， m << 1第k位没有被设置，所以按位或者可以被用来实现(m<<1) + (1<<k)的加法。 (((m<<1) | (1<<k)) << k)写成(2*m*(2**k) + 2**(2*k))移位和一个按位 – 或（后面减去获得new_diff ）。 也许还有更有效的方法来得到这个？ 无论如何，它比乘以m*m要好得多！ 与以上相比：

 >>> t('improved_i_sqrt(r(20000))', 'from __main__ import *', number = 5)/5. 0.10908999762373242 >>> all(improved_i_sqrt(n) == i_sqrt(n) for n in xrange(10**6)) True 

Answer 7

您的function在大input时失败：

 In [26]: isqrt((10**100+1)**2) ValueError: input was not a perfect square 

ActiveState网站上有一个配方，希望它更可靠，因为它只使用整数math。 它基于较早的StackOverflow问题： 编写自己的平方根函数

Answer 8

我发现这个线程前几天，重新写了nibot的解决scheme ，并通过减less一半的迭代次数，做一些其他小的性能调整，我能够提高性能的2.4倍：

 def isqrt(n): a = 0 # a is the current answer. r = 0 # r is the current remainder. for s in reversed(range(0, n.bit_length(), 2)): # Shift n by s bits. t = n >> s & 3 # t is the two next most significant bits of n. r = r << 2 | t # Increase the remainder as if no new bit is set. c = a << 2 | 1 # c is an intermediate value used for comparison. b = r >= c # b is the next bit in the remainder. if b: r -= c # b has been set, so reduce the remainder. a = a << 1 | b # Update the answer to include b. return (a, r) 

以下是timeit的结果：

 >>> timeit('isqrt(12345678901234567890)', setup='from __main__ import isqrt') 8.862877120962366 

然后，为了比较，我实现了最常用的整数平方根algorithm： 牛顿法 。 这个定义更加紧凑。

 def isqrt(n): x, y = n, n >> 1 while x > y: x, y = y, (y + n//y) >> 1 return (x, n - x*x) 

事实certificate，即使是长方根的优化版本也比牛顿的方法慢了约1.5倍。

 >>> timeit('isqrt(12345678901234567890)', setup='from __main__ import isqrt') 5.74083631898975 

所以，总而言之，如果您需要一个快速的纯Python整数平方根函数，那就不要再比上面提供的那个更进一步了。

编辑：我已经纠正了上述牛顿方法中的错误。 在我的机器上，它比user448810的解决scheme运行速度快了10％。

Answer 9

在计算机上不能精确地表示浮游物。 您可以testing所需的接近度设置epsilon到python的浮点精度内的一个小值。

 def isqrt(n): epsilon = .00000000001 i = int(n**.5 + 0.5) if abs(i**2 - n) < epsilon: return i raise ValueError('input was not a perfect square') 

Answer 10

我将这里给出的不同方法与一个循环进行了比较：

 for i in range (1000000): # 700 msec r=int(123456781234567**0.5+0.5) if r**2==123456781234567:rr=r else:rr=-1 

发现这个是最快的，不需要math导入。 很长时间可能会失败，但看看这个

 15241576832799734552675677489**0.5 = 123456781234567.0 

Answer 11

试试这个条件（不需要额外的计算）：

 def isqrt(n): i = math.sqrt(n) if i != int(i): raise ValueError('input was not a perfect square') return i 

如果你需要它返回一个int （而不是一个尾随零的float ），然后分配一个第二个variables或计算int(i)两次。