Week 3 - Solving second order ODEs

The aim of this program is to simulate the motion of a simple pendulum through animation and via a plot.

For this two programs are used, one to solve the ODE and the other one being the main program.

Since the input values are known, we have to declare it at the begining. The initial displacement is taken as zero and the initial velocity is 3.

The simulation time frame is 20 sec.

Using plot command, the plots for angular displacement and angular veocity are plotted.

In order to create the animation, we need a for loop for getting each frame of the image.

For capturing the frames, the movie command is used.

%Program for solving 2nd deg 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,100);

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

%Output plots
figure(1)
plot(t,results(:,1))
hold on
plot(t,results(:,2))
ylabel('Angular displacement & Angular velocity')
xlabel('time')
grid on

%Creating an animation of simple pendulum
figure(2)
ct = 1;
for i = 1:length(t_span)
    x0 = 0;
    y0 = 0;
    x1 = l*sind(results(i,1));
    y1 = -l*cosd(results(i,1));
    plot([x0 x1],[y0 y1],[-0.5 0.5],[0 0],'linewidth',3)
    hold on
    plot(x1,y1,'marker','o','markerSize',30,'markerFacecolor','g')
    axis([-1.5 1.5 -1.5 0.5])
    grid on
    pause(0.001)
    M(ct) = getframe(gcf);
    ct = ct+1;
    hold off
end
movie(M)
videofile = VideoWriter('ode_function.avi','Uncompressed AVI')
open(videofile)
writeVideo(videofile,M)
close(videofile)

%Program for solving ODE. 
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