Week 3 - Solving second order ODEs
Aim:
To write a program to simulate the transient behaviour of a simple pendulum and to create an animation of it's motion.
In the above equation,
g = gravity in m/s2,
L = length of the pendulum in m,
m = mass of the ball in kg,
b=damping coefficient.
(d^2 θ/dt^2) +(dθ/dt)*b/m+(g/l)*sinθ=0
θ=θ1
(d^2 θ1/dt^2 )+(dθ1/dt)*b/m+(g/l)*sinθ1=0
dθ1/dt= θ2
d^2 θ1/dt^2=dθ2/dt
θ2/dt+(b/m)*θ2+(g/l)*sinθ1=0
θ2/dt=-(b/m)*θ2-(g/l)*sinθ1
dθ1/dt=θ2
Code for simple pendulum
import math
import numpy as n1
import matplotlib.pyplot as plt
from scipy.integrate import odeint
def model(initial,t,b,l,g,m):
t1=initial[0]
t2=initial[1]
dt1=t2
dt2=(-b*t2/m)-(g*math.sin(t1)/l)
dt=[dt1,dt2]
return dt
b=0.05
l=1
g=9.81
m=1
initial=[0,3]
t=n1.linspace(0,20,200)
disvelo=odeint(model,initial,t,args=(b,l,g,m))
plt.plot(t,disvelo[:,0],'-',label='displacement')
plt.plot(t,disvelo[:,1],'--',label='velocity')
plt.legend()
plt.xlabel('time')
plt.ylabel('displacementvelocity')
plt.show()
y=disvelo[:,0]
ct=1
for ang in y :
x1=0;y1=0
x2=l*math.sin(ang)
y2=l*math.cos(ang)
file='test%05d.png'%ct
plt.figure()
plt.plot([-1.5,1.5],[0,0],linewidth=8)
plt.plot([x1,x2],[y1,y2],linewidth=3)
plt.plot(x2,y2,'o',markersize=15)
plt.xlim([-1.5,1.5])
plt.ylim([-1.5,1.5])
plt.savefig(file)
ct=ct+1
I have used odeint command to solve the 2 first order ODE equation in the function.The command have common syntax.It solve the 2 first order ODE equation and gives the results.the results are plotted.For animation of simple pendulum,we have separate the displacement value from the solution.For different displacement value the string position are saved as a image file.All images are stored in separate files.All images changed into animation using image magick software
Output graph