Solving 2nd Order ODEs

What is ODE?

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.

OBJECTIVE

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

CODE

clear all
close all
clc

b=0.2; %Damping coefficient
g=9.81; %gravity
l=1; %length
m=1; %mass

theta_0=[0;3] % 0 is the displacement and 3 is the velocity

t_span= linspace(0,10,100)

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

figure(1)
plot(t,results(:,1))
hold on
plot(t,results(:,2))
xlabel(\'Time ---->\')
ylabel(\'angular velocity ---->\')
legend(\'angular velocity\',\'angular acceleration\')


figure(2)
plot(results(:,1),results(:,2))
xlabel(\'angular velocity ---->\')
ylabel(\'angular acceleration ---->\')


ang_vel=results(:,1);
ct=1;
for i=1:length(ang_vel)
  ang_vel1=ang_vel(i)
  x0=0;
  y0=0;
  x1=1*sin(ang_vel1);
  y1=-1*cos(ang_vel1);
  
  figure(3)
  plot([-1 1],[0,0], \'linewidth\',2,\'color\',\'k\')
  hold on
  plot([0,0],[1 -1], \'linewidth\',1,\'linestyle\',\'--\')
  hold on
  plot([x0,x1],[y0 y1], \'linewidth\',1.5,\'color\',\'r\')
  hold on
  plot(x1,y1,\'o\',\'markers\',20,\'markerfacecolor\',\'y\')
  hold off
  axis([-1.5 1.5 -1.5 1.5])
  pause(0.00001)
  file_text=sprintf(\"pend%d.png\",ct)
  saveas (gca, file_text);
  ct=ct+1;
endfor

1st 3 lines set the environment ready

b,g,l,m are the predefined parameters

thete_0 is a variable that holds the displacement and the velocity

t_Span is a variable that takes 100 values between 0-10

[t, results] is a variable that stores the results obtained by the parameters once they pass through the function ode_func

ode45 uses t, theta ode_func t_span and theta_0 for solving the ODE
ode45 is a standard matlab ode solver

now the plot is made of the t and the 1st collum and then the 2nd column

where the x-axis is time and the y-axis is the angular velocity

ang_vel is a variable that stores the 1st column

ct is the counter

now the for loop below is used to fix one end of the string to be fixed and the other to move

a plot is made to fix a line, bob, and string

marker is used to make the bob

axis is fixed and the time for each plot is fixed

now the plots are being saved and after every save there would be an increment.

FUNCTION

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

the function saves the values of t, theta,b,g,l,m

the equation of 2nd order ODE is applied. where dtheta_d is a variable that stores the value of dtheta1_dt and dtheta2_dt

Results

https://www.youtube.com/watch?v=ThYyQDLxFN0