Solving second order ODE for pendulum motion

Objective - Write a program to simulate transient behaviour of simple pendulum and create animation.

Theory-

In Engineering, ODE is used to describe the transient behavior of a system. A simple example is a pendulum

The way the pendulum moves depends on the Newtons second law. When this law is written down, we get a second order Ordinary Differential Equation that describes the position of the \"ball\" w.r.t time.

Similarly, there are hundreds of ODEs that describe key phenomena and if you have the ability to solve these equations then you will be able to simulate any of theses systems.

Your objective is to write a program that solves the following ODE. This ODE represents the equation of motion of a simple pendulum with damping.

In the above equation,

g = gravity in m/s2,

L = length of the pendulum in m,

m = mass of the ball in kg,

b=damping coefficient.

Inputs-

L=1 metre,

m=1 kg,

b=0.05.

g=9.81 m/s2.

Simulate the motion between 0-20 sec,

for angular displacement=0,

angular velocity=3 rad/sec at time t=0.

Program-

A]Function-

Step 1) first we create function of ODE , this equation divide into 2 parts i.e. [dtheta1_dt=theta2]

[dtheta_2dt=-(b/m)*theta2 - (g/l)*sin(theta1)]

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

B]Main Program-

Step 1- Initially provide all inputs like bore, gravity , length ,mass and angular velocity is 3 and angular displacement is 20.

step 2- In step 2 use ode45 command to solve ode equation , and call the function of ODE.

step 3- from the above equation we plot the graph of angular displacement and angular veocity.

step 4 -ploting fig of pendulum we use the command for creating link x_end=l*sin(results(i,1)); and
y_end=-l*cos(results(i,1));

step 5- for ploting base the command is,

plot([x_base x_end],[y_base y_end],[-0.5 0.5],[0 0],\'linewidth\',3)

and for pendulum,

plot(x_end, y_end,\'Marker\',\'o\',\'MarkerSize\',30,\'MarkerFacecolor\',\'g\')

step 6- Now we fixing the axis x=[-1.5 , 1.5] and y=[-1.5 , 1.5].

and time interval between frame is 0.01 sec.

step 7- Now use videowriter command for save the animation file.

clear all
close all
clc

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

%initial codition
theta_0 =[0;3];

%time points
t_span= linspace(0,20,500);

%solve ode 
[t, results] = ode45(@(t,theta) ode_function(t,theta,b,g,l,m),t_span,theta_0);
figure
plot(t,results(:,1))
hold on
plot(t,results(:,2))
xlabel(\'time\')
ylabel(\'Angular displacement vs angular velocity\')
legend(\'Angular displacement\',\'angular velocity\')
hold off
grid on

%ploting pendulum
figure

ct=1;
for i=1:length(t_span)
    x_base = 0;
    y_base = 0;
    x_end=l*sin(results(i,1));
    y_end=-l*cos(results(i,1));
plot([x_base x_end],[y_base y_end],[-0.5 0.5],[0 0],\'linewidth\',3)
hold on
plot(x_end, y_end,\'Marker\',\'o\',\'MarkerSize\',30,\'MarkerFacecolor\',\'g\')
axis([-1.5 1.5 -1.5 0.5])
grid on
pause(0.01)
hold off

M(ct) = getframe(gcf);
ct=ct+1;
end

%video file
movie(M);
videofile = VideoWriter(\'ode_pendulum.avi\',\'Uncompressed AVI\');
open(videofile)
writeVideo(videofile,M)
close(videofile)

Plots-

Angular velocity and displacement vs time.

Pendulum-

youtube link-

error facing-

when program is run , showing this error . at line 19 i forget to insert (s) in result. for solving this error i use \'results\' instead of \'result\'.

Conclusion- The transient behaviour of pendulum motion with respect to time is successfully run by using ODE.