CODEWITHSUNDEEP
c++algorithmoopsquare-root
csciXAV_12
csciXAV_12

Reputation: 309

square root algorithm C++

I can not figure out why this algorithm enters an infinite loop if the number entered is over 12 digits. Can anyone see why it will never end? Thanks. I just updated the algorithm to use the fabs() function and still get an infinite loop. 

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;
}

Upvotes: 5

Views: 10317

Answers (5)

Tony Delroy
Tony Delroy

Reputation: 106246

I suggest waiting until you've got a stable answer, rather than fiddling with epsilon values:

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;
}

Upvotes: 4

denis-bu
denis-bu

Reputation: 3446

I'd say that when numbers are big enough you can't use absolute epsilon value because it doesn't fit into precision.

Try to use relative comparison instead. Consider following function to check if 2 doubles are equal:

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;
}

Upvotes: 0

ale64bit
ale64bit

Reputation: 6242

I believe you should use relative error for the termination, and not the absolute error. 

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

Otherwise it will take very long time (it's not an infinite loop) to compute square root of very long values.

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

Cheers!

EDIT: moreover, as pointed below in the comments, it is worthy to check if the guess was already guessed in order to avoid infinite loop with some specific corner cases.

Upvotes: 6

barak manos
barak manos

Reputation: 30146

This is not a direct answer to your question, but an alternative solution.

You can use 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 iterations should get you a good enough approximation, but you can also check the absolute difference.

Upvotes: 3

&#193;lvaro G&#243;mez
&#193;lvaro G&#243;mez

Reputation: 306

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

Heh, looks like you are not supposed to use abs() like that. There are some cases where it should stop but i won't since it's limited precision.

Instead use fabs()

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

Upvotes: 0

Related Questions