Infinite loop in babylonian method for square roots

Question

In some specific numbers my algorithm gets stuck. It never reaches the minimum approximation so we never get out of the while. I think I can either low my approximation requisites or use double for my numbers, but I'm trying to figure out other solutions.

I'm programming a babylonian algorithm to calculate roots. First I'm doing this in C and later I will do this in Assembly(University homework). When I try to find the root of numbers like 99999 the program iterates to infinity. I have already tried two different stop conditions, one of them I did exactly like this tutorial from geeks4geeks(the first one inside the site).

https://www.geeksforgeeks.org/square-root-of-a-perfect-square/

The second stop condition I tested was this:

while ((x*x - n) > e) {}

I tried something like this because it is more "relatable" to the method enunciation. The full code is showed below:

#include <stdio.h>
#include <math.h>

/*Returns the square root of n. Note that the function */

float squareRoot(float n)

{

    /*We are using n itself as initial approximation

   This can definitely be improved */

    float x = n;

    float y = 1;
    float e = 0.000001; /* e decides the accuracy level*/

    while ((x*x - n) > e) {

        x = (x + y) / 2;
        y = n / x;

//      if(prev_err == x-y){
//          printf("A aproximação por ponto flutuante alcançou o máximo possível para o caso\n");
//          return x;
//      }
//      prev_err = x-y;
        
        
    }


    return x;

}

 

/* Driver program to test above function*/

int main()

{

    int n;

    printf("Insira o número cuja raiz deseja calcular\n");

    scanf("%d", &n);


    printf("Square root of %d is %.8f\n", n, squareRoot(n));

    return 0;

}

Answer 1

float only can take 4 Size (bytes), so If you want to calculate the square root for number greater than 4 bytes You need to replace all float to double like this:

#include <stdio.h>
#include <math.h>

/*Returns the square root of n. Note that the function */

double squareRoot(double n)

{

    /*We are using n itself as initial approximation

   This can definitely be improved */

    double x = n;

    double y = 1;
    double e = 0.000001; /* e decides the accuracy level*/

    while ((x*x - n) > e) {

        x = (x + y) / 2;
        y = n / x;

//      if(prev_err == x-y){
//          printf("A aproximação por ponto flutuante alcançou o máximo possível para o caso\n");
//          return x;
//      }
//      prev_err = x-y;
        
        
    }


    return x;

}

 

/* Driver program to test above function*/

int main()

{

    int n;

    printf("Insira o número cuja raiz deseja calcular\n");

    scanf("%d", &n);


    printf("Square root of %d is %.8f\n", n, squareRoot(n));

    return 0;

}

if u want learn more about Primitive Data Types size in c: Primitive Data Types sizes

https://www.programiz.com/c-programming/c-data-types#:~:text=The%20size%20of%20float%20(single,data%20type)%20is%208%20bytes.

Answer 2

An absolute tolerance will never work. If n is large, x.x - n can remain large and the loop will never stop. If n is tiny, x.x - n can too quickly become small and the result will be quite inaccurate.

Your test has another big flaw: if x.x - n < e, the iterations will stop immediately if the LHS is negative, whatever its value.

---

The cure is to take the absolute value and a relative tolerance.

A better cure is to

• adapt the starting approximation to the magnitude of n (such as the nearest power of 4),

• use a fixed number of iterations (with a good starting approximation, 7 iterations are enough).