Solution of 2nd order ODE and animation of damped oscillations for a simple pendulum

The following code generates animation of pendulum motion with damped oscillations and an attached mass at its end. ode45 is used to solve the second order differential in time \'t\' and position \'theta\' which gives solution containing both position of the pendulum(theta) in one column and angular velocity of the mass in other column for various time. The string is plotted using the for loop and various points of the mass bob are obtained using relations x=l*cos(theta) and y=l*sin(theta) .

Code for plotting results and generating animation

%Motion of a simple pendulum with damping
clear all 
close all
clc
%given data(\'l\'in metre,\'m\'in kg,\'g\' in m/s^2)
l=1;
b=0.5;
g=9.81;
m=1;
x0=0;
y0=0;
ct=1;
%initial conditions
theta_i=[0;3];

%time inputs
t_span=linspace(0,20,500);

%solving
[t,results]=ode45(@(t,theta) pendulum_func(t,theta,l,b,g,m),t_span,theta_i);

%plotting results
figure(1)
plot(t,results)
xlabel(\'Time\')
ylabel(\'Results\')
grid on

%plotting string
    figure(2)
for i=1:length(t)
    THETA=results(i,1);
    x1=l*sin(THETA);
    y1=-l*cos(THETA);
    plot([x0 x1],[y0 y1],\'linewidth\',3)
    axis([-1 1 -1.5 1.5])
    
    %plotting rectangle
    rectangle(\'position\',[-0.5 0 1 0.5],\'linewidth\',5)
    
    %plotting circle
    viscircles([x1,y1],0.05,\'edgecolor\',\'r\',\'linewidth\',2)
    grid on
    
    M(ct)=getframe(gcf);
    ct=ct+1;
    pause(0.002)
end
movie(M)
video=VideoWriter(\'Motion_of_a_Pendulum.avi\',\'uncompressed AVI\');
open(video)
writeVideo(video,M)
close(video)

The function code contains the equations governing the motion of the pendulum.

Code for function:

function[dtheta_dt]= pendulum_func(t,theta,l,b,g,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

Animation file link:

https://drive.google.com/file/d/1dqLSryvZrX3jNdkbjqfEC-mZPLsjcQQf/view

The responses get damped down since b is quite substantial.For lower value of b,the oscillations of the pendulum keep on going for time t=20 seconds. It can be also observed from the plot that when the theta is at extreme position,angular velocity is zero and at mean position, angular velocity is maximum.