Solve Second Order ODE of Simple Pendulum and Animation of its motion

AIM: Write a program in MATLAB that solves the second order ODE of pendulum and animate the motion of pendulum.

INTRODUCTION: The Differential equation has many applications in various engineering applications and science. It is very important for engineers to be able to generate physical problems using mathematical calculations, and then solve these calculations to further the behavior of the systems involved.

How the pendulum moves depends on Newtons' second law. When this rule is written, we find the second-order Ordinary difference Equation that describes the position of the "ball" in relation to time.

GOVERNING EQUARIONS & GRAPHS:

`(d^2 theta)/(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

MATLAB have ODE45 command to solve the ODE. But for the second order ODE first we have to convert into first order ODE.

By converting this two equation in a matrix from,we get

Code:

clear all 
close all 
clc

                                                                    %inputs
b=0.05;                                               %damping co-efficient
l=1;                                           %Length of the Pendulum in m
m=1;                                                %mass of the ball in Kg
g=9.81;                                                   %gravity in m/s^2

                                                             %initial vaule
theta_0 =[0;3];
%initial anglar displacement is zero and initial angular velocity is 3 rad/s

                                                               %time points
t_span = linspace(0,20,600);               %total time for motion is 20 sec

                                                                 %solve ODE
[t, results]= ode45(@(t,theta) ode_func(t,theta,b,g,l,m),t_span,theta_0); 
                                    %ode45 to solve ODE with using ode_func

                                    
                                    
hold on
plot(t,results(:,1))
plot(t,results(:,2))

xlabel('Time')
ylabel('Anglar displacement & Angular velocity')

legend('Angular displacement','Angular velocity')
title('Transient behaviour of Simple Pendulum')
hold off

figure(2)
ct = 1
for i = 1:length(t_span)
    
    x0 = 0;
    y0 = 0;
    
    x1 = 1*sin(results(i,1));
    y1 = -1*cos(results(i,1));
    
   
    plot([-1 1],[0 0],'linewidth',10,'color','k')
    hold on
    plot([0 0],[0 -1.2],'linestyle','--','color','k')
    
    plot([x0 x1],[y0 y1],'linewidth',5,'color','r')
    plot(x1,y1,'marker','o','markerFacecolor','b','markersize',30)
    title('Simple Pendulum motion')
    
    axis([-1.5 1.5 -1.5 0.2])
    pause(0.0001)
   
    
    M(ct)= getframe(gcf);
    hold off
    ct=ct+1;
    
end

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

Code for sub-function:

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

Code explanation:

First, we plot the behavior of the pendulum for a moment rather than create images of a simple pendulum movement.

Step description of the pendulum behavior graph correction code:

First, we must provide the necessary input that measures (b), pendulum height (L), ball weight (m) and vitality (g).

We should give the theta value in the matrix form because it contains 2 values ​​namely angular velocity and angular velocity.

We have to allocate the estimated time and convert this space into an equivalent space with the use of the linspace command.

We should use the ode45 command to solve the first ODE with the help of the sub function (ode_func) by assigning the angular and angular velocity distribution at any time.

In implementation we solve our 2 equation by adding a given.

Now the plot is a graph of the angular displacement of the ball in relation to time and angular velocity in relation to time.

We plot these two graphs numerically so we use the myth command to easily understand the angular displacement and the angular velocity curve. Using the order title gives the graph title

The first of the pendulum's display distances in a different number. we set the sum.

We apply a loop to the change of display with displacement in relation to time.

Respond to a string connection at the end of a constant (x 0 , y 0 = (0,0)

We need to calculate the free end-of-wire coordinates (x 1 y 1 ) according to the number of time-dependent angles.

We have to set the pivot line the thread is on.

After that we have to arrange another line to show the ball balance.

Now edit the graph string.

We have to point the ball to the animation so we use the marker to edit the theme.

We must give the correct value for the X and Y axis.

To keep the rest in the ball we use the stop order.

We have to create a new M variable to keep the current mathematical frame.

We must provide an addendum to resist the command.

Then after the loop, we have to make a conversion movie that is used to save the frame. Therefore, the final images are constructed.

ERROR:

Many errors occurs due to typing mistakes, I manage to remove all those by reviewing the code multiple times.

RESULT:

Go through this link for animation:

https://youtu.be/0rGKbKK62Q4