Solving Simple Pendulum Motion ODE

Introduction:

In Engineering, ODE is used to describe the transient behavior of a system. A simple example is a pendulum.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.

This ODE represents the equation of motion of a simple pendulum with damping:

In the above equation,

g = gravity in m/s2,

L = length of the pendulum in m,

m = mass of the ball in kg,

b=damping coefficient.

We have taken the values as:

L=1 metre,

m=1 kg,

b=0.5.g=9.81 m/s2.

Now we have to simulate the mtion between 0-20 sec,

for angular displacement=0,angular velocity=3 rad/sec2 at time t=0.

Description:

Now to solve the second order ODE, we have to break the equation such that it becomes two equations of first order ODE. after that we have to create the 'Function file' which solves the two equations as our required ODE.

The MATLAB Function created for the this purpose looks like this:

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

It is clear from the function that our dtheta_dt is column vector of dtheta1_dt and dtheta2_dt.

Now the following programme is run which gives us Solved Second Order Ode and Video which represents the Pendulum motion in given conditions.

clear all
close all
clc
b =0.5;
g=9.81;
l=1;
m=1;
%input conditions
theta= [0;3]
%time
t_span=linspace(0,20,200);
%solve ODE
[t,results]=ode45(@(t,theta) ode_func(t,theta,b,g,l,m),t_span,theta);
x0 = 0;
y0 = 0;
THETA_1 = results(:,1);
ct= 1
for i = 1:length(t_span)
    r = THETA_1(i)
    x1 = 1*sin(r);
    y1 = -1*cos(r);
   
    plot([x0 x1],[y0 y1],'linewidth',1,'color','bla');
    viscircles([x1 y1],0.05)
    rectangle('Position',[-0.4 0 0.8 0.2],'facecolor','c')
    title('Simple Pendulum')
    line([0,0],[-5 0],'color','blue','linestyle','--')
    axis([-1 1 -1.4 0.2])
    M(ct) = getframe(gcf)
    ct = ct+1
 
    
end
movie(M)
videofile= VideoWriter('Pendulum_animation.avi','Uncompressed AVI')
open(videofile)
writeVideo(videofile,M)
close(videofile)

%displacement and velocity against time
close all
hold on 
plot(t,results)
plot(t,results)
grid on
xlabel('Time(in Sec)')
legend('displacement','velocity')

Above programme generates the graph which shows angular Displacement and angular velocity of Pendulum as shown:

YOUTUBE link of the video that represents SIMPLE PENDULUM MOTION:

https://youtu.be/lp4OBpLaIoY