Week 3 - Solving second order ODEs

AIM: Solving 2^nd order Ordinary Differential Equations.

Objectives of the Project:

1. Write a program that solves the ordinary differential equation of motion of a simple pendulum with damping.
2. Create an animation file that represents the position of a pendulum.

Intoduction to Differential Equation:

In mathematics, a differential equation is an equation that relates one or more functions and their derivatives. In applications, the functions generally represent physical quantities, the derivatives represent their rates of change, and the differential equation defines a relationship between the two.

In Engineering, ODE is used to describe the transient behavior of a system. An ordinary differential equation (ODE) is a differential equation containing one or more functions of one independent variables and the derivatives of those functions.

A second order differential equation is an equation involving the unknown function y, its derivatives y' and y'', and the variable x. We will only consider explicit differential equations of the form.

A typical approach to solving higher-order ordinary differential equations is to convert them to systems of first-order differential equations, and then solve those systems. The example uses Symbolic Math Toolbox to convert a second-order ODE to a system of first-order ODEs. Then it uses the PYTHON solver odeint to solve the system.

A second order linear differential equation to solve oscillating motion of pendulum. The way the pendulum moves depends on the Newton’s second law of motion. When this law is written down, we get a second ordinary Differential Equation that describes the position of ball w.r.t. time.

Consider a simple pendulum of mass m attached to string of length l fixed to a pivot at one of its end. And other end of string is displaced to an initial angle Θ and released, the pendulum will swing back and forth with periodic motion under gravity.

The Equation is given by:

`(d^2 θ)/dt^2 +b/m(dθ)/dt+g/l sin⁡θ=0`

Basic Fundamental equations are as follows:

Let’s assume, `θ=θ_1`

`(dθ_1)/dt=θ_2`

`(d^2 θ)/dt^2 =(d^2 θ_1)/dt^2 =d/dt ((dθ_1)/dt)=(dθ_2 )/dt`

`(dθ_2)/dt=-b/m θ_2-g/l sin⁡θ_1`

Here Inputs are as follows,

m=Mass of the ball;

l=Length of the string/pendulum;

b=damping coefficient;

g=Acceleration due to gravity;

Θ=Angle of inclination with horizontal.

Body of the Content:

• PYTHON is used to display the code of the program.

# To solve the ordinary differntial equation

# Importing Math library

import math
import matplotlib.pyplot as plt
import numpy as np 
from scipy.integrate import odeint

# Defining a function that returns dtheta/dt

def ode(theta,t,b,g,l,m):

	theta1=theta[0]
	theta2=theta[1]

	dtheta1_dt=theta2
	dtheta2_dt=-(b/m)*theta2 - (g/l)*math.sin(theta1)
	dtheta_dt=[dtheta1_dt,dtheta2_dt]
	return dtheta_dt

b=0.05
g=9.81
l=1
m=1

# Initial condition
theta_0=[0,3]

# time points
t=np.linspace(0,20,300)

# solve ODE
theta=odeint(ode,theta_0,t,args=(b,g,l,m))

#plot results
plt.plot(t,theta[:,0],'r',label='Displacement')
plt.plot(t,theta[:,1],'g--',label='Velocity')
plt.ylabel('Plot')
plt.xlabel('Time')
plt.legend()
plt.show()

# Animation
ct=1
for theta in theta[:,0]:
	x0=0
	y0=0
	x1=l*math.sin(theta)
	y1=-(l*math.cos(theta))

	filename=str(ct)+'.png'
	ct=ct+1

	plt.figure()
	plt.plot([-1.5,1.5],[0,0],linewidth=2)
	plt.plot([x0,x1],[y0,y1],linewidth=4)
	plt.plot(x1,y1,'o',markersize=15)
	plt.xlim([-1.5,1.5])
	plt.ylim([-1.5,1.5])
	plt.xlabel('X-axis')
	plt.ylabel('Y-axis')
	plt.savefig(filename)

Explanation of the above Programing code.

• Firstly, import math and import matplotlib command is used to perform math functions and plotting graph respectively.
• Also, import numpy, from scipy.integrate import odeint command are used to add two array and solving ordinary differential equations resp.
• Now, a function is defined named as ode where inputs are given such as theta,t,b,g,l and m. Under the function, the differential equations are defined which describes the position or displacement of the pendulum.
• Input parameters such as b, g ,l &m are defined and then the intial conditions as theta_0 is assigned in terms of displacement and velocity.
• A time array is also defined from 0 to 20 with 300 different points in between them.
• Then, the time and results (angular position and velocity) are calculated using odeint[which is an inbuilt function which integrates a system of non–stiff ordinary differential equations (non–stiff ODEs and DAEs) using second order Dormand-Prince method. This is a fourth–order accurate integrator. This solver requires six function evaluations per integration step]function which takes input parameters to calculate the output as the differentiation of variables w.r.t. time, initial conditions and time array.
• plot command is used to plot the graph of displacement and velocity. Also, plt.xlabel and plt.ylabel are used to give x and y label to axis.
• For creating an animation file, for loop is initiated for first and subsequent values of angular postion and velocity. Under this, initial value for x0 and yo are assigned and with the help of equation of x1 and y1,all values of results are calculated.
• The variable ct consider the count of each frame in the interval of substended angle.
• Now, plotting the position of fixed support and also pendulum string with the help of plt.plot command. Bob of pendulum is plotted with markersize’15’.
• The plt.xlim and plt.ylim commands are used to allot upper and lower limit to x and y axis for plotting the area under which pendulum is oscilating.
• To save the file, filename is created consist of string name with extension ‘.png’. plt.savefig command is used to save the obtained figure.
• The variable ct=ct+1 is used to give increment of 1 for each frame.
• After that, I have installed image magic software to convert all the images /plot into animation file .This can be done with the help of command prompt window(cmd.exe). First,I have find the path of directory where the plots/images are saved and later by using magick *.png output.gif command ,I have created an animation file successfully.

Results:

• Obtained values from command window:

• Plot of graph:

Result from Command Prompt window:

Result of Animation file:

The Simulation of motion of an oscillating pendulum is shown in below link:https://youtu.be/LYFUoPes2AU

Conclusion:

• The transient behavior of a pendulum system is described with the help of ordernary differential equation.
• Plot of graph of angular velocity and angular position vs time shows the visual representation of velocity and postion w.r.t. time.
• For the above simulation of an oscillating pendulum an animation file is created successfully.