SOLVING SECOND ORDER ODE's

AIM: TO SIMULATE THE TRANSIENT BEHAVIOUR OF A SIMPLE PENDULUM AND TO CREATE AN ANIMATION OF IT’S MOTION USING MATLAB

THEORY: In Engineering, ODE is used to describe the transient behaviour of a system. A simple example is a pendulum.

The way the pendulum moves depends on the Newtons second law. When this law is written down, we get a second order Ordinary Differential Equation that describes the position of the "ball" w.r.t time.

GOVERNING EQUATIONS:

• Ordinary Differential Equation (ODE) representing equation of motion of a simple pendulum with damping –

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

where,

g = gravity (m/s²)

L = length of the pendulum (m)

m = mass of the ball (kg)

b = damping coefficient

OBJECTIVES:

1. To study the transient behaviour of the pendulum system.
2. Animate the motion of the simple pendulum.

INPUT VALUES:

L = 1 m

m = 1 kg

b = 0.05

g = 9.81 m/s²

FUNCTION PROGRAM:

function [dtheta_dt] = ode_func(t,theta,b,g,L,m)
 
theta1 = theta(1);
theta2 = theta(2);
 
% Second order ODE broken down into two single order ODE for solving
dtheta1_dt = theta2;
dtheta2_dt = -(b/m)*theta2 - (g/L)*sin(theta1);
 
% The array consisting of two single order ODE's
dtheta_dt = [dtheta1_dt; dtheta2_dt];
end

MAIN PROGRAM:

clear all
close all
clc
 
% Given Parameters
% Length of pendulum
L=1;
% Mass of bob
m=0.1;
% Damping Co-efficient
b=0.05;
g=9.81;
 
% Initial condition
%[Displacement;Velocity]
theta_0 = [0;3];                                   

% Time points distributed over an array
t_span = linspace(0,20,200);
 
% Solve ODE
% Using inbuilt ODE solver
[t, results] = ode45(@(t,theta) ode_func(t,theta,b,g,L,m),t_span,theta_0);
 
% Plotting
figure(1)
 
% Plotting Displacement v/s Time
plot(t,results(:,1))
hold on
 
% Plotting Displacement v/s Time
plot(t,results(:,2))
 
title ('TRANSIENT BEHAVIOUR')
xlabel('TIME')
ylabel('PLOT (Displacement & Velocity)')
legend('Displacement','Velocity')
 
 
% Animating the motion of Simple Pendulum
figure(2)
ct = 1;
for i = 1:length(t_span)
    RESULTS = results(i);
    
    % Coordinates of start point
    x0 = 0;
    y0 = 0;
    
    % Coordinates of end point
    x1 = L*sind(RESULTS);
    y1 = -L*cosd(RESULTS);
    
    % Plotting
    % Plot representing the massless thread
    plot([x0 x1],[y0 y1],'linewidth',2,'color','g');
    hold on
    
    % Plot representing ceiling from which the pendulum hangs
    plot([-0.5 0.5], [0 0],'linewidth',5,'color','b');
    
    % Plot representing the bob
    plot(x1,y1,'marker','o','markersize',15,'markerfacecolor','r');
    
    title('SIMPLE PENDULUM MOTION')
    axis([-1.5 1.5 -1.5 0.5]);
    pause(0.00005);
    
    % Hiding the axes of the plot
    set(gca,'XTick',[], 'YTick', [])
    
    M(ct) = getframe(gcf);
    hold off
    ct=ct+1;
end
 
% Animating the motion of simple pendulum
movie(M)
videofile = VideoWriter('ode_main_vid','Uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)

EXPLANATION OF PROGRAM:

• Initially, writing a function program where a specified function is mentioned covering all the parameters and the variables. Remember that the name of the program should be saved with the function name.
• The given values are specified and the initial condition is imposed.
• Time span over which the pendulum’s transient behaviour has to be observed is also defined.
• Now, the ODE solver is used to calculate the results.
• Further, we plot the Displacement v/s Time and Velocity v/s Time graphs on the same figure to study the respective behaviours.
• After the plotting, animation of the pendulum motion is done using Movie Command and specifying the array over which the pendulum traverses.

TRANSIENT BEHAVIOUR:

ANIMATION: The link for the animation has been mentioned below.

Simple_Pendulum_Motion_Animation

CONCLUSION: Establishing a clear view of the transient behaviour of the simple pendulum by solving the second order ODE defining its motion. Also, animation of the pendulum motion is observed by converting the successive frames into a movie.