Solving second order ODEs (PYTHON)

AIM: PYTHON Program which will simulate the motion of pendulum by solving the given ODE (Ordinary Differential Equation).

THEORY:

A Pendulum is a body suspended from a fixed support so that it swings freely back and forth under the influence of gravity.

The way the 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.

The differential equation which represents the motion of a simple pendulum is

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

GOVERNING EQUATION:

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

where,

b is the damping co-efficient,

m is the mass of the ball,

g is acceleration due to gravity,

l is the length of the pendulum,

and θ is the angular displacement.

OBJECTIVE:

• Program to solve the given ODE which represents the equation of motion of a pendulum with damping.
• Simulate the motion of pendulum and create an animation.

INPUT VALUES:

l = 1 m

m = 1 kg

b = 0.05

`g = 9.81 m⁄s^2`

PROGRAM:

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

# Function that returns dz/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

b = 0.05
g= 9.81
l = 1
m = 1

# Initial condition
theta_0 = [0,3]

# Time points
t = np.linspace(0, 20, 300)

# Solve ODE 
theta = odeint(model, theta_0, t, args = (b,g,l,m))

# Plots
plt.plot(t, theta[:, 0], 'b-', label = r'$\frac{d\theta_1}{dt} = \theta2$')
plt.plot(t, theta[:, 1], 'r--', label = r'$\frac{d\theta_2}{dt} = -\frac{b}{m}\theta_2 -\frac{g}{L}sin\theta1$')
plt.xlabel('Time')
plt.ylabel('Plot')
plt.legend(loc = 'best')
plt.show()

#Animation using image stitching
# Defining counter variable
ct =1
for th in theta[:,0]:
	x_start = 0
	y_start = 0

	x_end = l*math.sin(th)
	y_end = -l*math.cos(th)

	# Naming the figures using the counter varible and saving them in '.png' format
	filename = 'test%05d.png' %ct
	ct = ct +1 

	# Plotting figures and saving them with particular filename
	plt.figure()


	# plotting the support
	plt.plot([-1,1],[0,0], linewidth =4, color = 'red')

	# plotting the string through which the ball is attached
	plt.plot([x_start, x_end][y_start,y_end], linewidth =2, color ='blue')

	#plotting the ball
	plt.plot(x_end,y_start, marker ='o', markersize =25, color = 'green')

	plt.xlim([-2, 2])
	plt.ylim([-1.5, 0.5])
	plt.savefig(filename)

PLOTS:

EXPALNATION OF THE CODE:

Initially, a function program is created where specified function is mentioned covering all the parameters and the variables.

The given values are specified and the initial condition are imposed.

The time period over which the pendulum’s transient behavior has to be observed is also defined.

Then, the ‘odeint’ command acts as a solver for the first-order differential equations extracted out of a single second order differential equation.

Animation of the pendulum motion is done using Imagemagik by stitching the images or plots obtained.

CONCLUSION: We have established and created the simulation of motion of simple pendulum by solving the given second order ODE.