Project

Breaking Ice with Air cushion Vehicle - Find minimum pressure with Newton-Raphson method using Python

GGurwinder Singh15 Jan 2020

Objective: -To determining the minimum cushion pressure needed to break a given thickness of ice using an air cushion vehicle by breaking down the equation according to Newton-Raphson method. Solve by using Newton Raphson method and Write a Python code: -

Finding minimum pressure by taking relaxation factor as 1 and h=0.6.
Print pressure, h values and no if iteration in a table format.
Find minimum relaxation factor using graph.

Theory: - produces successively better approximations to the roots (or zeroes) of a real- valued function. The most basic version starts with a single-variable function f defined for a real variable x, the function's derivative f ′, and an initial guess x0 for a root of f. If the function satisfies sufficient assumptions and the initial guess is close, then

Equation provided below is break down using Newton-Raphson method

Muller (“Ice Breaking with an Air Cushion Vehicle”) in Mathematical Modeling: Classroom Notes in Applied Mathematics, SIAM 1987) derived the equation.

Where,

p denotes the cushion pressure.
h denotes the thickness of the ice field.
r denotes the size of the air cushion.
σ denotes the tensile strength of the ice.
β is related to the width of the ice wedge

Inputs,

β = 0.5
r = 40 feet
σ = 150 pounds per square inch (psi). (convert to psi feet i.e. 150X144= 21600)
Construct a table with h in the first column, and p in the second column for h = 0.6, 1.2, 1.8, 2.4, 3.0, 3.6 and 4.2 feet

Coding: -

Import module – math, matplotlib, numpy as np and tabulate for creating tables.
For coding we have to define 2 function f(x) i.e. f and f’(x) i.e. fprime. f(x) is the equation given above and f’(x) is the derivative of f(x).
f(x) equation breaks in 4 terms and then create function f and similarly for fprime function.

So, Coding for Finding minimum pressure by taking relaxation factor as 1 and h=0.6.

After import module and define inputs, create f and fprime function and return values from them.
Take relaxation factor i.e. alpha = 1
Tolerance as tol = 1e-4
x_guess = 10 (x_guess – first guess value of pressure which we assume for the iteration process in N-R method).
Use new variable ‘it’ for storing the value of no of iterations.
Then, use while for calling function and do the N-R iteration for the given input.
Use append function for storing values in p and it arrays.
Print values.

import matplotlib.pyplot as plt 
import math
from tabulate import tabulate
import numpy as np 


#inputs
beta = 0.5
r = 40
sigma = 21600 #(150*144 for pound to psi feet)
h = 0.6
p = []

def f(p,h): # breaks equation into 4 parts
	a1 = pow(p,3)*(1-pow(beta,2))
	a2 = (0.4*h*pow(beta,2)-(sigma*pow(h,2)/pow(r,2)))*p**2
	a3 = pow(sigma,2)*pow(h,4)*p/(3*pow(r,4))
	a4 = pow((sigma*h*h/(3*r*r)),3)
	f1 = a1+a2+a3-a4
	return f1

def fprime(p,h): # breaks equation into 3 parts
	b1 = 3*pow(p,2)*(1-beta*beta)
	b2 = 2*p*((0.4*h*beta*beta)-(sigma*h*h/(r*r)))
	b3 = pow(sigma,2)*pow(h,4)/(3*pow(r,4))
	fprime1 = b1+b2+b3
	return fprime1

#Performing N-R iteration
alpha = 1  # relaxtion factor
tol = 1e-4
x_guess = 10
it=[] 
i = 1
while(abs(f(x_guess,h))>tol):
	x_guess = x_guess - alpha*(f(x_guess,h)/fprime(x_guess,h))
	i = i +1
	p.append(x_guess)
	it.append(i)
print(i,'|',x_guess)

Result: - The minimum pressure value is 4.23906 and it is obtained after 8(eight) iterations for the above equation using Newton Raphson method.

So, Coding for Print pressure, h values and no if iteration in a table format.

After import module and define inputs, create f and fprime function and return values from them.
Here, tabulate module uses for creating tables.
Take relaxation factor i.e. alpha = 1
Tolerance as tol = 1e-4
x_guess = 10 (x_guess – first guess value of pressure which we assume for the iteration process in N-R method).
Use new variable ‘it’ for storing the value of no of iterations.
For creating table, create header variable for table and give the specific name to the columns.
Here, the value of h is in array. So, for this for loop is used and then use while loop for calling function and do the N-R iteration for the given input.
Use append function for storing values in p and it arrays.
Use mydata = tuple(zip(it,h,p)) ; A Tuple is a collection of Python objects separated by commas. In some ways a tuple is similar to a list in terms of indexing, nested objects and repetition but a tuple is immutable unlike lists which are mutable. This command is basically used to arrange the values in order.
Then, use the print command to print the values in table format by use of mydata and header variable.

import matplotlib.pyplot as plt 
import math
from tabulate import tabulate
import numpy as np 


#inputs
beta = 0.5
r = 40
sigma = 21600 #(150*144 for pound to psi feet)
h = [0.6,1.2,1.8,2.4,3.0,3.6,4.2]
p = []

