网上找到的平方根源代码
Square root source code found on the net
我正在寻找一些计算平方根的算法并找到了这个源文件。我想尝试复制它,因为它看起来很简单,但我无法将它与某些已知算法(Newton、Babylon ...)联系起来。能告诉我名字吗?
int sqrt(int num) {
int op = num;
int res = 0;
int one = 1 << 30; // The second-to-top bit is set: 1L<<30 for long
// "one" starts at the highest power of four <= the argument.
while (one > op)
one >>= 2;
while (one != 0) {
if (op >= res + one) {
op -= res + one;
res += 2 * one;
}
res >>= 1;
one >>= 2;
}
return res;
}
就像 references, this is the classic "digit-by-digit" method done to calculate square root. Learned it in base 10小学时教的那样
OP 的代码也失败了 select 个数字。 sqrt(1073741824) --> -1
而不是预期的 32768。1073741824 == 0x40000000
。此外,它失败了大多数(所有?)值这个和更大的值。当然OP的sqrt(some_negative)
也是个问题
备选方案:也here
unsigned isqrt(unsigned num) {
unsigned res = 0;
// The second-to-top bit is set: 1 << 30 for 32 bits
// Needs work to run on unusual platforms where `unsigned` has padding or odd bit width.
unsigned bit = 1u << (sizeof(num) * CHAR_BIT - 2);
// "bit" starts at the highest power of four <= the argument.
while (bit > num)
bit >>= 2;
while (bit > 0) {
if (num >= res + bit) {
num -= res + bit;
res = (res >> 1) + bit; // Key dif between this and OP's code
} else {
res >>= 1;
}
bit >>= 2;
}
return res;
}
可移植性更新。需要4的最大次方。
#include <limits.h>
// greatest power of 4 <= a power-of-2 minus 1
#define POW4_LE_POW2M1(n) ( ((n)/2 + 1) >> ((n)%3==0) )
unsigned bit = POW4_LE_POW2M1(UINT_MAX);
我正在寻找一些计算平方根的算法并找到了这个源文件。我想尝试复制它,因为它看起来很简单,但我无法将它与某些已知算法(Newton、Babylon ...)联系起来。能告诉我名字吗?
int sqrt(int num) {
int op = num;
int res = 0;
int one = 1 << 30; // The second-to-top bit is set: 1L<<30 for long
// "one" starts at the highest power of four <= the argument.
while (one > op)
one >>= 2;
while (one != 0) {
if (op >= res + one) {
op -= res + one;
res += 2 * one;
}
res >>= 1;
one >>= 2;
}
return res;
}
就像
OP 的代码也失败了 select 个数字。 sqrt(1073741824) --> -1
而不是预期的 32768。1073741824 == 0x40000000
。此外,它失败了大多数(所有?)值这个和更大的值。当然OP的sqrt(some_negative)
也是个问题
备选方案:也here
unsigned isqrt(unsigned num) {
unsigned res = 0;
// The second-to-top bit is set: 1 << 30 for 32 bits
// Needs work to run on unusual platforms where `unsigned` has padding or odd bit width.
unsigned bit = 1u << (sizeof(num) * CHAR_BIT - 2);
// "bit" starts at the highest power of four <= the argument.
while (bit > num)
bit >>= 2;
while (bit > 0) {
if (num >= res + bit) {
num -= res + bit;
res = (res >> 1) + bit; // Key dif between this and OP's code
} else {
res >>= 1;
}
bit >>= 2;
}
return res;
}
可移植性更新。需要4的最大次方。
#include <limits.h>
// greatest power of 4 <= a power-of-2 minus 1
#define POW4_LE_POW2M1(n) ( ((n)/2 + 1) >> ((n)%3==0) )
unsigned bit = POW4_LE_POW2M1(UINT_MAX);