Week 3 - Solving second order ODEs

video animation link :https://drive.google.com/file/d/1FKuPLqRD60RoJDyM_9Ok2yPRcPgyoZiE/view?usp=sharing

explanation has been given in program for your convenience:

program 1 : solving ode

clear all
close all
clc

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

%initial conditions%
theta_0=[0,3]; %determination of initial conditions of pendulum%

%time values%
t_span=linspace(0,20,200); %the values of the time span for which the pendulum will oscillate% 

%solving ode% %the ouptupts of the ode ,the dependent and independent
%values ,the inputs the time span for which the ode is to be solved and the
%initial conditions are called%
[t,results]=ode45(@(t,theta) ode_function(t,theta,l,b,g,m),t_span,theta_0);
figure(1)
%plotting% %plot to show how the oscillations and damping occurs%
figure(1)
plot(t,results(:,1));
hold on
plot(t,results(:,2));
xlabel('time');
ylabel('velocity');
%creating animation to show working of pendulum%
ct=1; 
for i=1:length(t_span)
figure(2)
x0=0;%x coordinate for starting point of pendulum%
y0=0;%y coordinate for starting point of pendulum%

x1=l*sin(results(i,1));%x coordinate for ending point of pendulum%
y1=-l*cos(results(i,1));%y coordinate for ending point of pendulum%

plot([x0 x1],[y0 y1],'LineWidth',3);
hold on
plot(x1,y1,'o','markersize',15,'MarkerFaceColor','r');
axis([-2 2 -2 2])
hold off;
a(ct)=getframe(gcf);
ct=ct+1;
end
%avi creation%
movie(a)
videofile=VideoWriter('simple pendulum.avi');
open(videofile)
writeVideo(videofile,a)
close(videofile)

program 2 : ode function

function[dtheta_dt]=ode_function(t,theta,l,b,g,m); %creating function definition%
theta1=theta(1);%assigning theta1 value to theta%
theta2=theta(2);%assigning theta2 value to theta%
dtheta1_dt=theta2;%ode1%
dtheta2_dt=-(b/m)*theta2-(g/1)*sin(theta1);%ode2%
dtheta_dt=[dtheta1_dt;dtheta2_dt];%combining ode1 &ode 2%
end