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

The Newton Raphson method is an iterative technique that is used to find out the roots of a real valued function. It uses the idea that a continuous and differentiable function can be approximated by a straight-line tangent to it.

In this case Newton Raphson method is used to find out the cushion pressure required to cut a certain thickness of ice using the air cushion vehicle.

An Air cushion vehicle is similar to a Hovercraft which is used to break ice layers on the surface of rivers, lakes and even on the land in some countries .This vehicle uses air pressure to create cracks in the ice surface and then the crack propagates in the direction where there is the least resistance.

To create a crack in the ice surface minimum air pressure is required to be exerted by the vehicle. E.R.Muller presented an equation with few variables to calculate the pressure.This equation can be solved by Newton Raphson method.

Muller’s Equation

p3(1-β2) + (0.4hβ2 – σ h2/r2)p2 + (σ2 h4/3r4)p – (σh2 /3r2)


            p = cushion pressure

             β = width of the ice-wedge

             h = thickness of ice field

             σ = tensile strength of ice

             r = size of the air cushion

And the newton Raphson equation is given by,

              Xn+1 =Xn -f(x)/f’(x)


                    Xn+1 = solution of the nth iteration

To increase the accuracy of the obtained solution, a realaxed factor is multiplied to the ratio which gives us,

                 Xn+1 = Xn – a {f(x)-f’(x)}


  • Python program is shown below which solves the Muller’s equation by the Newton Raphson method.
  • All the variables except for the thickness of ice and pressure are assigned certain values.
import matplotlib.pyplot as plt 

def f(p,sigma,beta,r,h):
			return pow(p,3)*(1-pow(beta,2))+(0.4*H*pow(beta,2)-(sigma*pow(H,2)/pow(r,2)))*pow(p,2)+(pow(sigma,2)*pow(H,4)*p/(3*pow(r,4)))-pow((sigma*pow(H,2)/(3*pow(r,2))),2)

def fprime(p,sigma,beta,r,h):
			return 3*pow(p,2)*(1-pow(beta,2))+(0.4*H*pow(beta,2)-(sigma*pow(H,2)/pow(r,2)))*p*2+(pow(sigma,2)*pow(H,4)/(3*pow(r,4)))

sigma =150
beta = 0.5
r = 40
h = [0.6,1.2,1.8,2.4,3.0,3.6,4.2]
tol = 1e-14
iter = 1

#Peforming Newton raphson iteration

p_guess =120
P = []
print(" _______________________")
print("|Thickness |  Pressure |")
for i in range (0,len(h)):
		H = h[i]
		while (abs(f(p_guess,sigma,beta,r,H))>tol):
			p_guess = p_guess - alpha*(f(p_guess,sigma,beta,r,H)/fprime(p_guess,sigma,beta,r,H))
			iter = iter+1
			print("|%f  | %f |"%(H,P[i]))   
			p_guess = 120 

  • The below program shows a range of relaxation factors and optimum one is highlighted.
import matplotlib.pyplot as plt 
import math

def f(p,sigma,beta,r,h):
			return pow(p,3)*(1-pow(beta,2))+(0.4*h*pow(beta,2)-(sigma*pow(h,2)/pow(r,2)))*pow(p,2)+(pow(sigma,2)*pow(h,4)*p/(3*pow(r,4)))-pow((sigma*pow(h,2)/(3*pow(r,2))),2)

def fprime(p,sigma,beta,r,h):
			return 3*pow(p,2)*(1-pow(beta,2))+(0.4*h*pow(beta,2)-(sigma*pow(h,2)/pow(r,2)))*p*2+(pow(sigma,2)*pow(h,4)/(3*pow(r,4)))

sigma = 150
beta = 0.5
r = 40
h = 0.6
tol = 1e-2
#Peforming Newton Raphson Iterartion
alpha = []
i = 0.1
j = 0
		i = i+0.1
		j = j+1
		Iter =[]

		for i in range(0,len(alpha)):
				p_guess = 120
				ite =1 
						p_guess = p_guess-alpha[i]*(f(p_guess,sigma,beta,r,h)/fprime(p_guess,sigma,beta,r,h))
						ite = ite+1
						min_iter = min(Iter)
						for k in range(0,len(Iter)):
									pos = k
									print("The optimum relaxation factor is ",alpha[k])
									plt.legend(['Optimum Relaxation Factor'])
									plt.xlabel('Reduction Factor')
									plt.ylabel('Number of iterations')
									plt.title('No of Iterations v/s RF')

The above graph shows a decreasing trend in the number of iterations after 0, then the trend increases after a certain value of the relaxing factor. So, there must be minimum value for the number of iterations and its corresponding relaxation factor. This value is considered as the optimum relaxation factor because it reduces the calculation time in a calculating system.

Projects by Habishek Umapathy

AIM – To run a pipe flow simulation with an inlet Reynolds number of 100, 1000 and 10000. For each of these cases do the following Place the line probes at 95%, 90%, 85%of the pipe length. Compare the normalized velocity profile at each of these locations. N Read more

OBJECTIVE – We have to write the script which should take column numbers as the input and plot the respective columns as separate images. The plot labels should be extracted from the file. If there is a request for a plot between column 1(crank angle) and column 8 Read more

  OBJECTIVE – To write a program to perform curve fitting.   Curve Fitting is the process of constructing a curve or mathematical functions, which possess the closest proximity to the real series of data. By curve fitting, we can mathematically constr Read more

OBJECTIVE – To write a program to simulate the transient behaviour of a simple pendulum and to create an animation of its motion.   PROCEDURE – In Engineering, ODE is used to describe the transient behavior of a system. A simple example is a pendulum Read more

OBJECTIVE - To write a code for P_V diagram of an Otto cycle and calculate thermal efficiency




  • Initially state variables p1,t1,t3 and gamma a Read more

OBJECTIVE - 1. Write a function that extracts the 14 co-efficient and calculates the enthalpy, entropy and specific heats for all the species in the data file. Calculate the molecular weight of each species and display it in the command window. Plot the Cp, Ent Read more

OBJECTIVE - Write a code in MATLAB to optimise the stalagmite function and find the global maxima of the function. PROCEDURE- The genetic algorithm is a method for solving both constrained and unconstrained optimization problems that is based on natural sele Read more

OBJECTIVE – To write code to fit a linear and cubic polynomial for the Cp data and plot the linear and cubic fit curves along with the raw data points.   PROCEDURE – Curve Fitting is the process of constructing a curve or mathematical functions, whic Read more


The End