Simulating the Transient Behavior of a Simple Pendulum

Objective: To write a program in python to simulate the transient behaviour of a simple pendulum and to create an animation of its motion.

Steps:-

• After importing the required libraries in our code, we begin with defining a function named model which takes the input parameters of the pendulum, i.e., damping coefficient (b), acceleration due to gravity (g), length of the string (l), mass of the ball (m) and time (t) and then returns `(d)/(dt)`(theta), where theta is a `2xx1` matrix with the first row representing the pendulum\'s angle while the second row representing the angular velocity.
• We define the parameters b, g, l, m, the initial conditions of theta and the time for which the pendulum is to be evaluated.
• We apply odeint to solve the ODE returned by the function model.
• After getting the values of the pendulum\'s angle at different time, we use trigonometry to calculate the ball\'s position for different time instances.
• We plot the pendulum by joining the fixed position of the string to its moving end. We draw a small circle with radius of 0.1 unit to represent the ball.
• We save each figure in a folder maintaining proper sequence of filenames.
• Finally, we stitch all the images together to create an animation using the ImageMagick software.

Code:-

(The different steps have been explained as comments in the code.)

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

# defining a function that returns d/dt
def model(theta,t,b,g,l,m):

	theta1 = theta[0]
	theta2 = theta[1]
	# writing the equations of derivatives in the case of simple pendulum
	dtheta1_dt = theta2
	dtheta2_dt = -(b/m)*theta2 - (g/l)*math.sin(theta1)
	dtheta_dt = [dtheta1_dt,dtheta2_dt]
	return dtheta_dt

# The parameters of the pendulum

# Damping coefficient. Unit: kg/s
b = 0.05
# acceleration due to gravity. Unit: m/s^2
g = 9.81
# length of the pendulum. Unit: m
l = 1
# mass of the ball. Unit: kg
m = 1

# initial condition: first element denotes initial angle (in radians) and the second element denotes initial angular speed (in radian/second)
theta_0 = [0,3]

# time inputs: creating 201 elements between 0 and 20 seconds so that we get the matrix theta after each 0.1 second.
t = np.linspace(0,20,201)

# Now, solving the ODE using odeint 
theta = odeint(model,theta_0,t,args=(b,g,l,m))

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

# Cooridinate of the fixed point of pendulum: lets take it as origin
x0 = 0
y0 = 0

# count timer
ct = 1

# calculating the moving end of the pendulum after each 0.1 second
for i in range(0,200):
	# Coordinate of the moving end of the string
	x1 = l*math.sin(theta[i][0])
	y1 = -abs(l*math.cos(theta[i][0]))

	# creating unique file names for each position
	filename = \'test%05d.png\' %ct
	ct = ct + 1

	# plotting the string and the sphere
	# using subplot technique for creating the sphere of 0.1 m radius
	fig, ax = plt.subplots()
	plt.plot([x0,x1],[y0,y1])
	circle1 = plt.Circle((x1,y1),0.1,color=\'r\')
	plt.xlim([-1.25,1.25])
	plt.ylim([-2.5,0])
	ax.add_artist(circle1)
	# saving every image as a different file and maintaining the proper sequence of the files
	plt.savefig(filename)

Results:-

• We get the graph showing the pendulum\'s angle and its angular velocity at different time instances:-

(The blue continuous line represents the angle (in radian) made by the pendulum with the vertical axis while the red dotted line represents its angular velocity (in radian/second) with time.)

We can clearly see in the figure, that the angular velocity as well as the angle covered by the pendulum in each cycle is decreasing with time. This is due to the damping coefficient, b.

• After getting all the images and stitching them using ImageMagick software, we get the following GIF:-

This animation has also been uploaded on youtube: https://youtu.be/VMSnTmDquZ4

Errors Encountered While Programming:-

• Initially, I faced an error named \"SyntaxError: invalid syntax\", as can be seen in the following screenshot:-

I came to know that the error was because of the imbalance in brackets on line number 58. After correcting the brackets, the code ran smoothly.

• The second problem I faced was that the length of the pendulum was changing with the change in time. No error was being thrown on the command window while this was happening. The GIF showing this issue is as follows:-

I solved this error by setting the limit of the y-axis so that its scale does not change.

• While setting the limit for the y-axis after facing the precious error, I found out that the pendulum was now facing upwards instead of downwads. This was because I had set the limit to be [0,-2.5]. This can be seen in the following GIF:-

This error was happening because the limit should have been set to [-2.5,0] instead of [0,-2.5]. Correcting this solved this issue.

Screenshots of my Program:-

• Getting the graph showing the pendulum\'s angle and its angular velocity:-

• Code running successfully:-

• Getting all the images of the pendulum at different time instances:-