Solving a Pendulum's Second Order ODE

Here is a MATLAB program that solves a second order ODE which represents the equation of motion of a simple pendulum with damping along with the code to simulate the motion of a pendulum between 0-20 sec, for angular displacement=0, angular velocity=3 rad/sec at time t=0.

clear all
close all
clc

b=0.5;
m=1;
l=1;
g=9.81;

% Initial condition displacement and velocity
theta_0=[0;3];
%Time points
t_span=linspace(0,20,500);
% solve ODE
[t, results]=ode45(@(t, theta)ode_func(t,theta,b,g,l,m),t_span,theta_0);

figure(1)
plot(t,results(:,1))
hold on
plot(t,results(:,2))
ylabel('Plot')
xlabel('Time(sec)')
legend('Displacement(rad)','Velocity(rad/sec)')
hold off

%% Animation
figure(2)

O=[0 0];
axis(gca,'equal');
axis([-1 1.2 -1.5 .2]);
grid on
ct=1;
% Circle in Origin
O_circ= viscircles(O,0.01);
for i=1:length(t_span)
  
  %Mass Point
  P=l*[sin(results(i,1)) -cos(results(i,1))];
  %pendulum 
  pend= line([O(1) P(1)],[O(2) P(2)]);
  
  % Ball
  ball= viscircles(P,0.05);
  
  % Timme interval to update the plots
  pause(0.001);
  M(ct)=getframe(gcf);
  ct=ct+1;
  % Delete previous object if it is not the final loop
  if i<length(t_span)
  delete(pend);
  delete(ball);
  end
end

% Video making Commands
movie(M) %generates a movie with frames compiled from M
videofile = VideoWriter('Pendulum.AVI','Uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)

ODE Function

function [dtheta_dt]=ode_func(t,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

RESULTS:

• Plot

• Simulating the motion of the pendulum