编写自己的平方根函数

如何编写自己的函数来求整数的最准确的平方根?

在谷歌了它之后,我找到了 这个(从它的 原始链接归档) ,但是首先,我没有完全得到它,其次,它也是近似的。

假设平方根为最接近的整数(到实际的根)或浮点数。

180965 次浏览

当然是近似的,这就是浮点数的数学原理。

无论如何,标准的方法是使用 牛顿的方法。这和使用泰勒级数差不多,另一种方式立刻出现在我的脑海中。

一般来说,整数的平方根(例如2)可以近似为 只有(不是因为浮点数算法的问题,而是因为它们是无理数,无法精确计算)。

当然,有些近似值要比其他近似值好。我的意思是,当然,值1.732是一个更好的近似3的算术平方根,比1.7

你给出的那个链接的代码所使用的方法是采用第一个近似值,并用它来计算 好多了近似值。

这就是牛顿方法,你可以用每个新的近似值重复计算,它对你来说已经足够精确了。

事实上,必须的有一些方法来决定什么时候停止重复,否则它将永远运行。

通常,当近似值之间的差值为 少于时,您将停止。

编辑: 我认为没有比您已经找到的两个更简单的实现了。

让我指出一个非常有趣的方法来计算1/sqrt (x)的反平方根,这在游戏设计界是一个传奇,因为它快得令人难以置信。或者等等,看看下面的帖子:

Http://betterexplained.com/articles/understanding-quakes-fast-inverse-square-root/

PS: 我知道你只是想要平方根,但是优雅的地震克服了我的所有阻力:)

顺便说一下,上面提到的文章也讨论了无聊的牛顿-拉斐逊近似。

计算 X 的平方根的一个简单(但不是很快)的方法:

squareroot(x)
if x<0 then Error
a = 1
b = x
while (abs(a-b)>ErrorMargin)
a = (a+b)/2
b = x/a
endwhile
return a;

示例: squeroot (70000)

    a       b
1   70000
35001       2
17502       4
8753       8
4381      16
2199      32
1116      63
590     119
355     197
276     254
265     264

正如您所看到的,它为平方根定义了一个上边界和一个下边界,并缩小了边界,直到它的大小可以接受为止。

有更有效的方法,但这一个说明了过程,并很容易理解。

如果使用其他的整数,那么要注意把 Erroredge 设置为1,否则会产生一个无穷无尽的循环。

根据您的需要,可以使用简单的分而治之策略。它不会像其他方法那样收敛为 很快,但是对于新手来说,它可能更容易理解。另外,由于它是一个 O (log n)算法(每次迭代搜索空间减半) ,因此32位浮点数的最坏情况将是32次迭代。

假设你想要62.104的平方根。选择一个介于0和0之间的值,然后平方它。如果平方高于你的数字,你需要集中在数字少于中点。如果太低,就把注意力集中在高的地方。

使用真正的数学,您可以永远将搜索空间划分为两个(如果它没有一个有理平方根)。在现实中,计算机最终会用尽精度,你会有你的近似值。下面的 C 程序说明了这一点:

#include <stdio.h>
#include <stdlib.h>


int main (int argc, char *argv[]) {
float val, low, high, mid, oldmid, midsqr;
int step = 0;


// Get argument, force to non-negative.


if (argc < 2) {
printf ("Usage: sqrt <number>\n");
return 1;
}
val = fabs (atof (argv[1]));


// Set initial bounds and print heading.


low = 0;
high = mid = val;
oldmid = -1;


printf ("%4s  %10s  %10s  %10s  %10s  %10s    %s\n",
"Step", "Number", "Low", "High", "Mid", "Square", "Result");


// Keep going until accurate enough.


while (fabs(oldmid - mid) >= 0.00001) {
oldmid = mid;


// Get midpoint and see if we need lower or higher.


mid = (high + low) / 2;
midsqr = mid * mid;
printf ("%4d  %10.4f  %10.4f  %10.4f  %10.4f  %10.4f  ",
++step, val, low, high, mid, midsqr);
if (mid * mid > val) {
high = mid;
printf ("- too high\n");
} else {
low = mid;
printf ("- too low\n");
}
}


// Desired accuracy reached, print it.


printf ("sqrt(%.4f) = %.4f\n", val, mid);
return 0;
}

