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

AIM: To determine the minimum cushion pressure needed to break a given thickness of ice using air cushion vehicle by calculating using newton raphson method.

Case (i) relaxation factor 1

Case (ii) height (thickness of ice)= 0.6

INTRODUCTION:

The Newton-Raphson method (also known as Newton's method) is a way to quickly find a good approximation for the root of a real-valued function f(x) = 0f(x)=0. It uses the idea that a continuous and differentiable function can be approximated by a straight line tangent to it.

Suppose you need to find the root of a continuous, differentiable function f(x)f(x), and you know the root you are looking for is near the point x = x_0x=x0​. Then Newton's method tells us that a better approximation for the root is

`x_1 = x_0 - \frac{f(x_0)}{f'(x_0)}.x1​=x0​−f′(x0​)f(x0​)​`

This process may be repeated as many times as necessary to get the desired accuracy. In general, for any xx-value x_nxn​, the next value is given by

`x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)}.xn+1​=xn​−f′(xn​)f(xn​)​`

CASE 1: Alpha, Relaxation factor=1

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

#define muller's equation

def f(p,h):                                                  #b= beta s= sigma
    
    T1=pow(p,3)*(1-pow(b,2))                                
    T2=(0.4*h*pow(b,2)-(s*pow((h/r),2)))*pow(p,2)
    T3=((s**2)*(h**4)*p)/(3*(r**4))
    T4=pow((s*(h**2))/(3*(r**2)),3)
    F= T1 + T2 + T3 - T4
    return F
        
def fprime(p,h):
    T1=3*pow(p,2)*(1-pow(b,2))
    T2=2*p*(0.4*h*pow(b,2)-(s*pow((h/r),2)))
    T3=((s**2)*(h**4))/(3*(r**4))
    Fp= T1 + T2 + T3
    return Fp

#inputs

r= 40     # size of air cushion
b= 0.5     #width of ice 
s= 150*144     #tensile strength of ice 

height= [0.6,1.2,1.8,2.4,3,3.6,4.2]     #thickness of ice

p_guess= 10
tol= 1e-4
alpha= 1     # relaxation parameter
i=1     

#appending pressure values and no of iterations
pressure= []
iterations= []


for h in height :
    while(abs(f(p_guess,h))>tol):
        p_guess= p_guess - alpha*((f(p_guess,h))/(fprime(p_guess,h)))
        i= i + 1
    pressure.append(p_guess)
    iterations.append(i)
print(pressure)
print(iterations)

for i in range(0,7):
        if i==0 :
            print("height : pressure")
        print(height[i], ":" ,pressure[i])

headers= ['thickness','iterations','pressure']
table= tuple(zip(iterations,height,pressure))

#zip is a built-in function that takes two or more sequences
#and "zips" them into a list of tuples where each tuple contains one element rom each sequence.

print(tabulate(table, headers=headers, tablefmt="github"))

#plotting
plt.plot(height,pressure)
plt.xlabel('height in feets')
plt.ylabel('pressure in psf')
plt.show()

OUTPUT:

Values:

CASE 2: h= 0.6

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

#define muller's equation

def f(p,h):                                                  #b= beta s= sigma
    
    T1=pow(p,3)*(1-pow(b,2))                                
    T2=(0.4*h*pow(b,2)-(s*pow((h/r),2)))*pow(p,2)
    T3=((s**2)*(h**4)*p)/(3*(r**4))
    T4=pow((s*(h**2))/(3*(r**2)),3)
    F= T1 + T2 + T3 - T4
    return F
        
def fprime(p,h):
    T1=3*pow(p,2)*(1-pow(b,2))
    T2=2*p*(0.4*h*pow(b,2)-(s*pow((h/r),2)))
    T3=((s**2)*(h**4))/(3*(r**4))
    Fp= T1 + T2 + T3
    return Fp

#inputs

r= 40     # size of air cushion
b= 0.5     #width of ice 
s= 150*144     #tensile strength of ice 
h= 0.6    #thickness of ice

p_guess= 10
tol= 1e-4
alpha= np.linspace(0.01,1.99,200)
a=[]     #relaxation parameter
iterations= []

for j in alpha :
    i=1
    p_guess=10

    while(abs(f(p_guess,h))>tol):
        p_guess= p_guess - j*((f(p_guess,h))/(fprime(p_guess,h)))
        i= i + 1
        
    a.append(j)
    iterations.append(i)

print(iterations)
print(' The minimum value of iterations is:  ',min(iterations))

table= tuple(zip(iterations,a))
S= sorted(table)
print(S[0])
#zip is a built-in function that takes two or more sequences
#and "zips" them into a list of tuples where each tuple contains one element rom each sequence.


#plotting
plt.plot(a,iterations,color="blue",linewidth=5)
plt.xlabel('alpha',fontsize='large',fontweight='bold')
plt.ylabel('no. of iterations',fontsize='large',fontweight='bold')
plt.show()

OUTPUT:

Values:

CONCLUSION:

From the above result we can conclude that the minimum no of iterations required is 7 and the value of optimum relaxation factor is 1.0547