Department of Mathematics

LAB MANUAL

VI Semester B.Sc. Mathematics

MATDSCP6.2:
PRACTICAL’S ON NUMERICAL ANALYSIS
Contents
1 Program to find root of an equation using bisection and Regula-Falsi
methods 1

2 Program to find root of an equation using Newton-Raphson and Secant

methods 3

3 Program to solve system of algebraic equations using Gauss-elimination

method 5

4 Program to solve system of algebraic equation using Gauss-Jacobi method 5

5 Program to solve system of algebraic equation using Gauss-Seidel method 6

6 Program to evaluate integral using Trapezoidal rule 8

7 Program to evaluate integral using Weddle’s rule 9

8 Program to evaluate integral using Simpson’s one-third rule and Simp-

son’s three-eighth rule 10
1 Program to find root of an equation using bisection
and Regula-Falsi methods

1.1 Bisection Method

1.1.1 Program

# Find a root of the equation x ** 3 - 4 *x - 9 = 0 using Bisection Method .

def f ( x ) :
return x ** 3 - 4 * x - 9

a = float ( input ( " Enter the value of a = " ) )

b = float ( input ( " Enter the value of b = " ) )
if f ( a ) * f ( b ) > 0 . 0 :
print ( f " The root of f ( x ) does not lies in the interval [ { a } ,{ b } ] " )
else :
print ( f " The root of f ( x ) lies in the interval [ { a } ,{ b } ] " )
n = int ( input ( " Enter the number of iterations = " ) )
for i in range ( n + 1 ) :
x=(a+b)/2
print ( " Iteration : %d ,\ t x = % 0 . 6f \ t f ( x ) = % 0 . 6f " % (i ,x , f ( x ) ) )
if f ( a ) * f ( x ) < 0 :
b=x
else :
a=x
print ( f " The approximate root is { x } " )

1.1.2 Program

# Find a root of the equation cos ( x ) -x * e ** x = 0 using Bisection Method .

from sympy import *
def f ( x ) :
return cos ( x ) - x * exp ( x )

a = float ( input ( " Enter the value of a="))

b = float ( input ( " Enter the value of b="))
if f ( a ) * f ( b ) > 0 . 0 :
print ( f " The root of f ( x ) does not lies in the interval [ { a } ,{ b } ] " )
else :
print ( f " The root of f ( x ) lies in the interval [ { a } ,{ b } ] " )
n = int ( input ( " Enter the number of iterations = " ) )
for i in range ( n + 1 ) :
x=(a+b)/2

1
print ( " Iteration : %d ,\ t x = % 0 . 6f \ t f ( x ) = % 0 . 6f " % (i ,x , f ( x ) ) )
if f ( a ) * f ( x ) < 0 :
b=x
else :
a=x
print ( f " The approximate root is { x } " )

1.1.3 Exercise

1. Find the positive root of the equation x3 − 9x + 1 = 0 by bisection method.

2. Find the positive root of f (x) = x − cos x = 0 by bisection method.

3. Use the bisection method in four stages to find the real root of the equation
x log10 x − 102 = 0.

1.2 Regula-Falsi Method

1.2.1 Program

# Find a real root of the equation cos ( x ) = 3 *x - 1 using Regula - Falsi

Method .
from sympy import *
def f ( x ) :
return cos ( x ) - 3 * x + 1

a = float ( input ( " Enter the value of a = " ) )

b = float ( input ( " Enter the value of b = " ) )
if f ( a ) * f ( b ) > 0 . 0 :
print ( f " The root of f ( x ) does not lies in the interval [ { a } ,{ b } ] " )
else :
print ( f " The root of f ( x ) lies in the interval [ { a } ,{ b } ] " )
n = int ( input ( " Enter the number of iterations = " ) )
for i in range ( n + 1 ) :
x=(a*f(b)-b*f(a))/(f(b)-f(a))
print ( " Iteration : %d ,\ t x = % 0 . 6f \ t f ( x ) = % 0 . 6f " % (i ,x , f ( x ) ) )
if f ( a ) * f ( x ) < 0 :
b=x
else :
a=x
print ( f " The approximate root is { x } " )

2
1.2.2 Exercise

1. Find a real root of the equation x3 − 4x + 1 = 0 by Regula falsi method.

2. Compute the real root of x log10 x − 1.2 = 0 by the method of false position.

3. Use the method of false position to find the fourth root of 32.

4. Find the real root of the equation xex − 3 = 0 by Regula-Falsi method correct to
four decimal places.