这里有几个步骤,希望你能了解它是如何工作的。77:

pax> sqrt 77
Step      Number         Low        High         Mid      Square    Result
1     77.0000      0.0000     77.0000     38.5000   1482.2500  - too high
2     77.0000      0.0000     38.5000     19.2500    370.5625  - too high
3     77.0000      0.0000     19.2500      9.6250     92.6406  - too high
4     77.0000      0.0000      9.6250      4.8125     23.1602  - too low
5     77.0000      4.8125      9.6250      7.2188     52.1104  - too low
6     77.0000      7.2188      9.6250      8.4219     70.9280  - too low
7     77.0000      8.4219      9.6250      9.0234     81.4224  - too high
8     77.0000      8.4219      9.0234      8.7227     76.0847  - too low
9     77.0000      8.7227      9.0234      8.8730     78.7310  - too high
10     77.0000      8.7227      8.8730      8.7979     77.4022  - too high
11     77.0000      8.7227      8.7979      8.7603     76.7421  - too low
12     77.0000      8.7603      8.7979      8.7791     77.0718  - too high
13     77.0000      8.7603      8.7791      8.7697     76.9068  - too low
14     77.0000      8.7697      8.7791      8.7744     76.9893  - too low
15     77.0000      8.7744      8.7791      8.7767     77.0305  - too high
16     77.0000      8.7744      8.7767      8.7755     77.0099  - too high
17     77.0000      8.7744      8.7755      8.7749     76.9996  - too low
18     77.0000      8.7749      8.7755      8.7752     77.0047  - too high
19     77.0000      8.7749      8.7752      8.7751     77.0022  - too high
20     77.0000      8.7749      8.7751      8.7750     77.0009  - too high
21     77.0000      8.7749      8.7750      8.7750     77.0002  - too high
22     77.0000      8.7749      8.7750      8.7750     76.9999  - too low
23     77.0000      8.7750      8.7750      8.7750     77.0000  - too low
sqrt(77.0000) = 8.7750

62.104美元:

pax> sqrt 62.104
Step      Number         Low        High         Mid      Square    Result
1     62.1040      0.0000     62.1040     31.0520    964.2267  - too high
2     62.1040      0.0000     31.0520     15.5260    241.0567  - too high
3     62.1040      0.0000     15.5260      7.7630     60.2642  - too low
4     62.1040      7.7630     15.5260     11.6445    135.5944  - too high
5     62.1040      7.7630     11.6445      9.7037     94.1628  - too high
6     62.1040      7.7630      9.7037      8.7334     76.2718  - too high
7     62.1040      7.7630      8.7334      8.2482     68.0326  - too high
8     62.1040      7.7630      8.2482      8.0056     64.0895  - too high
9     62.1040      7.7630      8.0056      7.8843     62.1621  - too high
10     62.1040      7.7630      7.8843      7.8236     61.2095  - too low
11     62.1040      7.8236      7.8843      7.8540     61.6849  - too low
12     62.1040      7.8540      7.8843      7.8691     61.9233  - too low
13     62.1040      7.8691      7.8843      7.8767     62.0426  - too low
14     62.1040      7.8767      7.8843      7.8805     62.1024  - too low
15     62.1040      7.8805      7.8843      7.8824     62.1323  - too high
16     62.1040      7.8805      7.8824      7.8815     62.1173  - too high
17     62.1040      7.8805      7.8815      7.8810     62.1098  - too high
18     62.1040      7.8805      7.8810      7.8807     62.1061  - too high
19     62.1040      7.8805      7.8807      7.8806     62.1042  - too high
20     62.1040      7.8805      7.8806      7.8806     62.1033  - too low
21     62.1040      7.8806      7.8806      7.8806     62.1038  - too low
22     62.1040      7.8806      7.8806      7.8806     62.1040  - too high
23     62.1040      7.8806      7.8806      7.8806     62.1039  - too high
sqrt(62.1040) = 7.8806