def f(p,h): # breaks equation into 4 parts
	a1 = pow(p,3)*(1-pow(beta,2))
	a2 = (0.4*h*pow(beta,2)-(sigma*pow(h,2)/pow(r,2)))*p**2
	a3 = pow(sigma,2)*pow(h,4)*p/(3*pow(r,4))
	a4 = pow((sigma*h*h/(3*r*r)),3)
	f1 = a1+a2+a3-a4
	return f1

def fprime(p,h): # breaks equation into 3 parts
	b1 = 3*pow(p,2)*(1-beta*beta)
	b2 = 2*p*((0.4*h*beta*beta)-(sigma*h*h/(r*r)))
	b3 = pow(sigma,2)*pow(h,4)/(3*pow(r,4))
	fprime1 = b1+b2+b3
	return fprime1

#Performing N-R iteration
alpha = 1
tol = 1e-4
x_guess = 10
it=[]

headers = ["iteration","thickness","pressure"]

for i in h:
	iter = 1
	while(abs(f(x_guess,i))>tol):
		x_guess = x_guess - alpha*(f(x_guess,i)/fprime(x_guess,i))
		iter = iter + 1
	p.append(x_guess)
	it.append(iter)

mydata = tuple(zip(it,h,p))
print(tabulate(mydata,headers=headers,tablefmt='grid'))

Result: - The below result is the table we get after executing the code above code and show us the iteration, h-value(thickness) and the p-value(pressure) related data.

So, Coding for Find minimum relaxation factor using graph.

After import module and define inputs, create f and fprime function and return values from them.
Here, tabulate module uses for creating tables.
Take relaxation factor i.e. alpha = 1
Tolerance as tol = 1e-4
x_guess = 10 (x_guess – first guess value of pressure which we assume for the iteration process in N-R method).
Use new variable ‘it’ for storing the value of no of iterations.
For creating table, create header variable for table and give the specific name to the columns.
Here, the value of h is in array. So, for this for loop is used and then use while loop for calling function and do the N-R iteration for the given input.
Use append function for storing values in p and it arrays.
Use linspace command for alpha value.
Create two variables with empty array.
Then use for loop for different value of alpha find minimum pressure by N-P iteration.
Using append function save J and CT arrays.
Use tuple command to sort the values and the by using sorted command we can find the accurate minimum relaxation factor.
Print values.
Using plot command plot two graph one in between thickness and pressure and relaxation factor and iteration.
Using plt.show() command plot shows.

import matplotlib.pyplot as plt 
import math
from tabulate import tabulate
import numpy as np 

#inputs
beta = 0.5
r = 40
sigma = 21600 #(150*144 for pound to psi feet)
h = [0.6,1.2,1.8,2.4,3.0,3.6,4.2]
p = []

def f(p,h): # breaks equation into 4 parts
	a1 = pow(p,3)*(1-pow(beta,2))
	a2 = (0.4*h*pow(beta,2)-(sigma*pow(h,2)/pow(r,2)))*p**2
	a3 = pow(sigma,2)*pow(h,4)*p/(3*pow(r,4))
	a4 = pow((sigma*h*h/(3*r*r)),3)
	f1 = a1+a2+a3-a4
	return f1

def fprime(p,h): # breaks equation into 3 parts
	b1 = 3*pow(p,2)*(1-beta*beta)
	b2 = 2*p*((0.4*h*beta*beta)-(sigma*h*h/(r*r)))
	b3 = pow(sigma,2)*pow(h,4)/(3*pow(r,4))
	fprime1 = b1+b2+b3
	return fprime1

#Performing N-R iteration
alpha = 1
tol = 1e-4
x_guess = 10
it=[]

for i in h:
	iter = 1
	while(abs(f(x_guess,i))>tol):
		x_guess = x_guess - alpha*(f(x_guess,i)/fprime(x_guess,i))
		iter = iter + 1
	p.append(x_guess)
	it.append(iter)

Alpha = np.linspace(0.01,1.99,200)
J=[]
CT=[]
for j in Alpha:
	x_guess = 10
	ct = 1
	while(abs(f(x_guess,0.6))>tol):
		x_guess = x_guess - (j*f(x_guess,0.6)/fprime(x_guess,0.6))
		ct = ct + 1
	J.append(j)
	CT.append(ct)

W = tuple(zip(CT,Alpha))
w = sorted(W)

print(w[0])
plt.figure()
plt.xlabel("Thickness")
plt.ylabel("Pressue")
plt.plot(h,p)
plt.show()
plt.xlabel("Relaxation Factor")
plt.ylabel("Iteration")
plt.plot(J,CT)
plt.show()

Result: - The minimum relaxation factor for the given equation is 1.054 and at 7 iterations.

Output: -

Plots: - Relaxation factor and Iterations

Plots: - Pressure and Thickness

Conclusion: -

We learn how Newton-Raphson method work using python.
Learn how to create tables in python and plot results in graph.
Understand for and while loop.
Minimum pressure is 4.23906 and it is obtained after 8(eight) iterations at h = 0.6
Create table for pressure values for different given h-values and takin relaxation factor 1
Plot graph between thickness and pressure and relaxation factor and iteration
And also find the minimum relaxation factor.

Breaking Ice with Air cushion Vehicle - Find minimum pressure with Newton-Raphson method using Python

Ship work like this
in your own course.

Talk to a career counsellor

Ship work like thisin your own course.

Talk to a career counsellor

Ship work like this
in your own course.