Week 3 - Solving second order ODEs

Objective :-

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.

Explanation :-

What is ODE's :-

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.

Given Inputs :-

L=1 metre,

m=1 kg,

b=0.05.

g=9.81 m/s2.

Step 1 -

Convert second order ODE to first order

Lets

θ = θ1

so, dθ/dt = dθ1/dt = θ2 ..........1

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

Hence d^2θ/dt2 = dθ2/dt from Eq. 2

and dθ1/dt = θ2 from Ea. 1

So,

dθ2/dt+(b/m)*θ2+(g/l)*sinθ = 0

dθ2/dt = -(b/m)*θ2-(g/l)*sinθ .....3

Step 2 -

Create a function of ODE

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

Step 3 -

Create a simple pendulam program by using ode45 Function

- initial condition provided is angular displacement=0,angular velocity=3 rad/sec.

-Time interval 0-20 Sec. in 500 steps.

-Solving ODE equation by using ode45 function

clear all
close all
clc

b = 0.05;
g = 9.81;
l = 1;
m = 1;

% Initial condidition
theta_0 = [0;3];

% Time points
t_span = linspace(0,20,500);

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

% plot
plot(t,results(:,1));

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

Result graph -

In above graph blue line shows the angular displacement of pendulam w.r. to time as we see initial angular displacement is 0

and red line shows the angular Velocity as we see initial angular velocity is 3m/s and goes changing with time.

Final Program to create aniamtion of pendulum movement -

clear all
close all
clc

b = 0.05;
g = 9.81;
l = 1;
m = 1;

% Initial condidition
theta_0 = [0;3];

% Time points
t_span = linspace(0,20,100);

% slove ODE
[t_span, results] = ode45(@(t_span,theta) ode_Function(t_span,theta, b,g,l,m),t_span, theta_0);

% plot
figure(1)
hold on
plot(t_span,results(:,1));
plot(t_span,results(:,2));
xlabel('time')
ylabel('plot')
saveas(1,'Time_plot.png')

% Animation
X0 = 0;
Y0 = -0.5;
X1 = -2;
Y1 = -0.5;
X2 = 2;
Y2 = -0.5;
ct=1;
for i = 1:length(t_span)
    X = l*sin(results(i));
    Y = -l*cos(results(i));
    plot([X0,X1,X2],[Y0,Y1,Y2],'linewidth',8)
    hold on
    plot([X0,X],[Y0,Y],'linewidth',2)
    hold on
    plot(X,Y,'o','markersize',15,'markerfacecolor','b')
    axis([-3 3 -3 0])
    hold off
    pause(0.1);
    d(ct)= getframe(gcf);
    ct=ct+1;

end

movie(d)
videofile = VideoWriter("Simple_pendulum.avi",'compressed AVI');
open (videofile)
writeVideo(videofile,d)
close(videofile)

Error Faced during Programing :-

Unrecognised function or variable 'ode_func' error occurse due to Function not define so we define function to overcome this error.

Same error occurse after defining fuction for (t,theta) due to we used Function_ode to define function and in program we recall it as func_ode. so to overcome we define same name in program.

Saved file name with different than function name causing error occurse, corrected it by renaming as a function name.

Link of Google Drive for AVI File -

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