Solving Second Order ODEs using Matlab

Aim: to solve a second order ODE using Matlab and

to simulate the motion of pendulum using Matlab.

CODE 1: for function 'ODE_func'

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 2: for plotting time vs distance (amplitude) and for simulating pendulum motion

%solving second order ODEs 
%simulation of pendulum motion using Matlab

clear all
close all
clc

%input parameters
b=0.5;
g= 9.81;
l= 1;
m= 1;

%given initial conditons
theta_0= [0;3];

%time points
t_span= linspace(0,20,300);

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

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

%simulation of pendulum motion
%input coordinates

for i=1:length(results(:,1))
    res= results(i);
    
x0=0;
y0=0;

x1= l*sind(res);
y1= l*(-cosd(res));

figure(2)
plot([x0 x1], [y0 y1], 'linewidth',3,'color','g');
axis([-0.03 0.03 -1 1]);
M(ct)=getframe(gcf);
ct=ct+1;

end

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

Procedure:

1.the second order ODE used for this project (simulation of pendulum motion) is

`(d^2theta)/dt^2 + (b/m)((d theta)/dt) + (g/l)sintheta = 0`

where, b= damping value

m= mass (kg)

g= gravity (m/s^2)

l= length (m)

2. for solving a second order ODE, we split it into 2 first order ODEs and then solve these 2 first order ODEs to get the overall solution.

3. this is the procedure that we do in Matlab by using 'function' operator

4. function operator is used to split the second order ODE into 2 first order ODEs. the function operator is named as 'ODE_func'

5. then, we write the main code by giving input parameters for b, g, l and m.

6. here, b=0.5, g=9.81 m/s^2, l=1m, m=1 kg

7. then, we give the initial variable values. here, the initial time is 0 secs and velocity is 3 m/s.

8. the time range is given as 20 seconds with linspace command of 300 intermediate values.

9. then, we use ode45, which is an in-built ODE solver in Matlab, to give the ode_func in the main code.

10. counter variable i8s initiated as 1.

11. the plot between distance and time is written using plot and hold on command.

12. for the simulation of pendulum motion, we use a 'for' loop to give values of the coordinated of the position of end points of the pendulum.

13. the origin point is named as [x0, y0] and is given a value of (0,0).

14. a new variable 'res;' is given to find the other end point which oscillates with respect to time.

15. then, movie command is used to get the video file of the simulation. it is stored as 'uncompressed AVI' format.

16. the program is run and we successfully obtain the simulation of pendulum motion under given conditions.

ERRORS done and STEPS taken to correct them:

ERROR 1:

a typo was done here at y1. it was subsequently corrected.

ERROR 2:

the counter 'ct' variable was not initiated. this was corrected by giving ct=1 in the code.

RESULTS:

the following graph shows the relation between amplitude (distance) and time of the 2nd order ODE representing pendulum motion.

the following Video shows the simulation of the pendulum under given conditions.

https://drive.google.com/open?id=1dtfqeL23QsU-9216PjRH091v0Vqknrrx