Week 3 - Solving second order ODEs

Aim:-
To solve an ordinary differential equation (ODE) and use the calculated information
to animate the motion of the pendulum via MatLab.
Introduction:-
Pendulums have many applications and were utilized often before the digital age.
They are used in clocks and metronomes due to the regularity of their period, in
wrecking balls and playground swings, due to their simple way of building up and
keeping energy.
A simple pendulum consists of a point mass suspended on a string or wire that has
negligible mass. If the pendulum weight or bob is pulled to a relatively small angle
from the vertical and let go, it will swing back and forth at a regular period and
frequency. Figure 1 shows a simple pendulum model

The figure at the right shows an idealized pendulum, with a "massless" string or
rod of length L and a bob of mass m. The open circle shows the rest position of
the bob. When the bob is moved from its rest position and let go, it swings back
and forth. The time it takes the pendulum to swing from its farthest right position
to its farthest left position and back to its next farthest right position is
the period of the pendulum.
The primary forces acting on the bob are the gravitational force that makes it
move in the first place and the force exerted by the string to keep it moving
along a circular path. In addition, there may be a damping force from friction at
the pivot or air resistance or both. We will construct a model to describe how the
angle theta of the pendulum varies as a function of time t.
Let s(t) be the distance along the arc from the lowest point to the position of the
bob at time t, with displacement to the right considered positive. Let theta(t) be
the corresponding angle with respect to the vertical. The figure shows tangential
and radial components of the gravitational force on the pendulum bob. The
radial component is exactly balanced by the force exerted by the string, so the
only relevant force producing the motion is the tangential component of the
gravitational force. For the moment, we ignore the damping force, if any

An ordinary differential equation (ODE) is a differential equation containing one or
more functions of one independent variable and the derivatives of those functions.
The ODE equation and the calculations are shown below in Figure 2.

Code Explanation:-
The function code is shown in figure 3
% A function program which will be called in the main program to solve the
% ODE
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

A function is a group of statements that together perform a task. In MATLAB,
functions are defined in separate files. The name of the file and of the function
should be the same. Functions operate on variables within their own workspace,
which is also called the local workspace, separate from the workspace you access at
the MATLAB command prompt which is called the base workspace.
Function handles are simple to function codes that can be created without using the
regular function definition using the keywords ‘function’ and ‘end’. They do not
require any separate ‘.m’ files or not needed to be declared separately below/above
the main code. Function handles are directly declared in the main code itself and
also accessible from any point in that particular code, i.e they are localized to a
particular conditional. The declaration of a function handle is very similar to that of a
function declaration. It can be distinguished by the means of ‘@’ symbol before your
value/expression allocation. The syntax for creating a function handle is. A function
was defined in order to establish a relation between all the given input values, the
function code is based on the equation derived in figure 2.

Figure 4 shows the main code of the project.
% Objective: Solving second order ODEs of a penedulum
clear all
close all
clc
b = 0.05 % damping coefficient
g = 9.81 % gravity (m/s^2)
l = 1 % length of the spring (m)
m = 0.1 % mass of the bob (kg)
% initial condition
% angular displacement is 0 and angular velocity is 5
theta_0 = [0;5]
% time span
% Pendulum oscillates between 0‐20 seconds
t_span = linspace(0,20,600)
% solving ODE
[t,results] = ode45(@(t,theta) ode_func(t, theta, b,g,l,m), t_span, theta_0);
%Plotting
figure(1)
hold on
plot(t,results(:,1))
plot(t,results(:,2))
title('Damped Harmonic motion of a Pendulum')
legend('Angular displacement', 'Angular velocity')
ylabel('Plot')
xlabel('Time(seconds)')
% For loop
% Before starting the for loop, a counter value was defined
ct =1;
for i = 1:length(results(:,1))
x_start = 0
y_start = 0
x_end = l*sin(results(i,1))
y_end = ‐l*cos(results(i,1))
% plotting the animation
figure(2)
plot([‐1 1], [0 0], 'linewidth', 5)
hold on
plot([x_start x_end],[y_start y_end],'linewidth',3)
hold on
plot(x_end,y_end,'o','markersize',20,'markerFacecolor','y')
title('Animation of transient behaviour of a Pendulum')
axis([‐1.5 1.5 ‐1.5 1])
pause(0.03)
hold off
M(ct) = getframe(gcf)
ct = ct+1
end
% Animation of the pendulum
movie(M)
videofile=VideoWriter('Pendulum Animation.avi', 'uncompressed avi')
open(videofile)
writeVideo(videofile,M)
close(videofile)

Before writing the code, it's important to start with clear all, close all and clc
commands since
Clear all - removes data from the workspace, freeing up system memory
Close all- deletes the current figure or the specified figure
clc - clears all input and output from the Command Window display, giving you a
"clean screen
Next step is to define the inputs and variables. The damping coefficient, gravity,
length of the spring and mass of the bob have all been defined.
Next step is to define the initial conditions, the angular displacement and angular
velocity have been defined. Timespan has also been defined.
Next step is to solve the ODE, Ode45 was used, it is a function in Matlab that solves
ordinary or partial differential equations using a numerical method. The numerical
method in question here is the 4th order Runge Kutta method. The function ode45
takes 3 inputs. The 3 inputs are the function to be integrated, the time span and the
initial conditions. [t, results] means that the results that have been calculated with
respect to time.
After calculating the results, a graph has to be plotted, the graph contained two
plots. One for angular displacement with respect to time and one for angular
velocity with respect for time. The legend for the graph and the axis' labels have
been defined. Figure 5 shows the graph obtained

The second part of the code is to establish a for loop, A for loop is a repetition
control structure that allows you to efficiently write a loop that needs to execute a
specific number of times. a counter loop was added before the For loop. the start
and endpoints for x and y-axis have been defined. The loop is used where all the
values of timespan are assigned to a variable 'i' one at a time. The coordinates of
[ xo x1 y0 y1] and their corresponding relations are entered to form a 'string' of the
pendulum. Markersize command was used to form the bob of the pendulum.
A pause command was used for pausing the command for a certain time, for this
project a pause of 0.03 s was used. The M command was used to capture the
current frame every time the for loop runs. The final step of the loop is to add a +1 to
the ct variable; this will increment the frame every time the for loop run. End
command was used for the plot and for the for loop to stop both from running.
The final step of the code is to create a movie of the animation of the robotic arm, so
the movie command was used. The VideoWriter command was used to assign the
video format, for this project .avi was used as the format. The open command was
used to open the video. writevideo command was used to write the video with plot
frames. Finally, the close command was used to close the video file.
Conclusion:-
The x_start,y_start, x_end y_end values are shown below in figure 6
x_start =
0
y_start =
0
x_end =
‐0.0113
y_end =
‐0.9999