Scribd
Upload
25 views18 pages

Numerical Method

The document provides implementations of various numerical methods including Bisection, Regula Falsi, Newton Raphson, Secant, interpolation methods (Forward, Backward, Linear, Lagrange), and integration methods (Trapezoidal, Simpson's 1/3, Simpson's 3/8). Each method includes a function definition, input handling, and iterative calculations to find roots or perform interpolations and integrations. The document is structured with code snippets and explanations for each method.

Uploaded by

mandev Khatri
Copyright
© © All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd
25 views18 pages

Numerical Method

The document provides implementations of various numerical methods including Bisection, Regula Falsi, Newton Raphson, Secant, interpolation methods (Forward, Backward, Linear, Lagrange), and integration methods (Trapezoidal, Simpson's 1/3, Simpson's 3/8). Each method includes a function definition, input handling, and iterative calculations to find roots or perform interpolations and integrations. The document is structured with code snippets and explanations for each method.

Uploaded by

mandev Khatri
Copyright
© © All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd
You are on page 1/ 18

Bisection method

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

# Implementing Bisection Method


def bisection(x0,x1,e):
step = 1
print('\n\n*** BISECTION METHOD IMPLEMENTATION
***')
condition = True
while condition:
x2 = (x0 + x1)/2
print('Iteration-%d, x2 = %0.6f and f(x2) =
%0.6f' % (step, x2, f(x2)))

if f(x0) * f(x2) < 0:


x1 = x2
else:
x0 = x2

step = step + 1
condition = abs(f(x2)) > e

print('\nRequired Root is : %0.8f' % x2)

# Input Section
x0 = input('First Guess: ')
x1 = input('Second Guess: ')
e = input('Tolerable Error: ')
# Converting input to float
x0 = float(x0)
x1 = float(x1)
e = float(e)

#Note: You can combine above two section like this


# x0 = float(input('First Guess: '))
# x1 = float(input('Second Guess: '))
# e = float(input('Tolerable Error: '))

# Checking Correctness of initial guess values and


bisecting
if f(x0) * f(x1) > 0.0:
print('Given guess values do not bracket the
root.')
print('Try Again with different guess values.')
else:

bisection(x0,x1,e)

Regula falsi method


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

# Implementing False Position Method


def falsePosition(x0,x1,e):
step = 1
print('\n\n***
FALSE POSITION METHOD
IMPLEMENTATION ***')
condition = True
while condition:
x2 = x0 - (x1-x0) * f(x0)/( f(x1) - f(x0) )
print('Iteration-%d, x2 = %0.6f and f(x2) =
%0.6f' % (step, x2, f(x2)))

if f(x0) * f(x2) < 0:


x1 = x2
else:
x0 = x2

step = step + 1
condition = abs(f(x2)) > e

print('\nRequired root is: %0.8f' % x2)

# Input Section
x0 = input('First Guess: ')
x1 = input('Second Guess: ')
e = input('Tolerable Error: ')

# Converting input to float


x0 = float(x0)
x1 = float(x1)
e = float(e)

#Note: You can combine above two section like this


# x0 = float(input('First Guess: '))
# x1 = float(input('Second Guess: '))
# e = float(input('Tolerable Error: '))
# Checking Correctness of initial guess values and
false positioning
if f(x0) * f(x1) > 0.0:
print('Given guess values do not bracket the
root.')
print('Try Again with different guess values.')
else:

falsePosition(x0,x1,e)

Nr method
def f(x):
return x**3 - 5*x - 9

# Defining derivative of function


def g(x):
return 3*x**2 - 5

# Implementing Newton Raphson Method

def newtonRaphson(x0,e,N):
print('\n\n***
NEWTON RAPHSON METHOD
IMPLEMENTATION ***')
step = 1
flag = 1
condition = True
while condition:
if g(x0) == 0.0:
print('Divide by zero error!')
break

x1 = x0 - f(x0)/g(x0)
print('Iteration-%d, x1 = %0.6f and f(x1) =
%0.6f' % (step, x1, f(x1)))
x0 = x1
step = step + 1

if step > N:
flag = 0
break

condition = abs(f(x1)) > e

if flag==1:
print('\nRequired root is: %0.8f' % x1)
else:
print('\nNot Convergent.')

# Input Section
x0 = input('Enter Guess: ')
e = input('Tolerable Error: ')
N = input('Maximum Step: ')

# Converting x0 and e to float


x0 = float(x0)
e = float(e)

# Converting N to integer
N = int(N)

#Note: You can combine above three section like


this
# x0 = float(input('Enter Guess: '))
# e = float(input('Tolerable Error: '))
# N = int(input('Maximum Step: '))

# Starting Newton Raphson Method


newtonRaphson(x0,e,N)

Secant method

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

# Implementing Secant Method

def secant(x0,x1,e,N):
print('\n\n*** SECANT METHOD IMPLEMENTATION
***')
step = 1
condition = True
while condition:
if f(x0) == f(x1):
print('Divide by zero error!')
break

x2 = x0 - (x1-x0)*f(x0)/( f(x1) - f(x0) )


print('Iteration-%d, x2 = %0.6f and f(x2) =
%0.6f' % (step, x2, f(x2)))
x0 = x1
x1 = x2
step = step + 1

if step > N:
print('Not Convergent!')
break

condition = abs(f(x2)) > e


print('\n Required root is: %0.8f' % x2)

# Input Section
x0 = input('Enter First Guess: ')
x1 = input('Enter Second Guess: ')
e = input('Tolerable Error: ')
N = input('Maximum Step: ')

# Converting x0 and e to float


x0 = float(x0)
x1 = float(x1)
e = float(e)

# Converting N to integer
N = int(N)

#Note: You can combine above three section like


this
# x0 = float(input('Enter First Guess: '))
# x1 = float(input('Enter Second Guess: '))
# e = float(input('Tolerable Error: '))
# N = int(input('Maximum Step: '))

# Starting Secant Method


secant(x0,x1,e,N)
Forward interpolation
n = int(input('Enter number of data points: '))

# Making numpy array of n & n x n size and


initializing
# to zero for storing x and y value along with
differences of y
x = [Link]((n))
y = [Link]((n,n))

# Reading data points


print('Enter data for x and y: ')
for i in range(n):
x[i] = float(input( 'x['+str(i)+']='))
y[i][0] = float(input( 'y['+str(i)+']='))

# Generating forward difference table


for i in range(1,n):
for j in range(0,n-i):
y[j][i] = y[j+1][i-1] - y[j][i-1]

print('\nFORWARD DIFFERENCE TABLE\n');

for i in range(0,n):
print('%0.2f' %(x[i]), end='')
for j in range(0, n-i):
print('\t\t%0.2f' %(y[i][j]), end='')
print()
Backward interpolation
import numpy as np
import sys

# Reading number of unknowns


n = int(input('Enter number of data points: '))

# Making numpy array of n & n x n size and


initializing
# to zero for storing x and y value along with
differences of y
x = [Link]((n))
y = [Link]((n,n))

# Reading data points


print('Enter data for x and y: ')
for i in range(n):
x[i] = float(input( 'x['+str(i)+']='))
y[i][0] = float(input( 'y['+str(i)+']='))

# Generating backward difference table


for i in range(1,n):
for j in range(n-1,i-2,-1):
y[j][i] = y[j][i-1] - y[j-1][i-1]

print('\nBACKWARD DIFFERENCE TABLE\n');

for i in range(0,n):
print('%0.2f' %(x[i]), end='')
for j in range(0, i+1):
print('\t%0.2f' %(y[i][j]), end='')
print()

Linear interpolation
print('Enter first point:')
x0 = float(input('x0 = '))
y0 = float(input('y0 = '))

# Reading second point


print('Enter second point:')
x1 = float(input('x1 = '))
y1 = float(input('y1 = '))

# Reading calculation point


xp = float(input('Enter calculation point xp: '))

# Calculating interpolated value


yp = y0 + ((y1-y0)/(x1-x0)) * (xp - x0)

# Displaying result
print('Interpolated value at %0.4f is %0.4f'
%(xp,yp))

Lagrange interpolation

import numpy as np

# Reading number of unknowns


n = int(input('Enter number of data points: '))
# Making numpy array of n & n x n size and
initializing
# to zero for storing x and y value along with
differences of y
x = [Link]((n))
y = [Link]((n))

# Reading data points


print('Enter data for x and y: ')
for i in range(n):
x[i] = float(input( 'x['+str(i)+']='))
y[i] = float(input( 'y['+str(i)+']='))

# Reading interpolation point


xp = float(input('Enter interpolation point: '))

# Set interpolated value initially to zero


yp = 0

# Implementing Lagrange Interpolation


for i in range(n):

p = 1

for j in range(n):
if i != j:
p = p * (xp - x[j])/(x[i] - x[j])

yp = yp + p * y[i]
# Displaying output
print('Interpolated value at %.3f is %.3f.' % (xp,
yp))

Trapezoidal method
def f(x):
return 1/(1 + x**2)

# Implementing trapezoidal method


def trapezoidal(x0,xn,n):
# calculating step size
h = (xn - x0) / n

# Finding sum
integration = f(x0) + f(xn)

for i in range(1,n):
k = x0 + i*h
integration = integration + 2 * f(k)

# Finding final integration value


integration = integration * h/2

return integration

# Input section
lower_limit = float(input("Enter lower limit of
integration: "))
upper_limit = float(input("Enter upper limit of
integration: "))
sub_interval = int(input("Enter number of sub
intervals: "))

# Call trapezoidal() method and get result


result = trapezoidal(lower_limit, upper_limit,
sub_interval)
print("Integrationresult by Trapezoidal method is:
%0.6f" % (result) )

Simpson 1/3

def f(x):
return 1/(1 + x**2)

# Implementing Simpson's 1/3


def simpson13(x0,xn,n):
# calculating step size
h = (xn - x0) / n

# Finding sum
integration = f(x0) + f(xn)

for i in range(1,n):
k = x0 + i*h

if i%2 == 0:
integration = integration + 2 * f(k)
else:
integration = integration + 4 * f(k)

# Finding final integration value


integration = integration * h/3
return integration

# Input section
lower_limit = float(input("Enter lower limit of
integration: "))
upper_limit = float(input("Enter upper limit of
integration: "))
sub_interval = int(input("Enter number of sub
intervals: "))

# Call trapezoidal() method and get result


result = simpson13(lower_limit, upper_limit,
sub_interval)
print("Integrationresult by Simpson's 1/3 method
is: %0.6f" % (result) )

Simpson3/8

def f(x):
return 1/(1 + x**2)

# Implementing Simpson's 3/8


def simpson38(x0,xn,n):
# calculating step size
h = (xn - x0) / n

# Finding sum
integration = f(x0) + f(xn)

for i in range(1,n):
k = x0 + i*h
if i%3 == 0:
integration = integration + 2 * f(k)
else:
integration = integration + 3 * f(k)

# Finding final integration value


integration = integration * 3 * h / 8

return integration

# Input section
lower_limit = float(input("Enter lower limit of
integration: "))
upper_limit = float(input("Enter upper limit of
integration: "))
sub_interval = int(input("Enter number of sub
intervals: "))

# Call trapezoidal() method and get result


result = simpson38(lower_limit, upper_limit,
sub_interval)
print("Integrationresult by Simpson's 3/8 method
is: %0.6f" % (result) )
Eulers
def f(x,y):
return x+y

# or
# f = lambda x: x+y

# Euler method
def euler(x0,y0,xn,n):
# Calculating step size
h = (xn-x0)/n

print('\n-----------SOLUTION-----------')
print('------------------------------')
print('x0\ty0\tslope\tyn')
print('------------------------------')
for i in range(n):
slope = f(x0, y0)
yn = y0 + h * slope
print('%.4f\t%.4f\t%0.4f\t%.4f'%
(x0,y0,slope,yn) )
print('------------------------------')
y0 = yn
x0 = x0+h

print('\nAt x=%.4f, y=%.4f' %(xn,yn))

# Inputs
print('Enter initial conditions:')
x0 = float(input('x0 = '))
y0 = float(input('y0 = '))

print('Enter calculation point: ')


xn = float(input('xn = '))

print('Enter number of steps:')


step = int(input('Number of steps = '))

# Euler method call


euler(x0,y0,xn,step)

Rk4
def f(x,y):
return x+y

# or
# f = lambda x: x+y

# RK-4 method
def rk4(x0,y0,xn,n):

# Calculating step size


h = (xn-x0)/n

print('\n--------SOLUTION--------')
print('-------------------------')
print('x0\ty0\tyn')
print('-------------------------')
for i in range(n):
k1 = h * (f(x0, y0))
k2 = h * (f((x0+h/2), (y0+k1/2)))
k3 = h * (f((x0+h/2), (y0+k2/2)))
k4 = h * (f((x0+h), (y0+k3)))
k = (k1+2*k2+2*k3+k4)/6
yn = y0 + k
print('%.4f\t%.4f\t%.4f'% (x0,y0,yn) )
print('-------------------------')
y0 = yn
x0 = x0+h
print('\nAt x=%.4f, y=%.4f' %(xn,yn))

# Inputs
print('Enter initial conditions:')
x0 = float(input('x0 = '))
y0 = float(input('y0 = '))

print('Enter calculation point: ')


xn = float(input('xn = '))

print('Enter number of steps:')


step = int(input('Number of steps = '))

# RK4 method call


rk4(x0,y0,xn,step)

You might also like