APPLICATION OF NEWTON-RAPHSON METHOD TO OBTAIN THE MINIMUM PRESSURE REQUIRED TO BREAK THE ICE OF GIVEN THICKNESS
OBJECTIVE
To find the minimum pressure required to break the ice of given thickness using cushion vehicle.
INTRODUCTION
When the ice pack in one section of the river, it hinders the flow of water and causing floods above the packed ice. Self-powered air cushion vehicles can be used successfully to break it. An experiment with air cushion platforms tied to the front of ships near Thunder Bay, Ontario, Canada, has demonstrated that an air cushion platform can be a relatively efficient means of ice breaking compared to conventional displacement vessels with their large power and crew requirements. The Alexander Henry, a coast guard vessel, has ice-breaking capabilities of two knots in twelve inches of ice. When fitted with an air cushion platform across its bow, the vessel was able to move through ice seventeen inches thick at a speed of nine knots, through twenty-four-inch ice at five knots, and thirty-inch ice at two knots. When fitted with the air cushion platform the fuel consumption of the Alexander Henry was drastically reduced.
The question we wish to explore is how the cushion pressure of the air cushion vehicle affect its icebreaking capacities.
THEORY
At low speeds, a cushion of air forms under the ice and it breaks under its own weight, while at high speeds the ACV sets up a standing wave about half the craft length astern in the ice
cover, which moves with the speed of the craft, and breaks the ice at the crest of the wave.
In a low-speed ACV, when a weight(pressure) is placed on a sheet of ice bends and radial cracks develop. These cracks emanate from the position of the weight outwards. If the weight is sufficient, circumferential cracks will then appear. When the bending moment(weight) exceed those which can be borne by the ice, it breaks.
The relationship between bending moment and tensile strength is given by:
`sigma=(3pr^2)/h^2[1-beta^(2/3)(1-(0.4h)/p)^(1/3)]`
`sigma=\" Tensile strength of ice in Psf\"`
`\"p\"=\"The air pressure or cushion pressure\"`
`\"r\"=\"the radius of load or size of air cushion\"`
`\"h\"=\"The thickness of ice field\"`
`beta=sin(b_0/2)/(b_0/2) ,\" \" b_0 \" width of ice-wedge at unit distance\"`
For the minimum cushion pressure required to break a given thickness of ice, we explore p(r, h), and from (1) we obtain a cubic in p, namely
`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`
The most profitable way to solve for p is probably by some iterative procedure. Thus we use the Newton-Raphson method to obtain the root (pressure) of the equation.
Newton–Raphson method 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 x0 for a root of f. If the function satisfies sufficient assumptions and the initial guess is close, then
`x_1=x_0-f(x_0)/(f\'(x_0))`
is a better approximation of the root than x0. Geometrically, (x1, 0) is the intersection of the x-axis and the tangent of the graph of f at (x0, f (x0)): 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-f(x_n)/(f\'(x_n))`
until a sufficiently precise value is reached.
Relaxation factor is used to modify the step size. If it is greater than one then it is called over relaxation, it amplifies the step size. If the relaxation factor is less than one then it is called under relaxation, it reduces the step size.
CODE TO GET THE MINIMUM PRESSURE:
import math
import numpy as np
import matplotlib.pyplot as plt
b=0.5
sigma=150*144
r=40
def f(p,h):
first_term=pow(p,3)*(1-pow(b,2))
second_term=pow(p,2)*(0.4*h*pow(b,2)-sigma*pow(h,2)/pow(r,2))
third_term=(pow(sigma,2)*pow(h,4)/(3*pow(r,4)))*p
four_dash=sigma*pow(h,2)/(3*pow(r,2))
fourth_term=pow(four_dash,3)
return first_term+second_term+third_term-fourth_term
def f_d(p,h):
first_term=3*pow(p,2)*(1-pow(b,2))
second_term=2*p*(0.4*h*pow(b,2)-sigma*pow(h,2)/pow(r,2))
third_term=(pow(sigma,2)*pow(h,4)/(3*pow(r,4)))
return first_term+second_term+third_term
def pressure(h):
p_guess=80
al=1
#Newton-Raphson algorithm
while (abs(f(p_guess,h))>1e-5):
p_guess=p_guess-al*(f(p_guess,h)/f_d(p_guess,h))
return p_guess
min_press=[]
for i in range(7):
h= [0.6,1.2,1.8,2.4,3.0,3.6,4.2]
min_press.append(pressure(h[i]))
print(\"Thickness(Inches)\\tMinimum Pressure(Psf)\")
for i in range(7):
print(\'{:10.3f}\\t{:24.9f}\'.format(h[i],min_press[i]))
`f(p)=\"first term+second term+third term+fourth term\"`
`\"first term\"=p^3(1-beta^2)`
`\"second term\"=(0.4hbeta^2-(sigmah^2)/r^2)p^2`
`\"third term\"=(sigma^2h^4)/(3r^4)p`
`\"four dash\"=(sigmah^2)/(3r^2)`
`\"fourth term\"=\"four_dash\"^3`
`f_d(p)=\"first term+second term+third term\"`
`\"first term\"=3p^2(1-beta^2)`
`\"second term\"=(0.4hbeta^2-(sigmah^2)/r^2)2p`
`\"third term\"=(sigma^2h^4)/(3r^4)`
`p_(n+1)=p_n-f(p)/(\"f_d\"(p))`
We start the process with some arbitrary initial value p_guess. To improve the guess, function f(p) and f_d(p) gives the value of the function and the first derivative. The improved result is stored as a new p_guess. This process is continued until the value of the function f(p) is less than tolerance value.
We calculate the minimum pressure for the values of thickness and store it in \'min_press\'.
import math
import numpy as np
import matplotlib.pyplot as plt
def f(p):
b=0.5
sigma=150*144
h=0.6
r=40
first_term=pow(p,3)*(1-pow(b,2))
second_term=pow(p,2)*(0.4*h*pow(b,2)-sigma*pow(h,2)/pow(r,2))
third_term=(pow(sigma,2)*pow(h,4)/(3*pow(r,4)))*p
four_dash=sigma*pow(h,2)/(3*pow(r,2))
fourth_term=pow(four_dash,3)
return first_term+second_term+third_term-fourth_term
def f_d(p):
b=0.5
sigma=150*144
h=0.6
r=40
first_term=3*pow(p,2)*(1-pow(b,2))
second_term=2*p*(0.4*h*pow(b,2)-sigma*pow(h,2)/pow(r,2))
third_term=(pow(sigma,2)*pow(h,4)/(3*pow(r,4)))
return first_term+second_term+third_term
al=np.linspace(0.1,1.9,41)
sol=[]
no=[]
for n in range(41):
iteration=0
p_guess=10
while (abs(f(p_guess))>1e-9):
iteration=iteration+1
p_guess=p_guess-al[n]*(f(p_guess)/f_d(p_guess))
no.append(iteration)
sol.append(abs(f(p_guess)))
plt.figure()
plt.plot(al,sol)
z=[0]*41
plt.plot(al,z,\'r--\',label=\'zero line\')
plt.legend(loc=\'right\')
plt.xlabel(\'Relaxation Factor\')
plt.ylabel(\'|f(p)|\')
plt.show()
plt.plot(al,no)
plt.xlabel(\'Relaxation Factor\')
plt.ylabel(\'No of Iterations\')
plt.show()We vary the relaxation factor and calculate the value of the function at the root and the number of iterations.
RESULTS AND CONCLUSION
The figure shows the value of minimum pressure for the given thickness.
We can observe that the number of iterations required to reach a solution is minimum for relaxation factor=1. As we move on either side of 1 number of iteration increases.
From the above graph, the solution is accurate when relaxation factor=1.
Considering the number of iteration and accuracy of the solution, the optimum value of the relaxation factor is 1.Thus minimum pressure required for the thickness of 0.6 feet considering a relaxation factor of 1 is 4.239 Psf and number of iterations performed is 13.