49岁:

pax> sqrt 49
Step      Number         Low        High         Mid      Square    Result
1     49.0000      0.0000     49.0000     24.5000    600.2500  - too high
2     49.0000      0.0000     24.5000     12.2500    150.0625  - too high
3     49.0000      0.0000     12.2500      6.1250     37.5156  - too low
4     49.0000      6.1250     12.2500      9.1875     84.4102  - too high
5     49.0000      6.1250      9.1875      7.6562     58.6182  - too high
6     49.0000      6.1250      7.6562      6.8906     47.4807  - too low
7     49.0000      6.8906      7.6562      7.2734     52.9029  - too high
8     49.0000      6.8906      7.2734      7.0820     50.1552  - too high
9     49.0000      6.8906      7.0820      6.9863     48.8088  - too low
10     49.0000      6.9863      7.0820      7.0342     49.4797  - too high
11     49.0000      6.9863      7.0342      7.0103     49.1437  - too high
12     49.0000      6.9863      7.0103      6.9983     48.9761  - too low
13     49.0000      6.9983      7.0103      7.0043     49.0598  - too high
14     49.0000      6.9983      7.0043      7.0013     49.0179  - too high
15     49.0000      6.9983      7.0013      6.9998     48.9970  - too low
16     49.0000      6.9998      7.0013      7.0005     49.0075  - too high
17     49.0000      6.9998      7.0005      7.0002     49.0022  - too high
18     49.0000      6.9998      7.0002      7.0000     48.9996  - too low
19     49.0000      7.0000      7.0002      7.0001     49.0009  - too high
20     49.0000      7.0000      7.0001      7.0000     49.0003  - too high
21     49.0000      7.0000      7.0000      7.0000     49.0000  - too low
22     49.0000      7.0000      7.0000      7.0000     49.0001  - too high
23     49.0000      7.0000      7.0000      7.0000     49.0000  - too high
sqrt(49.0000) = 7.0000

我首先想到的是: 这是一个使用二进制搜索的好地方(灵感来自于这个伟大的 教程)

为了找到 vaule的平方根,我们在 (1..value)中搜索 number,其中预测器 我们选择的预测器是 number * number - value > 0.00001

double square_root_of(double value)
{
assert(value >= 1);
double lo = 1.0;
double hi = value;


while( hi - lo > 0.00001)
{
double mid = lo + (hi - lo) / 2 ;
std::cout << lo << "," << hi << "," << mid << std::endl;
if( mid * mid - value > 0.00001)    //this is the predictors we are using
{
hi = mid;
} else {
lo = mid;
}


}


return lo;
}

在 Python 中计算任意精度的平方根

#!/usr/bin/env python
import decimal


def sqrt(n):
assert n > 0
with decimal.localcontext() as ctx:
ctx.prec += 2 # increase precision to minimize round off error
x, prior = decimal.Decimal(n), None
while x != prior:
prior = x
x = (x + n/x) / 2 # quadratic convergence
return +x # round in a global context




decimal.getcontext().prec = 80 # desirable precision
r = sqrt(12345)
print r
print r == decimal.Decimal(12345).sqrt()

产出:

111.10805551354051124500443874307524148991137745969772997648567316178259031751676
True

下面计算 N > 0的楼层(sqrt (N)) :

x = 2^ceil(numbits(N)/2)
loop:
y = floor((x + floor(N/x))/2)
if y >= x
return x
x = y

