Week 3 - Solving second order ODEs
Green Color = Comments,
Blue Color = Program,
Black Color = Explanations,
Red Color = Instructions
# Youtube Link for the Program = https://youtu.be/RAy1n59IC_w
# Program
% A Program on Solving second order ODE with Pendulum Animation.
close all
To close all the plots in the figure window
clear all
To clear/delete all the saved data in the workspace, editor & command window
clc
To clear/delete all the temporary saved data in the current window
% Input
Comment to specify below are all the input data (taken from Goggle)
Below are the values of all the necessary data required to perform the program
Values such as
Damping Coefficient (b)
Gravity (g) in m/s^2
Length of Pendulum (l) in m
Mass of the Ball (m) in Kg
b = 0.05;
g = 9.81;
l = 1;
m = 0.2;
% Initial Condition of Angular displacement & velocity in m/s
Here 0.5 is Angular Displacement & 4 is velocity m/s
theta_0 = [0.5;4];
% Time Points
Here linspace command is used to specify that the motion of t_span will be from 0 to 40 with 600 values
t_span = linspace(0,40,600);
% Solving ODE
Solving the ODE with ode45 command and using the function program named "ode_func"
You can find the function program at the end of this Program
[t, results] = ode45(@(t,theta) ode_func(t, theta, b,g,l,m),t_span, theta_0);
After every value semi-colon (;) is used to not display the value in the command window
% Plotting Angular displacement, Velocity and time
This plot is named as figure (1)
figure(1)
plot(t,results(:,1))
hold on
grid on
plot(t,results(:,2))
xlabel('Time')
ylabel('Plot')
legend('Displacement','Velocity')
title('Displacement & Velocity plot')
% Animation Movie
Here S is used as a counter and given a value as 1
Here using a counter for plotting with a specific command
S = 1
Here in this program, multiple values are being taken also in a specific sequence without entering them manually, For that for command is being used for specifying the values i.e. i = from 1 to length of (results(:,1)).
for i = 1:length (results(:,1))
% Coordinates
Below are the values of the coordinates i.e. (x0 y0), (x1 y1)
x0 = 0
y0 = 0
x1 = l*sin(results(i,1));
y1 = -1*cos(results(i,1));
% Plotting
Comment for Plotting
Here plot command is used including the Coordinates, Linewidth, Markersize of the ball, and specifying the axis.
The pause & hold command is used to add an interval + hold the animation and make it smoother.
This plot is named as figure(2)
figure(2)
plot([-1.5 1.5],[0 0])
hold on
plot([x0 x1],[y0 y1],'linewidth',2,'Color','b')
hold on
axis([-1.5 1.5 -1.5 1.5])
plot (x1,y1,'O','Markersize',30,'Color','r')
Here naming the movie/animation/video as A of counter S as get frame for framing it in one.
A(S) = getframe(gcf);
hold off
This is to specify to add 1 after every plotting.
S = S+1;
end
% Creating animated movie from the above plot
Specifying movie named as A
movie(A)
Naming the video/animation/movie using VideoWriter command as second order ode Pendulum.m and specifying the format of the movie as Uncompressed AVI
videofile = VideoWriter('second order ode Pendulum.m','Uncompressed AVI');
specifying to open the program after saving it in a specific name and format
open(videofile)
this command is used to play the video/animation/movie named A in the program
writeVideo(videofile,A)
After successfully playing the video/animation/movie close command is used to close it.
close(videofile)
# Program for Defining function named ode_func
% Function program for 2nd order ODE
function [dtheta_dt] = ode_func(t,theta,b,g,l,m)
theta1 = theta(1)
theta2 = theta(2)
second order ode is breaked into 2 single order ode
dtheta1_dt = theta2;
dtheta2_dt = -(b/m)*theta2 - (g/l)*sin(theta1);
dtheta_dt = [dtheta1_dt; dtheta2_dt];
end
# Output
- You can find the animation/video/movie in the Youtube link mentioned at the top.
- The video shows the movement with respect to Angular displacement, Velocity, and Time of a Pendulum shown with the help of second order ODE.
- The values are taken randomly such as Damping Coefficient (b), Gravity (g) in m/s^2, Length of Pendulum (l) in m, Mass of the Ball (m) in Kg.
- As the Time increases (left to right), Angular Displacement & Velocity decreases (shown in figure(1)).
- This program shows that one can make animations and do calculations of ODE's in Matlab software with ease.
# Errors that were faced during the program
Line no. 53 (in main Program)
A(S) = getframe(gcf);
hold off
S = S+1;
The mistake is highlighted with yellow above
As the hold on command is used to hold every displacement in its place in the plot figure (2), accordingly hold must be set off too.
The error was that I didn't use the hold off command so that all the displacement was shown (Snap is attached)
The snap of the error is attached below named as Error 01