Simulating a 2nd order ODE system of a simple pendulum

Objective - To understand the motion of simple pendulum using MATLAB

Theory :

A simple pendulum consist of mass 'M' attached to one end of a light inextensible string of length 'l', and the other end of the string attached to a fixed point.

When the bob in equilibrium is subjected to an external force it displaces from its initial position and oscillates in a circular arc in the same vertical plane.

This motion can be described in mathematical form and a code can be written in MATLAB to get the actual representation of motion.

In Engineering, ODE is used to describe the transient behaviour 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 Diffrential Equation that describes th position of the bob w.r.t. time.

Mathematical Equation:

The following is the mathematical equation :

In the above equation,

g=9.81 , acceleration due to gravity in m/s^2,

L = 1 , length of the pendulum in m,

m = 1 , mass of the ball in kg,

b = 0.05 , damping coeffiicient.

PROCEDURE :

• In order to solve this, firstly the 2nd order ODE is converted into first order which can be then be used in function.
• Calling the function using the inbuilt ode45
• Inputting initial data, such as time, angular displacement and angular velocity along with parameters L,g,b and m.

The results are stored in 2 columns, in the first column angular displacement values and in 2nd column angular velocity values.

MATLAB CODE :

clear all
close all
clc

% Inputs

L = 1; % Length of the string
m = 1; % Mass of the bob
g = 9.81; % Acceleration due to gravity
b = 0.05; % Damping coefficient
ct = 1; % Counter variable

% Initial condition

dis_v = [0;3];

t_span = linspace(0,20,500); % Time

% ODE function calling

[t,output] = ode45(@(t,theta) ode_pendulum(t,theta,b,g,L,m),t_span,dis_v);

figure(1)
hold on
plot(t,output(:,1))
plot(t,output(:,2))
xlabel('X-axis')
ylabel('Y-axis')
grid on
title('variation of angular displacement and angular velocity w.r.t time','color','r')

theta2 = output(:,1);

for i = 1:length(theta2);

theta_x = theta2(i); % theta_x keeps changing w.r.t time

x0= 0;
y0= 0;

x1 = L*sin(theta_x);
y1 = -L*cos(theta_x);

% plotting the pendulum
figure(2)
plot([-0.5 0.5],[0,0],'color','b','linewidth',5);
hold on
plot([x0 x1],[y0 y1],'color','r','linewidth',2)
hold on
plot(x1,y1,'o','markers',15,'markerfacecolor','y')
hold off
grid on
axis([-1.5 1.5 -1.5 0.5]);
pause(0.002)
M(ct) = getframe(gcf);
ct = ct+1;

end

movie(M)
videofile = Videowriter('pendulum_challenge.avi','Uncompressed avi');
open(videofile)
writeVideo(videofile,M)
close(videofile)