Week 3 - Solving second order ODEs
Aim –
To simulate transient behaviour of simple pendulum and create animation of its motion.
Introduction –
Simple pendulum consists of bob attached to light and inextensible string attached to rigid point.
ODE are be used to describe the transient behaviour.
For simulating simple pendulum following ode is used.
Where,
g – Acceleration due to gravity in m/s2,
L – Length of the pendulum in m,
m- Mass of bob in kg,
b- Damping coefficient.
Time for simulation is 20s, for angular displacement is 0,angular velocity is 3 rad/sec at time t =0.
python code -
#programme to simulate transient behaviour of simple pendulum
import math
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import numpy as np
def model(theta,t,b,g,m,l):
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.812
l =1
m =0.1
theta_0 =[0,3]
t =np.linspace(0,20,400)
#solving second order ode
theta =odeint(model,theta_0,t,args=(b,g,l,m))
#plot results
plt.plot(t,theta[:,0],'g',label =r'displacement')
plt.plot(t,theta[:,1],'b',label =r'velocity')
plt.xlabel('plot')
plt.ylabel('time')
plt.legend()
plt.savefig('output.png')
plt.show()
#animation
ct =1
for theta in theta[:,0]:
x0 =0
y0 =0
x1=l*math.sin(theta)
y1=-l*math.cos(theta)
plt.figure()
plt.plot([-1,1],[0,0],"black",linewidth=5)
plt.plot([x0,x1],[y0,y1],"red")
plt.plot(x1,y1,'o',markersize = 8)
plt.xlim([-1,1])
plt.ylim([-1.5,1])
filename =str(ct)+'.png'
plt.savefig(filename)
ct =ct+1
Explaination of code -
importing libraries,defining the functions,descritising the sub vertical angle,solving ode,ploting results,saving images.
Output -
plot-
Animation-
Conclusion -
velocity and displacement amplitudes are decreasing due to damping medium.