Solving the Simple pendulum ODE and animating the motion of it using MATLAB

Aim: To write a program that simulates the transient behaviour of a simple pendulum and to create an animation of its motion.

Governing equations:

To implement the second-order equation of the simple pendulum, it is converted into two first-order ODEs.

`d^2/dt^2(theta)+b/m+g/l(sintheta)` where

b is the damping coefficient

g is the gravitational acceleration in m/sec^2

l is the length of the pendulum in metre

m is the mass of the bob in kg

Objectives:

1. To write a MATLAB program to implement the simple pendulum ODE.
2. To animate the motion of the pendulum.

MATLAB Program:

The below function is for the splitting up of the second-order simple pendulum ODE to two first-order ODEs.

Function Code:

function[dtheta_dt]=ode_function(t,theta,b,g,l,m);
theta1=theta(1);%The first term of theta vector
theta2=theta(2);%The second term of theta vector
dtheta1_dt=theta2;%when the second order ODE is split into two first order ODE
dtheta2_dt=-b/m*theta2-g/l*sin(theta1);
dtheta_dt=[dtheta1_dt;dtheta2_dt];%the array of the displacement and velocity 
end

MATLAB Code:

The below program is for the plotting displacement and velocity and for the animation of the motion of the simple pendulum.

Here ODE 45 which is the inbuilt solver configuration in MATLAB is used for solving the ODE.

Figure(1) is for the plot of displacement and velocity and Figure(2) is for the animation of the simple pendulum.

clear all 
close all
clc

b=0.05;%damping coefficient (no units)
g=9.81;%gravitational force in m/sec^2
l=1;%length of the pendulum in metres 
m=1;%mass of the ball in kg

%initial conditions
theta_0=[0,3];
 
%points of time
time_span=linspace(0,20,500);

%solving the ODE
[t,results]=ode45(@(t,theta)ode_function(t,theta,b,g,l,m),time_span,theta_0);
%ode45 is a programmed solver in MATLAB and the required parameters are
%given
 
figure(1)%this is required as the plots need to be in two different graphs
plot(t,results(:,1))%displacement with time plot
hold on%this command is used to plot multiple plots in the same graph
plot(t,results(:,2))%velocity with time plot
xlabel(\'time\')%to label the x-axis
ylabel(\'displacement and velocity\')%to label the y-axis

figure(2)%this is required as the plots needs to be in two different graphs
z=1;%count of the frame
for i=1:length(time_span)%loop inorder for the motion of simple pendulum
    %initial points 
     x0=0;
     y0=0;
     
     x1=1*sin(results(i,1));%final points
     y1=-1*cos(results(i,1));
      
     plot([-1,1],[0,0],\'linewidth\',2,\'color\',\'b\')%for the base line
      hold on%this command is used to plot multiple plots in the same graph
     
      plot([x0,x1],[y0,y1],\'linewidth\',2,\'color\',\'b\')%for the pendulum string
      
     
     plot(x1,y1,\'marker\',\'o\',\'markerfacecolor\',\'g\',\'markersize\',20)%for the bob
     
     axis([-2,2,-2,1])%limits of x-axes and y-axes
     
     
     pause(0.001)%each frame will pause for 0.001 seconds
     
     M(z)=getframe(gcf);%getframe command is used to get each frame and will be added to the array \'M\'
     hold off%for each frame to be clearer (or) overlapping of frames doesn\'t occur
     
     z=z+1;
end
movie(M)%to create a movie 
 videofile = VideoWriter(\'motion_of_simple_pendulum.avi\',\'uncompressed AVI\');%to specify the movie name and type of video
 open(videofile)%opens the videofile
 writeVideo(videofile,M)%to play the video
 close(videofile)%to close the videofile

The plot below is for the displacement and velocity of the simple pendulum with respect to change in time in seconds.

The link below is for the video of a simple pendulum in oscillation.

Youtube link:https://youtu.be/K-oofvptCgg

Errors Encoutered:

1. The error above is a typo error of \'markersize\'.
2. In the animation process, hold off command was not used initially. So, the video was not proper, then it was rectified
3. The bob was in the opposite axis i.e the bob was above the baseline initially and then when the cosine was made negative in the code, it was in the correct direction i.e below the baseline.

Conclusion:

Initially, there was a difficulty in doing the animation then I found the way of doing it. The oscillation of simple pendulum with the given parameters is shown in the video which is generated from MATLAB.