Simulating the transient behaviour of a Simple Pendulum using Python

AIM:

Simulate the transient behaviour of a simple pendulum using Python.

OBJECTIVES:

1. To write a program that solves the second order ODE which represents the equation of motion of a simple pendulum with damping.
2. To simulate the motion of the simple pendulum with damping and creating an animation from the plots obtained in the program.

GOVERNING EQUATIONS:

Second order ODE representing the equation of motion of a simple pendulum with damping:

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

In the above equation,

g = gravity in m/s^2

l = length of the pendulum in m

m = mass of the ball in kg

b = damping coefficient

THEORY:

A pendulum is a body suspended from a fixed support so that it swings freely back and forth under the influence of gravity. When a pendulum is displaced sideways from its resting, equilibrium position, it is subject to a restoring force due to gravity that will accelerate it back toward the equilibrium position. When released, the restoring force acting on the pendulum's mass causes it to oscillate about the equilibrium position, swinging back and forth.

Damping suggests that the angular amplitude decreases exponentially over time. This means that the amplitude of the pendulum will go on decreasing with time and eventually the pendulum will stop.

In Engineering, ODE is used to describe the transient behaviour of a system. A simple example is a pendulum. The way in which pendulum moves depends on the Newtons second law. When this law is written down, we get a second order Ordinary Differential Equation that describes the position of the "ball" w.r.t time.

In this project we will simulate the motion of the simple pendulum with damping by solving the second order ODE which represents the equation of motion of a simple pendulum with damping.

In this project we solve the second order ODE by considering the following:

`θ=θ_1`

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

`d/dx [[θ_1],[θ_2 ]]=[[θ_2],[-b/m*θ_2-g/l sin⁡θ_1 ]]`

This is a system of ODEs with two first order ODE in it. We use this to get our solution.

PROGRAM:

Programming language used is Python.

Code to simulate the transient behaviour of a simple pendulum and to create an animation of it's motion using second order Ordinary Differential Equation in Python.

# Program for solving second order ODEs and simulate its motion
import matplotlib.pyplot as plt 
import numpy as np 
from scipy.integrate import odeint 
 
 
def ode_function(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 

# Inputs  
# damping coefficient 
b = 0.05 
# acceleration due to gravity(m/s^2) 
g = 9.81   
# length of string (m)                                                                
l = 1   
# mass of bob (kg)                                                                   
m = 1       
 
# Initial Conditions 
theta_0 = [0,3]  
# time points 
t_span = np.linspace(0,20,100) 
 
# Solve ODE 
 
theta = odeint(ode_function,theta_0,t_span,args = (b,g,l,m)) 
 
# Plotting the graph 
plt.plot(t_span,theta[:,0],'b-',label=r'$\frac{d\theta1}{dt} = \theta2$') 
plt.plot(t_span,theta[:,1],'r--',label=r'$\frac{d\theta2}{dt} = -\frac{b}{m}\theta2 - \frac{g}{l}sin\theta1$') 
plt.xlabel('Time (second)') 
plt.ylabel('Plot') 
plt.title('Damped Harmonic motion of a Simple Pendulum') 
plt.legend(loc='best') 
plt.ylim([-3.5,3.5]) 
plt.xlim([0,20]) 
plt.grid() 
plt.show() 
 
# Iteration using for loop to create animation of motion of pendulum 
 
ct =1 
range = theta[:,0] 

for i in range:  
	x_start = 0                                                             
	y_start = 0  
	x_end =  l*math.sin(i)                                            
	y_end = -l*math.cos(i)  

	# Plotting for animation  
	plt.figure()    
	# Plotting the rigid support through which string of pendulum is attached  
	plt.plot([-0.5,0.5],[0,0],linewidth=8,color='y')  
	# Plotting String of pendulum  
	plt.plot([x_start,x_end],[y_start,y_end],linewidth=3)
	# Plotting Bob of pendulum  
	plt.plot(x_end,y_end,'o',markersize=20,color='r')  
	# Plotting dotted vertical line for refernce  
	plt.plot([0,0],[0,-1.5],'--',color='k')  
	plt.xlim([-1.25,1.25])  
	plt.ylim([-2,1])  
	plt.title('Transient behaviour of a Simple Pendulum')  
	filename = 'test%05d.png' % ct  
	plt.savefig(filename)  
	ct = ct + 1 

EXPLANATION OF THE PROGRAM:

1. Program started by importing modules required. We need math module for mathematical calculations, matplotlib module for plotting results, numpy module to evaluate time intervals, scipy odeint module to solve ODEs.
2. In order to solve ODEs, a function is considered which will return the value of dΘ/dt.

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

Above function consists of theta, time, damping co-effcient, length of pendulum, mass of bob, acceleration due to gravity. As explained in calculation part, Θ1 is displacement function and Θ2 is angular velocity function.

3. Input values are: l = 1m, m = 1kg, b = 0.05, g = 9.81 m/s^2.
4. Time interval are specified with help of numpy module. Time interval of 20 seconds with 200 step size.
5. The function defined in called back in order to calculate theta.

theta = odeint(ode_function,theta_0,t_span,args = (b,g,l,m))

Above function made use of odeint function of scipy module.

6. Result obtained from above function are plotted by using matplotlib module. Result obtained will have 2 columns and multiple rows. First columns are the values of displacement while that of in second column are the values of angular velocity.

therefore, dispacement results are plotted as;

plt.plot(t, theta[:,0]) = it implies theta should consider first column and all rows. It will give response of displacement curve w.r.t time with red color.

plt.plot(t, theta[:,1]) = it implies thata should consider second row and all rows. It will give response of velocity curve w.r.t time with green dashed line.

Legend function used will display annotations of curve at "best" suited location.

7. After plotting response of pendulum movement against time, next objective is to plot motion of simple pendulum.

8. After this position of the pendulum is plotted at particular time in the time span and captured frame by frame to make an animation.

9. Then we can use command line tool “Cygwin” and “imagemagick” program tool to combine all these positions to form an animation.

RESULTS:

1 Plot of angular displacement and angular velocity of the pendulum over its time span.

2 The plots made over the time span of the pendulum saved frame by frame.

3 The simulation of the transient behaviour of a simple pendulum is obtained in form of an animation, which can viewed by visiting the link mentioned below:

https://youtu.be/CNbPbsYLY3o

CONCLUSION:

1. Python can be used to solve several ODEs that describe key phenomena and is able to simulate any of these systems. Which can be used for better understanding of the phenomenon.
2. Simulation of the transient behaviour of a simple pendulum by solving the ODEs governing its motion ,using Python was done.