Week 3 - Solving second order ODEs

clear all
close all
clc

b=0.05;
g=9.81;
l=1;
m=1;

% initial condition
theta_0 = [0;3];

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

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

figure(1)
 % angular displacement
plot(t,results(:,1)) 
hold on 
 % angular velocity
plot(t,results(:,2)) 
ylabel('angular displacement & angualar velocity')
xlabel('time')
title('Simple pendulum')
grid on

% creating animation of simple pendulum

figure(2)
ct=1;

for i = 1:length(t_span)
    THETA = results(i,1);

x1 = l*sin(THETA);
y1 = -l*cos(THETA);

plot([-0.5,0.5],[0,0],'linewidth',3,'Color','r');
line([0 x1],[0 y1],'linewidth',3,'color','g');
hold on
plot(x1,y1,'o','Markers',20,'markerfacecolor','b');
grid on
axis([-2 2 -2 1.5]);
pause(0.0000000001)
hold off
M(ct)=getframe(gcf);
ct=ct+1;
end
movie(M)
videofile = VideoWriter('simple_pendulum.avi','uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)

ans 2nd order ode to represent equation for pendulum

given equation consist of b,g,l,m and sinθ parameters are given( V=3rad/sec and desplacement=0 at initial condition) and time span is 20sec

step1)as the equation is given for the problem and following parameters are also give which are considerate in the equation
as
b=0.05 damping coefficient
g=9.81 gravitational force
l=1 length of string
m=1 mass of pendulum

step2) equation as shown is consist of motion of simple pendulum
solving the equation it is splited it into single order equation one input parameter are as above

dtheta1_dt = theta2;
dtheta2_dt = -(b/m)*theta2 - (g/l)*sin(theta1);

step3) there are two ode for angular displacement and velocity and for solving ode we provide the initial condition
and assign it as theta_0
theta_0 = [0;3];.....(0 angular displacement for 3 initial angular velocity)

step4) for ode we give the time interval of integration in this program is 20sec,& t_span is interval integration assign to theta_0
t_span = linspace(0,20,500);.......linspace for interval time is ((20sec)) is equally divided in 500 intervals

step5) for solving ode we are using ode45 inbuilt command in matlab for solving ode is function of theta and time

step6) 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/1)*sin(theta1);
dtheta_dt = [dtheta1_dt; dtheta2_dt];
end

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/1)*sin(theta1);
dtheta_dt = [dtheta1_dt; dtheta2_dt];
end

step7)solve ODE by follwing expresssion
[t,results] = ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0);

step8) now plot time result to the angular displacement & the angular velocity
plot(t,results)
% angular displacement
plot(t,results(:,1))
hold on
% angular velocity
plot(t,results(:,2))
ylabel('angular displacement & angualar velocity')
xlabel('time')
title('Simple pendulum')
grid on

step9)
second plot will be animation of simple pendulum for, creating the aimation first counter is set to 1

ct=1;------------------------------------------------------------------------counter added

for i = 1:length(t_span)
THETA = results(i,1);

x1 = l*sin(THETA);
y1 = -l*cos(THETA);

plot([-0.5,0.5],[0,0],'linewidth',3,'Color','r');
line([0 x1],[0 y1],'linewidth',3,'color','g');
hold on
plot(x1,y1,'o','Markers',20,'markerfacecolor','b');
grid on
axis([-2 2 -2 1.5]);
pause(0.0000000001)
hold off
M(ct)=getframe(gcf);
ct=ct+1;
end
movie(M)
videofile = VideoWriter('simple_pendulum.avi','uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)

step10) to get the individual value of THETA which is equal to result of i,1 that means each time it will access
to the column one values and with changing i that is row number . and that theta is used to calculate the coordinate x1 &y1
i = 1:length(t_span)
THETA = results(i,1);

x1 = l*sin(THETA);
y1 = -l*cos(THETA);

step11) for plotting the line of pendulum each time we are using the x1 y1 coordinates and initial coordinate are given as 0,0
plot([-0.5,0.5],[0,0],'linewidth',3,'Color','r');
line([0 x1],[0 y1],'linewidth',3,'color','g');
hold on
plot(x1,y1,'o','Markers',20,'markerfacecolor','b');
grid on
axis([-2 2 -2 1.5]);

step12)for representing the ball of pendulum we using marker with 25 size and blue colour
plot(x1,y1,'o','Markers',20,'markerfacecolor','b');

step13) individual frames which is act is stored rin separate array and getframes command is used to collect the all frames and creating animation
and using movie (M) command animation of simple pendulum is created

movie(M)
videofile = VideoWriter('simple_pendulum.avi','uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)

step14) errors which is faced are counter not added,ode45making function,and during animation, got help from videos available on internet and skilllync videos

step15)following youtube link for animation of simple pendulum,

https://youtu.be/S6XPqqOLga0