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

AIM:-

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.

NEWTON RAPHSON METHOD:-

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.

f(p,h)=p3(1−β2)+(0.4hβ2−σh2r2)p2+σ2h43r4p−(σh23r2)3

f'(p,h)=3p2(1−β2)+(0.4hβ2−σh2r2)2p+σ2h43r4

In numerical analysis, Newton's method, also known as the Newton–Raphson method, named after Issac Newton and Joseph Raphson, is a root-finding algorithm that 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 x1=x0-{f(x0)f′(x0)}" role="presentation">x1=x0−{f(x0)f'(x0)}x1=x0-{f(x0)f′(x0)} 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 xn+1=xn-{f(xn)f′(xn)}" role="presentation">xn+1=xn−{f(xn)f'(xn)}xn+1=xn-{f(xn)f′(xn)} until a sufficiently precise value is reached. This algorithm is first in the class of Householder's method, succeeded by Halley's method. The method can also be extended to complex functions and to systems of equations.

CODE :-

# TO CALCULATE VALUES OF MINIMUM PRESSURE REQUIRED TO BREAK AN ICE WEDGE OF GIVEN THICKNESS

# USING NEWTON RAPHSON METHOD
import matplotlib.pyplot as plt         # Library for making plots 
             
# Given Inputs :-  

# Width of the ice wedge
beta = 0.5
# Size of the air cushion
r = 40
# Tensile Strength of the ice (pounds per square inch)
sigma = 150
# Tensile Strength of the ice (pounds per square feet)
Sigma = 150*144
# Thickness of the ice field
h = [0.6, 1.2, 1.8, 2.4, 3.0, 3.6, 4.2] 

# Defining function of the ice breaking equation
def function(p, h, beta, r, Sigma):
    
    term1 = pow(p, 3)*(1-pow(beta, 2))
    term2 = 0.4*h*pow(beta, 2)-(Sigma*pow((h/r), 2))
    term3 = pow(Sigma, 2)*pow((h/r), 4)*p
    term4 = Sigma*pow((h/r), 2)
    
    E = term1 + term2*pow(p, 2) + term3/3 - pow((term4/3), 3)
    
    return E

# Defining Derivative of the above ice breaking equation
def function_prime(p, h, beta, r, Sigma):
    
    term1 = 3*pow(p, 2)*(1-pow(beta, 2))
    term2 = 0.4*h*pow(beta, 2)-(Sigma*pow((h/r), 2))
    term3 = pow(Sigma, 2)*pow((h/r), 4)
    
    E1 = term1 + 2*p*term2 + term3/3
    
    return E1   
# USING NEWTON RAPHSON METHOD
# Guessing the pressure value
p1 = 30
# Considerign Relaxation factor (alpha)
alpha = 1
# Assigning Tolerance factor 
Tol = 1e-4
# Countor variables for incrementing the number of iterations
ct = 1
# Assigning empty array of cushion pressure applied(pounds per square feet)
c_p = []
# Assigning empty array of cushion pressure applied(pounds per square inch)
C_p = []
# Assigning empty array of the number of iterations
n_i = []
# Loop to calculate all minimum pressure at each value of thickness
for i in h:
    while(abs(function(p1, i, beta, r, Sigma)) > Tol):
        p1 = p1 - alpha*(function(p1, i, beta, r, Sigma)/function_prime(p1, i, beta, r, Sigma))
        ct = ct+1
        
    # Storing all cushion pressure for each values of the thickness of the ice
    c_p.append(p1)
    # Storing number of iteration used
    n_i.append(ct)
    
# Loop to calculate all minimum pressure at each value of thickness
for i in h:
    while(abs(function(p1, i, beta, r, sigma)) > Tol):
        p1 = p1 - alpha*(function(p1, i, beta, r, sigma)/function_prime(p1, i, beta, r, sigma))
        ct = ct+1
        
    # Storing all cushion pressure for each values of the thickness of the ice
    C_p.append(p1)    
    
# Plotting the minimum cushion pressure versus the isce wedge thickness
# Creating emptly figure window
plt.figure()
# Creating plot
plt.plot(h, c_p, color = 'b', marker = 'o', markerfacecolor = 'r')
# Giving title to our plot
plt.title(' Plotting the minimum cushion pressure versus the ice wedge thickness')
# Giving labels
plt.xlabel('Ice Wedge Thickness(feet)')
plt.ylabel('Minimum Cushion Pressure(pound per square feet)')
# Saving the plot
plt.savefig('Plot of ice wedge thickness vs minimum cushion pressure.png')
plt.show()
print('--------------------------------------------------------------------------------------------------------------------------------------')
print('Thickness of the Ice wedge','|','Minimum Cushion Pressure(psft)','|','Minimum Cushion Pressure(psi)','|','No. of Iterations')
print('--------------------------------------------------------------------------------------------------------------------------------------')
for j in range(0,7):
    print('              ',h[j],'                     ',c_p[j],'          ',C_p[j],'           ',n_i[j])
print('--------------------------------------------------------------------------------------------------------------------------------------')
    
 

STEPS:-