这是一个版本的牛顿的方法给出了克兰德尔 & Pomerance,“素数: 一个计算的视角”。你应该使用这个版本的原因是,那些知道自己在做什么的人已经证明,它会精确地收敛到平方根的底部,而且它很简单,所以出现实现错误的可能性很小。它也很快(尽管构造一个更快的算法是可能的——但是要正确地做到这一点要复杂得多)。适当实现的二进制搜索对于非常小的 N 来说可能更快,但是在那里您也可以使用一个查找表。

要舍入到 最近的整数,只需使用上面的算法计算 t = floor (sqrt (4N))。如果设置了 t 的最小有效位,则选择 x = (t + 1)/2; 否则选择 t/2。注意,这是对平局进行四舍五入; 您也可以通过查看余数是否为非零(即 t ^ 2 = = 4N)进行四舍五入(或四舍五入到偶数)。

请注意,您不需要使用浮点算术。事实上,你不应该。这个算法应该完全使用整数来实现(特别是,floor ()函数只是表示应该使用常规的整数除法)。

找到一篇关于 整数平方根的文章。

这是一个略有改进的版本:

unsigned long sqrt(unsigned long a){
int i;
unsigned long rem = 0;
unsigned long root = 0;
for (i = 0; i < 16; i++){
root <<= 1;
rem = (rem << 2) | (a >> 30);
a <<= 2;
if(root < rem){
root++;
rem -= root;
root++;
}
}
return root >> 1;
}

反过来,正如名字所说,但有时候“足够接近”就是“足够接近”; 无论如何,这是一个有趣的读法。

Quake3快速 InvSqrt ()的起源

我在学校学过一个算法,你可以用它来计算精确的平方根(如果平方根是一个无理数,那么它的精度可以高到任意程度)。它肯定比牛顿的算法慢,但它是精确的。 假设你想计算531.3025的平方根

首先,你要把从小数点开始的数字分成两组:
{5}{31}。{30}{25}
然后:
1)找到第一组最小的平方根或者等于第一组的实际平方根: sqrt ({5}) > = 2。这个平方根是你最终答案的第一个数字。让我们把我们已经找到的最后一个平方根的数字表示为 B。现在 B = 2。
2)接下来计算{5}和 B ^ 2之间的差: 5-4 = 1。
3)其后所有2位数字组均须符合以下条件:
将余数乘以100,然后将其加到第二组: 100 + 31 = 131。
找到根的 X.next 数字,这样131 > = ((B * 20) + X) * X. X = 3。43 * 3 = 129 < 131.现在 B = 23。另外,由于小数点左侧没有更多的2位数组,因此已经找到了最终根的所有整数位。
4)对{30}和{25}重复同样的操作。所以你有:
{30} : 131-129 = 2.2 * 100 + 30 = 230 > = (23 * 2 * 10 + X) * X-> X = 0-> B = 23.0
{25} : 230-0 = 230.230 * 100 + 25 = 23025.23025 > = (230 * 2 * 10 + X) * X-> X = 5-> B = 23.05
最终结果 = 23.05。
这种算法看起来很复杂,但是如果你在纸上使用你在学校里学过的“长除法”的同样符号,那么它会简单得多,只不过你不用除法,而是计算平方根。

这是 Facebook 等社交网站常问的一个面试问题。我不认为在采访中使用牛顿方法是个好主意。如果面试官问你牛顿方法的原理,而你并不真正理解它,那该怎么办?

我用 Java 提供了一个基于二进制搜索的解决方案,我相信每个人都能理解。

public int sqrt(int x) {


if(x < 0) return -1;
if(x == 0 || x == 1) return x;


int lowerbound = 1;
int upperbound = x;
int root = lowerbound + (upperbound - lowerbound)/2;


while(root > x/root || root+1 <= x/(root+1)){
if(root > x/root){
upperbound = root;
} else {
lowerbound = root;
}
root = lowerbound + (upperbound - lowerbound)/2;
}
return root;
}

你可以在这里测试我的代码: Leetcode: sqrt (x)

