Week 3 - Solving second order ODEs

clear all
close all
clc

% input parameters 
b = 0.05;
m = 0.1;
g = 9.81;
L =1;

% initial conditions 
initial_0 = [0;5];

% time points 
t_span = linspace(0,10,100);

% solving ODE 
[t, results] = ode45(@(t, theta) ode_func(t,theta,b,m,L,g),t_span, initial_0);

    figure(2);
for i = 1:length(results(:,1))
    x0 = 0
    y0 = 0
    x1= L*sin(results(i,1));
    y1= (-1)*L*cos(results(i,1));
    
    plot([-1 1], [0 0], 'linewidth', 4, 'color', 'g');
    axis ([-1 1 -2 0]);
    hold on 
    plot([x0 x1], [y0 y1], 'linewidth',2,'color', 'b');
    hold on 
    plot(x1, y1, '-o','markers', 20, 'markerfacecolor','r'); %plotting the bob using marker and facecolor
    s1 = 'Time= ' ;
    s2 = num2str(t(i));
    text(0,-0.5,[s1, s2]);
    M(i) = getframe(gcf);
    hold off
end 

movie(M)
videofile = VideoWriter('Pendulum.avi');
open(videofile);
writeVideo(videofile,M);
close(videofile);

figure(3)
plot(t, results(:, 1));
hold on;
plot(t, results(:, 2)); 

xlabel('time');
ylabel('plot');



%% ode_func 

function [dtheta_dt] = ode_func(t,theta,b,m,L,g)

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 



% Conclusion : The first code written solves the equation for two terms namely anugular dispalcement and velocity (2 columns of results).Also it gives the animation of the pendulum and saved as a movie in the same folder. For loop is used to capture each picture at time interval in between 0 to 10 seconds for 100 intervals. For second code, it is for defining of the function ode45 which includes two first order ODEs which are solved in that one.