SOLVING 2ND ORDER ODE REPRESENTING THE MOTION OF A SIMPLE PENDULUM USING MATLAB

AIM:

To write a MATLAB program to solve a ODE which describes the motion of the simple pendulum with damping and to animate the motion of the same.

The ODE that describes the motion of simple pendulum is

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

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.

The values considered are

L=1 metre,

m=1 kg,

b=0.05.

g=9.81 m/s2.

The time span is taken ass 0 seconds to 20 seconds. Also angular velocity at time = 0 is 3 rad/s.

In above equation, the symbol \'theta\' is the angular displacement.

CONVERSION OF 2ND ORDER ODE TO SYSTEM OF FIRST ORDER ODE\'S:

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

In the above equation, consider

`theta = theta1`

so,

`(d theta )/dt = (d theta1 )/dt = theta2` ...................................................................eqn (1)

Now the original eqn becomes,

`(d theta2)/dt = -(b/m)(theta2) - (g/ L)sin theta1` ....................................................eqn (2)

So we got two firt order equations eqn (1) and eqn (2).

MATLAB PROGRAM:

Those two equations are being written in matlab program to be solved.

The ODE function is a follows

function ydot = simple_pendulum(t,y,L,m,b,g)


y1=y(1);
y2=y(2);

dy1_dt=y2;

dy2_dt=(((-b*y2)/m) - (g*(sin(y1)))/L);

ydot = [dy1_dt ; dy2_dt];

end

NOTE: In the above snippet, instead of the symbol \'theta\' , the variable y is used.

The main program is

clear all
close all
clc

%Equation of simple pendulum
%(d^2)y/dt^2 + b*dy/m*dt + g*siny/L = 0

%g = gravity in m/s2,0
%L = length of the pendulum in m,
%m = mass of the ball in kg,
%b=damping coefficient,
%y = angular displacement in radians,
%t = time in seconds.

L=1;

m=1;

b=0.05;

g=9.82;

%Providing the number of time steps and initial conditions

t = linspace(0,20,500);
y0 = [0 3];

%Solving the ode

[t, results] = ode45(@(t,y)simple_pendulum(t,y,L,m,b,g), t,y0);


figure(1)
plot(t,results(:,1),\'r\')
hold on
plot(t,results(:,2), \'b\')

xlabel(\'Time span\')
ylabel(\'Angular spec\')

legend(\'Angular displacement\',\'Angular velocity\')

%pivot points
xo=0;
yo=0;
figure(2)
n=1;
for i = 1:500

    x1 = L*sin(results(i,1));
    y1 = -L*cos(results(i,1));
    hold off
    plot([-1,1],[0,0])
    hold on
    plot([xo x1],[yo y1])
    hold on 
    plot(x1,y1,\'marker\',\'o\',\'linewidth\',7)
    axis([-2 2 -2 2])
    pause(0.01)
    p(n)=getframe(gcf);
    n=n+1;
    
end

movie(p)
videofile=VideoWriter(\'Simple_pendulum.avi\',\'uncompressed AVI\');
open(videofile)
writeVideo(videofile,p)
close(videofile)

The timespan and initial conditions are being defined first. And then using the MATLAB inbuilt function called \'ode45\' is used to solve the two ode\'s which we derived previously. The function is called then with the appropriate syntax for ode45 function.

In the figure 1, the results are plotted

From the above plot, the decreasing fluctuation of the angular displacement and angular velocity shows the simple harmonic motion of the pendulum.

ANIMATING THE MOTION OF PENDULUM:

Then in figure2, the animation file is created.

The origin is kept as the pivot points.

xo=0;
yo=0;

x1,y1 are the coordinate points of the pendulum bob. The number of frames id defined in the variable n and it will be increasing after each time step. Thus a number of time step frames are compiled into a single movie file which is stored in variable p.

Then using the Videowriter function, the movie is written with its name and file type and saved and uploaded to youtube.

Thus the 2nd order ODE is converted to system of 1st order ODE\'s and is solved using MATLAB and the animation file representing the motion of Simple pendulum is being created.