Week 3 - Solving second order ODEs

% Usin Second order derivation to calculating pendulum motion
% Input
 b = 0.05; % Damping coefficient
 g = 9.81; % Acceleration due to gravity
 l = 1; % Length of the pendulum
 m = 0.1; % Mass of the ball

 % Initial condition
 theta_O = [0;3]; % theta is the initial condition for the pendulum
t_span = linspace(0,20,500); % t_span is the time variation of motion of the pendulum
 % Slove ODE
 [t,theta] = ode45(@(t,theta) ode_func(t,theta,b,g,l,m),t_span,theta_O); 
 figure(1) % Figure (1) shows the damping of the pendulum in curves with respective to time and plot
 hold on
 plot (t,theta(:,1),'Color','r')
 hold on
 plot (t,theta(:,2),'Color','b')
 xlabel('time');
 ylabel('plot');
 %displ = theta(:,1);
 %% Animation
 ct = 1; % ct is counter 
 for i = 1:length(t) % i is a veriable and for loop is started
x = l*sin(theta(i)); % x is used to determine the point in the x axis
y = -l*cos(theta(i)); % y is used to determine the point in the y axis
x0=0;
y0=0;
figure(2) % Figure (2) shows the motion of the pendulum in animation
plot([x0 x],[y0 y],'LineWidth',3)
hold on
plot(x,y,'.','MarkerSize',60)
xlabel('x axis');
ylabel('y axis');
axis([-1.5 1.5 -1.5 0.2])
pause(0.03)
hold off
Q(ct) = getframe(gcf); % Q is denoted as variable for getting the pendulum motion frame by frame
% Text = sprintf("pendulum%d.png,ct);
% Saves(1,Text)
ct = ct+1; % here the ct is incremented by 1 every time
 end % For loop is ended with a command end
movie(Q) % Movie command is used to create a movie of pendulum motion
videofile = VideoWriter('pendulum_animation.avi','Uncompressed avi'); 
% (Videofile command is used to create a videofile and videowriter is used
% to write a video name and its formate.)
open(videofile) % Open command is used to open the videofile 
writeVideo(videofile,Q) % Writre command is used to write the motion of the pendulum is a video formate
close(videofile) % Close command is used to close the video file and end the video.

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

ANIMATION MOTION OF PENDULUM :

Pendulum_animation