Numerical Methods

Newton-Raphson method

  Newton-Raphson is a technique to find the roots of algebraic and transcendental equations of the form `f(x)=0`. We can start by assuming that `x_0` is an approximate solution of \( f(x)=0 \) such that \(x_1=x_0+h\) is an exact solution. $$\therefore f(x_1)=0 $$ As, \(h\) is a very shall number, we can Taylor expand `f(x_1)` about `x_0` to get:

$$f(x_1)=f(x_0+h)=f(x_0)+\frac{h}{1!} f^{\prime}(x_0)+\frac{h^2}{2!} f^{\prime\prime}(x_0)+\cdots=0$$

As, `h` is very small we can neglect `h^2` and higher powers of `h` as well as higher derivatives of `f(x_0)` to obtain $$f(x_0)+hf^{\prime}(x_0) \cong 0$$ Or, if \(f^{\prime}(x_0)\neq 0\); $$ h=-\frac{f(x_0)}{f^{\prime}(x_0)}$$ Therefore, \(x_1=x_0-\frac{f(x_0)}{f^{\prime}(x_0)}\) is the first approximation to the root of the equation. Similarly, we can get the next approximate root (say `x_2`) by taking `x_1` as the trial solution. Hence, we can write $$x_2=x_1- \frac{f(x_1)}{f^{\prime}(x_1)}$$ In general, $$\boxed{x_{n+1}=x_n- \frac{f(x_n)}{f^{\prime}(x_n)},\; n=0,1,2,\cdots}$$ The above equation is the Newton-Raphson formula. Starting from an initial approximate root `x_0` we can get successive approximations `x_n` by applying this formula iteratively for `n` times.


Program in C
========================================================================= Program to find solutions to equations of the form f(x)=0 using Newton-Raphson Method LANGUAGE :: C Compiler: GNU GCC Compiler =========================================================================
#include <stdio.h> #include <stdlib.h> #include <math.h> double f(double x){ double y; y=(x*x*x)-(2*x)-5; return y; } double df(double x){ double y; y=3*(x*x)-2; return y; } int main() { int i,n; double x,x0,e; printf("Enter the maximum number of iterations:\n"); scanf("%d",&n); printf("\nEnter the trial solution:\n"); scanf("%lf",&x0); printf("\nEnter the tolerance:\n"); scanf("%lf",&e); printf("\nIterations:\n"); for(i=1;i<=n;i++){ x=x0-f(x0)/df(x0); printf("Iteration:%2d x=%10.6lf f(x)=%10.6lf \n",i,x,f(x)); if(fabs(x-x0)<=e) break; x0=x; } printf("\nTHE APPROXIMATE SOLUTION IS: %lf",x); return 0; }
OUTPUT of C program
Enter the maximum number of iterations: 
100

Enter the trial solution:
0

Enter the tolerance:
0.000001

Iterations:
Iteration: 1  x= -2.500000  f(x)=-15.625000
Iteration: 2  x= -1.567164  f(x)= -5.714632
Iteration: 3  x= -0.502592  f(x)= -4.121770
Iteration: 4  x= -3.820706  f(x)=-53.132488
Iteration: 5  x= -2.549393  f(x)=-16.470758
Iteration: 6  x= -1.608111  f(x)= -5.942390
Iteration: 7  x= -0.576100  f(x)= -4.039002
Iteration: 8  x= -4.597710  f(x)=-92.995258
Iteration: 9  x= -3.083543  f(x)=-28.151977
Iteration:10  x= -2.022194  f(x)= -9.224909
Iteration:11  x= -1.123764  f(x)= -4.171613
Iteration:12  x=  1.208652  f(x)= -5.651658
Iteration:13  x=  3.580790  f(x)= 33.751515
Iteration:14  x=  2.655233  f(x)=  8.409627
Iteration:15  x=  2.216106  f(x)=  1.451367
Iteration:16  x=  2.102125  f(x)=  0.084892
Iteration:17  x=  2.094584  f(x)=  0.000358
Iteration:18  x=  2.094551  f(x)=  0.000000
Iteration:19  x=  2.094551  f(x)= -0.000000

THE APPROXIMATE SOLUTION IS: 2.094551
 
Program in Fortran 95
===================================================================================
Program to find solutions to equations of the form f(x)=0 using                
Newton-Raphson Method                                                    
LANGUAGE :: FORTRAN 95  
Compiler :: GNU Fortran
===================================================================================
PROGRAM newtonRaphson
    INTEGER i,n
    REAL x,x0,e
    WRITE(*,*)'Enter the maximum number of iterations'
    READ*,n
    WRITE(*,*)'Enter the trial solution'
    READ*,x0
    WRITE(*,*)'Enter the tolerance'
    READ*,e
    WRITE(*,*)'Iterations:'
    DO i=1,n
        x=x0-f(x0)/df(x0)
        WRITE(*,10)i,x,f(x)
        IF (ABS(x-x0)<=e) EXIT
        x0=x
    END DO
10  FORMAT(/,1x,'Iteration',I2,':',4X,'x=',F12.9,4X,'f(x)=',F12.8)
    PRINT*,'THE APPROXIMATE SOLUTION IS',x
    STOP
 END PROGRAM

 REAL FUNCTION f(x)
 REAL x
 f=(x**3)-(2*x)-5
 RETURN
 END

REAL FUNCTION df(x)
REAL x
df=3*(x**2)-2
RETURN
END



OUTPUT of Fortran95 program
 Enter the maximum number of iterations
 100
 Enter the trial solution
 0
 Enter the tolerance
 0.000001
 Iterations:
 Iteration 1:    x=-2.500000000    f(x)=-15.62500000
 Iteration 2:    x=-1.567164183    f(x)= -5.71463251
 Iteration 3:    x=-0.502592385    f(x)= -4.12176943
 Iteration 4:    x=-3.820705891    f(x)=-53.13246536
 Iteration 5:    x=-2.549392939    f(x)=-16.47075081
 Iteration 6:    x=-1.608111024    f(x)= -5.94238663
 Iteration 7:    x=-0.576099694    f(x)= -4.03900290
 Iteration 8:    x=-4.597699642    f(x)=-92.99465179
 Iteration 9:    x=-3.083536386    f(x)=-28.15179825
 Iteration10:    x=-2.022189140    f(x)= -9.22485638
 Iteration11:    x=-1.123758674    f(x)= -4.17160273
 Iteration12:    x= 1.208699346    f(x)= -5.65154457
 Iteration13:    x= 3.580445528    f(x)= 33.73895264
 Iteration14:    x= 2.655045271    f(x)=  8.40602684
 Iteration15:    x= 2.216037750    f(x)=  1.45049477
 Iteration16:    x= 2.102116823    f(x)=  0.08480024
 Iteration17:    x= 2.094583511    f(x)=  0.00035763
 Iteration18:    x= 2.094551563    f(x)=  0.00000095
 Iteration19:    x= 2.094551563    f(x)=  0.00000095
 THE APPROXIMATE SOLUTION IS   2.09455156
 
Program in Scilab
======================================================================
 Program to find solutions to equations of the form f(x)=0 using                
Newton-Raphson Method                                                    										
 LANGUAGE  : SCILAB                                  	
=======================================================================

function y=f(x)
    y=(x^3)-(2*x)-5
endfunction

function y=df(x)
    y=3*(x^2)-2
endfunction

n = input("Enter the maximum number of iterations:")
x0 = input("Enter the trial solution:")
e = input("Enter the tolerance:")
mprintf("Iterations:\n")
mprintf("\n %d %f %f",n,x0,e)
for i=1:1:n
    x = x0 - f(x0)/df(x0)
    mprintf("\nIteration %2d   x=%10.6f   f(x)=%10.6f",i,x,f(x))
    if (abs(x-x0)<=e) then
        break
    end
    x0 = x
end
mprintf("\nThe approximate solution is %f",x)


OUTPUT of Scilab program
 Enter the maximum number of iterations:100

Enter the trial solution:0

Enter the tolerance:0.000001

Iterations:

Iteration  1   x= -2.500000   f(x)=-15.625000
Iteration  2   x= -1.567164   f(x)= -5.714632
Iteration  3   x= -0.502592   f(x)= -4.121770
Iteration  4   x= -3.820706   f(x)=-53.132488
Iteration  5   x= -2.549393   f(x)=-16.470758
Iteration  6   x= -1.608111   f(x)= -5.942390
Iteration  7   x= -0.576100   f(x)= -4.039002
Iteration  8   x= -4.597710   f(x)=-92.995258
Iteration  9   x= -3.083543   f(x)=-28.151977
Iteration 10   x= -2.022194   f(x)= -9.224909
Iteration 11   x= -1.123764   f(x)= -4.171613
Iteration 12   x=  1.208652   f(x)= -5.651658
Iteration 13   x=  3.580790   f(x)= 33.751515
Iteration 14   x=  2.655233   f(x)=  8.409627
Iteration 15   x=  2.216106   f(x)=  1.451367
Iteration 16   x=  2.102125   f(x)=  0.084892
Iteration 17   x=  2.094584   f(x)=  0.000358
Iteration 18   x=  2.094551   f(x)=  0.000000
Iteration 19   x=  2.094551   f(x)= -0.000000
The approximate solution is 2.094551