Writing a program in Matlab

Project report

Aim:

Writing a program in Matlab that will simulate the pendulum motion.

Objectives:

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

`(d^2 θ)/(dt^2 )+ b/m (dθ)/(dt)+ + g/L sin⁡θ=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.

Theory:

A pendulum is a weight suspended from a pivot so that it can swing freely. When a pendulum is displaced sideways from its resting, equilibrium position, it is subject to a restoring force due to gravity that will accelerate it back toward the equilibrium position. When released, the restoring force acting on the pendulum\'s mass causes it to oscillate about the equilibrium position, swinging back and forth. The time for one complete cycle, a left swing and a right swing, is called the period. The period depends on the length of the pendulum and also to a slight degree on the amplitude, the width of the pendulum\'s swing.

A simple pendulum consists of a mass m hanging from a string of length L and fixed at a pivot point P. When displaced to an initial angle and released, the pendulum will swing back and forth with periodic motion. By applying Newton\'s secont law for rotational systems, the equation of motion for the pendulum may be obtained

`τ=Iα= -mg sin⁡θL=mL^2 (d^2 θ)/(dt^2 ) ;`

After rearranging,

`(d^2 θ)/(dt^2 )+ g/L sin⁡θ = 0`

The second order differential equation of simple pendulum

`(d^2 θ)/(dt^2 )+ b/m (dθ)/dt+ + g/L sin⁡θ = 0 …………………………(A)

To solve the above second order differential equation in matlab it has to be converted into first order differential equation

Lets`θ= θ1`

`(dθ)/(dt)= (dθ1)/dt=θ2 ………………………….(1)`

`(d^2 θ)/(dt^2 )= (d^2 θ1 )/(dt^2 )= d/dt ( (dθ1)/dt)= d/dt( θ2) …………………………..(2)`

Substituting the above equations to (A)

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

`(dθ2)/dt+ b/m θ2+ g/L sin⁡θ1=0`

`(dθ1)/dt= θ2`

`( dθ2)/dt = -b/m θ2- g/L sin⁡θ1`

`d/dt [(θ1@θ2)]= [(θ2@-b/m θ2- g/L sin⁡θ1 )]`

Program:

Matlab program that will simulate the pendulum motion

clear all 
close all
clc

b=0.05;      % damping coefficient
g=9.81;      % gravity in m/s2
l=1;         % length of the pendulum in m
m=1;         %  mass of the ball in kg

%initial conditions
theta_0=[0;3];
%time_plots
t_span=linspace(0,20,500);

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


%Pendulum
count=1;
for i=1:length(results(:,1));
x0=0;`
y0=0;
x1=l*sin(results(i,1));
y1=-l*cos(results(i,1));
plot([-0.5,0.5],[0,0],\'linewidth\',3);
line([x0,x1],[y0,y1],\'linewidth\',3)
hold on
plot(x1,y1,\'o\',\'markers\',15,\'markerfacecolor\',\'r\');
grid on
axis([-1.5 1.5 -1.5 0.5]);
pause(0.01);
M(count)=getframe(gcf);
hold off
count=count+1;
end

movie(M)
videofile=VideoWriter(\'pendulum.avi\',\'Uncompressed AVI\')
open(videofile)
writeVideo(videofile,M)
close(videofile)

The function call:

The function is used for the differenciation of the given data

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

The output:

Step wise explaination of codes:

1. First the constants of the differential equations are defined, as they have given in the problem statement. The damping coefficient,accelaration due to gravity, length of the pendulum, mass of the ball is given and defined.
2. The initial conditions like angular displacement and angular velocity are assigned to the variable theta_0. The time span is given 0 to 20.
3. The solving differencial equation is by using the code :
   [t, results]=ode45(@(t,theta) ode_func(t,theta,b,g,l,m), t_span,theta_0);
   The ode45 is inbuilt differencial equation command is used to differenciate the function.
4. At the time function call is made to the function ode_func(t,theta,b,g,l,m), which will take the values of defined variables and calculates thecolumn matrics having two values.
5. By this the differenciation is done to the function.
6. The for loop is used to calculate the motion of pendulum at different time intervals and position of the pendulum is calculated using the formulas
   `x1=l*sin(results(i,1))`
   `y1=-l*cos(results(i,1))`
7. Markers command is used to creat the ball at the position of (x1,y1)
8. Getframe(gcf) to get all the images at every point and to assign it to variable M
9. Movie is created using the command movie.

Error occured and correction method followed

1. 
   
   The error occured line21 stating: Expression or statement is incorrect--possibly unbalanced (, {, or [.
   and is corrected by puting the colon (:)
   for i=1:length(results(:,1));
   
   
   
   
2. 
   
   Error is shown in the command window in the line of 16 , but after checking it is found to be all correct, so i checked the function ode_func( t,theta,b,g,l,m) there I found the error in the line 6 where i used colon(:) instead of using the semi-colon(;)
   dtheta_dt=[dtheta1_dt ;dtheta2_dt];
   
   
   
   
3. 
   
   Error occured in getting the correct animation of the pendulum motion as it is tracing the entire motional path.
   To correct this error hold off command is used.
   pause(0.01);
   hold off
   M(count)=getframe(gcf);
   
   

Result:

The required animational plot is obtained

https://youtu.be/XIyIsj2cK8g

Conclusion:

By this analysis the differential equation of the simple pendulum is studied. By altering different variables different results are obtained. By changing the damping factor, time taken to come to equilibrium position is noticed. The different commands used to creat the pendum animation is also understood. The use of differencial equations in matlab and their importance is also studied.