BREAKING ICE WITH AIR CUSHION VEHICLE - FIND MINIMUM PRESSURE WITH NEWTON RAPHSON METHOD USING PYTHON
1.AIM :
TO WRITE A PYTHON PROGRAM TO FIND MINIMUM PRESSURE REQUIRED TO BREAK ICE WITH AN AIR CUSHION VEHICLE USING NEWTON RAPHSON METHOD
2.OBJECTIVE :
- To write a python code to calculate the minimum cushion pressure needed to break a given thickness of ice using an Air Cushion Vehicle by using the Newton Raphson Method for finding the root.
- To calculate the optimum value of the relaxation factor for the thickness of ice to be 0.6 feet.
3.THEORY :
Newton–Raphson method
In numerical analysis, Newton's method, also known as the Newton–Raphson method, named after Isaac Newton and Joseph Raphson, is a root-finding algorithm which 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 `x_0` for a root of f. If the function satisfies sufficient assumptions and the initial guess is close, then
`x_{1}=x_{0}-{frac {f(x_{0})}{f'(x_{0})}}`
is a better approximation of the root than `x_0`. Geometrically, (`x_1`, 0) is the intersection of the x-axis and the tangent of the graph of f at (`x_0`, f (`x_0`)): that is, the improved guess is the unique root of the linear approximation at the initial point. The process is repeated as
`x_{n+1}=x_{n}-{frac {f(x_{n})}{f'(x_{n})}}`
until a sufficiently precise value is reached. This algorithm is first in the class of Householder's methods, succeeded by Halley's method. The method can also be extended to complex functions and to systems of equations.
Here in this work a python program is written to calculate the minimum cushion pressure needed to break a given thickness of ice using an Air Cushion Vehicle by using the Newton Raphson Method and further an optimum value of the relaxation factor is calculated for the thickness of ice to be 0.6 feet.
The program is divided in 2 cases, in case 1 the cushion pressure is calculated for a list of heights (0.6, 1.2, 1.8, 2.4, 3.0, 3.6, 4.2) and in case 2 a list of relaxation factor ranging between 0 to 2 is used for calculation and the most optimum relaxation factor value is cound.
Ice breaking with an Air Cushion Vehicle is governed by the formula:
`p^3(1-beta^2)+(0.4hbeta^2-(sigmah^2)/r^2)p^2+(sigma^2h^4)/(3r^4)p -((sigmah^2)/(3r^2))^3=0`
where p denotes the cushion pressure, h is the thickness of the ice field, σ the tensile strength of the ice, r is the size of the air cushion and β is related to the width of
the ice wedge.
Values taken in the program are β = 0.5, r = 40 feet and σ = 150 pounds per square inch (psi)
4.PYTHON PROGRAM :
import math
import numpy as np
import matplotlib.pyplot as plt
#inputs
sigma = 150*144
beta = 0.5
r = 40
#defining function for calculating f(p,h) (ice crusher function)
def f(p,h):
return (1-beta**2)*p**3 + (0.4*h*beta**2 - sigma*h**2/r**2)*p**2 + (sigma**2*h**4)*p/(3*r**4) - (sigma *h**2/(3*r**2))**3
#defining function for calculating derivative of f(p,h) function, f_(p,h)
def f_(p,h):
return (3*p**2*(1-beta**2) + 2*p*(0.4*h*beta**2 - sigma *h**2/r**2) + sigma**2*h**4/(3*r**4))
#set of heights inputs
h_range = [0.6,1.2,1.8,2.4,3.0,3.6,4.2]
# relaxation factor
alpha = 1
i = 1
#tolerance
tol = 1e-4
#initial guessed value of pressure
p_guess = 10
P = []
iterations = []
for h in h_range:
while (abs(f(p_guess,h))>tol):
p_guess = p_guess-alpha*f(p_guess,h)/f_(p_guess,h)
i = i + 1
iterations.append(i)
P.append(p_guess)
print(' pressure values for different heights ')
#printing the table
print('------------------------------------------------------------')
print('Height(ft) ','\t','Cushion pressure(psft)','\t','iterations')
print('------------------------------------------------------------')
for i in range(0,7):
print(' ',h_range[i],'\t\t\t\t','{:.5f}'.format(P[i]),'\t\t\t\t',iterations[i])
print('------------------------------------------------------------')
#plotting Pressure and height
plt.plot(h_range,P,color = 'blue')
plt.xlabel('Thickness of ICE (feet)')
plt.ylabel('Cushion pressure ' r'($ \frac{pounds}{sqft} $)')
plt.title('Cushion pressure v/s Ice Thickness')
plt.show()
#CASE2
#list of alphas from 0 to 2 for findig the optimum alpha
alpha_range = np.arange(0.02,2.0,0.04)
iterations =[]
i_min =10000
alpha_best = []
for alpha in alpha_range:
i=1
p_guess=10
while (abs(f(p_guess,0.6))>tol):
p_guess = p_guess-alpha*f(p_guess,0.6)/f_(p_guess,0.6)
i = i + 1
if (i_min > i):
i_min = i
alpha_opt = alpha
iterations.append(i)
print('==========================================================================')
print('optimum relaxation factor (ice thickness = 0.6 ft) ')
print('minimum number of iterations = ',i_min, 'and it occurs at alpha = ','{:.2f}'.format(alpha_opt))
#printing the table
print('--------------------------------------------------------------')
print('Relaxation factor',' ','iterations')
print('--------------------------------------------------------------')
for i in range(0,50):
print(' ','{:.2f}'.format(alpha_range[i]),'\t\t\t\t\t',iterations[i])
i=i+1
print('------------------------------------------------------------')
plt.plot(alpha_range,iterations)
plt.plot(alpha_opt,i_min,marker = '*',color = 'red',label = 'optimum relaxation factor value')
plt.legend()
plt.xlabel('relaxation factor')
plt.ylabel('iterations to converge')
plt.show()Explanation:
- Initially all the necessary libraries are imported and values to the inputs( β = 0.5, r = 40 feet and σ = 150 pounds per square inch (psi) or σ = 150*144 pounds per square feet) is assigned.
- Then, the governing equation and its derivative is defined as function.
`f(p,h)=p^3(1-beta^2)+(0.4hbeta^2-(sigmah^2)/r^2)p^2+(sigma^2h^4)/(3r^4)p -((sigmah^2)/(3r^2))^3`
`f'(p,h)=3p^2(1-beta^2)+(0.4hbeta^2-(sigmah^2)/r^2)2p+(sigma^2h^4)/(3r^4)`
- For case 1, a list of heights=[0.6,1.2,1.8,2.4,3.0,3.6,4.2] is defined, initial p_guess (guess value of pressure) value is taken to be 10 and relaxation factor is taken to be 1.
- Using for loop the values of cushion pressure corresponding to different heights is calculated, the formula `p_(guess) = p_(guess) - alpha(f(p_(guess),h)/(f'(p_(guess), h)))`is used to find the new p_guess after each iteration.
- Then the value of cushion pressure corresponding to each height is displayed.
- In case 2, the value of optimum relaxation factor is found for a height of 0.6.
- A list of relaxation factor(alpha) values ranging from 0 to 2 is defined.
- And the no of iterations needed to converge for each value is calculated using for loop and stored in an array.
- The alpha for which the no of iterations is minimum is the optimum relaxation factor.
- A graph is plotted to understand the variation of iterations with relaxation factor.
5.RESULTS :
Case 1 result-
Value of cushion pressure corresponding to height
This result can be plotted in a graph as shown below
As it can be seen from the trend, the cushion pressure increases with the increase in height (thickness of ice).
Case 2 results-
This result can be plotted in a graph as shown below
The optimum value of relaxation factor is 1.02 and it gives result in 8 iterations
6.CONCLUSIONS :
- The python program was written and executed successfully.
- The value of cushion pressure corresponding to the thickness of ice was calculated.
- The optimum value of relaxation factor is 1.02 and it gives result in 8 iterations
7.REFERENCES :