通知
有关Erlang
或C / C++
中的解决方案,请转到下面的试用版4。
维基百科文章
整数平方根
- "整数平方根"的定义可以在这里找到
计算平方根的方法
- 可以在这里找到一种实现"比特魔术">的算法
[试用1:使用库函数]
代码
isqrt(N) when erlang:is_integer(N), N >= 0 ->
erlang:trunc(math:sqrt(N)).
问题
这个实现使用了C库中的sqrt()
函数,因此它不适用于任意大的整数(注意,返回的结果与输入不匹配。正确的答案应该是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+
]
代码
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).
说明
此实施基于以下观察:
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
时,它似乎在我(非常快)的机器上永远运行。
[试用3:仅使用Bigint+1
、-1
和div 2
的二进制搜索]
代码
变体1(我最初的实现)
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(参考Vikram Bhat的答案,修改了上述代码,使边界改为Mid+1或Mid-1)
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 * 1111111111111111111111111111111111111111
),结果立即显示出来。然而,在我的机器上求解一百万(1000000)个61位数字(即从1000000000000000000000000000000000000000000000000000000000000
到1000000000000000000000000000000000000000000000000000001000000
)需要60秒(+-2秒)。我想做得更快。
[试验4:仅使用Bigint+
和div
的牛顿法]
代码
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++中的代码(出于您的兴趣)
递归变体
#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);
}
迭代变体
#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秒)解决了一百万(1000000)个61位数字,所以这比试用版3更快。它能更快吗?
关于我的机器
处理器:3.4 GHz Intel Core i7
内存:32 GB 1600 MHz DDR3
操作系统:Mac OS X版本10.9.1
相关问题
python 中的整数平方根
用户448810的答案使用"牛顿法">。我不确定用"整数除法"做除法是否可以。我稍后会尝试此更新。【更新(2015-01-11):这样做是可以的】
数学的答案是使用第三方Python包
gmpy
,这对我来说不是很好,因为我主要感兴趣的是在Erlang中只使用内置的工具来解决它。DSM的答案似乎很有趣。我真的不明白发生了什么,但似乎有"比特魔法">参与其中,所以它也不太适合我。
元整数平方根中的无限递归
- 这个问题是针对C++的,AraK(提问者)的算法看起来与上面的试验2的想法相同
下面这样的二进制搜索不需要浮点除法,只需要整数乘法(比牛顿方法慢):-
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)
迭代,因此即使是非常大的也不应该永远运行
Log10(目标)计算(如果需要):-
acc = target
log10 = 0;
while(acc>0) {
log10 = log10 + 1;
acc = acc/10;
}
注意:acc/10
是整数除法
编辑:-
有效界:-sqrt(n)的位数大约是n的一半,因此您可以通过high = 10^(log10(N)/2+1)
&;CCD_ 19得到更严格的约束,并且它应该提供2倍的速度。
评估界限:-
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);