Week 3 - Solving second order ODEs

The ODE function used for calculations:

ODE Function:

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 

Main Program:

close all
clear all
clc

b = 0.05;       %damping coefficient
g = 9.81;       %gravity
l = 1;          %length of pendulum
m = 1;          %mass of pendulum

%initial condition
theta_0 = [0;5];

%time points 
tspan = linspace(0,20,500);

%solve ODE
[t,results] = ode45(@(t,theta) ode_func(t,theta,b,g,l,m),tspan, theta_0);

%Plotting
figure(1)
plot(t,results(:,1))
hold on 
plot(t,results(:,2))
ylabel('Time')
xlabel('Velocity')
grid on
hold off

%displacement and velocity results-
displacement = results(:,1);
velocity = results(:,2);

%Considering (x0,y0) as the coordinates of the point of suspension- 
x0 = 1;
y0 = 1;
% for the increment of the frames 'ct' is used and set to value 1
ct = 1;

%for the trajectory of position of Bob of pendulum
for i = 1:length(displacement)
        x1 = 1 - l*sin(displacement(i));
        y1 = 1 - l*cos(displacement(i));

        figure(2)
        plot([0 2],[1 1],'LineWidth',15,'color','b')
        hold on
        plot([1 x1],[1 y1],[x1 x1],[y1 y1],'-O','markersize',25,'MarkerFaceColor','r','LineWidth',4)
        plot(x0,y0,'-O','markerSize',3,'markerfacecolor','b')
        %Taking (1,1) as origin and [0 2 -0.5 2] as the reference axis
        axis([0 2 -0.5 2])    
        
        %for changing the frame of the movie
        M(ct)= getframe(gcf); 
        %for incrementing the frame
        ct = ct+1;
        hold off
end

%for saving the animation-
movie(M)
videofile = VideoWriter('simple pendulum.avi','Uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)

Output:

Figure 1-

Figure 2-

https://drive.google.com/file/d/1qFSs10DevjIofdeizoECMKNSyR104jD-/view?usp=sharing