SIMULATION OF A SIMPLE PENDULUM BY SOLVING 2ND ORDER DIFFERENTIAL EQUATION USING MATLAB

AIM: To simulate the transient behaviour of a simple pendulum and to animate the oscillation using MATLAB.

OBJECTIVE:

1. To execute the program using 2nd order differential equation for calculating displacement and velocity with respect to time in MATLAB.

2. To animate the oscillation of a simple pendulum

3. To plot displacement and velocity with respect to time.

THEORY:

A simple pendulum consists of a large mass (bob) suspended by a light string from a rigid support . The length of the string is large compared to the dimensions of the bob. A stationary pendulum has its string in a vertical position with the bob hanging straight down. This is known as the pendulum's equilibrium position. Oscillations are set up in a pendulum by displacing the bob from its equilibrium positions and releasing it.

In this simple pendulum ,the oscillation comes to rest or equilibrium position after some time. this is because of the damping effect.

Damping is an influence within or upon oscillatory system that has the effect of reducing, restricting or preventing its oscillations. In a physical system, damping is produced by process that dissipate the energy stored in the oscillation.

In addition to this the same pendulum also describes the motion of spring, mass and damper sysytem.

2ND ORDER DIFFERENTIAL EQUATION:

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

here, assumed values are

b = damping co-efficient = 0.05

m = mass of bob = 1Kg

g = acceleration due to gravity = 9.81`m/s^2`

l = length = 1m

here we assume,

`theta=theta_1`then,

`(d theta)/dt=(d theta_1)/dt=theta_2`

`(d^2 theta)/dt^2=d/dt[(d theta)/dt]=(d theta_2)/dt`

by substituting these in 2nd order ODE equation, we get

`(d theta_2)/dt+b/m*theta_2+g/l *sin(theta_1)=0`

In order to solve this original 2nd order ODE equation, we have to break down this equation into series of 1st order ODE's.

1st ODE equation = `(d theta)/dt = theta_2`

2nd ODE equation = `(d theta_2)/dt= -b/m*theta_2-g/l *sin(theta_1)`

MAIN PROGRAM:

clear all
close all
clc

% inputs
b =0.6;
g = 9.81;
l = 1;
m = 1;

% initial condition
theta_0 = [0;3];

% time points
t_span = linspace(0,10,250);

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

figure(1)
plot(t,results(:,1))
hold on
plot(t,results(:,2))
ylabel('plot')
xlabel('time')
legend ('displacement','velocity')
hold off
grid on

% counter

ct = 1;
for i=1:length(results)
    THETA = results(i,1);
    x0=0;
    y0=0;
    x1=l*sin(THETA)
    y1=-l*cos(THETA)
    figure(2)
    plot(x1,y1,'o','markersize',20,'markerfacecolor','r')
    hold on
    plot([x0 x1],[y0 y1],'linewidth',3)
    plot([-.5 .5],[0 0],'linewidth',10)
    axis ([-1.5 1.5 -2 .5])
    pause(0.03)
    hold off
    M(ct)=getframe(gcf);
    ct=ct+1;
    
end

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

FUNCTION PROGRAM:

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

PLOT RESULTS:

SIMULATED VIDEO:

RESULTS AND CONCLUSION:

1. Dsplacement and Velocity are calculated and ploted as shown above.

2. Animation of simulated 2nd order ODE is done and shown above.

3. Program is executed by breaking 2nd order ODE into 1st order ODE using MATLAB.