2 Program to find root of an equation using Newton-

Raphson and Secant methods

2.1 Newton-Raphson Method

2.1.1 Program

# Find a real root of the equation 3x = cos ( x ) + 1 using Newton - Raphson

Method .
from sympy import *
def f ( x ) :
return 3 * x - cos ( x ) - 1

def df ( x ) :
return 3 + sin ( x )

a = float ( input ( " Enter the value of a = " ) )

b = float ( input ( " Enter the value of b = " ) )
if f ( a ) * f ( b ) > 0 . 0 :
print ( f " The root of f ( x ) does not lies in the interval [ { a } ,{ b } ] " )
else :
print ( f " The root of f ( x ) lies in the interval [ { a } ,{ b } ] " )
n = int ( input ( " Enter the number of iterations = " ) )
x0 = ( a + b ) / 2
for i in range ( n + 1 ) :
xn = x0 - f ( x0 ) / df ( x0 )
print ( " Iteration : %d ,\ t x = % 0 . 6f \ t f ( x ) = % 0 . 6f " % (i , xn , f ( xn ) ) )
x0 = xn
print ( f " The approximate root is { x } " )

3
2.1.2 Exercise

1. Use Newton-Raphson method to find a real root of the equation x3 − 2x − 5 = 0.

2. Applying Newton-Raphson method to find the real root of the equation x4 −x−10 =
0 which is near to x = 2.

3. Using Newton-Raphson method find a real root of the equation sin x = 1 − x (near
0.5).

4. Using Newton-Raphson method find a real root of the equation x2 + x = cos x (near
0.5).

2.2 Secant Method

2.2.1 Program

Find a real root of the equation cos x = 3x − 1 using secant Method.

# Find a real root of the equation cos ( x ) = 3 *x - 1 using secant Method .
from sympy import *
def f ( x ) :
return x ** 3 + x - 1
# return cos ( x ) -3 * x + 1

a = float ( input ( " Enter the value of a = " ) )