下面是一种用三角法求平方根的方法。虽然不是最快的算法,但很精确。代码是用 javascript 编写的:

var x = 5;
var sqrt = Math.cos(Math.asin((x-1)/(x+1)))*(x+1)/2;
alert(sqrt);
// Fastest way I found, an (extreme) C# unrolled version of:
// http://www.hackersdelight.org/hdcodetxt/isqrt.c.txt         (isqrt4)


// It's quite a lot of code, basically a binary search (the "if" statements)
// followed by an unrolled loop (the labels).
// Most important: it's fast, twice as fast as "Math.Sqrt".
// On my pc: Math.Sqrt ~35 ns, sqrt <16 ns (mean <14 ns)


private static uint sqrt(uint x)
{
uint y, z;
if (x < 1u << 16)
{
if (x < 1u << 08)
{
if (x < 1u << 04) return x < 1u << 02 ? x + 3u >> 2 : x + 15u >> 3;
else
{
if (x < 1u << 06)
{ y = 1u << 03; x -= 1u << 04; if (x >= 5u << 02) { x -= 5u << 02; y |= 1u << 02; } goto L0; }
else
{ y = 1u << 05; x -= 1u << 06; if (x >= 5u << 04) { x -= 5u << 04; y |= 1u << 04; } goto L1; }
}
}
else                                             // slower (on my pc): .... y = 3u << 04; } goto L1; }
{
if (x < 1u << 12)
{
if (x < 1u << 10)
{ y = 1u << 07; x -= 1u << 08; if (x >= 5u << 06) { x -= 5u << 06; y |= 1u << 06; } goto L2; }
else
{ y = 1u << 09; x -= 1u << 10; if (x >= 5u << 08) { x -= 5u << 08; y |= 1u << 08; } goto L3; }
}
else
{
if (x < 1u << 14)
{ y = 1u << 11; x -= 1u << 12; if (x >= 5u << 10) { x -= 5u << 10; y |= 1u << 10; } goto L4; }
else
{ y = 1u << 13; x -= 1u << 14; if (x >= 5u << 12) { x -= 5u << 12; y |= 1u << 12; } goto L5; }
}
}
}
else
{
if (x < 1u << 24)
{
if (x < 1u << 20)
{
if (x < 1u << 18)
{ y = 1u << 15; x -= 1u << 16; if (x >= 5u << 14) { x -= 5u << 14; y |= 1u << 14; } goto L6; }
else
{ y = 1u << 17; x -= 1u << 18; if (x >= 5u << 16) { x -= 5u << 16; y |= 1u << 16; } goto L7; }
}
else
{
if (x < 1u << 22)
{ y = 1u << 19; x -= 1u << 20; if (x >= 5u << 18) { x -= 5u << 18; y |= 1u << 18; } goto L8; }
else
{ y = 1u << 21; x -= 1u << 22; if (x >= 5u << 20) { x -= 5u << 20; y |= 1u << 20; } goto L9; }
}
}
else
{
if (x < 1u << 28)
{
if (x < 1u << 26)
{ y = 1u << 23; x -= 1u << 24; if (x >= 5u << 22) { x -= 5u << 22; y |= 1u << 22; } goto La; }
else
{ y = 1u << 25; x -= 1u << 26; if (x >= 5u << 24) { x -= 5u << 24; y |= 1u << 24; } goto Lb; }
}
else
{
if (x < 1u << 30)
{ y = 1u << 27; x -= 1u << 28; if (x >= 5u << 26) { x -= 5u << 26; y |= 1u << 26; } goto Lc; }
else
{ y = 1u << 29; x -= 1u << 30; if (x >= 5u << 28) { x -= 5u << 28; y |= 1u << 28; } }
}
}
}
z = y | 1u << 26; y /= 2; if (x >= z) { x -= z; y |= 1u << 26; }
Lc: z = y | 1u << 24; y /= 2; if (x >= z) { x -= z; y |= 1u << 24; }
Lb: z = y | 1u << 22; y /= 2; if (x >= z) { x -= z; y |= 1u << 22; }
La: z = y | 1u << 20; y /= 2; if (x >= z) { x -= z; y |= 1u << 20; }
L9: z = y | 1u << 18; y /= 2; if (x >= z) { x -= z; y |= 1u << 18; }
L8: z = y | 1u << 16; y /= 2; if (x >= z) { x -= z; y |= 1u << 16; }
L7: z = y | 1u << 14; y /= 2; if (x >= z) { x -= z; y |= 1u << 14; }
L6: z = y | 1u << 12; y /= 2; if (x >= z) { x -= z; y |= 1u << 12; }
L5: z = y | 1u << 10; y /= 2; if (x >= z) { x -= z; y |= 1u << 10; }
L4: z = y | 1u << 08; y /= 2; if (x >= z) { x -= z; y |= 1u << 08; }
L3: z = y | 1u << 06; y /= 2; if (x >= z) { x -= z; y |= 1u << 06; }
L2: z = y | 1u << 04; y /= 2; if (x >= z) { x -= z; y |= 1u << 04; }
L1: z = y | 1u << 02; y /= 2; if (x >= z) { x -= z; y |= 1u << 02; }
L0: return x > y ? y / 2 | 1u : y / 2;
}

