Week 3 - Solving second order ODEs

AIM-2nd ODE is used to describe the transient behavior of a system. A simple example is a pendulum

g = gravity in m/s2,

L = length of the pendulum in m,

m = mass of the ball in kg,

b=damping coefficient.

use,

L=1 metre,

m=1 kg,

b=0.05.

g=9.81 m/s2.

t_span= 0-20 sec,

angular displacement=0,

angular velocity=3 rad/sec at time t=0.

clear all
close all
clc

b=0.05
g=9.81
l=1
m=2
 % initial conditions
 theta_0 = [0;5];
 
 %time ponts
 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)
plot(t,results(:,1))
 hold on
 plot( t,results(:,2))
 hold off
 grid on
 ylabel('plot')
 xlabel('Time')

 legend('displacement(m)', 'velocity (m/s^2)');
 saveas(figure(1),'shm.jpg');
 theta= results(:,1);
 ct=1
 % animating the transient motion of the simple pendulam
 figure(2)
 
 for i=1:length(results(:,1))
     x0=0;
     y0=0;
     
     x1=1*sin(results(i,1));
     y1=-1*cos(results(i,1));
     
     plot([0,0],[-1.5 0.5])
     
     hold on
     plot([x0 x1],[y0 y1],'linewidth',2,'color','b')
     axis([-1.5 1.5 -1.5 0.5])
     plot(x1, y1, 'o','markersize',40,'color','r')
   
     pause(0.1)
     M(ct)=getframe(gcf)
     hold off
     ct=ct+1
   
 end
     
%making od1e_func videofile
movie(M)
videofile=VideoWriter('od1e_func.avi');
open(videofile)
writeVideo(videofile,M)
close(videofile)

first of all clear all close all & clc

this command used for clear all screen for preveious run command

t_span used for 0-20 sec in 500 times repeat

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

solving ode45 function ( t, theta) we call the function using parameter t,theta_0

and generate the graph angular velocity and angular displacement wrt time.

initilising the counter command ct

now use the for loop command 500 times values for pendulam

now to animation of video generate file

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

results

programe run sucessfully

output:-(link)

https://youtu.be/7IXe6nUf2u4