CODEWITHSUNDEEP
cinfinite-loopnewtons-method
user3306583
user3306583

Reputation: 119

Infinite loop in Newton Raphson method

I'm writing a program to calculate a value for lambda except it's getting stuck in the do loop. It doesn't seem to be updating the n value as I was hoping that by setting nMax = 100 it would at least finish quickly in case the other case (eps >= e_allow) never happened. Can anyone spot where the infinite loop is coming from?

Thanks a lot for your help. My program is: 

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

//Solution to Linear Dispersion Relationship for a Horizontal Bottom

double f ( double, double, double, double);
double df (double, double, double); 

int main(void) {


    float pi = 3.1415927;
    double g = 9.81;                //Acceleration due to gravity
    double omega;                  //Angular frequency of wave
    double kn, dk;             //Wave numbers
    double e_allow, eps;    //Tolerance used to find difference between kNew and kOld   
    int nMax, n;          // Setting up counter for iterations
    double lambda;                 //Value we need to solve for 
    double h, T;                   //Depth and Period. User needs to input these values.
    double fn, dfn;               //Declaring functions describing linear dispersion


    printf("Enter the value of the depth (m).\n");
    scanf("%lf", &h); 
    printf("\nEnter the value for the period (s).\n");
    scanf("%lf", &T);


    /* 
    Assuming shallow water conditions. Therefore, tanh(k*d) --> kd. kOld is calculated using 
    the relationship omega^2 = kold^2/(g*d) 
    */

    omega = (2.0 * pi)/T;
    printf("Angular velocity = %.5f\n", omega); 
    kn = sqrt((omega * omega)/(g * h));             
    printf("Wave number kn = %.5f\n", kn);

    e_allow = 1.0e-5;
    n = 0; 
    nMax = 100;

    //Following the Newton Raphson Method
    do {
        fn = f(kn, omega, g, h);
        dfn = df(kn, g, h);
        dk = -(fn/dfn);
        kn = kn + dk;
        eps = fabs(dk/kn);
        n = n + 1;
    } while (eps >= e_allow || n == nMax);


    //Calculating value of lambda using final kn value. 
    lambda = (2.0 * pi)/kn;


    //printf("\nfn = %.6f, dfn = %.6f, dx = %.6f, eps = %.6f\n", fn, dfn, dk, eps);


    printf("\nkn = %.5lf, Wavelength = %.5lf\n", kn, lambda); */

    //printf("fn = %.6lf, dfn = %.6lf, eps = %.6lf, dk = %.6lf, kn = %.6lf\n", fn, dfn, eps, dk, kn);



    return 0;

}


double f (double kn, double omega, double g, double h) {
    return ((omega*omega) - (g * kn * tanh(kn * h)));
}

double df (double kn, double g, double h) {
    return ((( g * kn * h )/( 1 + ( h * h * kn * kn ))) - g * tanh( kn * h));
}

Upvotes: 0

Views: 475

Answers (1)

Lutz Lehmann
Lutz Lehmann

Reputation: 25992

Your derivative is partially wrong. Somehow you plugged in the derivative formula for the arctan function.

The derivative of tanh(x) is 1/square(cosh(x))=1-square(tanh(x)).

And keep track of the sign, the minus sign does not get lost in the first term. The correct sign alone makes the code work, without correcting the rest of the derivative.

Upvotes: 1

Related Questions