Simulation of transient behaviour of a simple pendulum and its animation

AIM : Program to simulate the transient behaviour of a simple pendulum and to create an animation of it's motion.

THEORY : A pendulum is a weight suspended from a pivot so that it can swing freely.[1] When a pendulum is displaced sideways from its resting, equilibrium position, it is subject to a restoring force due to gravity that will accelerate it back toward the equilibrium position. When released, the restoring force acting on the pendulum's mass causes it to oscillate about the equilibrium position, swinging back and forth.

The 2nd order differential equation (with damping) is given as:

where g = gravity in m/`s^2`, L = length of the pendulum in m, m = mass of the ball in kg, b=damping coefficient , `theta`= Angular Displacemnt, `d(theta)` = Angular Velocity,`d^2(theta)`) = Angular Acceleration.

MATLAB CODE :

1) Code to simulate transient behaviour of simple pendulum :

%PROGRAM TO SIMULATE THE TRANSIENT BEHAVIOR OF A SIMPLE PENDULUM AND THE
%ANIMATION OF ITS MOTION
clear all
close all
clc

% INPUTS
b = 0.05;                         % damping coefficient
g = 9.81;                         % gravity constant
l = 1;                            % length of the pendulum from the fixed point
m = 1;                            % mass of the ball
theta_0 = [0;3];                  % initial conditions ( angular_displacement , velocity)             
t_span = linspace(0,20,350);      % time plots of 10 seconds
[t,results] = ode45(@(t,theta)ode_function(t,theta,b,g,l,m),t_span,theta_0); % 2nd ODE of pendulum is being called from a function

% GRAPH PLOT
figure(1)                         
plot(t,results(:,1),'color','r')
hold on 
plot(t,results(:,2),'color','b')  
xlabel('time')
ylabel('plot')
axis([0 10 -5 5])

% INPUT POINTS OF PENDULUM
x0 = 0;                             % fixed points 
y0 = 0;
no_of_frames = 1;                   % let initial number of frames be 1


for i = 1:length(t_span)            % for loop to plot at each an every point
    x1 =  1*sin(results(i,1));      % coordinates of bob on x-axis  
    y1 = -(l*cos(results(i,1)));    % coordinates of bob on y-axis
    
   % ANIMATION PLOT
   figure(2)    
   plot([-1,1],[0,0],'linewidth',2,'color','r')
   hold on 
   plot([x0 x1],[y0 y1],'linewidth',3,'color','b')
   hold off
   hold on
   plot(x1,y1,'o','markersize',20,'markerfacecolor','g','markeredgecolor','g')
   hold off
   axis([-1 1 -1.5 0.5])
   grid on
   pause(0.001)
   M(no_of_frames) = getframe(gcf);   
   no_of_frames = no_of_frames + 1;
end
    
% ANIMATION OF PENDULUM
movie(M)
v = VideoWriter('PENDULUM.avi','Uncompressed AVI');
open(v)
writeVideo(v,M)
close(v)

2) Code for calling the function :

Here , we convert 2nd ODE to 1st ODE

function [dtheta_dt] = ode_function(t,theta,b,g,l,m) 

theta1 =theta(1)
theta2 = theta(2)
dtheta1_dt=theta2;
dtheta2_dt = -(b/m)*theta2 - (g/1)*sin(theta1);
dtheta_dt = [dtheta1_dt;dtheta2_dt];

end

OUTPUT :

ANIMATION LINK : https://www.youtube.com/watch?v=OhICSOltb_Y&feature=youtu.be

ERRORS :

corrected error : add figure(2) before plotting .