利用内置函数计算一个数的平方根

# include"iostream.h"
# include"conio.h"
# include"math.h"
void main()
{
clrscr();
float x;
cout<<"Enter the Number";
cin>>x;


float squreroot(float);
float z=squareroot(x);
cout<<z;




float squareroot(int x)
{




float s;
s = pow(x,.5)
return(s);
}

一个简单的解决方案,可以处理浮点平方根和任意精度使用二进制搜索

用红宝石编码

include Math


def sqroot_precision num, precision
upper   = num
lower   = 0
middle  = (upper + lower)/2.0


while true do
diff = middle**2 - num


return middle if diff.abs <= precision


if diff > 0
upper = middle
else diff < 0
lower = middle
end


middle = (upper + lower)/2.0
end
end


puts sqroot_precision 232.3, 0.0000000001

假设我们想找到2的算术平方根而你找到了 估计是1.5。 我们说 a = 2,x = 1.5。为了计算更好的估计,我们将 a 除以 x。这给出了一个新的值 y = 1.333333。 然而,我们不能仅仅把这作为我们的下一个估计(为什么不呢?).我们需要把它和之前的估计值平均起来。所以我们下一个估计,xx 将是(x + y)/2,或者1.416666。

Double squareRoot(Double a, Double epsilon) {
Double x = 0d;
Double y = a;
Double xx = 0d;


// Make sure both x and y != 0.
while ((x != 0d || y != 0d) && y - x > epsilon) {
xx = (x + y) / 2;


if (xx * xx >= a) {
y = xx;
} else {
x = xx;
}
}


return xx;
}

Epsilon 决定了近似值的精确度。函数应该返回满足 abs (x * x-a) < epsilon 的第一个近似值 x,其中 abs (x)是 x 的绝对值。

square_root(2, 1e-6)
Output: 1.4142141342163086

使用二进制搜索

public class FindSqrt {


public static void main(String[] strings) {


int num = 10000;
System.out.println(sqrt(num, 0, num));
}


private static int sqrt(int num, int min, int max) {
int middle = (min + max) / 2;
int x = middle * middle;
if (x == num) {
return middle;
} else if (x < num) {
return sqrt(num, middle, max);
} else {
return sqrt(num, min, middle);
}
}
}

已经有很多答案了,这是我的,这是最简单的一段代码(对我来说) ,这是它的 算法

Python 2.7中的代码:

from __future__ import division
val = 81
x = 10
def sqr(data,x):
temp = x - ( (x**2 - data)/(2*x))
if temp == x:
print temp
return
else:
x = temp
return sqr(data,x)
#x =temp
#sqr(data,x)
sqr(val,x)