Functions and floating point comparison

Question

#include<stdio.h>
#include<stdlib.h>
#define abs(a) ((a)>0 ? a: -a)
#define eps_sqrt 0.00000000000001
#define it 100

float sqrt(float x)
/*The Square Root Function using the Newton's Method*/
{
    int it_sqrt=0;
    float a_sqrt = x/2;
    while ((abs((a_sqrt*a_sqrt)-(x))>=eps_sqrt) && (2.0*a_sqrt != 0) && (it_sqrt<=it))
    {
        a_sqrt = a_sqrt - ((a_sqrt*a_sqrt)-(x)/(2.0*a_sqrt));
        it_sqrt++;
    }
    return a_sqrt;
}

int main()
{
    printf("%.5f\n", sqrt(5));
    system ("pause");
}

i tried using the Newton's iteration method to find the square root on Python and it worked, perfectly well. I'm new on C and I don't understand why this function didn't work for me. Whenever I run it, it returns "-1.#INF0A" Any help will be appreciated.

---

Edit: I tried changin the eps to 0.000001 and it also didn't work.

Answer 1

double mysqrt(double x){
    double eps=pow(10,-10);
    double x0 = 0.0;
    double x1 = x/2.0;
    while(fabs(x1 - x0)>eps){
        x0 = x1;
        x1 = x0 + (x - x0*x0)/x0/ 2.0;
    }
    return x1;
}

macro expansion
abs((a_sqrt*a_sqrt)-(x))
expansion (((a_sqrt*a_sqrt)-(x))>0 ? (a_sqrt*a_sqrt)-(x): -(a_sqrt*a_sqrt)-(x))
NG: -(a_sqrt*a_sqrt)-(x)

abs((a_sqrt*a_sqrt- x))
expansion (((a_sqrt*a_sqrt- x))>0 ? (a_sqrt*a_sqrt- x): -(a_sqrt*a_sqrt- x))

rewrite
#define abs(a) ((a)>0 ? a: -a)
to
#define abs(a) ((a)>0 ? a: -(a))

Answer 2

Changing this line:

                a_sqrt = a_sqrt - ((a_sqrt*a_sqrt)-(x)/(2.0*a_sqrt));

to

                a_sqrt = a_sqrt - ((a_sqrt*a_sqrt - x)/(2.0*a_sqrt));

works for me.

Answer 3

Try to use a bigger epsilon, maybe python uses doubles instead of floats.

Answer 4

This is one of the rare cases where using double actually makes sense. Note that the precision of float is significantly lower than eps_sqrt:

[mic@mic-nb tmp]$ cat tmp2.c
#include <stdio.h>
#include <math.h>

int main() {
    double a = sqrtl(2.0);
    printf("%1.20f\n", a - (float) a);
}
[mic@mic-nb tmp]$ gcc tmp2.c; ./a.out
0.00000002420323430563
vs. your value of:
0.00000000000001

So your program will, in most cases, never terminate.