CREATING MOTION OF SIMPLE PENDULUM WITH DAMPING USING MATLAB ODE SOLVER

Pendulum is a body suspended from a fixed support so that it swings freely back and forth under the influence of gravity. 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 mathematics of pendulums is in general quite complicated. Simplifying assumptions can be made, which in the case of a simple pendulum allow the equations of motion to be solved analytically for small-angle oscillations.
Damped Simple Harmonic Motion
When the motion of an oscillator reduces due to an external force, the oscillator and its motion are damped. These periodic motions of gradually decreasing amplitude are damped simple harmonic motion. An example of a damped simple harmonic motion is a simple pendulum.
In the damped simple harmonic motion, the energy of the oscillator dissipates continuously. But for a small damping, the oscillations remain approximately periodic. The forces which dissipate the energy are generally frictional forces.

(d²θ/dt²) +(b/m)(dθ/dt) + (g/l)sinθ =0

It is the equation for the damped simple harmonic motion.

Where

b=damping coefficient of the string material

m=mass of the bob

g=gravity

l=length of the string

this equation is solved by converting in to the matrix as shown below

d/dt[θ1 θ2]=[θ2 (-b/m)θ2-(g/l)sinθ]

the MATLAB Program for the animation of pendulum as follow

2^nd Order Differential function solver:

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

The main program for the animation as follows:

%Cleaning the screen, command line and workspace, closing all plots
clear all
close all
clc
 
%input parameters
l=1;
g=9.81;
m=1;
b=0.05;
 
%initial conditions
theta_0=[0;3]; % initial displacement=0 and the angular velocity=3rad/sec
 
%time
t_span=linspace(0,20,500); %time span ranging between 0-20sec and divided in to 500 parts
 
%calculating the 2nd order ODE with the help of funtion
[t,results]=ode45(@(t,theta) ode_func(t,theta,b,g,l,m),t_span,theta_0);
 
%creating the loop for the oscillations of pendulum
 
ct=1;
for i=1:length(results(:,1)); %as the results contains both time and angular velocity
        x0=0;
        y0=0;
        x1=l*sin(results(i,1)); %assigning the direction of oscillation of the bob
        y1=-l*cos(results(i,1));
        %creating the regid support
        plot([-1,1],[0,0],\'linewidth\',4,\'color\',\'k\');
        hold on;    %erase the previous stored data dynamically
        
        %creating the rope
        plot([x0,x1],[y0,y1],\'linewidth\',2,\'color\',\'b\');
        
        %creating the shape and size of the bob using marker and facecolor command
        plot(x1,y1,\'-o\',\'markers\',20,\'markerfacecolor\',\'r\'); 
        grid on;
        axis([-1.5,1.5,-1.5,0.5]);  %assaigning the plot area
        pause(0.0003);
        hold off;
        q(ct)=getframe(gcf); %capturing the each frame for loop ct
        ct=ct+1;
end
movie(q)
videofile=VideoWriter(\'pendulum_oscilation.avi\');
open(videofile)
writeVideo(videofile,q)
close(videofile)

The YouTube link for the animation as follows:

https://youtu.be/G2qyqTbKgfo