Week 3 - Solving second order ODEs

Aim: To simulate the simple pendulum motion

Objective:

• To write a program that solves ordinary differential equation
• Create an animation showing the motion of pendulum

Problem Statement:

Equation,

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

where,

b = damping coefficient

m = mass of the ball, kg

g = gravity, m/s²

l = length of the pendulum, m

• Simple pendulum motion depends on Newton’s second law, which says ‘Total force being applied to a moving body represents the motion of the body’
• The key factor in this governing equation is Torque, and it comes from gravity, since no additional force is applied to pendulum
• Here the governing equation is of second order ODE, to solve this differential equation we use “odeint” function from “scipy” module
• But for “odeint” to work, we must break original second order ODE into two first order ODE’s
• “odeint” solves the equation of motion of a simple pendulum with damping, and provides the results

Solution Procedure:

• Import necessary modules such as “numpy”, “matplotlib”, “math”, and odeint from “scipy”
• Define a function to solve two first order equations of motion and return the solution values
• Define given inputs such as damping coefficient(b), gravity(g), length of pendulum(l), mass of the ball(m) and Initial conditions such as of angular displacement and angular velocity
• Then we define time span, which gives the range to integrate the differential equation with provided initial conditions
• Then we use “odeint” function to solve the differential equation, in which we make a function call of the function named “model” and obtain the solution results which are angular displacement and angular velocity

• Next, we plot the solution results, showing the wave curves of angular displacement and angular velocity against time
• And then we plot another series of figures using the ‘for loop’ for the purpose of animation of pendulum motion
• Then “magick” tool is used inside command prompt to make animation out of the generated figures from python

Code:

# Solving 2nd order ODE 

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

# function that returns dz/dt
def model(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    # Damping coefficient
g = 9.81    # Gravity in m/s^2
l = 1       # Length of the pendulum in m
m = 1       # Mass of the ball in kg

# initial conditions
theta_0 = [0,3]

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

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

# plot results
plt.plot(t,theta[:,0],'b-',label=r'$\frac{d\theta_1}{dt}=\theta2$')
plt.plot(t,theta[:,1],'r--',label=r'$\frac{d\theta2}{dt}=-\frac{b}{m}\theta_2-\frac{g}{L}in\theta_1$')
plt.xlabel('time')
plt.ylabel('Plot')
plt.legend(loc='best')
plt.show()


# Pendulum simulation
# Angular displacements(theta_d)
theta_d = []
for displacements in theta[:,0]:
	theta_d.append(displacements)

# Counter,ct
ct = 1
for ang_rad in theta_d:
	x0 = 0
	y0 = 0

	x1 = l*math.sin(ang_rad)
	y1 =-l*math.cos(ang_rad)

	filename = str(ct) + '.png'
	ct = ct+1
	plt.figure()
	plt.plot([-0.5,0.5],[0,0])
	plt.plot([x0,x1],[y0,y1])
	plt.plot(x1,y1,'-o')
	plt.xlim(-1.5,1.5)
	plt.ylim(-1.5,1)
	plt.title('pendulum simulation')
	plt.savefig(filename)
		

Results:

• Angular displacement and Angular Velocity against time plot
  

• Animation of the simple pendulum motion

The above video shows the animation of simple pendulum motion which was created by solving the governing equation of motion which is second order ODE, and the function “odeint” is used to solve it and the solutions were plotted and animation was made using the figures. We can see that the pendulum is slowing down with time as the damping effect is considered in the program, which means after some time the pendulum will stop oscillating.

Youtube link:- https://youtu.be/tj9Ni9HGwJk

Discussion: The graph of Angular velocity and Angular displacement against time has been plotted, and we can see the damping effect as the values of angular displacement and angular velocity keeps decreasing. We have also created a animation showing the simple pendulum motion, from which we can observe that oscillation is decreasing and the back and forth motion slowing down over time. In this program we learnt solving the second order ODE using “odeint” function and creating an animation using “magick” tool in command prompt.