Week 3 - Solving second order ODEs

close all 
clc

% Main inputs 
b=0.05
g=9.81
l=1
m=1
 % Simulate the motion between 0-20 sec, for angular displacement=0,angular velocity=3 rad/sec at time t=0.

% Initial condition are given 
theta_0 = [0;4];

% Tiime points allocated 
t_span = linspace(0,30,200);

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

plot(t,results (:,1))
hold on
plot(t,results(:,2))

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

ct=1;



for i = 1:length(t)
    
    % Initial possition of the pendulum 
    
    x0 = 0;
    y0 = 0;
    
    %Final possition of the pendulum 
    
    x1 =  l*sin(results(i,1));
    y1 = -1*cos(results(i,1));
    

    figure(2)
    % Base line are located 
    plot([x0 x1],[y0 y1],[-1.5,1.5],[0,0],'linewidth',3);
    hold on;
    % ploting of pendulum ball
    plot(x1,y1,'d','markers',30,'markerfacecolor','g');
    plot(x1,y1,'+');
    grid on
    
    axis([-1.2 2 -1.4 1]);
    pause(0.01);
    hold off;
    M(ct) = getframe(gcf);
    ct= ct+1;
   
end








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

    
    
    
    
function [dtheta_dt] = ode_func(~,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