Simple Pendulum by Solving Second Order ODEs

AIM: To ceate an aniamtion of simple pendulum by solving second oder ODE function for a time inteval of 0-20secs with angular velocity of 3m/s at t=0.

INTRODUCTION:

Simple pendulum equation

`(d^2 theta)/(dt^2 )+b/m (d*theta)/(dt)+g/L sin⁡theta=0` (*)

where,

g = gravity in m/s2,

L = length of the pendulum in m,

m = mass of the ball in kg,

b = damping coefficient.

For solving this ODE function in MATLAB, we have to convert the above second order ODE into two first order ODE which will then be specified to two column vectors. These two columns vectors would be used in the function 'ode45' that is predefined in MATLAB

Assuming,

`theta= theta_1`

we get,

`(d*theta)/(dt)= (d*theta_1)/dt= theta_2` Equation (1)

`(d^2theta)/(dt^2)=(d^2 theta_1)/(dt^2)=d/dt((d theta_1)/dt)=d/dt(theta_2)` Equation (2)

Substituting these two equations in the main equation we get,

`(d theta_2)/(dt)+b/m theta_2+g/L sin theta_1 =0` Equation (3)

Thus, the second order ODE is conveted into first oder which can now be used in the function.

PROCEDURE:

• Write the function for the solved first order ODE which would be used as calling function within the pre-existing ode45 function.
• Specify the time interval using linspace command.
• Specify the inital condition for the simple pendulum, angular velocity = 3m/s at(displacement=0,t=0).
• Specify the input parameters L g m b.
• The results are stored in two columns denoting displacement in column 1 and angular velocity in column 2 which is the output assigned from the ode45 function.
• Using Plot command we can generate the graph seperately for both the parameters with respect to time.
• These graphs can be made to appear in the same window using a subplot command.
• In subplot three animation of the simple pendulum is constructed with the help of for loop function
• The result of angular velocity is assgined to a variable called angle which is then used to determine the length of the counter variable in for loop
• The fixed point for one end of the pendulum is specified in co-ordinates x1 and y1
• The moving end of the pendulum is specified in co-ordinates x2 and y2 using a mathematical expression such the movement is along the specified angle toa certain distance that is generated by the for loop function.
• These co-ordinates are plotted and the ends of the pendulum is displayed using a marker command.
• The axis are used to place the pendulum animation in a graph and pause command determines the speed of the animation.
• These frames ae captured using a get frame command and stored in a variable called M(ct)
• Then we can take the output of these frames from M using the movie command and videowriter command.
• We can also specify the required video format and specify the file name.

CODE 1: ODE Function

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

CODE 2: Main Function

clear all
close all
clc

% Simple pendulum equation:
% (d^2(theta)/dt^2)+((b/m)*(d(theta)/dt))+((g/L)*sin(theta))=0

%Input
b=0.
g=9.81
l=1
m=1

%Initial Conditions,
%Time (0-20secs)
t_span= linspace(0,20,500) 
%Displacement=0, Angular Velocity= 3 (at t=0)
theta_0= [0,3];

% %ODE
[t,result]= ode45(@(t,theta) pendulum_func(theta, b,g,l,m),t_span,theta_0);

%Plotting Angular Displacement and Angular Velocity
subplot(3,1,1)
plot(t,result(:,1),'color','g')
title('Simple Pendulum')
xlabel('Time [s]')
ylabel({'Angular Displacement';'[rad]'})
legend ('Time','Angular Displacement')

subplot(3,1,2)
plot(t,result(:,2),'color','b')
title('Simple Pendulum')
xlabel('Time [s]')
ylabel({'Angular Velocity';'[rad/s]'})
legend ('Time','Angular Velocity')

%Animation of simple pendulum
subplot(3,1,3)
angle=result(:,2);
ct=1;
for i=1:length(result(:,2))
     % Co-ordinates of fixed points
   x1=10;
   y1=20;
     % Co-ordinates of second points
   x2=10-15*sind(angle(i))
   y2=20-15*cosd(angle(i))
   
   plot([x1 x2],[y1 y2],'-o','markerfacecolor','r')
   axis([0 20 0 20])
   pause(0.03);
   M(ct)=getframe(gcf);
   ct=ct+1;
end

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

OUTPUT:

The video of the output is specified in this youtube link,

https://youtu.be/7oylWM2o_1g

ERRORS:

The movement of the pendulum was not steady, this was because the angles specified in x2 and y2 cordinates. It should be specified as sind and cosd.

This error occurred because the paramter l was entered as L in one point of the calling function. Thus the parameters used did not match with the other places it had been used which caused the error.

RESULT:

Thus, the the ODE is solved and the animation for a simple pendulum is generated.