simulation of Kinematics of simple pendulum using ODE in MATLAB

MATLAB FUNCTION CODE FOR FIRST ODE OF SIMPLE PENDULUM:-

function [dtheta_dt] = f(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

MATLAB CODE FOR KINEMATICS OF SIMPLE PENDULUM USING ODE:-

For time interval ; 0 to 100 seconds

% Kinematics of simple pendulum %

clear all
close all
clc 

m = 1;  % mass of ball (kg)
b = 0.05; % damping coefficient
g = 9.81;  % acceleration due to gravity (m/s^2)
l = 1; % length of rod (m)


%% solving ode of the pendulum


% initial angle condition of the pendulum

init_theta = [0;pi];


% solve ode. the time interval is from 0 to 100 seconds

[t,theta] = ode45(@(t,theta) f(t,theta,b,g,l,m),[0,100],init_theta);

%% Plot for Angular displacement vs time

figure(1)
grid on
plot(t,theta(:,1),'color','r');
pause(0.001)
hold on
plot(t,theta(:,2),'color','b');
pause(0.001)
xlabel('Time(sec)')
ylabel('Angular Displacement(rad) & Angular Speed(rad/s)')
legend('Angular displacement(rad)','Angular Speed(rad/s)')
title('Angular Displacement & Angular Speed vs Time')


%% Animation and Movie creation for simple pendulum motion

figure(2)
ct = 1;   % counter
for i = 1:length(theta);
    THETA = theta(i);

        clf;    % to clear the previous plot
        
        
        % coordinates of the pendulum
        
        % origin point & fixed end point of pendulum rod
        x0 = 0;
        y0 = 0;
        
        % point for the location of the ball & free end of the pendulum rod
        
        x1 = l*sin(THETA);
        y1 = -l*cos(THETA);
        
        
        % plot for pendulum
        line([x0 x1],[y0 y1],'linewidth',2,'color','k');
        hold on
        plot(x1,y1,'ro','markers',15,'MarkerFacecolor','r');
        title('Kinematics of pendulum with Damping effect')
        grid on
        axis(gca,'equal');
        axis([-1.5 1.5 -1.5 0.5])
        pause(0.001)
        M(ct) = getframe(gcf);
        ct = ct+1; 
end

% movie
movie(M)
videofile = VideoWriter('Kinematics of pendulum with Damping effect','MPEG-4');
open(videofile)
writeVideo(videofile,M)
close(videofile)

Results:-

1. Plot for angular displacement and angular velocity vs Time

From the above plot

At time t=0,

Angular displacement = 0 rad,

Angular Velocity = 3.14 rad/s

2. Animation of the simple pendulum with damping effect