Week 3 - Solving second order ODEs

Aim : To create a code for simple pendulum by using ode.

Object : to create code for simple pendulum and create animition video for pendulum motion and displacement and find velocity.

Explanation of simple pendulum and ode

Simple pendulum :

A heavy point mass suspended by a weightless, perfectly flexible and inextensible string, from a perfectly rigid support is called “Simple pendulum”.

The tension in the string depends both on mg and centripetal force

mv^2/ l

The acceleration of the simple pendulum will be given by,

a=−gsinθ

The time period of the simple pendulum is,

T=2*π/ω

The value of the theta (θ) and the rate of theta (θ) at which it is changing for a simple pendulum can be determined or solved by using the set of “ordinary differential equations”.

Ordinary Differential Equations (ODE) :

An ordinary differential equation is an equation containing an unknown function of one real or complex variable x. The unknown function is generally represented by a variable which therefore depends on x.

The examples for the ode are,

dy/dx=F(x)

dy/dx=F(y/x)

d^2y/dx^2=F(y)

In engineering, ODE is used to describe the transient behaviour of a system. A example is a simple pendulum. The following is the second order ode of a simple pendulum. The way the pendulum moves depends on the newtons second law .

Newtons second law of motion states that the acceleration of the object is directly proportional to the net force acting on it, and inversely proportional to its mass. When this law is written in differential equations we will get the ordinary differential equation.

Second order differential equation for simple pendulum :

The position(theta) of the bob with respect to time of a simple pendulum can be determined by using the second order differential equations using the mat lab. The value theta can be solved by using the second order ode and this can be achieved by breaking the higher order ode i.e., second order ode into no of first order ode.

The second order ode for simple pendulum will be as follows,

Where,

b = damping coefficient

m = mass of the bob kg,

g = acceleration due to gravitym/s^2,

l = length of the string m.

θ=θ1

dθ/dt=dθ1/dt=θ2−eq1

d^2θ/dt^2=d^2θ1/dt^2=d/dt(dθ1/dt)=d(θ2)−eq2

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

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

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

dθ1/dt=θ2

dθ^2/dt=−b/m*θ2−g/l*sinθ1

The simple pendulum problem can be solved by using the second order ode functions and by solving these we will get the theta1 an theta2 values of the simple pendulum.

Description

here is given the challenges is some inputs that is value of mass of ball, penedulum length etc is below

lengthof the string l = 1m,

mass of the bob, m = 1kg,

damping coefficient,b = 0.05,

acceleration due to gravity, g = 9.81.

so get initial inputs write code according to it and mention input in code after that

define initial displacement & velocity
given data - simple motion between 0-20sec for angular displacement =0,
& angular velocity = 3 rad/sec at time t =0. so write code according it
theta_0 = [0;3]

after that select time here time is mention in challanges

time input (given)
time span is 0 to 20 sec and get value 500 intervals so write code as below
t_span = linspace(0,20,500)

so here is 0 to 20 sec and value get at 500 interval in between 0 to 20 sec

after that we write differintial equation also create differintial equation function due to clearly mention above explanation to slove ode for simple pendulum and get theta 1 & 2 value so 1 st write here function differintial equation here

differential equation
[t,results]= ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0) ... ode solver one degree equation.

after that write code for plot damping and oscillation value

plotting the damping & the oscillation values.
figure(1)
plot(t,results(:,1))
hold on
plot(t,results(:,2))

so in figure one show a graph for value of damping and oscillation

set the counter to the one
ct = 1

executing for loop from results of theta from linspace
for i = 1:length(results(:,1))
after that we select initial coordinates for pendulum or fixed point of pendulum
x0 = 0
y0 = 0

here is always set at 0 after that write code for find x and y position value differnt angle theta so use foemula as below for verical and horizontal is below
formula for next posion x1 & y1 ( vertical and horizontal )
x1 = l*sin(results(i,1)) for vertical
y1 = -l*cos(results(i,1)) for horizontal

y is put in negative sign because pendulum downward position

after that write code for plotting of pendulum movement & moving and adding bob given below and show in figure 2
figure(2)
plot([-1 1],[0 0],'linewidth',6,'color','r') in get line width and color also
hold on
line ([x0 x1],[y0 y1],'linewidth',4,'color','b')

plot(x1,y1,'v','markers',20,'markerfacecolor','r')
hold off
grid on
axis([-1.5 1.5 -1.5 1.5]) defining axis of graph
pause(0.005) select pause value so animation video run fast so select 0.005 sec

after that write code create animation of pendulum movement given below
M(ct) = getframe(gcf);
ct = ct+1; % incremeting counter
end

animation video save as avi file so write below code

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

details code :

clear all
close all
clc

% program on one degree equation for simple pendulum 

% given inputs 
l = 1; % length of pendulum in meter 
m = 1; % mass of ball in kg
b = 0.05; % damping coefficients
g = 9.81; % gravity in m/s^2

% define initial displacement & velocity 
% given data - simple motion between 0-20sec for angular displacement =0,
% angular velocity = 3 rad/sec at time t =0.
theta_0 = [0;3];

% time input (given)
% time span is 0 to 20 sec and get value 500 intervals.
t_span = linspace(0,20,500);

% differential equation 
[t,results]= ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0); % ode solver one degree equation

% plotting the damping & the oscillation values.
figure(1)
plot(t,results(:,1));
hold on 
plot(t,results(:,2));

% set the counter to the one 
ct = 1;
for i = 1:length(results(:,1))  % executing for loop from results of theta from linspace 
    % initial coordinates for pendulum or fixed point of pendulum
    x0 = 0;
    y0 = 0;
    % formula for next posion x1 & y1 ( vertical and horizontal )
    x1 = l*sin(results(i,1)); % for vertical 
    y1 = -l*cos(results(i,1)); % for horizontal 


% plotting of pendulum movement & moving and adding bob 
figure(2)
plot([-1 1],[0 0],'linewidth',6,'color','r');
hold on 
line ([x0 x1],[y0 y1],'linewidth',4,'color','b');

plot(x1,y1,'v','markers',20,'markerfacecolor','r');
hold off 
grid on
axis([-1.5 1.5 -1.5 1.5]); % defining axis of graph 
pause(0.005);

% create animation of pendulum movement 
M(ct) = getframe(gcf);
ct = ct+1; % incremeting counter 
end

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

For solving the ode we are providing the formulas or equations by using a function so write function code and get value theta 1 and 2

function code in below

% ode_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

Results :

in this challange simulate the pendulum that is create a animation video

link of video :

https://drive.google.com/file/d/1q6Ey_5JzVm2z6dqbkac2dqFh-qHn-tod/view?usp=sharing

and get value of theta 1 &2 and also get graph for velocity and displacement

velocity and displacement graph

Error get during the create the code

1.invalid expression. when calling a function or indexing a variable, use parentheser. otherwise, check for mismatched delimiters. in line 27

this error occurce due to not completed the bracket '()' i put bracket and complete them so error is solve

screen shot of error

2.using 'pendulum' as both the name of variable and the name of scripts is not supported

this error occurce due to i save videofile as name of script name file so this error show so i change video file name and slove error

screen shot of error

3. invalid expression. when calling a function or indexing a variable, use parentheser. otherwise, check for mismatched delimiters. in line 23

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

in this write wrong i write function for theta so after that the ode45 @ write the t and theta and ode45 use to slove differiential equation and also provide backet after @ put all think in i get following code

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

and excute them and solve error

screen shot of error