b = float ( input ( " Enter the value of b = " ) )
if f ( a ) * f ( b ) > 0 . 0 :
print ( f " The root of f ( x ) does not lies in the interval [ { a } ,{ b } ] " )
else :
print ( f " The root of f ( x ) lies in the interval [ { a } ,{ b } ] " )
i=1
while True :
# for i in range (1 , n ) :
x=(a*f(b)-b*f(a))/(f(b)-f(a))
c=f(a)*f(b)
print ( " Iteration : %d ,\ t x = % 0 . 6f \ t f ( x ) = % 0 . 6f " % (i ,x , f ( x ) ) )
if ( c = = 0 ) or ( abs ( a - x ) < 0 . 00001 ) :
print ( f " The approximate root is { round (x , 6 ) } after { i }
iterations " )
break ;
else :
i+=1
a=b
b=x

4
3 Program to solve system of algebraic equations us-
ing Gauss-elimination method

3.1 Program
Gauss elimination method to solve 10x+2y+z = 9; x+10y−z = −22; −2x+3y+10z = 22.
# Gauss elimination method to solve 10x + 2y + z = 9 ; x + 10y - z = - 22 ; - 2x + 3y + 10z =
22
import numpy as np
a = np . array ( [ [ 10 ,2 , 1 ] ,[1 , 10 , - 1 ] ,[ -2 ,3 , 10 ] ] )
b = np . array ( [9 , - 22 , 22 ] )
n = len ( b )
# Elimination Phase
for k in range (0 , n - 1 ) :
for i in range ( k +1 , n ) :
if a [i , k ] ! = 0 . 0 :
lam = a [i , k ] / a [k , k ]
a [i , k + 1 : n ] = a [i , k + 1 : n ] - lam * a [k , k + 1 : n ]
b [ i ] = b [ i ] - lam * b [ k ]
# Back substitution
for k in range ( n -1 , -1 , - 1 ) :
b [ k ] = ( b [ k ] - np . dot ( a [k , k + 1 : n ] ,b [ k + 1 : n ] ) ) / a [k , k ]
print ( " The solution of the given system is " ,b )

3.2 Exercise
Solve the following system of equations by Gauss elimination method.

1. 5x + 2y + z = 12; x + 4y + 2z = 15; x + 2y + 5z = 0

2. 10x − 2y + z = 12; x + 9y − z = 10; 2x − y + 11z = 20

3. 10x + y + z = 12; 2x + 10x + z = 13; 2x + 2y + 10z = 14

4 Program to solve system of algebraic equation us-

ing Gauss-Jacobi method

4.1 Program
Solve by Gauss-Jacobi method 10x + 2y + z = 9; x + 10y − z = −22; −2x + 3y + 10z = 22.

5
# Solve by Gauss - Jacobi method 10x + 2y + z = 9 ; x + 10y - z = - 22 ; - 2x + 3y + 10z = 22
f1 = lambda x ,y , z : ( 9 - 2 * y - z ) / 10
f2 = lambda x ,y , z : ( - 22 - x + z ) / 10
f3 = lambda x ,y , z : ( 22 + 2 * x - 3 * y ) / 10
# Intitial conditions
x0 , y0 , z0 =0 ,0 , 0
itr = 1
e = float ( input ( " Enter the tolerable error : " ) )
condition = True
print ( " \ n \ itr \ t x \ t y \ t z \ n " )
while condition :
x1 = f1 ( x0 , y0 , z0 )
y1 = f2 ( x0 , y0 , z0 )
z1 = f3 ( x0 , y0 , z0 )
print ( " % d \ t % 0 . 4f \ t % 0 . 4f \ t % 0 . 4f \ n " % ( itr , x1 , y1 , z1 ) )
e1 = abs ( x0 - x1 ) ; e2 = abs ( y0 - y1 ) ; e3 = abs ( z0 - z1 )
itr + = 1
x0 = x1 ; y0 = y1 ; z0 = z1
condition = e1 > e and e2 > e and e3 > e
print ( " \ n Solution : x = % 0 . 4f , y = % 0 . 4f , z = % 0 . 4f \ n " % ( x1 , y1 , z1 ) )

4.2 Exercise
Solve the following system of equations by Gauss-Jacobi method.

1. 5x + 2y + z = 12; x + 4y + 2z = 15; x + 2y + 5z = 0

2. 10x − 2y + z = 12; x + 9y − z = 10; 2x − y + 11z = 20

3. 10x + y + z = 12; 2x + 10x + z = 13; 2x + 2y + 10z = 14

5 Program to solve system of algebraic equation us-

ing Gauss-Seidel method
5.1 Program
# Solve by Gauss - Seidel method 10x + 2y + z = 9 ; x + 10y - z = - 22 ; - 2x + 3y + 10z = 22
f1 = lambda x ,y , z : ( 9 - 2 * y - z ) / 10
f2 = lambda x ,y , z : ( - 22 - x + z ) / 10
f3 = lambda x ,y , z : ( 22 + 2 * x - 3 * y ) / 10
# Intitial conditions
x0 , y0 , z0 =0 ,0 , 0

6
itr = 1
e = float ( input ( " Enter the tolerable error : " ) )
condition = True
print ( " \ n \ itr \ t x \ t y \ t z \ n " )
while condition :
x1 = f1 ( x0 , y0 , z0 )
y1 = f2 ( x1 , y0 , z0 )
z1 = f3 ( x1 , y1 , z0 )
print ( " % d \ t % 0 . 4f \ t % 0 . 4f \ t % 0 . 4f \ n " % ( itr , x1 , y1 , z1 ) )
e1 = abs ( x0 - x1 ) ; e2 = abs ( y0 - y1 ) ; e3 = abs ( z0 - z1 )
itr + = 1
x0 = x1 ; y0 = y1 ; z0 = z1
condition = e1 > e and e2 > e and e3 > e
print ( " \ n Solution : x = % 0 . 4f , y = % 0 . 4f , z = % 0 . 4f \ n " % ( x1 , y1 , z1 ) )

5.2 Exercise
Solve the following system of equations by Gauss-Seidel method.

1. 5x − 2y + z = −4; x + 6y − 2z = −1; 3x + y + 5z = 13

2. 9x + 2y + 4z = 20; x + 10y + 4z = 6; 2x − 4y + 10z = −15

3. 10x + y + z = 12; x + 10y + z = 12; x + y + 10z = 12

7
6 Program to evaluate integral using Trapezoidal rule
6.1 Program
R6 1
Evaluate dx taking n = 6 using trapezoidal rule.
0 1 + x2

# Evaluate int_0 ^ 6 ( 1 /( 1 + x ^ 2 ) ) dx taking n = 6 using trapezoidal rule

x0 = int ( input ( " Lower limit of the interval : " ) )
xn = int ( input ( " Upper limit of the interval : " ) )
n = int ( input ( " Number of sub intervals : " ) )
def f ( x ) :
return ( 1 / ( 1 + x ** 2 ) )
h = ( xn - x0 ) / n
sum1 = f ( x0 ) + f ( xn )
for i in range (1 , n ) :
sum1 = sum1 + 2 * f ( x0 + i * h )
print ( " Estimated value of given integration : " ,h / 2 * sum1 )

6.2 Exercise
R5
1. Evaluate log10 xdx taking n = 6.
1

R3 1
2. Evaluate dx taking n = 6.
1 1+x
R6 1
3. Evaluate dx taking n = 6.
0 1 + x2

8
7 Program to evaluate integral using Weddle’s rule
7.1 Program
R6 1
Evaluate dx taking n = 6 using Weddle’s rule.
0 1 + x2

# Evaluate int_0 ^ 6 ( 1 /( 1 + x ^ 2 ) ) dx taking n = 6 using Weddle 's rule

x0 = int ( input ( " Lower limit of the interval : " ) )
xn = int ( input ( " Upper limit of the interval : " ) )
n = int ( input ( " Number of sub intervals : " ) )
def f ( x ) :
return ( 1 / ( 1 + x ** 2 ) )
h = ( xn - x0 ) / n
sum1 = f ( x0 ) + f ( xn )
for i in range (1 ,n , 6 ) :
sum1 = sum1 + ( 5 * f ( x0 + i * h ) + f ( x0 + ( i + 1 ) * h ) + 6 * f ( x0 + ( i + 2 ) * h ) + f ( x0 + ( i + 3 ) * h ) +
5 * f ( x0 + ( i + 4 ) * h ) )
print ( " Estimated value of given integration : " ,3 * h / 10 * sum1 )

7.2 Exercise
R5
1. Evaluate log10 xdx taking n = 6.
1

R3 1
2. Evaluate dx taking n = 6.
1 1+x
R6 1
3. Evaluate dx taking n = 6.
0 1 + x2

9
8 Program to evaluate integral using Simpson’s one-
third rule and Simpson’s three-eighth rule

8.1 Simpson’s one-third rule

R6 1
Evaluate dx taking n = 6 using Simpson’s one-third rule.
0 1 + x2

8.1.1 Program

# Evaluate int_0 ^ 6 ( 1 /( 1 + x ^ 2 ) ) dx taking n = 6 using Simpson 's one - third

rule
x0 = int ( input ( " Lower limit of the interval : " ) )
xn = int ( input ( " Upper limit of the interval : " ) )
n = int ( input ( " Number of sub intervals : " ) )
def f ( x ) :
return ( 1 / ( 1 + x ** 2 ) )
h = ( xn - x0 ) / n
sum1 = f ( x0 ) + f ( xn )
for i in range (1 , n ) :
if i % 2 = = 0 :
sum1 = sum1 + 2 * f ( x0 + i * h )
else :
sum1 = sum1 + 4 * f ( x0 + i * h )
print ( " Estimated value of given integration : " ,h / 3 * sum1 )

8.1.2 Exercise
R5
1. Evaluate log10 xdx taking n = 6.
1

R3 1
2. Evaluate dx taking n = 6.
1 1+x

R6 1
3. Evaluate dx taking n = 6.
0 1 + x2

8.2 Simpson’s three-eighth rule

R6 1
Evaluate dx taking n = 6 using Simpson’s three-eighth rule.
0 1 + x2

8.2.1 Program

10
# Evaluate int_0 ^ 6 ( 1 /( 1 + x ^ 2 ) ) dx taking n = 6 using Simpson 's three - eighth
rule
x0 = int ( input ( " Lower limit of the interval : " ) )
xn = int ( input ( " Upper limit of the interval : " ) )
n = int ( input ( " Number of sub intervals : " ) )
def f ( x ) :
return ( 1 / ( 1 + x ** 2 ) )
h = ( xn - x0 ) / n
sum1 = f ( x0 ) + f ( xn )
for i in range (1 , n ) :
if i % 3 = = 0 :
sum1 = sum1 + 2 * f ( x0 + i * h )
else :
sum1 = sum1 + 3 * f ( x0 + i * h )
print ( " Estimated value of given integration : " ,3 * h / 8 * sum1 )

8.2.2 Exercise
R5
1. Evaluate log10 xdx taking n = 6.
1

R3 1
2. Evaluate dx taking n = 6.
1 1+x
R6 1
3. Evaluate dx taking n = 6.
0 1 + x2

11