1. After importing the necessary libraries, we will define the functions of an polynomial and its derivatives. and return it to the fuction.
2. Define the value of inputs: β = 0.5, r = 40 feet and sigma= 150 pounds per square inch (psi) or Sigma= 150*144 pounds per square feet(psft).
3. Now, we will use newton raphson method to calculate desired output and assign all the inputs.
4. Considering initial pressure value was taken as 30 and also assign relaxation factor and tolerance factor and also assign the empty array of cushion pressure and no. of iterations.
5. A for loop was used to get the values of cushion pressure corresponding to different height , used to find the new p1 after each iteration. Alpha was assumed to 1 in this case and then incrementing countor variable to count all Iterations
6. Now we will make one more loop same as above to calculate cushion pressure for the pounds per square inch.
7. Plotting of desired output and printing our desired output table.

ERRORS:-

The above error is developed due import error in the tabulate library and it can corrected by the installing tabulate module.

OUTPUT:-

Result Table:-

Graph:-

Conclusion:-

We can see that above result depict that we can an minimum cushion pressure applied to the ice wedge for breaking purpose by considering constant relaxation factor is at the 11th iteration which is nearly zero.

So from these results we can conclude that the minimum pressure will depend on the minimium thickness of the ice field.And for a unity relaxation factor it seems that the deterministic values are varying gradually and its under-relaxing.

Now we will find Optimum relaxation factor for this case by plotting number of Iterations and Relaxation factor:-

CODE :-

# TO CALCULATE VALUES OF OPTIMUM RELAXATION FACTOR 
# USING NEWTON RAPHSON METHOD
import matplotlib.pyplot as plt   # Library for making plots                
import numpy as np
# Given Inputs :- 
# Width of the ice wedge
beta = 0.5
# Size of the air cushion
r = 40
# Tensile Strength of the ice (pounds per square feet)
sigma = 150*144
# Thickness of the ice field
h = 0.6

# Defining function of the ice breaking equation
def function(p, h, beta, r, sigma):
    
    term1 = pow(p, 3)*(1-pow(beta, 2))
    term2 = 0.4*h*pow(beta, 2)-(sigma*pow((h/r), 2))
    term3 = pow(sigma, 2)*pow((h/r), 4)*p
    term4 = sigma*pow((h/r), 2)
    
    E = term1 + term2*pow(p, 2) + term3/3 - pow((term4/3), 3)
    
    return E

# Defining Derivative of the above ice breaking equation
def function_prime(p, h, beta, r, sigma):
    
    term1 = 3*pow(p, 2)*(1-pow(beta, 2))
    term2 = 0.4*h*pow(beta, 2)-(sigma*pow((h/r), 2))
    term3 = pow(sigma, 2)*pow((h/r), 4)
    
    E1 = term1 + 2*p*term2 + term3/3
    
    return E1   
 
# USING NEWTON RAPHSON METHOD
# Guessing the pressure value
p1 = 30
# Considerig Relaxation factor (alpha)
alpha = np.linspace(0.1,1.9,19)
# Assigning Tolerance factor 
Tol = 1e-4
# Countor variables for incrementing the number of iterations
ct = 1
# Assigning empty array of cushion pressure applied
c_p = []
# Assigning empty array of the number of iterations
n_i = []
 

# Loop to calculate all minimum pressure at each value of thickness
for i in range(len(alpha)):
    while(abs(function(p1, i, beta, r, sigma)) > Tol):
        p1 = p1 - alpha[i]*(function(p1, i, beta, r, sigma)/function_prime(p1, i, beta, r, sigma))
        ct = ct+1
        
    # Storing all cushion pressure for each values of the thickness of the ice
    c_p.append(p1)
    # Storing number of iteration used
    n_i.append(ct)
    
    
# Plotting the minimum cushion pressure versus the isce wedge thickness
# Creating emptly figure window
plt.figure()
# Creating plot
plt.plot(alpha, n_i)
# Giving labels
plt.xlabel('Relaxation factor')
plt.ylabel('No. of Iteration occurs')
plt.show()
 

STEPS:-

1. After importing the necessary libraries, we will define the functions of an polynomial and its derivatives. and return it to the fuctions.
2. Define the value of inputs: β = 0.5, r = 40 feet and sigma= 150*144 pounds per square feet(psft).
3. Now, we will use newton raphson method to calculate desired output and assign all the inputs.
4. Considering initial pressure value was taken as 30 and also assign relaxation factor and tolerance factor and also assign the empty array of cushion pressure and no. of iterations and alpha range was inputed with the help of n numpy module using linspace command
5. A for loop runs to the length of the alpha was used to get the values of cushion pressure corresponding to different Relaxation factor , used to find the new p1 after each iteration. heigh of the iced wedge was assumed to 0.6 in this case and then incrementing countor variable to count all Iterations.
6. Plotting of desired output.

RESULT:-

The above graph describes the relationship between the relaxation factor relative to its corresponding number of iterations.

The curve depicted in the graph is an exponentially decreasing curve. From the graph we can understand that the no. of iterations is very high for lower values of the relaxation factor,therefore it is acting as an under-relaxation entity(,i.e, it is taking very small interval values for the iterations), but as the value of the relaxation factor increases there is a rapid drop in the no. of iterations(,i.e the values are taken with larger interval gap) and as the relaxation factor reaches an approximate value of above than 1 there is a slight increase again in the no. of iterations,this is indicating the overthrow nature of the iteration depending on the alpha value.

Therefore for this case the optimum value for the relaxation factor is 1 at the no. of iteration be 8.