SOLVING SECOND ORDER ODE USING PYTHON
1.AIM :
WRITE A PROGRAM TO SOLVE THE ODE REPRESENTING THE MOTION OF SIMPLE PENDULUM USING PYTHON AND ANIMATE THE MOTION
2.THEORY AND GOVERNING EQUATIONS :
A simple pendulum consists of a bob of mass m hanging from a string of length l and fixed at a pivot point. When displaced to an initial angle and released, the pendulum will swing back and forth in simple harmonic motion. By applying Newton's secont law for rotational systems, the ODE representing the equation of motion of a simple pendulum with damping is given by
`(d^ 2 θ)/ (d t^ 2) + b/ m ⋅ (d θ)/ (d t) + g/ l ⋅ sin θ = 0`
where
g = gravity in m/s^2
l = length of the pendulum in m
m = mass of the ball in kg
b=damping coefficient
3.OBJECTIVES :
- To solve the second order ODE using python.
- To plot the angular velocity v/s time and angular displacement v/s time graph.
- To animate the motion of the simple pendulum.
4.PYTHON PROGRAM :
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import math
#function that returns d(theta)/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
#damping coefficient
b=0.05
#gravity
g=9.81
#length of pendulum
l=1
#mass
m=1
#initial condition (at time t=0)
theta_0 = [0,3]
#time points(total time of 0-20sec divided into 150 parts)
t=np.linspace(0,20,150)
#solve ODE
theta = odeint(model,theta_0,t,args = (b,g,l,m))
#plot results
#plotting angular displacement
plt.plot(t,theta[:,0],'b-',label=r'Angular displacement')
#plotting angular velocity
plt.plot(t,theta[:,1],'r--',label=r'Angular velocity')
plt.ylabel('plot')
plt.xlabel('time')
plt.legend(loc='best')
plt.show()
#creating animation
ct = 1
for theta in theta[:,0]:
x0 = 0
y0 = 0 # x0,y0 fixed end of pendulum
x1 = l*math.sin(theta)
y1 = -l*math.cos(theta)
plt.figure()
plt.plot([-0.25,0.25],[0,0],"black",linewidth=2)
plt.plot([x0,x1],[y0,y1],"red")
plt.plot(x1,y1,'o',markersize=20)
plt.xlim([-1,1])
plt.ylim([-1.5,1])
ct = ct +1
filename = 'test%05d.png' % ct
plt.savefig(filename)
EXPLANATION:-
- Initially the required modules are imported using 'import' command.
- Then a function is defined which simplifies the ODE.
- Next, the input for damping coefficient, length, mass, initial displacement and velocity is given.
- The timespan over which the moiton is to be studied is set using linspace function.
- Then 'odeint' command is used to solve the ODE, which returns angular displacement and angular velocity values.
- Angular displacement and Angular velocity is plotted with respect to time.
- Using a for loop, the position of pendulum is calculated at every angle and the plot is saved as image.
- All the saved images are stiched together and converted into a video.
5.RESULTS :
The damping effect can be seen from reducing amplitude of the plot.
Simple pendulum motion animation -
6.CONCLUSION :
- The second order ODE was solved.
- The graph for angular velocity v/s time and angular displacement v/s was time obtained.
- The animation for the movement of simple pendulum was obtained.