Implementing Square Root function in C not Working-CodePudding

double _sqrt(double n)
{
    double sqrt = 0, c = 0, prec = 1;

    for (;; sqrt  = prec) //increments sqrt per precision
    {
        c = sqrt * sqrt;
        if (c == n)
        {
            return sqrt;
        }
        if (c > n) // if square is greater then..
        {
            sqrt -= prec; // decrement squareroot by previous precision
            prec *= 0.1; // increase precision eg. 1, 0.1, 0.01, 0.001 .....INF
        }
    }
}

This is the square root function that I've done. It works for some numbers but for others its just blank, doesn't return a thing. Where am I getting this wrong?

CodePudding user response：

Your code has a few problems in it, the first one being that your code may infinitely loop as you try to have an infinite accuracy for a (possibly) irrational number. Although doubles do not have an infinite accuracy, I certainly wouldn't recommend trying to evaluate that function to that high of a degree of accuracy. The solution to this is to add some sort of 'precision' or 'epsilon' value (as mentioned in the comments below). Or even better to break out of the loop when increment becomes too small (as @Pete Becker suggested):

double _sqrt(const double& n)
{
    const double precision = 0.000001;
    double increment = 1.0, sqrt = 0.0;

    while (true)
    {
        if (increment <= precision)
            break;
        else if (sqrt * sqrt > n)
        {
            sqrt -= increment;
            increment *= 0.1;
        }
        else 
            sqrt  = increment;
    }

    return sqrt;
}

You can of course do this with the standard library, using std::sqrt, but I assume you are doing this for fun.

If you are interested in other square root algorithms, there are some scary looking ones on Wikipedia here.

One that I personally like is Newton's Method, but this is for finding the root's of functions in general. An application of this for the square root function that I copied and modified from here is:

double _sqrt(const double& n)
{
    const double precision = 0.00001;
    double sqrt, x = n;

    while (true)
    {
        sqrt = 0.5 * (x   (n / x));

        if (std::abs(sqrt - x) < precision)
            break;

        x = sqrt;
    }

    return sqrt;
}

CodePudding user response：

What you are trying to do is similar to binary search. Mine main concern is with prec = 1; and sqrt = prec because once the sqrt is really big number the expression sqrt = prec will not change the sqrt in any way due to rounding and your for loop get stuck...

so what you should do is realize the starting prec value and how many times the precision should increase by prec *= 0.1;

for numbers |x|>1 we now that magnitude of sqrt(x) has half of integer bits that are needed for int(x) so we can derive the prec from that

y = log2(|x|)
prec = y*0.1

the log2 does not need to be real log2 you can simply take exponent and divide it by 2 and construct y as:

y = 2^((exponent_extracted_from_x>>1) 1)

for (|x|<=1) start with for example: prec = 0.5

Now the ending condition you are missing is to stop once sqrt stops changing for example:

...
for (sqrt0=sqrt-prec; sqrt!=sqrt0; sqrt0=sqrt,sqrt =prec)
 {
 ... your original code ...
 }
return sqrt;
}

Also I do not see any sign handling (you know in case x<0)

for more info see: