Week 3 - Solving second order ODEs

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


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

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


#initial condition
theta_0 = [0,5]

t = np.linspace(0,10,150)

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


#plot
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}\theta2-\frac{g}{l}sin\theta1$')
plt.xlabel('time')
plt.ylabel('plot')
plt.legend(loc = 'best')
plt.show()


ct = 1

x1 =[]
y1 = []
for i in theta[:,0]:
	x0 = 0
	y0 = 0

	x1 = (l*math.sin(i))
	y1 = (-l*math.cos(i))


	#plot   
	    
	plt.figure(1)
	plt.plot([-1,1],[0,0],linewidth = 2)
	  
	#allocting co-ordinates for fixed support
	plt.plot(0,0,'o',Markersize = 10,color = 'b',markerfacecolor='b') 
	   
	#allocting co-ordinates for pendulum strin
	plt.plot([0,x1],[0,y1],linewidth = 2,color='b') 
	plt.axis([-1,1.5,-1.5,1])
	plt.plot(x1,y1,'o',Markersize =20,markerfacecolor='r')
	filename = 'pendulum%d05.png'%ct 
	ct =ct+1;
	plt.savefig(filename)
	plt.clf()




explanation

1. for solving ode in python we use the 'odeint' command from the 'scipy' library 
   with help of plot command plot the graph for 'ode'.
2. to create a pendulum animation we use for loop,marksize,markerfacecolor,linewidth,clf 
   commmads.