Week 3 - Solving second order ODEs

AIM:Write a program in Matlab\Octave that will simulate the pendulum motion.

Theory:A pendulum is a device that is found in wall clocks. It consists of a weight (bob) suspended from a pivot by a string or a very light rod so that it can swing freely. When displaced to an initial angle and released, the pendulum will swing back and forth with a periodic motion. By applying Newton’s second law of motion for rotational systems, the equation of motion for the pendulum may be obtained.

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 theses systems.

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

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.

Write a program in Matlab\Octave that will simulate the pendulum motion, just like the one shown in the start of this challenge.

use,

L=1 metre,

m=1 kg,

b=0.05.

g=9.81 m/s2.

Here we are simulating the motion between 0-20 sec, for angular displacement=0,angular velocity=3 rad/sec at time t=0.

code:

ODE Function code:

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 explaination:

1) In this code the Second order ODE is Simplified into first order ODE, by introducing the terms theta1 and theta2.

2)Here the function declare ode function for t, theta, b, g,l,m.

3)Then this theta1 and theta2 are diffrentiated and dtheta_dt gives a column Matrix.

Main code:

clear all
close all
clc

b=0.05;
g=9.81;
l=1;
m=1;
%initial condition
theta_0=[0;3];

%timepoints 
t_span=linspace(0,20,200)

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

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

ct=1;
for i=1:length(result(:,1))                    
    x0=0;
    y0=0;

x1=l*sin(result(i,1));
y1=-l*cos(result(i,1));

figure(2)
%plotting 
plot([x0 x1],[y0 y1],'linewidth',3)
axis([-2 2 -2 2])
hold on
plot([-1 1],[0 0],'linewidth',4)

plot(x1,y1,'.','MarkerSize',30)

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

end
    movie(M)
    videofile=VideoWriter('ode_function_for_loop.avi','uncompressed AVI');
    open(videofile)
    writeVideo(videofile,M)
    close(videofile)

code explaination:

1)First in the Editor window the parameters for b,g,l,m are given

2)Here the theta_0 represents angular displacement i.e is 0 and angular velocity i.e 3 rad/sec.

3)In the t_span 200 values are taken between 0 to 20 time period .

4)Now Using ODE 45 Command ODE function is given.

5)Now using plot Command angular displacement and angular velocity is plotted.

6)Now we are using for loop command as 'i' i.e for angular displacement.

7)For plotting the bob we are using point i.e '.' and increasing size by MarkerSize upto 30.

8)Using the Movie Command the plot are saved in uncompressed avi.

Plot:

Animation of simple pendulum:

https://drive.google.com/file/d/1vJADDxdJafq5lZbAiAFEX2GyL6NNSwbI/view?usp=sharing

Result:

here ,we simulated simple pendulum using ODE function.