平方根算法 C++

square root algorithm C++

我不明白为什么如果输入的数字超过 12 位,这个算法会进入无限循环。谁能明白为什么它永远不会结束?谢谢。我刚刚更新了算法以使用 fabs() 函数,但仍然得到一个无限循环。

double squareroot(double x)

{ /* computes the square root of x */

/* make sure x is not negative .. no math crimes allowed! */
assert( x >= 0 );
if (x==0) return 0;

/* the sqrt must be between xhi and xlo */
double xhi = x;
double xlo = 0;
double guess = x/2;

/* We stop when guess*guess-x is very small */

while (abs(guess*guess-x) > 0.00001 )
{
    if (guess*guess > x){
        xhi = guess;
    }

    else {
        xlo = guess;
    }

    guess = (xhi + xlo)/2;
}
return guess;
}

http://floating-point-gui.de/errors/comparison/

嘿,看来你不应该那样使用 abs() 。在某些情况下它应该停止,但我不会,因为它的精度有限。

改为使用 fabs()

http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm

我相信你应该使用相对误差来终止,而不是绝对误差。

while (abs((guess*guess-x) / guess) > 0.00001)

否则计算非常长的值的平方根需要很长时间(不是无限循环)。

http://en.wikipedia.org/wiki/Approximation_error

干杯!

编辑: 此外,正如下面评论中所指出的,值得检查 guess 是否已经被猜到,以避免某些特定的无限循环边角案例。

我会说,当数字足够大时,您不能使用绝对 epsilon 值,因为它不符合精度。

尝试改用相对比较。考虑以下函数来检查 2 个双打是否相等:

bool are_eql(double n1, double n2, double abs_e, double rel_e)
{
  // also add check that n1 and n2 are not NaN
  double diff = abs(n1 - n2);
  if (diff < abs_e)
    return true;
  double largest = max(abs(n1), abs(n2));
  if (diff < largest * rel_e)
    return true;
  return false;
}

这不是您问题的直接答案,而是替代解决方案。

您可以使用 Newton's method for finding roots:

assert(x >= 0);
if (x == 0)
    return 0;

double guess = x;
for (int i=0; i<NUM_OF_ITERATIONS; i++)
    guess -= (guess*guess-x)/(2*guess);
return guess;

24 次迭代应该可以获得足够好的近似值,但您也可以检查绝对差异。

我建议等到你得到一个稳定的答案,而不是摆弄 epsilon 值:

double squareroot(double x)
{
    if (x < 1) return 1.0 / squareroot(x);  // MSalter's general solution

    double xhi = x;
    double xlo = 0;
    double guess = x/2;

    while (guess * guess != x)
    {
        if (guess * guess > x)
            xhi = guess;
        else
            xlo = guess;

        double new_guess = (xhi + xlo) / 2;
        if (new_guess == guess)
            break; // not getting closer
        guess = new_guess;
    }
    return guess;
}