一种计算任意大整数的整数平方根（isqrt）的有效算法 [英] An efficient algorithm to calculate the integer square root (isqrt) of arbitrarily large integers
问题描述
对于 Erlang 或 C / C ++ ，请转到下面的试用版4 。
For a solution in Erlang or C / C++, go to Trial 4 below.
- 可以在这里找到整数平方根的定义
- 可以在这里找到bit magic的算法
- An algorithm that does "bit magic" could be found here
isqrt(N) when erlang:is_integer(N), N >= 0 ->
erlang:trunc(math:sqrt(N)).
问题
此实现使用 sqrt（）函数，因此它不适用于任意大的整数（请注意，返回的结果与输入不符，正确答案应为 12345678901234567890 ）：
Problem
This implementation uses the sqrt() function from the C library, so it does not work with arbitrarily large integers (Note that the returned result does not match the input. The correct answer should be 12345678901234567890):
Erlang R16B03 (erts-5.10.4) [source] [64-bit] [smp:8:8] [async-threads:10] [hipe] [kernel-poll:false]
Eshell V5.10.4 (abort with ^G)
1> erlang:trunc(math:sqrt(12345678901234567890 * 12345678901234567890)).
12345678901234567168
2>
[试用2：使用Bigint <$ c $
[ Trial 2 : Using Bigint + Only ]
Code
isqrt2(N) when erlang:is_integer(N), N >= 0 ->
isqrt2(N, 0, 3, 0).
isqrt2(N, I, _, Result) when I >= N ->
Result;
isqrt2(N, I, Times, Result) ->
isqrt2(N, I + Times, Times + 2, Result + 1).
描述
这个实现是基于以下观察结果：
Description
This implementation is based on the following observation:
isqrt(0) = 0 # <--- One 0
isqrt(1) = 1 # <-+
isqrt(2) = 1 # |- Three 1's
isqrt(3) = 1 # <-+
isqrt(4) = 2 # <-+
isqrt(5) = 2 # |
isqrt(6) = 2 # |- Five 2's
isqrt(7) = 2 # |
isqrt(8) = 2 # <-+
isqrt(9) = 3 # <-+
isqrt(10) = 3 # |
isqrt(11) = 3 # |
isqrt(12) = 3 # |- Seven 3's
isqrt(13) = 3 # |
isqrt(14) = 3 # |
isqrt(15) = 3 # <-+
isqrt(16) = 4 # <--- Nine 4's
...
问题
此实现涉及只有 bigint添加，所以我预计它运行得很快。但是，当我使用 1111111111111111111111111111111111111111 * 1111111111111111111111111111111111111111 ，它似乎在我的（非常快）的机器上永远运行。
Problem
This implementation involves only bigint additions so I expected it to run fast. However, when I fed it with 1111111111111111111111111111111111111111 * 1111111111111111111111111111111111111111, it seems to run forever on my (very fast) machine.
isqrt3(N) when erlang:is_integer(N), N >= 0 ->
isqrt3(N, 1, N).
isqrt3(_N, Low, High) when High =:= Low + 1 ->
Low;
isqrt3(N, Low, High) ->
Mid = (Low + High) div 2,
MidSqr = Mid * Mid,
if
%% This also catches N = 0 or 1
MidSqr =:= N ->
Mid;
MidSqr < N ->
isqrt3(N, Mid, High);
MidSqr > N ->
isqrt3(N, Low, Mid)
end.
变式2（修改上面的代码，以便边界与Mid + 1或Mid-1相反，参考 Vikram Bhat的回答）
Variant 2 (modified above code so that the boundaries go with Mid+1 or Mid-1 instead, with reference to the answer by Vikram Bhat)
isqrt3a(N) when erlang:is_integer(N), N >= 0 ->
isqrt3a(N, 1, N).
isqrt3a(N, Low, High) when Low >= High ->
HighSqr = High * High,
if
HighSqr > N ->
High - 1;
HighSqr =< N ->
High
end;
isqrt3a(N, Low, High) ->
Mid = (Low + High) div 2,
MidSqr = Mid * Mid,
if
%% This also catches N = 0 or 1
MidSqr =:= N ->
Mid;
MidSqr < N ->
isqrt3a(N, Mid + 1, High);
MidSqr > N ->
isqrt3a(N, Low, Mid - 1)
end.
问题
现在解决了79数字（即 1111111111111111111111111111111111111111 11111111111111111111111111111111111111 ）以减轻速度，结果立即显示。但是，我的机器需要60秒（+ - 2秒）来解决一百万（100万）61位数字（即从$ code> 1000000000000000000000000000000000000000000000000000000000000 到 1000000000000000000000000000000000000000000000000000001000000 ）。我想做的更快。
Problem
Now it solves the 79-digit number (namely 1111111111111111111111111111111111111111 * 1111111111111111111111111111111111111111) in lightening speed, the result is shown immediately. However, it takes 60 seconds (+- 2 seconds) on my machine to solve one million (1,000,000) 61-digit numbers (namely, from 1000000000000000000000000000000000000000000000000000000000000 to 1000000000000000000000000000000000000000000000000000001000000). I would like to do it even faster.
isqrt4(0) -> 0;
isqrt4(N) when erlang:is_integer(N), N >= 0 ->
isqrt4(N, N).
isqrt4(N, Xk) ->
Xk1 = (Xk + N div Xk) div 2,
if
Xk1 >= Xk ->
Xk;
Xk1 < Xk ->
isqrt4(N, Xk1)
end.
C / C ++中的代码（为了您的兴趣）
递归变体
Code in C / C++ (for your interest)
Recursive variant
#include <stdint.h>
uint32_t isqrt_impl(
uint64_t const n,
uint64_t const xk)
{
uint64_t const xk1 = (xk + n / xk) / 2;
return (xk1 >= xk) ? xk : isqrt_impl(n, xk1);
}
uint32_t isqrt(uint64_t const n)
{
if (n == 0) return 0;
if (n == 18446744073709551615ULL) return 4294967295U;
return isqrt_impl(n, n);
}
迭代变体
Iterative variant
#include <stdint.h>
uint32_t isqrt_iterative(uint64_t const n)
{
uint64_t xk = n;
if (n == 0) return 0;
if (n == 18446744073709551615ULL) return 4294967295U;
do
{
uint64_t const xk1 = (xk + n / xk) / 2;
if (xk1 >= xk)
{
return xk;
}
else
{
xk = xk1;
}
} while (1);
}
问题
Erlang代码在我的机器上在40秒（+ - 1秒）内解决了一百万（100万）61位数字，所以这比试用3 更快。可以更快吗？
Problem
The Erlang code solves one million (1,000,000) 61-digit numbers in 40 seconds (+- 1 second) on my machine, so this is faster than Trial 3. Can it go even faster?
处理器： 3.4 GHz Intel Core i7
Processor : 3.4 GHz Intel Core i7
内存： 32 GB 1600 MHz DDR3
Memory : 32 GB 1600 MHz DDR3
操作系统： Mac OS X版本10.9.1
OS : Mac OS X Version 10.9.1
-
user448810的回答使用牛顿法。我不知道是否使用整数除法进行分割是否可行。我会稍后再来这个更新。 [UPDATE（2015-01-11）：可以这样做）
The answer by user448810 uses "Newton's Method". I'm not sure whether doing the division using "integer division" is okay or not. I'll try this later as an update. [UPDATE (2015-01-11): It is okay to do so]
数学回答涉及使用第三方Python包 gmpy ，这对我来说不是很有利，因为我主要感兴趣
The answer by math involves using a 3rd party Python package gmpy, which is not very favourable to me, since I'm primarily interested in solving it in Erlang with only builtin facilities.
回答者： DSM 似乎很有趣。我真的不明白发生了什么，但是似乎这个位魔法涉及到，所以也不适合我。
The answer by DSM seems interesting. I don't really understand what is going on, but it seems that "bit magic" is involved there, and so it's not quite suitable for me too.
- 这个问题是针对C ++的，而AraK（提问者）的算法看起来和以前一样，
推荐答案
如下二进制搜索如下，不需要浮动划分整数乘法（比牛顿法慢）： -
How about binary search like following doesn't need floating divisions only integer multiplications (Slower than newtons method) :-
low = 1;
/* More efficient bound
high = pow(10,log10(target)/2+1);
*/
high = target
while(low<high) {
mid = (low+high)/2;
currsq = mid*mid;
if(currsq==target) {
return(mid);
}
if(currsq<target) {
if((mid+1)*(mid+1)>target) {
return(mid);
}
low = mid+1;
}
else {
high = mid-1;
}
}
这适用于 O（logN）迭代，所以不应该永远运行甚至非常大的数字
This works for O(logN) iterations so should not run forever for even very large numbers
Log10（target）如果需要， -
acc = target
log10 = 0;
while(acc>0) {
log10 = log10 + 1;
acc = acc/10;
}
注意： acc / 10 是整数除法
编辑： -
有效绑定： - sqrt（n）的大小是n的一半，所以你可以通过 high = 10 ^（log10（N）/ 2 + 1 ）&& low = 10 ^（log10（N）/ 2-1）得到更紧的绑定，应该提供2倍的速度。
Efficient bound :- The sqrt(n) has about half the number of digits as n so you can pass high = 10^(log10(N)/2+1) && low = 10^(log10(N)/2-1) to get tighter bound and it should provide 2 times speed up.
评估界限： -
bound = 1;
acc = N;
count = 0;
while(acc>0) {
acc = acc/10;
if(count%2==0) {
bound = bound*10;
}
count++;
}
high = bound*10;
low = bound/10;
isqrt(N,low,high);
这篇关于一种计算任意大整数的整数平方根（isqrt）的有效算法的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
