Week 3 - Solving second order ODEs

Objective: Program to solve 2nd Order ODE which represents the equation of motion of a simple Pendulum with damping

In Above Equation,

• g=Gravity in m/s2
• L=Lenght of Pendulum in m
• m=Mass of the ball in Kg
• b=Damping Coefficient

As per question,

L=1m; m=1Kg; g=9.81m/s2; b=0.05

Now, We have to Simulate the motion between 0-20 sec, for angular displacement=0,angular velocity=3 rad/sec at time t=0.

% Program to stimulate the Pendulum motion
clear all
close all
clc
b=0.05;
g=9.81;
l=1;
m=1;

%Initial Condition
theta_0 = [0;3];

%Time Points
t_span = linspace(0,10,500);

function[dtheta_dt] = ode_func(~, theta, b,g,l,m)
             theta1=theta(1);
             theta2=theta(2);
             dtheta1_dt=theta2;
             dtheta2_dt=-(b/m)*theta2-(g/l)*sin(theta1);
             dtheta_dt=[dtheta1_dt;dtheta2_dt];
         end

%SolveODE
[t, results] = ode45(@(t,theta) t_span , theta_0, ode_func(t, theta, b,g,l,m));
plot(t,results(:,1))
Hold On
plot(t,results(:,2))
ylabel('Plot')
xlabel('Time')

% for Animation
figure(1)
for i=1:length(t)
i;
clf
plot([-1.5,1.5],[0,0],'b-','LineWidth',10)
hold on
plot([0,0],[0,-1.5],'b--')
hold on
plot(x,y,'r--')
hold on
plot([0,x(i)],[0,y(i)],'k-','LineWidth',3)
plot(x(i),y(i),'r-o','markersize',30,'markerfacecolor','r')
plot(x(i),y(i),'r--')
axis([-1.5 1.5 -1.5 0])
 
%pause(0.1)
%if i~=length(t)
%clf
%end
movieVector(i)=getframe;
end

myWriter=VideoWriter('Pendulum');
%myWriter.FrameRate=10; % Frame per second
open(myWriter);
writeVideo(myWriter,movieVector);
close(myWriter);