Week 3 - Solving second order ODEs

Aim : A program to solve ODE (ordinary diffrential equation)

Objective : Writing a program that solves the following ODE. This ODE represents the equation of motion of a simple pendulum with damping

Explanation :

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.

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.

defining about second order ODE

The second-order ODE expresses the fact that the moving mass has both the kinetic and potential energies, Motion of an inertial mass under gravitational field.

explanation about Pendulum

A simple pendulum consists of a ball (point-mass) m hanging from a (massless) string of length L and fixed at a pivot point P. when dispaced to an initial angle and released , the pendulum will swing back and forth with periodic motion. By applying Newton's law for rotational systems, the equation of motion for the pendulum may be obtained.

Inertia force = m*d^theta/dt^2

Damping force = b*dtheta/dt

Gravity force = g/l *sin(theta)

Solving the ODE of a pendulum by using inputs and equations

Inputs :

Damping coefficient (b) = 0.05

Acceleration due to gravity (g) = 9.81 m/s2

Length of the string (l) = 1 metre

Mass of bob (m) = 1 kg

Initial displacement = 0

Initial velocity = 5

Time.(t_span) = 0 to 20 seconds

Simulaton

The motion between 0-20 sec, for angular displacement=0,angular velocity=3 rad/sec at time t=0.

solving the diffrential equation using ODE45 which is inbuilt in MATLAB

[t,r] = ode45(@(t,theta)ode_func(t,theta,b,g,l,m),tspan,theta_0)

where ,

t and r is time output variables

t,theta is definition variables

b,g,l,m are input parameters

now we are going to talk about diffrential equation

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

assume θ=θ1

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

if θ2=dθ1/dt

then

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

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

dθ2/dt = -b/m*θ2-g/l*sinθ1 this is equaion 2

dθ1/dt = θ2 this is eqution 1

equations for displacement and velocity

sinθ = x1-x0/l

x1 = lsinθ

cosθ = y1-y0/l

y1 = lcosθ1

we got all equations now we will code

Code :

clear
close
clc

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

% initial condition angular dispalcement and angular velocity
theta_0 = [0;3];

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

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

% plotting 
figure(1)

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

% Inputs 
x0 = 0;
y0 = 0;
disp_value = results(:,1);
ct = 1;
for i = 1:length(disp_value)
    x1 = l*sin(disp_value(i));
    y1 = -l*cos(disp_value(i));
    figure(2)
    plot([-2,2],[0,0])
    hold on
    plot([x0,x1],[y0,y1])
    plot(x1,y1,'.','Markersize',50)
    hold off
    axis([-2 2 -2 2])
    M(ct) = getframe(gcf);
    ct = ct+1;
   

end

movie(M)
videofile = VideoWriter('second_order_ode.avi','Uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)

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

Result :

link : https://youtu.be/az1L5nx6l-U

Got errors while using :

1) showing error due to inputed undefinied name

got solved after using proper name input (disp_value)

2) showing pendulum improper

got solved using y1 = -l*cos(disp_value(i))

3) Got an error not used M in writeVideo(videofile)

afte got solved using writeVideo(videofile,M)

Conclusion : successfully got ploted for angular displacement and angular velocity

Animated for pendulum with given inputs