Simulation of a Simple Pendulum

b = 0.5 Damping factor (kg/s)

m = 1 Mass (kg)

g = 9.81 Acceleration due to gravity (m/s^2)

L = 1 Length of Pendulum (m)

Angular displacement (rad)

Angular velocity (rad/sec)

simple_pendulum_run

% This script simulates simple pendulum motion

% Clearing
close all;
clear;
clc;

% Inputs
b = 0.5; % Damping factor (kg/s)
m = 1;   % Mass (kg)
g = 9.81; % Acceleration due to gravity (m/s^2)
L = 1;    % Length of Pendulum

% Initial condition (At t = 0)

theta0 = [0 ; 3]; % theta0 = [Angle of pendulum axis from vertical(rad) ; Angular velocity(rad/s)] 

% Simulation Time(seconds)
  % Simulation runs from o sec to 10 sec
  % with time step of 0.02 sec
t_span = linspace(0,20,500);

% Solving ODE for simple pendulum

[t, result] = ode45(@(t, theta) simple_pendulum_ode(theta, b, m, g, L), t_span, theta0);
Ang_Disp = result(:,1);
Ang_Vel  = result(:,2);

figure(1);
plot(t,Ang_Disp);
hold on;
plot(t, Ang_Vel);
hold off;
xlabel('Time in seconds');
ylabel('Angular displacement and Angular velocity');
legend('Angular displacement', 'Angular velocity');
% Animation

figure(2);

ct = 1;
for i = 1:size(t)
    figure(2);
    plot([0 L * sin(Ang_Disp(i))], [1.5 (1.5-L * cos(Ang_Disp(i)))], '-g', 'LineWidth', 1);
    hold on;
    plot(L * sin(Ang_Disp(i)), 1.5-L * cos(Ang_Disp(i)), 'Ob', 'MarkerSize', 20);
    hold off;
    axis([-1 1 0 1.5]);
    M(ct) = getframe(gcf);
    ct = ct + 1;
end
hold off;
movie(M);
videofile = VideoWriter('simple_pendulum.avi','Uncompressed AVI');
open(videofile);
writeVideo(videofile,M);
close(videofile);



simple_pendulum_ode

function [ dtheta_dt ] = simple_pendulum_ode(theta, b, m, g, L)
%SIMPLE_PENDULUM 1st order ODE set 
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

Simulation video:-