Week 3 - Solving second order ODEs

Objective: The aim is to obtain the animation plot of the pendulum by programming using Matlab.

Code:

clear all
close all
clc

%inputs
b = 0.05;
g = 9.81;
l = 1;
m = 0.1;

% initial condition
theta_0 = [0;5];

% time points
t_span = linspace(0,10,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))
xlabel('time')
ylabel('angular displacement & angular velocity')
title('variation of angular displacement & angular velocity w.r.t time')

theta0 = results(:,1);
ct = 1;

for i=1:length(theta0)
    THETA = theta0(i);
    x0 = 0;
    y0 = 0;
    x1 = l*sin(THETA);
    y1 = -l*cos(THETA);
    
    figure(2)
    plot([x0 x1],[y0 y1],'linewidth',3,'color','r')
    hold on
    plot([-0.5 0.5],[0 0],'linewidth',3,'color','b')
    hold on
    plot(x1,y1,'.','markersize',50)
    hold off
    axis([-1 1 -1 1])
    pause(0.1)
    M(ct) = getframe(gcf);
    ct = ct+1;
end

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

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

Explanation: Here, I provided the necessary input of pendulum motion requires such as mass of the bob, damping coefficient, string length, gravity, time span & theta for the mathematical equation of “d²Ѳ/dt² + (b/m)(dѲ/dt) + (g/l)sinѲ” which I’ve taken to solve second order ODE using programming for getting the pendulum motion animation plot. Usually, when the external force applied to the bob it’s motion is in the circular with respect to all the axis. To determine the variable of angular displacement & angular velocity with respect to time I used the Matlab inbuilt syntax “ode45” by giving input functions. Inside function, I converted 2^nd order ODE into 1^st order ODE by coding formula for getting the plot seems great. Then, for the pendulum animation plot result I used “for loop” command and theta as result of angular displacement as an input. Then, I provided the coordinates & contour inside the command for the looping purpose which changes with respect to time. Finally, I used movie syntax for saving the pendulum animation result as a video. Generally, the function script must be at the end of the program. This is what I get as an error. Other than this, all others are finding simple errors like missing letters as per assigned in input.

Result: https://youtu.be/sDDI7RnIBhE