Week 3 - Solving second order ODEs

Animation of Simple Pendulum

clear all
close all
clc


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

%initial condition
theta_0 = [0;5];

%time span
t_span = linspace(0,10,500);

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

plot(t,results(:,1))
hold on
plot(t,results(:,2))
ylabel('Plot')
xlabel('time')

The above code will create a graph of a simple damped pendulum.


%%%%%%% the ODE function %%%%%%%%%%%%

function [dtheta_dt] = ode_fun(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

The above code is to define the function which can be called in the simple pendulum graph

%%%%%%%%%%% simple pendulum animation %%%%%%%%%

clear all
close all
clc


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

%initial condition
theta_0 = [0;5];

%time span
t_span = linspace(0,10,500);

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

%Animation
%Origin
x(1)=0;
y(1)=0;


 

 %for loop
 ct =1;
 for i=1:length(results)
     x(2)= l*sin(results(i));
     y(2)= -l*cos(results(i));
     plot(x,y,'k-');
     hold on;
     plot(x(2),y(2),'k.','MarkerSize',35);
     hold off;
     grid on;
     axis([-1.5 1.5 -1.5 1.5]);
     pause(0.01);
     
     M(ct) = getframe(gcf);
     ct = ct+1;
   
 end
    movie(M) 
    videofile = VideoWriter('Simple_Pendulum','Uncompressed AVI');
    open(videofile)
    writeVideo(videofile,M)
    close(videofile)