Week 3 - Solving second order ODEs

AIM : To solve the ODE which represents the equation of motion of a simple pendulum with samping.

OBJECTIVES : Write a program that solves the following ODE. This ODE represents the equation of motion of a simple pendulum with damping. And also write a program to simulate the transient behavior of a simple pendulum and to create an animation of its motion.

Theory: - In Engineering, ODE is used to describe the transient behavior 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. Similarly, there are hundreds of ODES that describe key phenomena and if you have the ability to solve these equations then you will be able to simulate any of these systems.

GOVERNING EQUATIONS:

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 (ode) that describes the position of the "ball" w.r.t time.

The 2nd order ode below represents the equation of motion of a simple pendulum with damping

where, b-damping coefficient.

m = mass of the ball in kg,

g= gravity in m

L= length of the pendulum in m.

Since the given function is in the second-order differential equation we will convert it into first order to solve the ODE.

Inputs

• L=1 metre

• m=1 kg

• b=0.05

• g=9.81 m/s^2

• Time span= 0-20 sec

• angular displacement=0

• angular velocity=3 rad/sec2 at time t=0

CODE

clear all 
close all
clc

b=0.02;
g=9.81;
l=1;
m=0.1;

% initial condition
theta_0 = [0;5];

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

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

plot(t,results(:,1))
hold on
plot(t,results(:,2))
ylabel('plot')
xlabel('time')

% Iteration using for loop to create animation of simple pendulum

ct=1;
theta_1 = results(:,1);

for i=1:length(results (:,1))

% initial co-ordinates

x_start= 0;

y_start= 0;

theta_2 = theta_1(i);

% final co-ordinates

x_end = 1*sin(theta_2);

y_end = -1*cos(theta_2);

% ploting for animation 
figure(2)
% base of pendulum
plot([-1,1],[0,0],'linewidth',4,'color', 'black')
hold on 
% string of pendulum
line([x_start x_end],[y_start y_end], 'linewidth',6, 'color', 'R')
hold on

% bob of pendulum

plot(x_end,y_end, 'marker', 'o', 'markersize', 28, 'markerfacecolor', 'Y')

axis([-2 2 -2 2])

pause(0.03)

hold off

M(ct)=getframe(gcf);
ct=ct+1;

end

% creating movie from captured frames
movie(M)
videofile = VideoWriter('simple_pendulum.avi','uncompressed AVI');
open(videofile);
writeVideo (videofile,M);
close(videofile);

CODE FOR FUNCTION PROGRAM

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/1)*sin(theta1);
dtheta_dt = [dtheta1_dt;dtheta2_dt];
end

Detail explanation each step of the program:

1. In the first step Input and initial condition of a simple pendulum are fed.

2. Then the program is written to solve the ODE in which Matlab inbuilt function (ode45) is used with constant theta, t, and another function is used, which is defined separately in function program.

3. In the next step result of solved ODE is plotted w.r.t time using plot command.

4. Now iteration of co-ordinate is fed using for loop to create animation.

5. Then plot command is used to plot pendulum string and bob.

6. In the next step, the move command is used to create an animated video of the oscillating pendulum.

Output:

The output obtained is graphical representation velocity curve and displacement curve w.r.t Time. Angular Displacement and Angular Velocity plots were generated for a time range of 20 secs.

Result ;

Conclusion ;

Conclusion:

The x_start,y_start, x_end y_end values a

x_start =

0

y_start =

0

x_end =

-0.7779

y_end =

-0.6283

ANIMATION VIDEO LINK

Here is the link for the video of the animation that I got for this project.

https://youtu.be/t09utcGCMiw