MAT 275 Laboratory 5

The Mass-Spring System

In this laboratory we will examine harmonic oscillation. We will model the motion of a mass-
spring system with differential equations.
Our objectives are as follows:
1. Determine the effect of parameters on the solutions of differential equations.
2. Determine the behavior of the mass-spring system from the graph of the solution.
3. Determine the effect of the parameters on the behavior of the mass-spring.
The primary MATLAB command used is the ode45 function.

Mass-Spring System without Damping

The motion of a mass suspended to a vertical spring can be described as follows. When the
spring is not loaded it has length `0 (situation (a)). When a mass m is attached to its lower
end it has length ` (situation (b)). From the first principle of mechanics we then obtain
mg + −k(` − `0 ) = 0. (1)
|{z} | {z }
downward weight force upward tension force

The term g measures the gravitational acceleration (g ≈ 9.8m/s2 ≈ 32f t/s2 ). The quantity k is
a spring constant measuring its stiffness. We now pull downwards on the mass by an amount y
and let the mass go (situation (c)). We expect the mass to oscillate around the position y = 0.
The second principle of mechanics yields
d2 (` + y) d2 y
mg + −k(` + y − `0 ) = m 2
, i.e., m + ky = 0 (2)
|{z} | {z } | dt dt2
weight
{z }
upward tension force
acceleration of mass

using (1). This is a second-order ODE.

THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED 1
2021 v12 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
Equation (2) is rewritten
d2 y
+ ω02 y = 0 (3)
dt2
where ω02 = k/m. Equation (3) models simple harmonic motion. A numerical solution in the
interval 0 ≤ t ≤ 15 with initial conditions y(0) = 0.8 meter and y 0 (0) = −0.3 , (i.e. the mass is
initially stretched downward 0.8 meters and released with an initial upward velocity of 30 cm/s;
see setting (c) in figure), m = 4, k = 9, is obtained by first reducing the ODE to first-order
ODEs (see Laboratory 4).
Let v = y 0 . Then v 0 = y 00 = −ω02 y = − 94 y. Also, v(0) = y 0 (0) = −0.3. The MATLAB
program LAB05ex1 implements the problem.

LAB05ex1.m
clear all ; % clear all variables
m = 4; % mass [ kg ]
k = 9; % spring constant [ N / m ]
omega0 = sqrt ( k / m );
y0 =0.8; v0 = -0.3; % initial conditions
[t , Y ] = ode45 (@ f ,[0 ,15] ,[ y0 , v0 ] ,[] , omega0 ); % solve for 0 <t <15
y = Y (: ,1); v = Y (: ,2); % retrieve y , v from Y
figure (1); plot (t ,y , 'ro - ' ,t ,v , 'b + - ' ); % time series for y and v
grid on ; axis tight ;
% ---------------------------------------------------
function dYdt = f (t ,Y , omega0 ); % function defining the DE
y = Y (1); v = Y (2);
dYdt =[ v ; - omega0 ^2* y ];
end

Note that the parameter ω0 was passed as an argument to ode45 rather than set to its value
ω0 = 32 directly in the function f. The advantage is that its value can easily be changed in the
driver part of the program rather than in the function, for example when multiple plots with
different values of ω0 need to be compared in a single MATLAB figure window.

NOTE: Use the provided Live Script template for the following exercises. Just like for all the
other lab reports, unlesss otherwise specified, include in your lab report all M-files, figures,
MATLAB input commands, the corresponding output, and the answers to the questions.

1. From the graph in Fig. 1 answer the following questions.

(a) Which curve represents y = y(t)? How do you know?

(b) What is the period of the motion? Answer this question first graphically (by reading
the period from the graph) and then analytically (by finding the period using ω0 ).
(c) We say that the mass comes to rest if, after a certain time, the position of the mass
remains within an arbitrary small distance from the equilibrium position. Will the
mass ever come to rest? Why?
(d) What is the amplitude of the oscillations for y?

THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED 2
2021 v12 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
 Figure 1: Harmonic motion

(e) What is the maximum velocity (in magnitude) attained by the mass, and when is
it attained? Make sure you give all the t-values at which the velocity is maximal and
the corresponding maximum value. The t-values can be determined by magnifying
the MATLAB figure using the magnify button , and by using the periodicity of
the velocity function.
How far is the mass from the equilibrium when the maximum velocity is attained?
(f) How does the size of the mass m and the stiffness k of the spring affect the motion?
Support your answer first with a theoretical analysis on how ω0 (and therefore the
period of the oscillation) is related to m and k, and then graphically by running
LAB05ex1.m first with m = 10 and k = 9, and then with m = 4 and k = 24 . Include
the two corresponding graphs.

2. The energy of the mass-spring system is given by the sum of the kinetic energy and the
potential energy. In the absence of damping, the energy is conserved.

(a) Add commands to LAB05ex1 to compute and plot the quantity E = 21 mv 2 + 12 ky 2

as a function of time. What do you observe? (pay close attention to the y-axis scale
and, if necessary, use ylim to get a better graph). Include at least one plot. Does
the graph confirm the fact that the energy is conserved?
dE
(b) Show analytically that dt = 0. (Note that this proves that the energy is constant).
(c) Add commands to LAB05ex1 to plot v vs y (phase plot). Include the plot. Does the
curve ever get close to the origin? Why or why not? What does that mean for the
mass-spring system?

Mass-Spring System with Damping

When the movement of the mass is damped due to viscous effects (e.g., the mass moves in a
cylinder containing oil, situation (d)), an additional term proportional to the velocity must be

THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED 3
2021 v12 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
added. The resulting equation becomes

d2 y dy d2 y dy
m 2
+ c + ky = 0 or 2
+ 2p + ω02 y = 0 (4)
dt dt dt dt
c
by setting p = 2m . The program LAB05ex1 is updated by modifying the function f. Here we
set c = 4 .

LAB05ex2.m
clear all ; % clear all variables
m = 4; % mass [ kg ]
k = 9; % spring constant [ N / m ]
c = 4; % friction coefficient [ Ns / m ]
omega0 = sqrt ( k / m ); p = c /(2* m );
y0 =0.8; v0 = -0.3; % initial conditions
[t , Y ] = ode45 (@ f ,[0 ,15] ,[ y0 , v0 ] ,[] , omega0 , p ); % solve for 0 <t <15
y = Y (: ,1); v = Y (: ,2); % retrieve y , v from Y
figure (1); plot (t ,y , 'ro - ' ,t ,v , 'b + - ' ); % time series for y and v
grid on ; axis tight ;
% ---------------------------------------------------
function dYdt = f (t ,Y , omega0 , p ); % function defining the DE
y = Y (1); v = Y (2);
dYdt =[ v ; ?? ]; % fill - in dv / dt
end

Figure 2: Damped harmonic motion

3. Fill in LAB05ex2.m to reproduce Fig. 2 and then answer the following questions.

(a) For what minimal time t1 will the mass-spring system satisfy |y(t)| < 0.06 for all
t > t1 ? You can answer the question either by magnifying the MATLAB figure using
the magnify button , or use the following MATLAB commands (explain):

THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED 4
2021 v12 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
 for i =1: length ( y )
M ( i )= max ( abs ( y ( i : end )));
end
i = find (M <0.06); i = i (1)
disp ([ '| y | <0.06 for t > t1 with ' num2str ( t (i -1)) ' < t1 < ' num2str ( t ( i ))])
(b) What is the largest (in magnitude) velocity attained by the mass, and when is it
attained? Answer by using the magnify button.
(c) How does the size of c affect the motion? To support your answer, run the file
LAB05ex2.m for c = 6, c = 12, and c = 30. Include the corresponding graphs with a
title indicating the value of c used.
(d) Determine analytically the smallest (critical) value of c such that no oscillation ap-
pears in the solution.

4. (a) Modify LAB05ex2.m to compute and plot the quantity E = 12 mv 2 + 12 ky 2 as a function

of time. What do you observe? Is the energy conserved in this case? Include the
plot.
dE dE
(b) Show analytically that dt < 0 for c > 0 while dt > 0 for c < 0.
(c) Modify LAB05ex2.m to plot v vs y (phase plot). Include the plot. Comment on the
behavior of the curve in the context of the motion of the spring. Does the graph
ever get close to the origin? Why or why not?

THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED 5
2021 v12 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University