Open Vibrations

Open Vibrations
Austin R.J. Downey

Department of Mechanical Engineering
Department of Civil and Environmental Engineering
University of South Carolina, Columbia SC, USA
Laura Micheli
Department of Civil and Environmental Engineering
University of South Carolina, Columbia SC, USA
December 24, 2023

Open Vibrations CONTENTS
Contents
Preface iv
Cover Art . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iv
License . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
Source Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
I Foundational Concepts 1
1 Basic Concepts in Vibrations 3
1.1 Single Degree-of-Freedom Systems . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.1.1 Spring-Mass Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.1.2 Linear Springs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.2 Equivalent Stiffness . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.2.1 Equivalent Stiffness of Structural Systems . . . . . . . . . . . . . . . . . . 8
1.2.2 Springs in Series and Parallel . . . . . . . . . . . . . . . . . . . . . . . . 11
1.3 Equation of Motion for an Oscillating System . . . . . . . . . . . . . . . . . . . . 13
2 Free Vibrations 18
2.1 Mathematical Modeling of Free Vibration . . . . . . . . . . . . . . . . . . . . . . 18
2.1.1 Solve for the Natural Frequency (ωn ) of the System . . . . . . . . . . . . 19
2.1.2 Solve for Initial Phase (φ ) of the System . . . . . . . . . . . . . . . . . . 20
2.1.3 Solve for Amplitude (A) of the System . . . . . . . . . . . . . . . . . . . 21
2.1.4 Response for Simple Harmonic Motion . . . . . . . . . . . . . . . . . . . 22
2.1.5 Special Considerations for No Initial Velocity (v0 = 0) . . . . . . . . . . . 22
2.2 General Solution for Vibrating Systems . . . . . . . . . . . . . . . . . . . . . . . 25
2.2.1 Formulating the General Solution for a 1-DOF Spring-Mass System . . . . 28
2.2.2 Solution of 1-DOF System in Three Forms . . . . . . . . . . . . . . . . . 31
2.3 Damping . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
2.3.1 Modeling Vibrating Systems with Damping . . . . . . . . . . . . . . . . . 34
2.3.2 Modeling Underdamped Motion . . . . . . . . . . . . . . . . . . . . . . . 35
2.3.3 Modeling Overdamped Motion . . . . . . . . . . . . . . . . . . . . . . . . 42
2.3.4 Modeling Critically Damped Motion . . . . . . . . . . . . . . . . . . . . 43
2.3.5 Standard Form of the EOM . . . . . . . . . . . . . . . . . . . . . . . . . . 44
2.4 Logarithmic Decrement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
3 Forced Vibrations 52
3.1 Harmonic Excitations of Undamped Systems . . . . . . . . . . . . . . . . . . . . 53
3.2 Harmonic Resonance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
3.3 Harmonic Excitations of Underdamped Systems . . . . . . . . . . . . . . . . . . . 61
3.4 Frequency Response of Underdamped Systems . . . . . . . . . . . . . . . . . . . 66
3.5 Base Excitation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
3.5.1 Displacement Transmissibility Solution for Base Excitation . . . . . . . . 74
3.5.2 Force Transmissibility Solution for Base Excitation . . . . . . . . . . . . . 78
i
3.6 Numerical Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83
4 Transfer Function Approach 86

4.1 Transfer Function Method (Generic) . . . . . . . . . . . . . . . . . . . . . . . . . 86
4.2 Transfer Function Method for Solving Vibrating Systems . . . . . . . . . . . . . . 89
4.2.1 Free Vibration for Undamped Systems . . . . . . . . . . . . . . . . . . . . 89
4.2.2 Impulse Response Function . . . . . . . . . . . . . . . . . . . . . . . . . 90
4.2.3 Unit Step function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
4.2.4 Undamped Spring-mass System . . . . . . . . . . . . . . . . . . . . . . . 93
4.2.5 Under Damped Spring-mass System . . . . . . . . . . . . . . . . . . . . . 93
4.3 System Response to Arbitrary Inputs . . . . . . . . . . . . . . . . . . . . . . . . . 96
4.4 Transfer Function for Response to Random Inputs . . . . . . . . . . . . . . . . . . 100
4.4.1 Defining the Transfer Function H(s) . . . . . . . . . . . . . . . . . . . . . 100
4.4.2 Transfer Function Method (Steady-state Solution) . . . . . . . . . . . . . . 101
4.4.3 Response to Random Inputs . . . . . . . . . . . . . . . . . . . . . . . . . 105
5 Multiple Degree-of-freedom Systems 113

5.1 General Discussion on Mode Shapes . . . . . . . . . . . . . . . . . . . . . . . . . 114
5.2 Modeling Undamped Two Degree of Freedom Systems . . . . . . . . . . . . . . . 116
5.2.1 Solution for the Two-Degree-of-Freedom System . . . . . . . . . . . . . . 117
5.2.2 Deriving the Solution for the Two-Degree-of-Freedom System . . . . . . . 117
5.3 Explicit method for Solving 2-DOF Systems . . . . . . . . . . . . . . . . . . . . . 127
5.4 Eigenvalue-based Solution for Natural Frequencies and Mode Shapes . . . . . . . 130
5.4.1 Deriving the Eigenvalue-based Solution . . . . . . . . . . . . . . . . . . . 131
5.5 Transfer-function Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137
5.6 Multiple Degrees of Freedom . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139
5.7 Modal Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142
5.8 Numerical Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 149
II Applied Topics 153

6 Vibration Control 154
6.1 Vibration Nomograph . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 154
6.2 Vibration Isolation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157
6.3 Vibration Absorption . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 158
6.3.1 Vibration Absorption for Undamped Systems . . . . . . . . . . . . . . . . 158
6.3.2 Vibration Absorption for Damped Systems . . . . . . . . . . . . . . . . . 161
6.4 Active Vibration Suppression . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 162
6.4.1 Metrics for Vibration Control . . . . . . . . . . . . . . . . . . . . . . . . 162
6.4.2 Position-Derivative (PD) Control . . . . . . . . . . . . . . . . . . . . . . . 163
6.4.3 Proportional-Integral-Derivative (PID) Control . . . . . . . . . . . . . . . 165
ii
7 Experimental Vibrations 167

7.1 Sensing and Data Acquisition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 168
7.1.1 Accelerometers . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 168
7.1.2 Data Acquisition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 170
7.2 Controlled Force Excitation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 171
7.2.1 Modal Hammers . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 171
7.2.2 Electrodynamic Shakers . . . . . . . . . . . . . . . . . . . . . . . . . . . 172
7.3 Digital Signal Processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172
7.3.1 Sampling and Quantization . . . . . . . . . . . . . . . . . . . . . . . . . . 173
7.3.2 Aliasing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 174
7.3.3 Time-Frequency Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 176
8 Structural Dynamics 177

8.1 Single-story frame . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 177
8.2 Duhamel’s Integral . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 180
8.3 Two-story frame . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 182
iii
Open Vibrations
Preface
This open-source text is designed to offer a complete introduction to the field of vibrations, specifi-
cally tailored for undergraduate students. It covers the fundamental principles of vibrations, includ-
ing single and multi-degree of freedom systems, transfer function approaches, vibration control,
along with measurement and instrumentation. Each chapter includes examples and case studies to
reinforce the concepts presented. With its simple and clear explanations and practical approach,
this textbook serves as an essential resource for undergraduate students studying vibrations in en-
gineering disciplines.
Cover Art
The B-57 Canberra is an American-built copy of the British English Electric Canberra which first
flew in 1949. During initial high-speed flight testing, excessive vibrations were measured on the
canopy and a small fairing was added behind the canopy to reduce the aerodynamic load on the
canopy and thereby reduce the vibrations. Overall, airframes flight testing is said to have gone
very smoothly.
The B-57 was initially a twin-engined tactical bomber and reconnaissance aircraft but over the
years, various versions were produced or modified from the original stock. These include the WB-
57F, a specialized strategic reconnaissance version developed for the U.S. Air Force that is still
flown by NASA for scientific missions. Of note, in 2011 NASA determined that they needed a
third WB-57F to support their mission and an additional WB-57 (s/n 63-13298) was removed from
the Air Forces Boneyard in Tucson Arizona after 40 years of storage and returned to operational
status. As of 2022, three airframes are still flying for NASA.
The airframe on the cover is the SN 52-1516 and is an EB-57B “Night Intruder” which is an
electronic countermeasure (ECM) version of the B-57 and is on static display at the Air Force
Armament Museum at Eglin Air Force Base, Florida.
iv
Open Vibrations License
License
This work is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License[cc-
by-sa 4.0]. More information on the Attribution-ShareAlike 4.0 International (CC BY-SA 4.0)
license can be found at https://creativecommons.org/licenses/by-sa/4.0/. Unless oth-
erwise denoted, all text, figures, diagrams, and photos used in this work are the sole property of
the authors and are released under CC BY-SA 4.0 both in part and in whole. Reworks and redistri-
butions of this work that fall within the CC BY-SA 4.0 licenses are encouraged.
Source Code
The source code for this text is available at https://github.com/austindowney/Open-Vibrations.
1
Open Vibrations
Part I
Foundational Concepts
The Arthur Ravenel Jr. Bridge spans the Cooper River outside of Charleston South Carolina
(USA) with a main span of 471 m (1,546 feet) and uses dampers on the center cables to mitigate
wind and traffic-induced vibrations in the structure. It is the third longest among cable-stayed
bridges in the Western Hemisphere (2023).
2
Open Vibrations
1 Basic Concepts in Vibrations

The study of vibrations, within the broader field of classical mechanics, is the investigation of
oscillations that occur about an equilibrium point. Vibrations, both desired and undesired, are
present in all mechanical systems and can be helpful (e.g. a soil sieve, rotary sander) or destructive
(e.g. an aircraft frame in resonance). The oscillations that form a vibrating system may be periodic
(e.g., pendulum) or random (e.g. turbulence in an airplane), or a combination of the two.
Vibrations impact our daily lives in a variety of ways, from the sound made by banjo strings
that vibrates between 140 and 400 Hz to the 4-6 Hz vibration felt by a passenger in a car seat. The
consideration of vibrations, and their associated mathematical modeling, are important factors in
the design of mechanical systems. In the material that follows, the fundamental theories of vibra-
tion are presented and modeled using fundamental physical principles such as Newton’s three laws
of motion. These models and analyzed using the mathematical tools of calculus and differential
equations.
Vibration Case Study 1.1 Why study vibrations? One day, it could save your life! The
British Aircraft Corporation TSR-2 was a strike and reconnaissance aircraft developed dur-
ing the Cold War by the British Aircraft Corporation (BAC), for the Royal Air Force (RAF).
During the second flight test of air-frame XR219, vibration from one of the plane’s fuel
pump causes a vibration at the resonant frequency of the human eyeball which caused a mo-
mentary loss of vision. Test pilot Roland Beamont throttled back one engine to restore his
full vision. Roland Beamont was an expert in vibrations, having led the vibration program of
the Hawker Typhoon during WW II. In this role, he fit vibrographs to airplanes to determine
the effectiveness of propeller balancing and the testing of seats with vibration isolators.
Figure 1.1: The only BAC TSR-2 prototype to fly, picture taken in 1966 at what is now BAE
Warton Lancashire.a
a RuthAS, CC BY 3.0 <https://creativecommons.org/licenses/by/3.0>, via Wikimedia Commons
3
Open Vibrations 1.1 Single Degree-of-Freedom Systems
Review 1.1 Newton’s three laws of motion:
1. In an inertial frame of reference, an object either remains at rest or continues to move

at a constant velocity unless acted upon by a force.
2. In an inertial reference frame, the vector sum of the forces F on an object is equal to
the mass m of that object multiplied by the acceleration of the object: F = ma. (It is
assumed here that the mass m is constant)
3. When one body exerts a force on a second body, the second body simultaneously
exerts a force equal in magnitude and opposite in direction on the first body.
1.1 Single Degree-of-Freedom Systems

In its simplest form, the phenomenon of vibration is the exchange of energy between potential and
kinetic energy. Therefore, a vibrating system must have a component that stores potential energy.
This component must also be capable of releasing the energy as kinetic energy. This kinetic energy
is stored in the movement of a mass where the measure of this movement is the velocity of the
system and the continuous interchange between potential and kinetic energy is the vibration of the
system. The simplest vibrating systems can be modeled as a single-degree-of-freedom (1-DOF)
system. In a 1-DOF system, one variable can describe the motion of a system. Potential examples
of 1-DOF systems include:
1. yo-yo
2. pogo stick
3. door swinging on axis
4. throttle (gas pedal)
Variables often used for describing 1-DOF systems are x(t), y(t), z(t), and θ (t). Examples of
1-DOF systems are presented in figure 1.2 where the assumption of small displacements is made.
Note: we will often drop the “(t)” for simplicity in this material.
4
Figure 1.2: Examples of single degree of freedom (DOF) systems showing: (a) a vertical spring-
mass system; (b) a simple pendulum; and (c) a rotational spring-mass system.
1.1.1 Spring-Mass Model
“All models are wrong, but some are useful ”

George E.P. Box (1919 - 2013)
Newtonian physics describes the motion of particles in terms of displacement x, velocity ẋ,
and acceleration ẍ vectors. Moreover, Newton’s second law of motion says that the change in the
velocity of mass in motion is a product of the force acting on the mass. A simple way to express
this phenomenon is through a spring-mass model as presented in figure 1.3. These spring-mass
models neglect the mass of the spring and concentrate all the mass of the system into a single
point. Note that in this case the force vector and mass-acceleration vectors lie on the same axis and
as such are collinear. Therefore, these vectors can be easily treated as scalers simplifying the math
used in the modeling of the system.
Figure 1.3: A single-degree-of-freedom (1-DOF) spring-mass model showing: (a) annotated

schematic of a mass-spring system; and (b) the equivalent free-body diagram represented as a
point-mass system.
5
1.1.2 Linear Springs

Springs are mechanical devices that store energy, moreover, an ideal spring is a theoretical repre-
sentation of this mechanical device that is massless and responds with a linear increase in force for
a unit increase in displacement (i.e. F = kx). For simplicity, the springs in the spring-mass models
considered in this text are always assumed to be ideal linear springs. A graphical representation of
the idealized linear spring is presented in figure 1.4 where a unit force F applied to the free end of
the spring results in a unit displacement x of the spring. The resulting mathematical relationship,
F = kx, is known as Hooke’s Law. Nonlinear springs add considerable complexity to the modeling
of spring-mass systems, therefore, these are not considered in this introductory work.
Figure 1.4: Force-displacement plot for a linear spring.
Vibration Case Study 1.2 Why study vibrations? Because vibrations form an integral part
of how we interact with our world and as such, are an important consideration in product.
For example, vibrations in the automotive industry fall within a field of expertise termed
Noise Vibration and Harshness (NVH). NVH is important because, within a single company,
different levels of NVH will be desired for different market segments and products.
With a proper understanding of NVH, engineers can design cars that can adapt to their
environment or desired use case. Consider the 2019 VW Golf GTI shown in figure 1.5(a)
equipped with a dynamic suspension system where the driver can select between ‘comfort’
‘normal’ and ‘sport’ suspension options. To investigate the effect of these suspension set-
tings, an engineer can install an accelerometer (a sensor used for measuring acceleration) as
shown in figure 1.5(b). An important consideration in measuring acceleration is where and
how to mount the accelerometer. Here, the accelerometer is mounted in the cup holder to
measure the vertical acceleration in the center of the car.
6
Figure 1.5: VW Golf GTI with three suspension modes, showing: (a) the car, and; (b) the
accelerometer and data acquisition system used for measuring vibrations.
Figure 1.6 shows the measured acceleration in the frequency domain for the three sus-
pension modes during 5 minutes of interstate driving. While we will delve into the technical
of power spectral density, for now, consider the area under the curve to be representative of
the amount of energy in the measured vibrations for each suspension setting.
Figure 1.6: Response spectrum of the vibrational energy in a VW GTI measured using the
experimental setup shown in figure 1.5(b).
The sport mode is by far the suspension mode with the firmest ride and the highest
amount of measured vibration energy. While a stiff ride is beneficial during spirited driv-
ing on a track, the associated NVH level is tiring during prolonged driving. However, the
comfort mode adds a considerable amount of damping to the suspension, resulting in a ride
quality that is much more amenable to everyday driving. An engineer, using their knowledge
of vibrations, could develop systems that enable a single product (such as an automobile or
an airplane) to function well in multiple use cases; thereby increasing its usefulness and
7
Open Vibrations 1.2 Equivalent Stiffness
marketability.
1.2 Equivalent Stiffness

The generalized concept of stiffness can be directly related to mechanical systems and structural
components through Hooke’s law.
Review 1.2 Hooke’s Law states that the force (F) needed to extend or compress a spring
by some distance x scales linearly with respect to that distance. This law can be extended to
tensional stress of a uniform and elastic bar where the length, area, and Young’s modulus of
the bar are represented by l, A, and E, respectively. Knowing the tensile stress in the bar:
F
σ= (1.1)
A
and the definition of strain:
∆l
ε= (1.2)
l
Hooke’s law can be expanded to represent a uniform and elastic bar:
σ = Eε (1.3)
It follows that the change in length ∆l can be expressed as:

Fl
∆l = εl = (1.4)
AE
Note: Hooke’s law is often expressed using the convention that F is the restoring force
exerted by the spring on the applied force at the free end. Defining the stiffness and dis-
placement as k = AEl and ∆l = x, respectively. The equation for Hooke’s Law becomes:
F = −kx (1.5)
since the direction of the restoring force is opposite the spring displacement.
1.2.1 Equivalent Stiffness of Structural Systems

For a rod with a uniform cross-section, a direct representation of the system can be developed
as expressed in figure 1.7 where the vibration along the axis of the rod is to be considered. The
stiffness of the rod, k, is a measure of the resistance offered by an elastic body to deformation.
8
Figure 1.7: Equivalency between a vertical bar with a mass attached to the bottom and a spring-
mass model of the system.
For this 1-DOF system, the equation of a spring can be rearranged such that the stiffness can
be defined as:
F
k= (1.6)
x
The stiffness of the spring can be more closely related to material properties of the bar A, E, and l
considering that Hooke’s Law for the uniform tension on a bar can be expressed as:
σ = Eε (1.7)
This expression can be expanded into the form:

F x
=E (1.8)
A l
F
rearranging the terms and recalling the expression k = x leads to:
EA
k= (1.9)
l
In a similar fashion, we can also solve the equivalent system for a mass at the end of a cantilever
beam.
Figure 1.8: Equivalency between a cantilever beam and a spring mass system.
9
From engineering mechanics, we can compute the deflection at the point of a beam δ with a
point load P. This expression is typically expressed as:
Pl 3
δ= (1.10)
3EI
If we transform this equation into our variable system by exchanging P for F and δ for x. There-
after, the point load is replaced with the equivalent force F generated by the mass and the pull of
gravity(mg). As before, knowing that the stiffness of the system can be expressed as k = F/x we
can show that:
3EI
k= 3 (1.11)
l
Example 1.1 Considering the rod diagrammed below; calculate an equivalent spring constant
for the rod using the length of the rod l, its area A, and Young’s modulus E for a compressive
force F that compresses the rod a distance x. Additionally, is a linear spring a useful model
for a rod under compression? What if the rod is under tension?
Figure 1.9: Compressed cantilever rod.
Solution:
The rod shortens by a distance x under the axial force F, this can be related to the equation
of a linear spring F = kx by recalling from solid mechanics that the elongation (or shortening)
of a rod is expressed as
x σ Fl
x = l = εl = l = (1.12)
l E AE
where ε = xl is the strain value and σ = F/A is the stress induced in the rod. Combining
this expression with the equation of a linear spring yields:
F AE
k= = (1.13)
x l
As per the usefulness of the linear spring to represent an axial rod under compression
or tension, this would be application-specific but could generally be considered an excellent
first-order approximation.
10
1.2.2 Springs in Series and Parallel

In many cases, it becomes necessary to model a mechanical system as a set of springs (e.g., a
composite material, a table with multiple legs). For these systems, or for systems with more than
one spring acting on a body, equivalent stiffness can be calculated as:
Figure 1.10: Equations for calculating the equivalent stiffness of two springs (k1 and k2 ); (a) in
series; and (b) in parallel.
These are derived considering the displacement δ of the systems. For two springs in series:
Figure 1.11: Two springs k1 and k2 combined in series.
where the total displacement is

δac = δab + δbc (1.14)
Using the equation for stiffness k = F/δ , this converts to:
F F F
= + (1.15)
kac k1 k2
As F is the same throughout the system, we can cancel out F. Solving for the equivalent stiffness
yields:
1
kac = 1 1
(1.16)
k1 + k2
Similarly for a system of springs in parallel:
Figure 1.12: Two springs k1 and k2 combined in parallel.
11
The displacement in both springs is the same, so the total displacement is
δab = δ1 = δ2 = δ (1.17)
The forces in the direction of spring elongation sum to zero, therefore:
Fab = F1 + F2 (1.18)
Substituting the displacement and stiffness into the force equation yields:
δ kab = δ k1 + δ k2 (1.19)
this simplifies to:

kab = k1 + k2 (1.20)
Example 1.2 Calculate the equivalent stiffness of the following system:
Figure 1.13: Equivalent stiffness for springs in series and parallel.
The springs are combined as shown, using the equations defined before. Now, considering
that the displacement (δ ) of the top spring, and the bottom spring are the same we can state
the total stiffness k, which is the summation of the two. Therefore,
12
Open Vibrations 1.3 Equation of Motion for an Oscillating System
Figure 1.14: A spring mass system simplified down form springs in series and parallel..
1
where the final addition, (k1 + k2 ) + (k5 + 1 1 ) is applied at two springs in parallel as each
k3 + k4
spring is connected between the mass and the fixity. Rearranging this new expression to get a
common denominator:
(k1 + k2 + k5 )(k3 + k4 ) + k3 k4
k= (1.21)
k3 + k4
1.3 Equation of Motion for an Oscillating System

An Equation of Motion (EOM) is an equation that provides a basis for modeling a vibrating system
about its equilibrium point and relates the transfer of the potential energy from the spring to the
kinetic energy mass. In developing the EOM we assume that any surfaces are frictionless and
as such, no energy is extracted from the vibrating system. Referencing the 1-DOF system in
figure 1.15(a), and assuming the mass only moves in the x direction, the only force acting on the
mass in the x direction is the force that results from the elongation of the spring as annotated in
figure 1.15(b). Therefore, the sum of forces along the x axis must equal the mass (m) times the
acceleration of the mass (aẋ).
Figure 1.15: A spring-mass model of a 1-DOF system showing: (a) a schematic of the system; (b)
a free-body diagram of the system at its initial position.
Considering that positive displacements are to the right, the standard form of the equation of
motion for an undamped system without any excitation is expressed as:
s1 ẍ + s2 x = 0 (1.22)
13
where s1 and s2 are constants to be determined for the specific system. A systematic approach
to obtaining the free-body diagram (FBD) of a system under vibration can be expressed in three
steps:
1. Draw a free-body diagram (FBD) at the system’s equilibrium and displaced position (without
a displacing force).
2. Apply Newton’s second law to both FBDs ( equilibrium and displaced).
3. Combine the equations to write the EOM in standard form with the forcing component on
the right-hand side. For free vibration, the forcing component is 0.
Solving these three steps for 1-DOF system presented in figure 1.15 results in the EOM:
mẍ + kx = 0 (1.23)
Review 1.3 A second-order linear homogeneous differential equation has the form:
aẍ + bẋ + cx = 0 (1.24)

The EOM for a 1-DOF system under a free vibration is a second-order differential equation
due to acceleration (ẍ) being the second derivative of displacement (x) and homogeneous
as the forcing function (right-hand side of the equations) is zero. In EOM’s current form,
a = m, b = 0, and c = k. In future work, b will account for damping in the vibrating system.
Example 1.3 Considering the system:
Figure 1.16: A 1 DOF spring mass system with movement in the horizontal direction
Step-1 Define the direction of displacement, and draw the FBD for the equilibrium and
displaced state.
14
equilibrium state displaced state
Figure 1.17: Equivalent forces for a 1 DOF spring mass system with movement in the hori-
zontal direction
The equation for the equilibrium state is:
+
→ ∑ Fx = 0 (1.25)
and in the displaced state:

+
→ ∑ Fx = −kx (1.26)
This equation does not equal zero as the FBD does not account for the restoring force.
Step-2 Apply Newton’s second law (we want to store energy in the kinetic state) of motion to
the sum of forces for the displaced position we get:
+ ∑ Fx = −kx
ma = mẍ = → (1.27)
mẍ = −kx (1.28)

Step-3 Rearrange in the Equation to construct an EOM:
mẍ + kx = 0 (1.29)
Example 1.4 Some systems will have an initial displacement, as the system will oscillate
around this position we need to define the EOM about this position. Considering the system:
Figure 1.18: A 1 DOF spring mass system with movement in the vertical direction
15
Step-1 Define the direction of displacement, and draw the FBD for the equilibrium and dis-
placed state.
equilibrium state displaced state
Figure 1.19: Equivalent forces for a 1 DOF spring mass system with movement in the vertical
direction
+↓ ∑ Fx = mg − kδ = 0 (1.30)

+↓ ∑ Fx = mg − k(δ + x) (1.31)
This equation does not equal zero as the FBD does not account for the restoring force.
Step-2 Apply Newton’s second law (we want to store energy in the kinetic state) of motion to
the sum of forces for the displaced position we get:
mẍ = +↓ ∑ Fx = mg − kδ − kx (1.32)
We can than use the information from the equilibrium state to cancel out some terms, this
becomes:
mẍ = −kx (1.33)
Step-3 Rearrange in the Equation to construct an EOM:
mẍ + kx = 0 (1.34)
Vibration Case Study 1.3 Why study vibrations? One day, it could save your job! For a
project to be successful they need to be completed on time and within budget. Consider
the Ling-Temco-Vought (LTV) XC-142 which was a tilt-wing experimental aircraft devel-
oped in the ’60s for the US military and later turned over to NASA. During testing, the
cross-linked driveshaft produced excessive vibration and noise which resulted in a high pi-
lot workload. In general, the aircraft’s cross-linked driveshaft was the main technical issue
that caused the military to lose interest in the project.
16
Figure 1.20: A Ling-Temco-Vought XC-142A tested at the NASA Langley Research Center
in 1969. a
a NASA, Photograph published in Winds of Change, 75th Anniversary NASA publication, by James
Schultz, Public domain, via Wikimedia Commons
17
Open Vibrations
2 Free Vibrations
Vibrations (i.e. the exchange of potential and kinetic energy) require oscillatory motion that may
repeat itself regularly or irregularly. A motion that is repeated at time intervals is called periodic
motion. If this motion has a single frequency and amplitude it is called simple harmonic motion
are represents the most basic form of oscillatory motion as depicted in figure 2.1. For a 1-DOF sys-
tem, simple harmonic motion is defined as a periodic motion where the restoring force is directly
proportional to the displacement and acts in the direction opposite to that of displacement.
Figure 2.1: Oscillatory motion for a single degree of freedom system showing (a) periodic motion;
and (b) simple harmonic motion.
Given the nature of simple harmonic motion, constant amplitude, and frequency, the wave
starting at the origin O can be modeled at a point on the end of a vector with length A rotating at a
constant angular velocity ωn where the angle from the origin of the vector is φ , defined as φ = ωt.
Where ω is the lowercase Greek letter Omega and φ is the lowercase Coptic letter phi. This is
similar to a Greek phi (ϕ) and either can be used in this context. The subscript n on ω denotes
that this frequency relates to the natural frequency of the system, the only frequency in simple
harmonic motion. A visualization of the harmonic motion obtained from projecting the point on
the edge of a vector onto the ωnt space is presented in figure 2.2.
Figure 2.2: Harmonic motion represented at the projection of a point on the end of a vector moving
on a circle. Note the axis ωnt.
2.1 Mathematical Modeling of Free Vibration

The Development of a mathematical model for a system under free vibration would enable the
practitioner to predict, or model, the vibrating system of interest. Therefore, considering the fol-
18
Open Vibrations 2.1 Mathematical Modeling of Free Vibration
lowing system,
Figure 2.3: 1-DOF spring-mass system.
can be modeled expressed with the following EOM
mẍ(t) + kx(t) = 0 (2.1)
it becomes prudent to solve this homogeneous ordinary differential equation (ODE) to obtain a
model of the vibrating system. The simplest method for solving an ODE is to propose a solution
based on observations of a vibrating physical systema . Figure 2.4 reports and annotates the key
components from an observation of a vibrating system.
Figure 2.4: Summary of the temporal response for a 1-DOF system.
where x0 and v0 are the is the displacement and velocity at t=0 (i.e. the initial displacement).
A mathematical expression can now be formulated to represent the observed simple harmonic
motion. This expression can be based on the projection of a point on a vector (transposed into the
time domain) or assembled from constituent parts as done in what follows. Solving for a location
x, at a time t; x(t), the various characteristics of the expression can be identified:
2.1.1 Solve for the Natural Frequency (ωn ) of the System

• System oscillates → a sin function models this
• System oscillates at different speed → use a parameter to adjust ωn in rad/s.

a Calledan ansatz solution. Oxford Languages definition: noun, MATHEMATICS, an assumption about the form
of an unknown function which is made in order to facilitate solution of an equation or other problem.
19
• Systems have different amplitudes → use a parameter to adjust A in meters.
• System has different starting points → use a parameter to adjust φ in rad.
Using these four constituent components, an equation can be proposed:
x(t) = Asin(ωnt + φ ) (2.2)
Take the derivative to get velocity:
ẋ(t) = Aωn cos(ωnt + φ ) (2.3)
Take the derivative again to get acceleration:
ẍ(t) = −Aωn2 sin(ωnt + φ ) (2.4)
Substituting x and ẍ into the EOM for the considered 1-DOF system (mẍ(t) + kx(t) = 0) yields:
m − Aωn2 sin(ωnt + φ ) + k Asin(ωnt + φ ) = 0

(2.5)
Thereafter, dividing both sides by Asin(ωnt + φ ) results in the expression:
−mωn2 + k = 0 (2.6)
This expression can be rearranged into the more useful standard form:
r
k
ωn = (2.7)
m
Equation 2.7 represents a solution to the EOM presented in equation 2.1. This solution is not in
the form of an ODE so, therefore, we can experientially prove that this is the correct solution.
For example, we could build a system with known mass and stiffness and measure the natural
frequency of the system. Equation 2.7 equation leads to:
2π
T= (2.8)
ωn
where T is the period of oscillations and
ωn
fn = (2.9)
2π
where fn is the frequency of the oscillations.
2.1.2 Solve for Initial Phase (φ ) of the System

The EOM is a second-order ODE so there needs to exist two initial conditions (constants) to solve
it. For the systems under consideration, the displacement (x) and velocity (ẋ or v) at t = 0 are the
initial conditions. For simplicity, these are written as
x(0) = x0 (2.10)
20
ẋ(0) = v(0) = v0 (2.11)

Setting the equation to its initial state t = 0, the equations for displacement and velocity can be
simplified to:
x(0) = x0 = Asin(ωn 0 + φ ) = Asin(φ ) (2.12)
ẋ(0) = v0 = Aωn cos(ωn 0 + φ ) = Aωn cos(φ ) (2.13)
Thereafter, mathematical meanings for φ and A can be derived. To do this, φ can be solved for by
rearranging equations 2.12 and 2.13 for A:
x0
A= (2.14)
sin(φ )
and:
v0
A= (2.15)
ωn cos(φ )
Setting these two equations equal to each other cancels out A and creates:
x0 ωn v0
= (2.16)
sin(φ ) cos(φ )
therefore:
x0 ωn sin(φ )
= (2.17)
v0 cos(φ )
finally:
−1 x0 ωn
φ = tan (2.18)
v0
2.1.3 Solve for Amplitude (A) of the System

The amplitude of the vibrating system (A) is solved for in a similar manner to φ where the expres-
sions for x and ẋ are solved for at t = 0 and rearranged as to isolate φ . This operation results in the
equations:
x0
sin(φ ) = (2.19)
A
and:
v0
cos(φ ) = (2.20)
ωn A
From these equations a value for φ can be obtained knowing that sin(φ )2 +cos(φ )2 = 1. Therefore:
2 2
x0 v0
+ =1 (2.21)
A ωn A
ωn
multiplying each expression by 1 (also expressed as ωn ), gives the equation:
2 2 2
ωn x0 v0
+1 = 1×1 (2.22)
ωn A ωn A
21
which becomes: 2 2
ωn x0 v0
+ =1 (2.23)
ωn A ωn A
Further simplification is obtained by multiplying each side by (ωn A)2 to obtain:
ωn2 x02 + v20 = A2 ωn2 (2.24)
Solving for A, this expression rearranges to:

q s
ωn2 x02 + v20 2
v0
A= 2
= x0 + (2.25)
ωn ωn
2.1.4 Response for Simple Harmonic Motion

The time-varying displacement of a 1-DOF vibrating system under free response is expressed by
the equation x(t) = Asin(ωnt + φ ). Substituting in the expressions for A and φ results in:
q
ωn2 x02 + v20
x ω
0 n
!
x(t) = sin ωnt + tan−1 (2.26)
ωn v0
This equation provides a mathematical solution that relates displacement of the mass to the initial
conditions x0 and v0 . The solution is considered a free-response because no input is applied after
t=0. The relationship between the initial conditions (x0 and v0 ) and the amplitude and phase of the
response can be expressed using the Pythagorean theorem, a2 + b2 = c2 , as annotated in figure 2.5.
Figure 2.5: Trigonometric relationship between the initial conditions (x0 and v0 ), amplitude A, and
phase φ for free vibration of a 1-DOF system.
2.1.5 Special Considerations for No Initial Velocity (v0 = 0)

Upon close inspection of the temporal solution in equation 2.26, it becomes evident that any system
without initial velocity (i.e. v0 = 0) results in an undefined number for (x0 ωn )/v0 . A solution to
22
this challenge lies in the fact that the limit of tan−1 (x) approaches −π/2 at −∞ and π/2 at ∞, as
depicted in figure 2.6. Therefore, the solution at −∞ and ∞ is undefined, resulting in the expression:

x0 ωn π
= ± , when v0 = 0 (2.27)
v0 2
This step is applied in IEEE floating-point arithmetic (IEEE 754) and results in either π/2 or
±π/2 depending on the rounding format used. From the practitioner’s side, it becomes important
to recognize the situation v0 = 0 and correct this value as needed.
Figure 2.6: Response of tan−1 (or arctan) for x=-60 to 60 showing that the tan−1 is undefined as x
approaches −∞ and ∞.
Example 2.1 A vehicle wheel, tire, and suspension can be modeled as an SDOF spring and
mass as depicted below: The mass of the wheel and tire is measured to be 300 kg and its
frequency of oscillation is observed to be 10 rad/sec. What is the stiffness of the wheel as-
sembly?
23
Figure 2.7: Modeling of a vehicle wheel, tire, and suspension showing: (a) Graphical repre-
sentation; and (b) a spring-mass model.
Solution:
Considering: r
k
ωn = (2.28)
m
therefore, k = mωn2 = (300 kg)(10 rad/s)2 = 30 KN/m. Note: radians are a dimensionless
quantity and as such the units of mωn2 become kg · m where the unit value m
s2 m m is added such that
kg·m 1 N
the stiffness of the spring can be expressed as s2 · m = m .
Example 2.2 Consider the following 1-DOF system, where k = 857.8 N/m and m = 49.2 ×
10−3 kg, and calculate the natural frequency in rad/s and Hz. Also, find the period of oscilla-
tions and the maximum displacement if the spring is initially displaced 10 mm with no initial
velocity.
Figure 2.8: 1-DOF spring-mass system.
Solution:
24
Open Vibrations 2.2 General Solution for Vibrating Systems
r r
k 857.8
ωn = = = 132 rad/sec (2.29)
m 49.2 × 10−3
In Hz, this is:
ωn
fn = = 21 Hz (2.30)
2π
The period is:
2π
T== 0.0476 s (2.31)
ωn
The maximum displacement will happen when sin(ωnt + φ )=
√ 0, therefore, the value of A is
n 0 0 ω 2 x2 +v2
the maximum displacement. For an undamped system, A = ωn ,
√
q
ωn2 x02 + v20 1322 0.012 + 02
A= = = 0.01 m (2.32)
ωn 132
2.2 General Solution for Vibrating Systems

The EOM for a vibrating system has many solutions and can be expressed in various forms in-
cluding a general solution. These forms offer different mathematical approaches to solve the same
1-DOF spring-mass system and relate to each other through Euler’s equations.
Review 2.1 Vibration analysis uses complex numbers to solve the EOM’s differential equa-
tion. In this text the imaginary number is termed j (sometimes referred to as i): such that:
√
j= −1 (2.33)
and:
j2 = −1 (2.34)
A general complex number, x, can be expressed as:
x = a+bj (2.35)
here, a is referred to as the real number and b is the imaginary part of the number x. Such
complex numbers can be represented in the complex plane, also called an Argand plot. The
absolute value or modules is defined as |x| presented on the complex plot.
25
Figure 2.9: A conjugate pair of numbers (A and A∗) represented on the complex plane.
A and A∗ prime are complex conjugate pairs. In mathematics, the complex conjugate
of a complex number is the number with an equal real part and an imaginary part equal in
magnitude but opposite in sign. In other words, a conjugate pair is a + b j and a − b j.
Definition - con·ju·gate (adjective): Coupled, connected, or re-

lated.
Review 2.2 Euler’s (pronounced oy-ler) formula, is a mathematical formula in complex

analysis that establishes the fundamental relationship between the trigonometric functions
and the complex exponential function. Euler’s formula states that for any real number x,
e jψ = cos(ψ) + jsin(ψ) (2.36)

√
where j = −1. This equation can also be expressed as:
e− jψ = cos(ψ) − jsin(ψ) (2.37)
This can be expressed in terms of polar coordinates as:
26
Figure 2.10: Euler’s formula illustrated on the unit circle in the complex plane.
Figure 2.11: A Soviet Union stamp from 1957 with a Portrait of Leonhard Euler who worked
in various branches of the Imperial Russian Academy of Sciences and Imperial court during
his lifetime a .
Euler’s formula is named after the Swiss engineer and mathematician Leonhard Euler
(1707-1783), who among other things popularized the use of the Greek letter π to denote
the ratio of a circle’s circumference to its diameter, wast the first to use the expression f (x)
to denote a function, and correctly defining the base of the natural logarithm e; which is
now known as Euler’s number. While Euler developed “Euler’s formula” in 1748, it was
27
not used to describe points in a complex for another 50 years when the Danish-Norwegian
mathematician and cartographer Caspar Wessel presented to the Danish Academy in 1797
b.
a Post
of the USSR, Public domain, via Wikimedia Commons
b Whittaker,
Edmund Taylor, and George Neville Watson. A course of modern analysis: an introduction to
the general theory of infinite processes and of analytic functions; with an acount of the principal transcendental
functions. University Press, 1927.
2.2.1 Formulating the General Solution for a 1-DOF Spring-Mass System

We can also solve the following EOM as an elementary differential equation:
mẍ + kx = 0 (2.38)
in a more analytical manner using the theory of elementary differential equations. To do this the
form:
x(t) = aeλt (2.39)
is assumed, where a and t are nonzero constants that need to be determined. Using successive
differentiation, the assumed solution becomes:
ẋ(t) = λ aeλt (2.40)
and
ẍ(t) = λ 2 aeλt (2.41)
therefore, mẍ(t) + kx(t) = 0 becomes:
mλ 2 aeλt + kaeλt = 0 (2.42)
Next, the above expressions is divide by aeλt to obtain the characteristic equation:
mλ 2 + k = 0 (2.43)
This can be done because aeλt is never zero, therefore, the expression is never divided by zero.
The quadratic formula gives us:
r r
k k
λ =± − =± j = ±ωn j (2.44)
m m
q
remember that ωn = mk . Notice that the ± tells us there are two solutions to this problem. So,
putting λ back into the assumed solution results in two solutions (one positive, one negative):
x(t) = a1 e+ωn jt (2.45)
and
x(t) = a2 e−ωn jt (2.46)
As these solutions only consider, and are only valid for, linear systems; the sum of the solutions is
also a solution. This simplification results in:
x(t) = a1 e+ωn jt + a2 e−ωn jt (2.47)
where a1 and a2 are constants of integration that scale the unit Euler’s vector. The positive and
negative values in the exponent indicate that the terms are a conjugate pair.
28
Example 2.3 Show that x(t) = a1 e+ωn jt + a2 e−ωn jt is equal to Asin(ωnt + φ ).
Solution:
This equation was derived using Euler’s formula and it can be shown that this equation is
equivalent to the Asin(ωn + φ ). To recover the previously assumed solution, the knowledge
that a1 and a2 are complex congregate pairs and as such the magnitude can be expressed as
a1 = a2 is leveraged. Using Euler’s polar notation, a1 and a2 can be expressed as
a1 = a2 = ae jψ (2.48)
where a and ψ are real numbers, the equation becomes:
x(t) = ae j(ωnt+ψ) + ae− j(ωnt+ψ) (2.49)
this becomes:
x(t) = a(e j(ωnt+ψ) + e− j(ωnt+ψ) ) (2.50)
Remembering Euler’s equations from before, this becomes:

x(t) = a cos(ωnt + ψ) + jsin(ωnt + ψ) + cos(ωnt + ψ) − jsin(ωnt + ψ) (2.51)
combining the “cos” terms and canceling out the “sin” terms this becomes:
x(t) = 2a · cos(ωnt + ψ) (2.52)
This is equivalent to x(t) = Asin(ωnt + φ ) considering that A = 2a and knowing sin(φ ) =

cos(φ + ψ). To expand, this is because the sin and cos are only differentiated by a phase shift.
Next, a general solution for the EOM is obtained. Using the previous solution:
x(t) = a1 e+ωn jt + a2 e−ωn jt (2.53)
we can expand this into the form:

x(t) = a1 cos(ωnt) + jsin(ωnt) + a2 cos(ωnt) − jsin(ωnt) (2.54)
using trigonometric functions. This equates to:
x(t) = (a1 + a2 ) · cos(ωnt) + (a1 − a2 ) j · sin(ωnt) (2.55)
As x(t) is always real, A1 and A2 can be defined as:
A1 = (a1 + a2 ) (2.56)
and
A2 = (a1 − a2 ) j (2.57)
Lastly, the general solution is written as:
x(t) = A1 cos(ωnt) + A2 sin(ωnt) (2.58)
29
This is the general solution for the EOM (mẍ + kx = 0) of the considered oscillating system where
A1 and A2 are defined as: q
A= A21 + A22 (2.59)
and
A1
−1
φ = tan (2.60)
A2
These are obtained from a trigonometric relationship, similar to that used before:
Figure 2.12: Trigonometric relationship between the initial conditions, amplitude, and phase, for
free vibration of a 1-DOF system expressed with: (a) variables for initial conditions; and (b)
generic variables A1 and A2 .
again, A and φ are: q v !2

ωn2 x02 + v20
u
u v0
A= = tx02 + (2.61)
ωn ωn

−1 x0 ωn
φ = tan (2.62)
v0
Example 2.4 Using the general solution:
Calculate the values of A1 and A2 in terms of their initial conditions x0 and v0 .

Solution:
Knowing the following for x and ẋ:
ẋ(t) = −A1 ωn sin(ωnt) + A2 ωn cos(ωnt) (2.65)

Now apply the initial conditions, x(0) = 0 and v(0) = 0, this yields:
x(0) = x0 = A1 (2.66)
30
Open Vibrations 2.3 Damping
ẋ(0) = v0 = A2 ωn (2.67)
Solving for A1 and A2 shows us:
v0
A1 = x0 , and A2 = (2.68)
ωn
thus:
v0
x(t) = x0 cos(ωnt) + sin(ωnt) (2.69)
ωn
2.2.2 Solution of 1-DOF System in Three Forms

Form one, for mẍ + kx = 0 subject to nonzero initial conditions can be written as:
x(t) = a1 e+ωn jt + a2 e−ωn jt (2.70)
where a1 and a2 are complex terms. Form two is:
x(t) = Asin(ωnt + φ ) (2.71)
while form three is:
where A, φ , A1 , and A2 , are all real-valued constants. Each set of constants can be related to each
other by:

−1 A1
q
2
A = A1 + A2 2 φ = tan (2.73)
A2
A1 = (a1 + a2 ) A2 = (a1 − a2 ) j (2.74)
A1 − A2 j A1 + A2 j
a1 = a2 = (2.75)
2 2
Which follow from trigonometric identities and Euler’s formula.
2.3 Damping
The response of a spring-mass system predicts that a system will oscillate indefinitely. However,
we know that this is not true from observing real-world solutions. So based on real-world obser-
vations and mathematical conveniences, we need to add a term that will remove “energy” from
the system with time. To do this the idea of the ideal dashpot is introduced. A linear dashpot is
diagrammed in figure 2.13 and is a mechanical device that resists motion via viscous friction and
therefore converts the mechanical energy of the system into thermal energy that is dissipated.
Figure 2.13: Schematic of a liner dashpot showing the damping force ( fc ) acting in the opposite
direction of the displacement (x).
31
Just as spring forms a physical model of the cause vibration, through its storage and release
of energy, a dashpot (sometimes called a damper) forms a physical model for dissipating energy.
Dashpots create a resisting or damping force that acts opposite to the direction of travel (as anno-
tated in figure 2.13) and is proportional to the velocity. Therefore, the damping forces fc can be
computed as:
fc = cẋ (2.76)
the constant c, called the damping coefficient, has the units of kg/s. Dashpots are a mathematical
representation of viscous dampers installed in automobiles, aircraft, structures, and other mechan-
ical devices. However, all systems have inherent damping not just systems with physical dampers.
The spring-mass system can be used as a representation of real-world systems with inherent damp-
ing as demonstrated by the rubber engine mount depicted in figure 2.14.
Figure 2.14: Modeling of a rubber engine mount as an spring-dashpot-mass model showing (a) the
rubber engine mount; (b) idealized model of the rubber month; and (c) the FBD of the idealized
model
Depending on the amount of damping present in a system, the temporal response of the system
will represent itself in various ways, as represented in figure 2.15. To reiterate, an undamped case
will oscillate around the equilibrium and does not decay. If a limited amount of damping is present
in a system it will oscillate around the equilibrium and slowly decay with time to the equilibrium
position, this is termed underdamped. If an excessive amount of damping is present, the system
will not oscillate but decay directly to the equilibrium position, this is termed the overdamped
case. Lastly, there exists a special case that results in the system converging as quickly as possible
to the equilibrium position without oscillations; this case is termed the critically damped case.
Furthermore, the amount of damping required to obtain a critically damped system is the damping
value between the underdamped and overdamped cases for a specific system. To recap, the key
types of damping are:
• Undamped - Oscillates around the equilibrium and does not decay.
• Underdamped - Oscillates around the equilibrium and slowly decays and is the most com-
mon case.
• Overdamped - Does not pass the equilibrium position and is a simple decay with no oscil-
lation.
• Critically damped - provides the quickest approach to zero amplitude for a damped oscil-
lator.
32
Figure 2.15: Temporal responses for the three types of damping: underdamped, overdamped, and
critically damped.
Vibration Case Study 2.1 Dampers are used to extract energy from systems in an effort to
reduce their vibrations. The Author Ravel Junior Bridge in Charleston South Carolina is
a cable-stayed bridge over the Cooper River with a main span of 471 m (1,546 feet). The
bridge uses two dampers connected to the cables in the middle of the bridge to dissipate
excess energy in the cables that would otherwise cause unwanted vibrations from wind,
traffic, and seismic activity. The inclusion of these dampers on the bridges is a proactive
measure to enhance the bridge’s performance, safety, and user experience by reducing the
effects of vibrations.
Figure 2.16: Dampers installed on the center cables of the Author Ravel Junior Bridge in
Charleston South Carolina, showing: (a) the main span of the bridge taken from the watera
with the damper structure annotated; and (b) close up of the structure that holds the damper.
a originalun-annotated image by bbatsell, CC BY-SA 2.5 <https://creativecommons.org/licenses/by-
sa/2.5>, via Wikimedia Commons
33
2.3.1 Modeling Vibrating Systems with Damping

The spring-mass system of Chapter 1 can be expanded to a spring-dashpot-mass system that con-
siders the damping component of the system. A mathematical model of the spring-dashpot-mass
system can be developed for the case present in figure 2.17. Using the FBD for the system, it can
conclude that the EOM for this system:
Figure 2.17: Spring-dashpot-mass model showing: (a) a schematic of the system; and (b) the FBD
of the system.
is:
mẍ(t) = − fc − fk (2.77)
Rearranging into standard form and concerting forces into parameters c and k results in:
mẍ(t) + cẋ(t) + kx(t) = 0 (2.78)
This system is subject to the same initial conditions as before, x(0) = x0 and ẋ(0) = v0 . Again,
choosing to model it this way for convinces, so let’s solve it in a similar manner to the EOM
without damping. Again, assume the solution:
x(t) = aeλt (2.79)
here, a and t are nonzeros constants that need to be determined. Using successive differentiation,
we get:
ẋ(t) = λ aeλt (2.80)
and
ẍ(t) = λ 2 aeλt (2.81)
therefore, mẍ + cẋ + kx = 0 becomes:
mλ 2 aeλt + cλ aeλt + kaeλt = 0 (2.82)
Now we divide by aeλt to obtain the characteristic equation:
mλ 2 + cλ + k = 0 (2.83)
We can do this because aeλt is never zero, therefore, we never divide by zero. The quadratic
formula gives us: √
−c ± c2 − 4km −c 1 p 2
λ1,2 = = ± c − 4km (2.84)
2m 2m 2m
Some key points from this equation:
34
• The ± tells us there are two solutions to this problem
• if c2 − 4km < 0, system is Underdamped, solutions are complex conjugate pairs
• if c2 − 4km = 0, system is critically damped, solutions are equal negative real numbers
• if c2 − 4km > 0, system is Overdamped, solutions are distinct negative real numbers
From this, we can see that c2 − 4km = 0 is a special value, let us define a value for c that will give
us this critical damping number. We will call it the critical damping coefficient (ccr ). So setting
the equation as:
c2cr − 4km = 0 (2.85)
giving us:
c2cr = 4km (2.86)
next, we can derive the function:
√ m √
√
ccr = 2 km = 2 √ km = 2mωn (2.87)
m
q
remember that ωn = mk for an undamped system. Next, we generate a non-dimensional number
(ζ ), pronounced ‘zeta’ that will allow us to distinguish between different types of damping. ζ is
called the critical damping ratio.
c c c
ζ= = √ = (2.88)
ccr 2 km 2mωn
Now if we put the ζ back into the characteristic equation and resolve using the quadratic equation
we get: q
λ1,2 = −ζ ωn ± ωn ζ2 −1 (2.89)
From this equation, it becomes clear that ζ determines whether the roots are complex or real, this,
in turn, determines the nature of the response of the structure. Listing the possible responses we
get: For each damping case, we will have a different solution to the problem.
damping case critical damping ratio radicand solutions

under damped 0<ζ <1 c2 − 4km < 0 complex conjugate pairs
critically damped ζ =1 c2 − 4km = 0 equal negative real numbers
over damped 1<ζ c2 − 4km > 0 distinct negative real numbers
2.3.2 Modeling Underdamped Motion

In the case that 0 < ζ < 1, a√complex conjugate pair of roots are the solutions to the characteristic
equation after pulling out a −1:
q
λ1 = −ζ ωn + ωn 1 − ζ 2 j (2.90)
35
and: q
λ2 = −ζ ωn − ωn 1 − ζ 2 j (2.91)
Where the j is pulled out because:
q q q
1 − ζ 2 j = (1 − ζ 2 )(−1) = ζ 2 − 1 (2.92)
Next, let us “arbitrarily” define: q

ωd = ωn 1 − ζ 2 (2.93)
where ωd is the damped natural frequency. Therefore, the equations become:
λ1 = −ζ ωn + ωd j (2.94)
and:
λ2 = −ζ ωn − ωd j (2.95)
Again, we have two solutions to a linear problem, so we can combine these into one solution and
insert λ into the assumed solution aeλt to obtain:
x(t) = a1 e−ζ ωnt+ωdt j + a2 e−ζ ωnt−ωdt j (2.96)
where a1 and a2 are complex valued constants. This can now be simplified into:
x(t) = e−ζ ωnt (a1 eωdt j + a2 e−ωdt j ) (2.97)
Using Euler’s equations, (same as before) and choosing:
A1 = (a1 − a2 ) j (2.98)
and
A2 = (a1 + a2 ) (2.99)
Note that the A1 and A2 defined here are the reverse of those defined in Eq. 2.73. This is done to
allow the general form to be in the same format as before, however, assuming the same A1 and A2
would not change the final solution expressed below. The general form of this solution is then:
x(t) = e−ζ ωnt A1 sin(ωdt) + A2 cos(ωdt)

(2.100)
Recall that for undamped 1-DOF systems we showed
x(t) = Asin(ωnt + φ ) = A1 sin(ωnt) + A2 cos(ωnt) (2.101)
As e−ζ ωnt accounts for the damping, our current solution becomes:
x(t) = Ae−ζ ωnt sin(ωdt + φ ) (2.102)
Now that we have x and ẋ, we can solve for the boundary conditions x0 and v0 by setting t = 0,
we get:
x(0) = x0 = Asin(φ ) (2.103)
36
and taking the directive of x(t) using the product rule (fg)’= f’g+fg’, we get:
ẋ(t) = −ζ ωn Ae−ζ ωnt sin(ωdt + φ ) + Ae−ζ ωnt ωd cos(ωdt + φ ) (2.104)
ẋ(0) = v0 = −ζ ωn Asin(φ ) + Aωd cos(φ ) (2.105)
a simplification can be made to the prior equation by letting A = x0 /sin(φ ). This gives us the
equation:
x0 x0
ẋ(0) = v0 = −ζ ωn sin(φ ) + ωd cos(φ ) (2.106)
sin(φ ) sin(φ )
that can be simplified to:
ẋ(0) = v0 = −ζ ωn x0 + x0 ωd cot(φ ) (2.107)
The above equation related v0 to φ using terms that are known for a giving system (ζ , ωn , x0 , and ωd ).
Therefore, this expression can be used to solve for φ :
v0 + ζ ωn x0
cot(φ ) = (2.108)
x0 ωd
and as tan(φ ) = 1/cot(φ ): !
x0 ωd
φ = tan−1 (2.109)
v0 + ζ ωn x0
Thereafter, we can solve for A considering the fact that we sent A = x0 /sin(φ ). Using the trigono-
metric relationship between expressed in equation 2.108 and visualized in figure 2.18:
Figure 2.18: Trigonometric relationship between the initial conditions (x0 and v0 ), amplitude A,
and phase φ for underdamped motion of a 1-DOF system.
we show that sin(φ ) can be expressed as:

x0 ωd
sin(φ ) = p (2.110)
(v0 + ζ ωn x0 )2 + (x0 ωd )2
and applying A = x0 /sin(φ ) we get:
v !2
p u
(v0 + ζ ωn x0 )2 + (x0 ωd )2 u v0 + ζ ωn x0
A= = tx02 + (2.111)
ωd ωd
37
Finally, collecting all of our important equations:

√
• Critical damping coefficient: ccr = 2 km = 2mωn
c √c c
• Damping ratio: ζ = ccr = = 2mωn
2 km
p
• Damped natural frequency: ωd = ωn 1−ζ2
• Solution for underdamped system: x(t) = Ae−ζ ωnt sin(ωdt + φ ), where:

p !
(v0 + ζ ωn x0 )2 + (x0 ωd )2 x0 ωd
A= φ = tan−1
ωd v0 + ζ ωn x0
Figure 2.19: Four example responses for an underdamped 1-DOF system (ζ = 0.142) with various
initial conditions.
Example 2.5 Consider the following 1-DOF system, where k = 857.8 N/m, c = 7.8 kg/s, and
m = 49.2 × 10−3 kg, calculate the damped frequency in rad/s and Hz. What damping case is
this system?
38
Figure 2.20: 1-DOF spring-dashpot-mass system.
Solution:
Calculate the undamped frequency:
r r
k 857.8
ωn = = = 132 rad/s (2.112)
m 49.2 × 10−3
The systems critical damping value:
√ p
ccr = 2 km = 2 k = 857.8 · 49.2 × 10−3 = 12.993 kg/s (2.113)
And the critical damping ratio:

c 7.8
ζ= = = 0.600 (2.114)
ccr 12.993
This can also be expressed as 60% damped, this is an underdamped system, and the system
will oscillate. Now we can calculate the damped frequency:
q p
ωd = ωn 1 − ζ 2 = ωn 1 − 0.6002 = 105.6 rad/s (2.115)
Therefore, the system oscillates at 105.6 rad/sec or 16.81 Hz
Example 2.6 For a damped one DOF system where m, c, and k are known to be m = 1 kg,
c = 2 kg/s, and k = 10 N/m. Calculate the value of ζ and ωn . Is the system overdamped,
underdamped, or critically damped?
Solution:
The natural frequency is calculated as
r r
k 10
ωn = = = 3.16 rad/s (2.116)
m 1
39
The damping can be calculated as:

c 2 1
ζ= = q = √ = 0.316 (2.117)
2ωn m 2 10 10
1 (1)
So the damped natural frequency is equal to:

s 2
√

1
q
2
ωd = ωn 1 − ζ = 10 1 − √ = 3 rad/s (2.118)
10
As 0 < ζ < 1 the system is underdamped.
Example 2.7 Figure 2.21 shows an industrial device consisting of a mass isolated from its
fixtures by two rubber dampers and an offset spring with an angle α. Provide an estimate of
the system’s damped natural frequency in the vertical direction. Assume the rubber dampers
add damping and only negligible stiffness to the system and that the spring is long enough
such that the angles remain constant.
Figure 2.21: Industrial device (mass) connected to a fixed point with a rubber damper and
spring at an angle.
Solution:
First and foremost, we need to develop a mass-spring-dashpot representation of the sys-
tem. This is presented in what follows:
40
Figure 2.22: Mass-spring-dashpot representation of the industrial system represented if figure

2.21.
where the damping in the vertical direction provided by the rubber damper is modeled as a
dashpot in the vertical direction. As we only want an estimate of the frequency, the assumption
that the is small and as such α of the displaced state is equal α of the equilibrium state. This
leads to the FBD for the equilibrium and displaced states:
equilibrium position displaced position “x”
+↓ ∑ Fx = mg − kδ cos(α) = 0
+↓ ∑ Fx = mg − cẋ − kcos(α)(∆l + δ )
Applying Newton’s second law and combining these equations yields:
mẍ + cẋ + k∆lcos(α) = 0
41
Looking at the triangles formed by the dashpot and spring it can be shown that:
cos(α) = h/l = x/∆l
As we assumed the displacement is small and α remains unchanged. Therefore the prior
equation becomes:
x
mẍ + cẋ + k∆l = 0
∆l
This simplifies to the “normal” EOM for a 1-DOF system:
mẍ + cẋ + kx = 0
Therefore, once the values for the system are measured the system’s damped natural frequency
in the vertical direction can be estimated as
q
ωd = ωn 1 − ζ 2
Vibration Case Study 2.2 Epoxy-based damping is used in the automotive field to reduce
vehicle noise and vibration harshness (NVH). This passive form of damping is essential for
increasing comfort in modern light-wight vehicles.
Figure 2.23: Experimental modal analysis of an automotive (Jaguar) body in white, typically
done to reduce vehicle noise and vibration harshness.a
a Cjp24, CC BY-SA 3.0 <https://creativecommons.org/licenses/by-sa/3.0>, via Wikimedia Commons.
2.3.3 Modeling Overdamped Motion

In the case of overdamped systems, 1 < ζ , the solutions for λ are distinct real roots that are written
as: q
λ1 = −ζ ωn − ωn ζ 2 − 1 (2.119)
42
and: q
λ2 = −ζ ωn + ωn ζ 2 − 1 (2.120)
The solution for the EOM using the assumed solution then becomes:
√ √
2 2
x(t) = e−ζ ωnt (a1 e−ωnt ζ −1 + a2 e+ωnt ζ −1 ) (2.121)
This equation represents a non-oscillating response of the system. Again, a1 and a2 are solved for
using known boundary conditions x0 and v0 such that:
p
2
−v0 + − ζ + ζ − 1 ωn x0
a1 = p (2.122)
2ωn ζ 2 − 1
p
2
v0 + ζ + ζ − 1 ωn x0
a2 = p (2.123)
2ωn ζ 2 − 1
Typical responses for an overdamped system with various initial conditions are shown below in
figure 2.24.
Figure 2.24: Four example responses for an overdamped 1-DOF system (ζ = 2.371) with various
initial conditions.
2.3.4 Modeling Critically Damped Motion

In the case of critically damped systems, ζ = 1, the solutions for λ will be equal negative real
numbers, therefore from before:
q
λ1,2 = −ζ ωn ± ωn ζ 2 − 1 (2.124)
We get:
λ1 = λ2 = −ωn (2.125)
43
Because both solutions (a1 and a2 ) are the same, we multiply the second solution by t so the
solution for a critically damped system is in the same form as before. The solution for the EOM
using the assumed solution then becomes:
x(t) = a1 e−ωnt + a2te−ωnt (2.126)
This simplifies into:

x(t) = (a1 + a2t)e−ωnt (2.127)
This equation represents a non-oscillating response of the system. Again, a1 and a2 are solved for
using known boundary conditions x0 and v0 such that:
a1 = x0 (2.128)
a2 = v0 + ωn x0 (2.129)
2.3.5 Standard Form of the EOM

The EOM for a damped 1-DOF system is written in a “standard form” in which the effect of the
damping ratio and natural frequencies are more obvious. To get to the standard form, the normal
form of the EOM:
mẍ + cẋ + kx = 0 (2.130)
is divided by what the constant terms associated with the acceleration term. In this example, this
is m. Dividing every term by m yields:
c k
ẍ + ẋ + x = 0 (2.131)
m m
Numerical manipulations can be undertaken to get the coefficients of the velocity and displacement
terms into coefficients that more clearly express the characteristics of the vibrating system:
ẍ + 2ζ ωn ẋ + ωn2 x = 0 (2.132)
Example 2.8 An engine valve assembly is depicted in figure 2.25 where J is the inertia caused
by the right-hand side of the rocker arm. Derive an analytical solution for the natural frequency
of the rocker arm. Use the assumptions sin(θ ) = θ and cos(θ ) = 1.
44
Figure 2.25: Rocker arm assembly of an internal combustion engine showing: (a) a diagram
of the system and; (b) the FBD of the system.
Solution:
Taking the sum of the moments about O and considering the inertia caused by the right-
hand side of the rocker arm, J, the FBDs can be written as:
equilibrium position displaced position “x”
+ ∑ Mo = mgl − kl 2 δ = 0
↶
+ ∑ Mo = mgl − kl 2 δ − kl 2 θ − cl 2 θ̇ = 0
↶
Applying Newton’s second law and combining these equations yields:
(J + ml 2 )θ̈ + cl 2 θ̇ + kl 2 θ = 0 (2.133)
45
Therefore, the standard form of the EOM is:
cl 2 kl 2
θ̈ + θ̇ + θ =0 (2.134)
J + ml 2 J + ml 2
Results in the following analytical solution for the natural frequency:
s
kl 2
ωn = rad/s (2.135)
J + ml 2
Vibration Case Study 2.3 On August 17th , 2009 Turbine 2 of the hydroelectric power sta-
tion of the Sayano-Shushenskaya Dam near Sayanogorsk in Russia failed catastrophically.
The failure flooded the turbine hall and collapsed the ceiling. Killing 75 people, many of
whom were in the turbine hall to celebrate the anniversary of the plant’s general director.
Turbines of the type used at Sayano-Shushenskaya are designed to have high efficiency but
a very narrow working band. When they operate outside the designed working band, they
vibrate due to the pulsation of water flow and water strokes. These vibrations degrade the
turbine over time.
Turbine 2 had experienced excessive vibrations for a long time, ever since its installation
in 1979. Through the early ’80s several issues were fixed, along with substantial repairs in
2000 and 2005. In July 2009 the turbine again exceeded the allowed vibration specification
but stayed in operation. Over the years, the operating staff simply came to accept the higher
level of vibration. The final government report stated that the accident was caused by turbine
vibrations which led to fatigue damage in a turbine mount.
Figure 2.26: Sayano Shushenskaya’s turbine hall before the accident where turbine 2 (the
turbine that failed) is in the foreground of the image. a
46
Open Vibrations 2.4 Logarithmic Decrement
Figure 2.27: Sayano Shushenskaya’s turbine hall after the accident. b
a 4044415 Russian: Andrey Korzun English: Andrey Korzun, CC BY-SA 3.0

<https://creativecommons.org/licenses/by-sa/3.0>, via Wikimedia Commons
b The original source of this photograph is believed to be Jaffaa, a user of forums.drom.ru, who uploaded it
on August 17, 2009, a few hours after the accident. This image is a faithful digitization of a unique historic
image, and the copyright for it is most likely held by the person who created the image or the agency employing
the person. It is believed that the use of this image may qualify as fair use under the copyright law of the United
States.
2.4 Logarithmic Decrement

For a vibrating system, the mass (m) and stiffness (k) can be measured using scales and static
deflection tests. However, the damping coefficient (c) is a more difficult quantity to determine.
From k and m we can compute the natural frequency (ωn ) and the critical damping coefficient
(ccr ). Therefore, knowing that the critical damping ratio (ζ ) is defined as:
c c c
ζ= = √ = (2.136)
ccr 2 km 2mωn
if we calculate ζ , we can obtain c for the system of interest. This is made possible because ccr can
be calculated from k and m. Observing the temporal response for the underdamped system,
47
Figure 2.28: Measuring the peak displacement points in an experimental system with decay caused
by damping.
we mark three points of maximum amplitude, x1 , x2 , and x3 that happen at t1 , t2 , and t3 , respectively.
Considering displacement values for the first two points x1 and x2 , separated by a complete period
(T ). Knowing that one cycle is 2π, the time period for this complete cycle is given by:
2π 2π
t2 − t1 = = p (2.137)
ωd ωn 1 − ζ 2
where ωd is the damped natural frequency. This is the time period (T ) of damped oscillations. If
we derive an equation for the values of the peaks, also called the envelope of maximum values, we
get:
xpeaks = Ae−ζ ωnt (2.138)
Knowing that the system is underdamped, A can be solved for using the initial conditions x0 and
v0 , therefore: p
(v0 + ζ ωn x0 )2 + (x0 ωd )2
A= (2.139)
ωd
In terms of t1 and t2 , we can express the displacement at these times as:
x1 = Ae−ζ ωnt1 (2.140)
and
x2 = Ae−ζ ωnt2 (2.141)
therefore:
x1 e−ζ ωnt1
= −ζ ω t = eζ ωn (t2 −t1 ) (2.142)
x2 e n2
48
2π √2π
However, from before we know that t2 − t1 = ωd = . Therefore, we can express this last
ωn 1−ζ 2
equation as:
x1 √2πζ
=e 1−ζ 2 (2.143)
x2
Next, we take the natural log of both sides to get the logarithmic decrement, denoted by δ :

x1 x(t1 ) 2πζ
δ = ln = ln =p (2.144)
x2 x(t1 + T ) 1−ζ2
This shows us that the ratio of any two successive amplitudes for an underdamped system, vibrating
freely, is constant and is a function of the damping only. Sometimes, in experiments, it is more
convenient/accurate to measure the amplitudes after say “n” peaks rather than two successive peaks
(because if the damping is very small, the difference between the successive peaks may not be
significant). The logarithmic decrement can then be given by the equation

1 x1 1 x(t1 ) 2πζ
δ = ln = ln =p (2.145)
n xn+1 n x(t1 + nT ) 1−ζ2
Once we use the experimental data to obtain δ , and knowing that:
2πζ
δ=p (2.146)
1−ζ2
we can calculate the value of ζ :

δ
ζ=√ (2.147)
4π 2 + δ 2
Therefore, having ζ we can solve for the coefficient of damping, c, as:
√
c = ζ 2 km (2.148)
Example 2.9 Calculate the damping coefficient for the system with the measured amplitude
as expressed below given that m = 3 kg and k = 43 N/m. Use t1 = 1 sec, and tn+1 = t4 = 6 sec.
Use the peaks as marked in figure 2.29.
49
Figure 2.29: Response from an experimental system with noise.
Solution:
First, from the plot we can determine that x1 = −9.5 mm and x4 = −1.8 mm where n = 3.
Thereafter, we can solve for δ :

1 x1 1 −9.5
δ = ln = ln = 0.554 (2.149)
3 x4 3 −1.8
Next, we can calculate ζ , as:
δ 0.554
ζ=√ =√ = 0.0879 (2.150)
4π 2 + δ 2 4π 2 + 0.5542
And lastly: √ √
c = ζ 2 km = 0.0879 · 2 43 · 3 = 2.0 kg/s (2.151)
Example 2.10 The free-response of a 1000-kg automobile with stiffness of k = 400,000 N/m
is observed to be underdamped. Modeling the automobile as a single-degree-of-freedom os-
cillation in the vertical direction, as annotated in figure 2.7, determine the damping coefficient
if the displacement at t1 is measured to be 2 cm and 0.22 cm at t2 .
Solution:
Knowing x1 = 2 cm and x2 = 0.22 cm and t2 = T + t1 , therefore:
x1 2
δ = ln = ln = 2.207 (2.152)
x2 0.22
50
and:
δ 2.207
ζ= √ = √ = 0.331 (2.153)
4π 2 + δ 2 4π 2 + 2.2072
therefore, we can obtain the damping coefficient as
√ p
c = 2ζ km = 2(0.331) 400, 000 · 1, 000 = 13, 256 kg/s (2.154)
51
Open Vibrations
3 Forced Vibrations
Mechanical systems are subjected to external loading. For example, a piston in an engine when
forced up and down by a crankshaft or a seat in an airplane may vibrate due to the movement of
the jet engines transmitted through the aircraft structure. In real-world situations, structures are
subjected to complex loading that are hard to measure or not fully understood.
Vibration Case Study 3.1 Tall mast light poles are excited by a wind excitation and re-
spond across their entire frequency domain. Consider the light pole located in the state of
Kansas in the United States shown figure 3.1. The structure responds more at some frequen-
cies than other frequencies, as dictated by the structure’s geometry and material properties.
Studying how structures responded to forced inputs allows for a better design of the struc-
ture.
Figure 3.1: Tall mast light pole in the central United States showing: (a) the light mast; (b)
the measured temporal response of the light pole, and; (c) the frequency domain response of
the light pole. Light pole data provided by Jian lia and discussed in detail in Shaheen et al.
b.
a Jianli, CC BY-SA 4.0, Light pole data

b Shaheen, Mona, et al. “Wind-induced vibration monitoring of high mast illumination poles.” Sensors and
Smart Structures Technologies for Civil, Mechanical, and Aerospace Systems 2022. Vol. 12046. SPIE, 2022.
52
Open Vibrations 3.1 Harmonic Excitations of Undamped Systems
3.1 Harmonic Excitations of Undamped Systems

Investigating a single-degree of freedom system for a harmonic input is useful as it can be solved
mathematically with straightforward techniques. Consider the system:
Figure 3.2: 1-DOF system with an external force (F(t)) applied, showing: (a) the system configu-
ration; and (b) the free body diagram
where F(t) is the external force applied to the mass. For simplicity, let us consider a harmonic
excitation for F(t) such that:
F(t) = F0 cos(ωt) (3.1)
note that here, ω has no subscript and is the frequency in rad/sec of the driving force. ω is often
called the input frequency, driving frequency, or forcing frequency. F0 represents the magnitude of
the applied force. Building the EOM for the system in figure 3.2 yields:
mẍ(t) + kx(t) = F0 cos(ωt) (3.2)
For convenience, we drop the “(t)” to make the writing easier. Then, we convert the EOM to the
standard form by dividing the equation by m:
ẍ + ωn2 x = f0 cos(ωt) (3.3)
where:
F0
f0 = (3.4)
m
The EOM in this form is a second-order, linear nonhomogeneous differential equation. It is non-
homogeneous because there are no terms related to x on the right-hand side of the equation. One
way to solve such an ODE is to recall that the solution for a nonhomogeneous equation is the sum
of the homogeneous and particular solutions.
x = xh + x p (3.5)
again, noting that this is a temporal solution where “(t)” is implied. First, knowing that the solution
is the sum of two parts: 1) oscillations caused by the spring/mass system; and 2) vibrations caused
by the forcing function. The oscillations caused by the spring/mass system will form the homoge-
neous while the vibrations caused by the forcing function will form the particular solution. As we
know the solution for oscillations caused by the spring/mass system from our prior investigation
of unforced systems we set the equation for the homogeneous solution to be:
xh = Asin(ωnt + φ ) (3.6)
53
Next, we will denote the particular solution as x p . x p can be determined by assuming that it is in
form of the forcing function, therefore:
f0 cos(ωt) (3.7)
becomes:
x p = Xcos(ωt) (3.8)
where, x p is the particular solution and X is the amplitude of the forced response. Our total solution
for the harmonic excitations of undamped systems now becomes:
x(t) = Asin(ωnt + φ ) + Xcos(ωt) (3.9)
This approach, of assuming that x p = Xcos(ωt), in order to determine the particular solution is
called the method of undetermined coefficients. To calculate X, first we take the equations for
x p and ẍ p :
x p = Xcos(ωt) (3.10)
ẍ p = −ω 2 Xcos(ωt) (3.11)
and substituting these into the equation of motion in standard form yields:
−ω 2 Xcos(ωt) + ωn2 Xcos(ωt) = f0 cos(ωt) (3.12)
As long as cos(ωt) ̸= 0, solving for X yields:

f0
X= (3.13)
ωn2 − ω 2
Therefore, as long as ωn ̸= ω, the particular solution can take the form:
f0
xp = cos(ωt) (3.14)
2
ωn − ω 2
This then expands to the total form:
f0
x(t) = Asin(ωnt + φ ) + cos(ωt) (3.15)
ωn − ω 2
2
Expanding this to the general form for the homogeneous solution obtains the equation:
f0
x(t) = A1 sin(ωnt) + A2 cos(ωnt) + cos(ωt) (3.16)
ωn − ω 2
2
As before, we need to determine the values for the coefficients A1 and A2 by enforcing the initial
conditions x0 and v0 . Setting the time to zero (t = 0) and solving the initial displacement leads to:
f0
x(0) = x0 = A2 + (3.17)
ωn − ω 2
2
or:
f0
A2 = x0 − (3.18)
ωn − ω 2
2
54
again, solving the equation in terms of velocity:

f0
ẋ(t) = A1 ωn cos(ωnt) − A2 ωn sin(ωnt) − ω sin(ωt) (3.19)
ωn2 − ω 2
and solving for the initial velocity at t = 0:
ẋ(0) = v0 = A1 ωn (3.20)
or:
v0
A1 = (3.21)
ωn
Therefore, combining the equations we get:
v f0 f0
0
x(t) = sin(ωnt) + x0 − 2 2
cos(ωnt) + 2 cos(ωt) (3.22)
ωn ωn − ω ωn − ω 2
As before, we can relate A1 and A2 to each other through the basic trigonometric identities. This
yields,
x(t) = Asin(ωnt + φ ) + Xcos(ωt) (3.23)
s
v0 2
A= + (x0 − X)2 (3.24)
ωn

−1 ωn (x0 − X)
φ = tan (3.25)
v0
f0
X= 2 (3.26)
ωn − ω 2
Example 3.1 For the 1-DOF system:
Figure 3.3: 1-DOF spring-mass system subjected to an external force F(t).
with k = 10 N/m, m = 2.5 kg, ω = 4 rad/sec, F0 = 0.1 kN, x0 = 1 mm, and v0 = 0 mm/s plot
the temporal responses of the system considering the free-vibration case and the excited case.
Plot these on a single plot to compare the responses.
Solution:
The free-vibration response can be plotted using the expression:
v0
ωn
55
Open Vibrations 3.2 Harmonic Resonance
while the force vibration is expressed using:

v f0 f0
0
x(t) = sin(ωnt) + x0 − 2 2
cos(ωnt) + 2 cos(ωt) (3.28)
ωn ωn − ω ωn − ω 2
These temporal responses are plotted as (Note that the forcing function uses the axis on the
right):
Figure 3.4: Comparison of the temporal response for a 1-DOF system; expressing how the
forcing function changes the vibrational temporal response of the system.
3.2 Harmonic Resonance

Recall that our solution from before assumed that ωn ̸= ω, however, if ωn = ω then the system will
develop the phenomenon of resonance. Mathematically, this means the amplitude of the vibrations
becomes unbounded. The prior choice of Xcos(ωt) for the particular solution fails as it is also a
solution for a homogeneous equation. Therefore, a new particular solution is needed for the case
where ωn = ω. This new particular solution can be written as:
x p (t) = tXsin(ωt) (3.29)
Substituting this into the EOM of the system in standard form equation (from Boyce and DiPrima
(1997)) and solving for X yields:
f0
x p (t) = tsin(ωt) (3.30)
2ω
thus, the total solution can now be written as:
f0
x(t) = A1 sin(ωt) + A2 cos(ωt) + tsin(ωt) (3.31)
2ω
56
Note that ωn = ω, therefore, the frequencies are all in terms of the driving frequency ω. Again,
evaluating the solution at t = 0 for the initial conditions x0 and v0 yields:
v f0
0
x(t) = sin(ωt) + x0 cos(ωt) + tsin(ωt) (3.32)
ω 2ω
Where the first two terms account for the oscillations while the third term accounts for the con-
tinued increase of the maximum amplitude. The following plot shows the forced response of a
spring-mass system driven harmonically at its natural frequency.
Figure 3.5: Temporal response of a system in resonance showing the enveloped maximum ampli-
tude of displacement.
Example 3.2 Compute solutions for the homogeneous and particular solution separately, then
compute the total response of a spring-mass system with the following values: k = 1000 N/m,
m = 10 kg, subject to a harmonic force of magnitude F0 = 100 N and frequency of 8.162 rad/s,
and initial conditions given by x0 = 0 m and v0 = 0 m/s. Plot the response.
Solution: p
First, make sure that the system is not in resonance. Calculating that ωn = 1000/10 = 10
57
shows us that ωn ̸= ω. Next knowing that f0 = Fo /m = 10 we can find the homogeneous and
particular solutions as:
xh (t) = Asin(ωnt + φ ) (3.33)
x p (t) = Xcos(ωt) (3.34)
also:
x(t) = xh (t) + x p (t) (3.35)
where: s 2
v0
A= + (x0 − X)2 = (3.36)
ωn

−1 ωn (x0 − X)
φ = tan (3.37)
v0
f0
X= (3.38)
2
ωn − ω 2
This leads to the following results.
Figure 3.7: Temporal response for example problem where the envelope of the total solution
is a “beat” with a period of approximately 6 seconds.
Example 3.3 Considering the following system, write the equation of motion and calculate
the response assuming a) that the system is initially at rest, and b) that the system has an initial
displacement of 0.005 m. Use k = 2000 N/m, m = 100 kg, F(t) = 10sin(10t) N.
58
Solution:
The equation of motion is
mẍ + kx = 10sin(10t) (3.39)
or in standard form:
ẍ + ωn2 x = f0 sin(ωt) (3.40)
Note that the forcing function is in terms of sin, not cos as before, so we will have to resolve for
the constants A1 and A2 . Again, setting the particular solution to x p = Xsin(ωt) and solving
for X as before yields:
f0
x(t) = A1 sin(ωnt) + A2 cos(ωnt) + sin(ωt) (3.41)
2
ωn − ω 2
Now we can solve for A1 and A2 by setting the initial conditions x0 and v0 to t = 0. First,
setting t = 0 in the equation for x(t) yields:
A2 = x0 (3.42)
Then, a function for the velocity of the system is obtained:

f0
ẋ(t) = v0 = A1 ωn cos(ωnt) − A2 ωn sin(ωnt) + ω cos(ωt) (3.43)
ωn − ω 2
2
This allows us to obtain:

v0 ω f0
A1 = − · 2 (3.44)
ωn ωn ωn − ω 2
at t = 0. These lead to the full equation for the general solution:
v ω f0 f0
0
x(t) = − · 2 2
sin(ωnt) + x0 cos(ωnt) + 2 sin(ωt) (3.45)
ωn ωn ωn − ω ωn − ω 2
Also, knowing:
k √
r
ωn = = 20 rad/sec = 4.472 rad/sec (3.46)
m
59
and
F0 F0
fo = = = 0.1 N/kg (3.47)
m m
a) using the initial conditions x0 = 0 m and v0 = 0 m/s and the general expression obtained
above: 10 0.1 √ 0.1
x(t) = 0 − √ · 2
sin( 20t) + 0 + sin(10t) (3.48)
20 20 − 10 20 − 102
b) using the initial conditions x0 = 0.005 m and v0 = 0 m/s and the general expression obtained
above:
10 0.1 √ √ 0.1
x(t) = 0 − √ · 2
sin( 20t) + 0.05cos( 20t) + sin(10t) (3.49)
20 20 − 10 20 − 102
Figure 3.9: Temporal response for example problem.
Vibration Case Study 3.2 The Millennium Bridge is a pedestrian suspension bridge in
London over the River Thames. The supporting cables of the bridge are abnormally low
and rest below the deck level, giving a very shallow profile. This was required by London’s
protected Vistas which necessitates a clear line of view from Alexandra Palace to Saint
Paul’s Cathedral; as well as behind Saint Paul’s Cathedral where the bridge sits.
When opened on 10 June 2000, 2,000 pedestrians at 1.5 people per square meter used the
bridge. The bridge started to rock in the lateral direction at frequencies of between 0.5 Hz
and 1.1 Hz with accelerations up to 0.25 gn , this caused people on the bridge to try and
brace themselves by moving their body mass in sync with the bridge’s movement. This bio-
dynamic coupling created a forced lateral vibration in the bridge that would persist when
sufficient people were on the bridge.
To mitigate the vibrations, 37 dampers of 7 different types were installed to control the
lateral modes, with some also controlling vertical and torsional modes. After the installation
60
Open Vibrations 3.3 Harmonic Excitations of Underdamped Systems
of dampers, peak measured accelerations from 0.25 gn to 0.006 gn and no observable bio-
dynamic feedback occurred. In total, this retrofit took almost 2 years and added an extra
£5 million to the initial £18.2 million cost of the bridge.
Figure 3.10: View of Millennium Bridge in London UKa .

a David Martin / Under the Millennium Bridge / CC BY-SA 2.0
3.3 Harmonic Excitations of Underdamped Systems

Consider the system:
Figure 3.11: Damped 1-DOF system with an external force (F(t)) applied, showing: (a) the system
configuration; and (b) the free body diagram
Again, for simplicity, let us consider a harmonic excitation for F(t) such that:
F(t) = F0 cos(ωt) (3.50)
Building the EOM for the above system results in:
mẍ(t) + cẋ(t) + kx(t) = F0 cos(ωt) (3.51)
For convinces we can convert this to the standard form:
ẍ(t) + 2ζ ωn ẋ(t) + ωn2 x(t) = f0 cos(ωt) (3.52)
61
again, where:
F0
f0 = (3.53)
m
Recall that one way to solve such an equation is to obtain the sum of the homogeneous and partic-
ular solutions.
x(t) = xh (t) + x p (t) (3.54)
However, now that we have damping force to consider, our particular solution will have to consider
this damping. Therefore:
x p (t) = Xcos(ωt − φ p ) (3.55)
where φ p represents the phase shift. Note: φ p is represented in other texts as θ , θ p , or even just φ
but we will use φ p throughout the remainder of this text. Again, the phase shift is expected because
of the effect of the damping force. Now, our total equation is:
x(t) = Ae−ζ ωnt sin(ωd t + φ ) + Xcos(ωt − φ p ) (3.56)
We can use the method of undetermined coefficients to obtain X and φ p for the particular solution.
First, considering that we write the particular solution in the equivalent form:
x p (t) = Xcos(ωt − φ p ) = As cos(ωt) + Bs sin(ωt) (3.57)
Taking the derivative of the assumed forms of the particular solution yields:
x p (t) = As cos(ωt) + Bs sin(ωt) (3.58)
ẋ p (t) = −ωAs sin(ωt) + ωBs cos(ωt) (3.59)

ẍ p (t) = −ω 2 As cos(ωt) − ω 2 Bs sin(ωt) (3.60)
Recall that the homogeneous and particular solutions are each solutions on their own, therefore,
the EOM can be used to describe just the particular solution. Substituting x p . ẋ p , and ẍ p for x. ẋ,
and ẍ in the EOM in standard form:
ẍ + 2ζ ωn ẋ + ωn2 x = f0 cos(ωt) (3.61)
yields:
− ω 2 As cos(ωt) − ω 2 Bs sin(ωt) + 2ζ ωn − ωAs sin(ωt) + ωBs cos(ωt) +

(3.62)
ωn2 As cos(ωt) + Bs sin(ωt) = f0 cos(ωt)

and rearranging in terms of sin(ωt) and cos(ωt) yields:
(−ω 2 As + 2ζ ωn ωBs + ωn2 As − f0 )cos(ωt)+ (3.63)
(−ω 2 Bs − 2ζ ωn ωAs + ωn2 Bs )sin(ωt) = 0

From this expression, it is clear that there are two special moments in time where cos(ωt) and
sin(ωt) equal zero. First, considering that t = π/(2ω) results in cos(ωt)=0, sin(ωt)=1 and the
equation simplifies to:
(−2ζ ωn ω)As + (ωn2 − ω 2 )Bs = 0 (3.64)
62
Additionally, at t = 0, sin(ωt)=0 and cos(ωt)=1. Therefore, the equation yields
(ωn2 − ω 2 )As + (2ζ ωn ω)Bs = f0 (3.65)
We can solve two equations for two unknowns. Writing the two linear equations as the singular
matrix equation yields:
2
ωn − ω 2 2ζ ωn ω As

f
2 2 = 0 (3.66)
−2ζ ωn ω ωn − ω Bs 0

A
This can be solved by computing this system of equations for s . This gives us:
Bs
(ωn2 − ω 2 ) f0
As = 2 (3.67)
(ωn − ω 2 )2 + (2ζ ωn ω)2
2ζ ωn ω f0
Bs = (3.68)
2
(ωn − ω 2 )2 + (2ζ ωn ω)2
From trigonometric relationships we can see that,
q
X = A2s + B2s (3.69)

−1 Bs
φ p = tan (3.70)
As
We can now derive values for our particular solution x p :
f0
X=p (3.71)
(ωn2 − ω 2 )2 + (2ζ ωn ω)2

−1 2ζ ωn ω
φ p = tan (3.72)
ωn2 − ω 2
Now we can build a solution for the particular equation (x p ), therefore, the total solution becomes:
x(t) = xh (t) + x p (t) (3.73)

Note for larger values of t, the homogeneous solution approaches zero resulting in the particular
solution becoming the total solution. Therefore, the particular solution is sometimes called the
steady-state response, and the homogeneous solution is called the transient response. Solving for
the constants A and φ using boundary conditions (x0 = 0 and v0 = 0) results in a total solution
expressed as:
x0 − Xcos(φ p )
A= (3.75)
sin(φ )

−

−1 ω d x 0 Xcos(φ p )
φ = tan (3.76)
v0 + x0 − Xcos(φ p ) ζ ωn − ωXsin(φ p )
63
Finally, assembling all the terms:
x0 − Xcos(φ p )
A= (3.78)
sin(φ )

−1 ωd (x0 − Xcos(φ p ))
φ = tan (3.79)
v0 + (x0 − Xcos(φ p ))ζ ωn − ωXsin(φ p )
f0
X=p (3.80)
(ωn2 − ω 2 )2 + (2ζ ωn ω)2

−1 2ζ ωn ω
φ p = tan (3.81)
ωn2 − ω 2
Example 3.4 Consider the damped 1-DOF system below, and plot the total, steady state, and
transient responses for the following system configurations with no initial conditions. For
each configuration, comment on the temporal response and how it differs from the response
of the previous configuration.
a) k = 100 N/m, m = 10 kg, c = 10 kg/s, F0 = 1 N, and ω = 8.162.
b) k = 100 N/m, m = 10 kg, c = 10 kg/s, F0 = 3 N, and ω = 8.162.
c) k = 100 N/m, m = 10 kg, c = 10 kg/s, F0 = 3 N, and ω = 3.162.
Figure 3.12: Damped 1-DOF system with an external force (F(t)) applied, showing: (a) the
system configuration; and (b) the free body diagram
Solution:
The total response for the damped 1-DOF system subjected to an external force is modeled
using equations 3.77 through 3.81 while the transient response consists of the first half of
equation 3.77 and the steady state response consists of the second half of equation 3.77.
Solution a): Therefore, plotting the temporal responses for the configuration yields:
64
Figure 3.13: Temporal responses for a underdamped system with k = 100 N/m, m = 10 kg,
c = 10 kg/s, F0 = 1 N, and ω = 8.162.
Solution b): Configuration b increases the forcing function F0 to 3 N. This results in a similar
response to configuration a but with a linearly scaled amplitude:
c = 10 kg/s, F0 = 3 N, and ω = 8.162.
Solution c): Now, using ω = 3.162 rad/sec we put the system into resonance as ω = ωn .
However, unlike the undamped system, the amplitude of the displacement is not unbounded
as the damper absorbs energy from the system. Therefore, after about 7 seconds the system
65
Open Vibrations 3.4 Frequency Response of Underdamped Systems
enters an equilibrium state where any additional increase in amplitude caused by the system
entering into resonance is canceled out by the damping in the system as demonstrated in the
plot below:
c = 10 kg/s, F0 = 3 N, and ω = 3.162.
3.4 Frequency Response of Underdamped Systems

From equations 3.77 through 3.81 and the figures in example 3.4 we can see that for larger values
of t the transient response dies out while only the steady-state response controls the displacement
of the total response. This is always true if the system has any significant damping. Therefore, it is
often prudent to ignore the transient part and focus only on the steady-state response. Considering
the equation for the particular solution:
x p (t) = Xcos(ωt − φ p ) (3.82)
and knowing the values for X and φ p :
f0
X=p (3.83)
(ωn2 − ω 2 )2 + (2ζ ωn ω)2

−1 2ζ ωn ω
φ p = tan (3.84)
ωn2 − ω 2
We want to find a way to plot the responses of the system only in terms of the system’s natural
and driving frequencies, and its damping. First, we define a frequency ratio as the dimensionless
quantity
ω
r= (3.85)
ωn
66
Another common way to express r is β . Next, Recall that:

0 F
f0 m
X=p =p (3.86)
(ωn2 − ω 2 )2 + (2ζ ωn ω)2 (ωn2 − ω 2 )2 + (2ζ ωn ω)2
If we factor out ωn2 from the denominator and substitute in ωn2 = k/m and r = ω/ωn , we get:
F0 F0
m k
X= q = p (3.87)
ωn2 1 − ( ωωn )2
2
+ (2ζ ωωn )2 (1 − r )2 + (2ζ r)2
2
this becomes:
Xk Xωn2 1
= =p (3.88)
F0 f0 (1 − r2 )2 + (2ζ r)2
in a similar fashion, if we manipulate the equation for φ p we can get φ p in terms of r:

−1 2ζ r
φ p = tan (3.89)
1 − r2
If we solve for a few key values of r we can get the following data points. On the board, we can
solve for a few different frequency responses for a few different damping coefficients.
frequency ratio (r)

0 0.25 0.5 0.75 1 1.25 1.5 1.75 2.0
ζ = 0.1 1.00 1.07 1.32 2.16 5.00 1.62 0.78 0.48 0.33
ζ = 0.25 1.00 1.06 1.27 1.74 2.00 1.19 0.69 0.45 0.32
ζ = 0.5 1.00 1.03 1.11 1.15 1.00 0.73 0.51 0.37 0.28
ζ = 0.7 1.00 1.00 0.97 0.88 0.71 0.54 0.41 0.31 0.24
If we plot the values of the normalized amplitude vs r we obtain figure 3.16 where it can be
seen that the normalized amplitude is a function of damping in the system. However, it should be
noted that damping is only effective around resonance,
√ as below and above resonance, all damping
cases converge on similar values. Note that ζ ≥ 1/ 2 is the changeover point from where the max
normalized displacement is at r = 0 vs around resonance.
67
Figure 3.16: Normalized amplitude response for frequency ratio (r= (ω/ωn )) from 0 to 2 for a
variety of critical damping ratios.
And again, if we plot the values of the phase vs r we get figure 3.17. Note that all systems pass
through 90◦ at resonance. This means that when a system is under resonance, the position of the
system will lag the input force by 90◦ . This phase lag is also called quadrature as the system lags
the input by 90◦ at resonance.
68
Figure 3.17: Phase response for frequency ratio (r) from 0 to 2 for a variety of critical damping
ratios.
note that the dashed blue line is there because the phase values after π/2 need to be adjusted to
obtain a continuous plot. An astute observer would notice that the maximum amplitude is not
at ω = ωn . While resonance is defined as ω = ωn , this does not define the point of maximum
displacement of the steady-state response. Let us solve for the frequency ratio with the maximum
displacement. This will happen when
!
d Xk
=0 (3.90)
dr F0
We can show that: !

1 d
p =0 (3.91)
(1 − r2 )2 + (2ζ r)2 dr
when q
ωp √
rpeak = 1 − 2ζ 2 = , ζ < 1/ 2 (3.92)
ωn
√ √
however, this is only true for underdamped systems in which ζ < 1/ 2. If ζ ≥ 1/ 2 then the
value is imaginary and the peak value is at r = 0. In these cases, the maximum displacement
69
is a function of only ωn . ω p represents the driving frequency that corresponds to the maximum
amplitude ( Xk
F0 ) and is called the peak frequency, and can be calculated as:
q √
ω p = ωn rpeak = ωn 1 − 2ζ 2 , ζ < 1/ 2 (3.93)
Example 3.5 Consider the simple spring-mass system,
Figure 3.18: Damped 1-DOF spring-mass system subjected to an external force F(t).
where ωn = 132 rad/sec and ζ = 0.0085. Calculate the displacements of the steady-state
response for ω=132 and 125 rad/sec. In both cases, use f0 = 10 N/kg.
Solution: From before, we know the solution for the displacement of the particular solution
for ω=132 rad/sec is:
f0 10
X=p = = 0.034 m (3.94)
(ωn2 − ω 2 )2 + (2ζ ωn ω)2 2(0.0085)(132)2
while for ω=125 rad/sec X is:

f0 10
X=p =p = 0.005 m (3.95)
(ωn2 − ω 2 )2 + (2ζ ωn ω)2 (1799)2 + (280.5)2
Therefore, a slight change in the driving frequency (about 5%) results in an 85% change in
the amplitude of the steady-state response.
Example 3.6 The steady-state response for an engineered system must not surpass 1 cm, if
the system can be modeled as the spring and mass system below, what value of c must be
used?
70
Open Vibrations 3.5 Base Excitation
Use k = 2000 N/m, m = 100 kg, F(t) = 20 cos(6.3t) N.

Solution:
The steady state solution is:
x p (t) = Xcos(ωt − φ p ) (3.96)
knowing the amplitude is controlled by X:

f0
X=p (3.97)
(ωn2 − ω 2 )2 + (2ζ ωn ω)2
and recalling from the EOM in standard form that 2ζ ωn = c/m we can obtain:
f0
X=q (3.98)
(ωn2 − ω 2 )2 + ( mc ω)2
rearranging for c gives:

s s
2 2
f02 ωn2 − ω 2 F02 2
ωn2 − ω 2
c=m − = − m (3.99)
ω 2X 2 ω2 ω 2X 2 ω2
Therefore, if we set X = 0.01 m we can solve the above equation to yield c = 55.7 kg/s.
3.5 Base Excitation

Often, loading is not applied directly to the mass, but rather the mass of the system is excited when
the base of the mount that it is attached to is excited. This is called base excitation or sometimes
support motion. Examples of base excitation, or where base excitation is considered, include:
• machines on rubber mounts
• automobiles excited by the road
• building under earthquake loading
• hospital equipment
Consider the following system base excited system
71
Figure 3.20: Damped 1-DOF spring-mass system subjected to a displacement controlled base
excitation showing the FBDs for the equilibrium and displaced positions.
where x is the displacement of the mass and y is the displacement of the base. Note that
we consider positive upward here. The EOM can be constructed the same as before, but now
considering that the displacement of the springs and damper is x − y. In the equilibrium state,
where a positive x is up and the base displaces down:
+↑ ∑ Fx = −kδ − mg = 0 (3.100)
Note that these are both negative because the base displacing down “pulls” the mass down with
a force kδ (i.e. if you hold the mass and let the base “fall”). Conversely, the equation for the
displaced state is:
+↑ ∑ Fx = −k(δ + x − y) − mg − c(ẋ − ẏ) (3.101)
Apply Newton’s second law about the mass (mẍ) of motion to the sum of forces for the displaced
position we get:
+↑ ∑ Fx = mẍ = −kδ − kx + ky − mg − cẋ + cẏ (3.102)
applying the equation −kδ − mg = 0, and rearrange into the EOM yields:
mẍ + cẋ + kx = cẏ + ky (3.103)
Now as before we assume an input for the base excitation. For simplicity we assume:
y(t) = Y sin(ωbt) (3.104)
Taking the derivative of the assumed input yields:
ẏ(t) = Y ωb cos(ωbt) (3.105)
where Y is the amplitude and ωb is the frequency of the base excitation. Adding these terms into
our EOM yields:
mẍ + cẋ + kx = cY ωb cos(ωbt) + kY sin(ωbt) (3.106)
72
We can get this in standard form if we divide by m and apply the equations for the critical damping
ratio and natural frequency:
ẍ + 2ζ ωn ẋ + ωn2 x = 2ζ ωn ωbY cos(ωbt) + ωn2Y sin(ωbt) (3.107)
This equation can be related to a spring-mass-damper system with two harmonic inputs, one cos,
and one sin as shown below:
ẍ + 2ζ ωn ẋ + ωn2 x = Ccos(ωbt) + Dsin(ωbt) (3.108)
where C and D are arbitrary coefficients.
Vibration Case Study 3.3 Earthquakes are a classic and devastating example of base exci-
tation. On August 24th 2016 an earthquake hit Central Italy approximately 75 km (47 mi)
southeast of the city of Perugia. 299 people were killed and the town of Amatrice was heav-
ily damaged. A close look at the town center of Amatrice post-event, as shown in figure 3.21
shows that the town’s bell tower is still standing when the shorter residential buildings have
collapsed. A simplified explanation for the robustness of the bell tower can be found in the
fact the tall and slender bell tower has a natural frequency lower than that of the excitation
force of the earthquake. In comparison, the shorter and stiffer residential structures tend to
have a higher natural frequency that more closely aligns with the excitation frequency of the
earthquake, thereby resulting in these structures being excited closer to resonance.
Figure 3.21: The town center of Amatrice Italy after the August 24th 2016 earthquake that
measured 6.2 on the moment magnitude scale; note that the bell tower (lower natural fre-
quency) is still standing while shorter stiffer structures (higher natural frequency) have suf-
fered extensive damage. a
73
The architectural and cultural importance of bell towers leads to considerable efforts to
protect and preserve these historic structures, in addition to ensuring their safety to protect
the public post-event. During the August 24th earthquake, a team at the University of Perugia
was actively monitoring the bell tower at Basilica di San Pietro in the city of Perugia with
the intention of tracking the tower’s dynamics through time to better understand the tower’s
state; thereby enabling better preservation of the tower. Figure 3.22(a) shows the bell tower,
while figure 3.22(b) shows a sensor placed within the tower. Lastly, figure 3.22(c) shows
Italian researcher Nicola Cavalagli inspecting the data recorded from the accelerometer on
the bell tower on the morning of August 24th . A visual inspection of the monument the day
after the event did not result in the identification of damage. However, by comparing the
vibration signal from before and after the event, researchers were able to detect anomalies
in the tower’s structural behavior through statistical analysis of the vibration data.b This
statistical data is then matched with a finite element model of the system tower to infer
likely locations of damage.
Figure 3.22: Bell tower at Basilica di San Pietro, showing: (a) the bell tower, (b) a sensor
in the bell tower, and (c) data collected during the Central Italy earthquake of August 24,
2016. c
a Image cropped from original photo by Leggi il Firenzepost, CC BY 3.0
<https://creativecommons.org/licenses/by/3.0>, via Wikimedia Commons
b Giordano, P. F., Ubertini, F., Cavalagli, N., Kita, A., & Masciotta, M. G. (2020). Four years of structural
health monitoring of the San Pietro bell tower in Perugia, Italy: two years before the earthquake versus two
years after. International Journal of Masonry Research and Innovation, 5(4), 445-467.
c Austin R.J. Downey, CC BY-SA 3.0 <https://creativecommons.org/licenses/by/3.0>
3.5.1 Displacement Transmissibility Solution for Base Excitation

The steady-state solution is often more important than the transient solution when designing sys-
tems for continuous use. The particular solution for the base excited system annotated in figure
74
3.20 with the EOM presented in equation 3.108 can be expressed as x p (t). To solve for this expres-
sion we will use the linearity of the system and solve for a solution that is the sum of two particular
solutions. Resulting in:
(1) (2)
x p (t) = x p (t) + x p (t) (3.109)
Recall that the steady state solution for a harmonically excited spring-mass-damper can be
expressed as x p (t) = Xcos(ωt − φ p ), as denoted in equation 3.55. For the base excitation problem,
we will convert this expression to x p (t) = Xcos(ωbt − φ1 ). Therefore, for a base excited problem,
the forcing function can be expressed as the sum of particular solutions:
(1) (2)
Ccos(ωbt) + Dsin(ωbt) = x p = x p + x p (3.110)
where we dropped the (t) term from the expression for simplicity in writing and:
(1)
x p = X (1) cos(ωbt − φ1 ) (3.111)
(2)
x p = X (2) sin(ωbt − φ1 ) (3.112)
(1) (2)
Note that x p uses a cos term while x p uses a sin term. Both solutions use φ1 as the damping term
as the phase angle is independent of the excitation amplitude and the sin and cos terms account for
the difference in phase.
(1) (1)
For x p , we use the method of undetermined coefficients to obtain a solution for x p = X (1) cos(ωbt −
φ1 ). This can be as simple as setting 2ζ ωn ωbY equal to f0 from equation 3.71 that defines X for
underdamped systems. Again, 2ζ ωn ωbY comes from the EOM in standard form as presented in
equation 3.107. We can do this because both terms can be considered a “driving force”. This
results in the equation:
(1) 2ζ ωn ωbY
xp = q cos(ωbt − φ1 ) (3.113)
(ωn2 − ωb2 )2 + (2ζ ωn ωb )2
where:
−1 2ζ ωn ωb
φ1 = tan (3.114)
ωn2 − ωb2
(2)
Next, the particular solution associated with x p = X (2) sin(ωbt − φ1 ) can be obtained using the
same method of undetermined coefficients and setting f0 from equation 3.71 to the driving force
(2)
for x p in equation 3.107, ωn2 . This results in:
(2) ωn2Y
xp = q sin(ωbt − φ1 ) (3.115)
2 2 2
(ωn − ωb ) + (2ζ ωn ωb )2
As both equation 3.113 and 3.115 have the same argument (ωbt − φ1 ), these can be added as:
(1) (2)
xp = xp + xp (3.116)
75
to obtain: s
ωn2 + (2ζ ωb )2
x p = ωnY cos(ωbt − φ1 − φ2 ) (3.117)
(ωn2 − ωb2 )2 + (2ζ ωn ωb )2
and:
−1 ωn
φ2 = tan (3.118)
2ζ ωb
where φ2 is added to account for the cos and sin terms being combined. Again, the (t) has been
dropped for simplicity.
As before, if we want to investigate how a frequency input will affect the response (frequency
response) we can substitute substitute
ωb
r= (3.119)
ωn
into the temporal response to obtain:
s
1 + (2ζ r)2
X =Y (3.120)
(1 − r2 )2 + (2ζ r)2
Next, if we divide by Y we can obtain a normalized expression for the displacement:

s
X 1 + (2ζ r)2
= (3.121)
Y (1 − r2 )2 + (2ζ r)2
Plotting this for several critical damping ratios:
Figure 3.23: Displacement transmissibility for an underdamped 1-DOF system.
76
Around resonance, the maximum amount √ of displacement is transmitted to the mass. Addi-
tionally, the above plot shows that at r = 2 the displacement transmissibility X/Y √ is 1. Note
the “flip” where overdamped systems have a greater response to excitations after r = 2 than do
underdamped systems.
Example 3.7 A very common example of base motion is the SDOF model of a vehicle wheel
driving over a “rough” road as shown below. For this, let’s consider a generic modern sports
sedan that we can diagram as below
Figure 3.24: A 1-DOF “car” traveling over an uneven road.
where k = 300,000 N/m, m = 1600 kg, c = 15,000 kg/s, the period of road roughness = 3 m,
and the height of road roughness = 0.01 m. What is the deflection experience by the car at v =
50 km/h?
Solution:
The road is applying a base excitation that can be approximated as
Y = 0.005 m (3.122)
! !
1000 m 1 hours
v m/s = 50 km/hr = 13.888 m/sec (3.123)
1km 3600 s
! ! !
13.88 m 1 cycle 2π rad
ωb = = rad/s = 29.08 rad/s (3.124)
s 3m cycle
Therefore, the sinusoidal for the base excitation is then:
y(t) = (0.005)sin(29.08t) (3.125)
77
Next, we can calculate the natural frequency:

r r
k 300, 000
ωn = = = 13.69 rad/s (3.126)
m 1600
Therefore:
ωb
r= = 2.124 (3.127)
ωn
and:
c 15, 000
ζ= √ = √ = 0.342 (3.128)
2 km 2 1600 · 300, 000
Then it can be found that the maximum deflection of the car is:
s s
1 + (2ζ r)2 1 + (2 · 0.3423 · 2.124)2
X =Y = Y
(1 − r2 )2 + (2ζ r)2 (1 − 2.1242 )2 + (2 · 0.3423 · 2.124)2 (3.129)
= 0.0023 m
3.5.2 Force Transmissibility Solution for Base Excitation

For some systems, such as those with weak connections, the force transmitted to the mass is more
important than the displacement of the mass. The force transmitted to the mass is the sum of the
forces applied by the spring and damper. From the FBD above,
F(t) = k(x − y) + c(ẋ − ẏ) (3.130)
where this force is counteracted by the inertial force of the mass:
F(t) = −mẍ(t) (3.131)
Only considering the steady state we found that

s
ωn2 + (2ζ ωb )2
x p (t) = ωnY cos(ωbt − φ1 − φ2 ) (3.132)
(ωn2 − ωb2 )2 + (2ζ ωn ωb )2
if we differentiate this twice, to obtain ẍ(t) and combine this with F(t) = −mẍ(t) we get:
s
ωn2 + (2ζ ωb )2
F(t) = mωb2 ωnY cos(ωbt − φ1 − φ2 ) (3.133)
(ωn2 − ωb2 )2 + (2ζ ωn ωb )2
where the negative sign F(t) = −mẍ(t) as the force transmitted to the mass is both positive and
negative and we are solving for the amplitude of the transmitted force. Again applying:
ωb
r= (3.134)
ωn
this becomes:
F(t) = FT cos(ωbt − φ1 − φ2 ) (3.135)
78
where FT is the magnitude of the transmitted force and is

s
1 + (2ζ r)2
FT = kY r2 (3.136)
(1 − r2 )2 + (2ζ r)2
Again, this can be converted to force transmissibility to provide a normalized response such that:
s
FT 1 + (2ζ r)2
= r2 (3.137)
kY (1 − r2 )2 + (2ζ r)2
Plotting this for several critical damping ratios:
Figure 3.25: Force transmissibility for an underdamped 1-DOF system.

√ √ FT
Again, note the key location r = 2. At r = 2 the force transmitted to the system is 2 kY .
However,
√ also note that the normalized force does not necessarily fall off for r values greater than
r = 2.
Vibration Case Study 3.4 The Convair F2Y Sea Dart was a prototype seaplane fighter de-
veloped by the United States Navy in the early 1950s to enable sea-based jet fighters. One
key technical issue with the aircraft’s development was the violent forces induced into the
plane when the hydro-skis contacted the uneven surfaces of the water. Furthermore, adding
damping to the skies proved to be changed as the damping required changed significantly as
a function of the hydro-skis contact with the water. Significant work went into the skies and
shock-absorbing struts, which helped to improve the situation but it was never fully repaired.
79
Figure 3.26: The Convair F2Y Sea Dart, showing: a) XF2Y-1 Sea Dart (BuNo 135762) dur-
ing landing. This airframe disintegrated in mid-air over San Diego Bay, California (USA)
during a demonstration flight on November 4th, 1954 killing test pilot Charles E. Richbourg
after the airframe limitations were exceeda , and b) the hydro-skis undergoing extensive test-
ing on a pantograph mounted on a speed boat to study the forces transmitted to the airframe
from the hydro-skisb .
a Public Domain U.S. Navy National Museum of Naval Aviation photo No. 1996.253.7213.010
a Image from “The Impossible Takes Longer”, a film by Convair about Sea Dart development. The copy-
right of the image is unknown but may be held by the successor entities of Convair. It is believed that the use
of this image qualifies as fair use under the copyright law of the United States.
Example 3.8 For the system given below and excited at the base, should the system be excited
above or below the natural frequency if the transmitted force is the design limitation? Consider
the under-damped case with ζ = 0.1, and the over-damped case with ζ = 2 conditions.
Figure 3.27: Force transmissibility for an underdamped 1-DOF system.
80
Solution:
We can plot the transmissibility of both the force and displacement onto one plot. For
ζ = 0.1
Figure 3.28: Force and displacement transmissibility for the considered base excited system
with ζ = 0.1.
it is clear that to minimize the force, the system should be driven with a frequency below the
natural frequency. Next for ζ = 2:
Figure 3.29: Force and displacement transmissibility for the considered base excited system
with ζ = 2.
81
it can be seen that the same rationale applies. Therefore, for both ζ = 0.1 and ζ = 2 the
system should be excited below the natural frequency.
Example 3.9 A single-story building is subjected to a harmonic ground motion, ÿ(t) = Acos(ωbt).
a) Find the steady-state solution for the structure. b) If a damper was added between the base
and the floor, and r = 2, what would be the ideal critical damping coefficient to ensure the
safety of the building? (Think of safety as limiting displacement and transmitted force.)
Figure 3.30: A 1-DOF latterly excited system that represents a 1-story building.
Solution (a):
For simplicity, we can rearrange the system as what follows:
Figure 3.31: A base excited 1-DOF spring-mass system.
solving for the EOM yields:

mẍ + kx = ky (3.138)
Notice that this is the same as the EOM for a damped 1-DOF system if c = 0.
mẍ + cẋ + kx = +cẏ + ky → mẍ + kx = ky (3.139)
Therefore, we can use the solution:

s
ωn2 + (2ζ ωb )2
x p (t) = ωnY cos(ωbt − φ1 − φ2 ) (3.140)
(ωn2 − ωb2 )2 + (2ζ ωn ωb )2
82
Open Vibrations 3.6 Numerical Methods
where:
2ζ ωn ωb
−1
φ1 = tan (3.141)
ωn2 − ωb2

−1 ωn
φ2 = tan (3.142)
2ζ ωb
Now we have, or can easily get, values for ωn , ωb , and ζ . However, we do not have an
expression for Y . We can extract the displacement (and therefore the Y ) from the acceleration
as:
ÿ(t) = Acos(ωt) (3.143)
A
ẏ(t) = sin(ωt) +C1 (3.144)
ω
A
y(t) = − cos(ωt) +C1t +C2 (3.145)
ω2
Resulting in
A
Y =− (3.146)
ω2
Solution (b): From the plots we solved for before, we can see that we want a critical damping
coefficient that is as low as possible. This means any damping added to the system will
decrease its safety. This may seem counter-intuitive, but this is because we are attempting to
drive the structure at a frequency higher than its natural frequency, something that does not
commonalty happen. Typically excitations for a structure are well below its natural frequency.
3.6 Numerical Methods

Numerical methods can be used to solve the response of a system subjected to forced vibrations.
While not the most computationally efficient method, the EOM is an ODE that can be solved
directly while considering the initial directions to obtain the response of the system.
Example 3.10 Directly Solving the Ordinary Differential Equation
Using the EOM for the system in figure 3.32 solve for its temporal response by directly
solving the ODE for a system initially at rest with m = 1 kg, c = 0.2, k = 2.0, and F(t) =
1/2 sin(2πt).
Solution:
83
In MATLAB, ode45 is a versatile ODE solver and is one of the first solvers you should
try for most problems. The solver is setup as [t,y] = ode45(odefun,tspan,y0), where
tspan = [t0 tf], integrates the system of differential equations y’=f(t,y) from t0 to tf
with initial conditions y0. Each row in the solution array y corresponds to a value returned in
column vector t. The ODE is re-organized as
ẍ = ( ft − cẋ − kx)/m (3.147)
for the ode45 solver. Listing 4 reports the code needed to solve the time response of the
system shown in figure 3.32.
Listing 1: MATLAB code for solving the EOM through time.
% Time span for simulation
tspan = [0 , 10]; % Start time and end time
% Initial conditions [x , x ']

in it ia l_c on di tio ns = [0 , 0];
% Use ode45 to solve the system of ODEs

[t , y ] = ode45 (@ equation_of_motion , tspan , ini ti al _co nd it ion s ) ;
% Extract displacement and velocity

x = y (: , 1) ;
x_dot = y (: , 2) ;
The code in listing 4 needs to be combined with the functions in listing 5 and plotting code
to obtain the results shown in figure 3.33.
Listing 2: Functions called from the main code in listing 4.
% Equation of motion for the system
function dydt = e qua ti on _of _m ot ion (t , y )
% Mass , damping coefficient , and spring constant
m = 1.0; % Mass
c = 0.2; % Damping coefficient
k = 2.0; % Spring constant
% Unpack the state variables

x = y (1) ;
x_dot = y (2) ;
% Define the force excitation function f ( t )

f_t = f o r c e _ e x c i t a t i o n _ f u n c t i o n ( t ) ;
% Equation of motion
x_dotdot = ( f_t - c * x_dot - k * x ) / m ;
% Pack the derivatives into the output vector dydt

dydt = [ x_dot ; x_dotdot ];
end
% Force excitation function f ( t ) for a sinusoidal force excitation
84
function f_t = f o r c e _ e x c i t a t i o n _ f u n c t i o n ( t )
f_t = 0.5 * sin (2 * pi * t ) ;
end
Figure 3.33: Displacement response of the 1-DOF system in in figure 3.32.
85
Open Vibrations
4 Transfer Function Approach

Thus far, this text has only considered forced vibrations for 1-DOF systems excited with forcing
functions that can be easily expressed using either sin or cos examples. Therefore, the previously
developed solutions are only acceptable for systems with known and simple excitations. This
chapter will introduce the concept of transfer functions for solving vibration-related problems.
The transfer function, in particular the Laplace transfer function, is an important tool in the study
of vibrations as it allows the practitioner to solve for the temporal response of a system for a variety
of inputs using a single approach. Examples of force excitation that can be calculated include using
this method include:
• sinusoidal
• base excitation
• impulse
• arbitrary input
4.1 Transfer Function Method (Generic)

Consider the following system
Figure 4.1: Generic system H subjected to an input F and its corresponding output X.
where F is the input, H is the system, and X is the output from the system. This formulation is
called the transfer-function approach and is commonly used for the formulation and solution of
dynamic problems in the control literature. It can also be used for solving various forced-vibration
problems including those from complex or stochastic inputs.
Review 4.1 Laplace transforms, or more broadly integral transforms, are a procedure for
integrating the time (t) dependence of a function into a function of position or space (s). By
transforming the whole differential equation from the time domain into a lower-order func-
tion of space the problem becomes easier to solve as the function can often be manipulated
algebraically.
86
Open Vibrations 4.1 Transfer Function Method (Generic)
Figure 4.2: Portrait of Pierre-Simon Laplace by Johann Ernst Heinsius (1775).a
The Laplace transform is named after mathematician and astronomer Pierre-Simon Laplace
(23 March 1749 - 5 March 1827 ). Pierre-Simon Laplace was one of the greatest scientists
of all time and is often considered the French Newton. He taught Napoleon at the École
Militaire in 1784, became a count of the empire in 1806, and a marquis in 1817 after the
restoration of the monarchy. He is credited with advancements in engineering, mathematics,
statistics, physics, astronomy, and philosophy; however, maybe his greatest achievement is
not only surviving but benefiting from the change from the Ancien Régime → Bonaparte →
Bourbon Restoration.
Of interest to this class is the Laplace transform (L [ ]) of the function f (t), expressed
as L [ f (t)]. Here, a Laplace transform is used as a method of solving the differential equa-
tions of motion by reducing the computation needed to that of integration and algebraic
manipulation.
The definition of the Laplace transform of the function f (t) is:
Z ∞
L [ f (t)] = F(s) = f (t)e−st dt (4.1)
0
where s represents a variable in the complex plane (also called the s-plane) and f (t) = 0 for
all values of t < 0. Here, the s is a complex value. Lastly, the term F(s) is a generic term
that represents the input to a system. As this class needs the derivatives of the base function,
we will calculate these next:
∞ d[ f (t)]
Z Z ∞
L f˙(t) = f˙(t)e−st dt = e−st

dt (4.2)
0 0 dt
87
Open Vibrations 4.1 Transfer Function Method (Generic)
integration by parts yields:

∞
Z ∞
L f˙(t) = e−st f (t) e−st f (t)dt

+s (4.3)
0 0
Astutely, it can be noticed that the second term s 0∞ e−st f (t)dt is the input to the system
R
F(s). With a little rearranging, this becomes:
L f˙(t) = sF(s) − f (0)

(4.4)
Taking the derivative of again yields:
L f¨(t) = s2 F(s) − s f (0) − f˙(0)

(4.5)
A few key points of the Laplace transforms are:
• The domain of the problem changes from the real number line (t) to the complex plane
(s-plane).
• The integration of the Laplace transform changes differentiation into multiplication.
• The transform procedure is linear. Therefore, the transform of the linear combination
of two transforms is the same as the linear transformation of these functions.
• To move from the time domain to the complex number plane we typically use tables
of pre-solved integral.
• The function x(t) can be obtained by taking the inverse Laplace transform defined as
x(t) = L [X(s)]−1
The Laplace transform can be calculated in symbolic form. In particular interest to this
text is the Laplace form of the system input F(s) and output X(s). To expand the symbolic
form of the Laplace transform for the system inputs are and for system outputs:
L [ f (t)] = F(s) (4.6)
L f˙(t) = sF(s) − f (0)

(4.7)
L f¨(t) = s2 F(s) − s f (0) − f˙(0)

(4.8)
here, f (0) and f˙(0) are the initial values of the function f (t). Furthermore, the system
outputs are:
L [x(t)] = X(s) (4.9)
88
Open Vibrations 4.2 Transfer Function Method for Solving Vibrating Systems
L [ẋ(t)] = sX(s) − x(0) (4.10)

L [ẍ(t)] = s2 X(s) − sx(0) − ẋ(0) (4.11)
here, x(0) and ẋ(0) are the initial values of the function x(t).
a Johann Ernst Heinsius, CC BY-SA 4.0 <https://creativecommons.org/licenses/by-sa/4.0>, via Wikimedia
Commons
4.2 Transfer Function Method for Solving Vibrating Systems

As mentioned in the introduction to this chapter, a variety of systems can be solved for using the
transfer function method. The procedure for using the Laplace transform to solve equations of
motion expressed as an inhomogeneous ordinary differential equation is:
1. Take the Laplace transform of both sides of the EOM while treating the time derivatives
symbolically.
2. Solve for X(s) in the obtained equation.
3. Apply the inverse transform x(t) = L [X(s)]−1
4.2.1 Free Vibration for Undamped Systems

Consider the undamped single-DOF system:
Figure 4.3: A spring-mass model of a 1-DOF system.
The EOM for this system is a homogeneous differential equation because the right-hand side is
equal to zero:
mẍ(t) + kx(t) = 0 (4.12)
Here we will leave the “(t)” for clarity to differentiate the time domain solution from Laplace
solution “(s)” in the s-plane, as discussed in review 4.1. The EOM can be rewritten in standard
form as:
ẍ(t) + ωn2 x(t) = 0 (4.13)
where the initial conditions at t = 0 are x(0) = x0 and ẋ(0) = v0 . Taking the Laplace transforms, in
symbplic form using equations 4.9 - 4.11, of both sides of the EOM yields:
2
s X(s) − sx0 − v0 + ωn2 X(s) = 0

(4.14)
89
using equations 4.9 and 4.11 from section 4.1. Solving for the output of the system X(s) yields:
sx0 + v0
X(s) = (4.15)
s2 + ωn2
We can expand this form of X(s) to obtain equations listed in our Laplace Transform table:
sx0 v0 ωn
X(s) = + · (4.16)
s + ωn2
2 s + ωn2
2 ωn
This becomes:
s v0 ωn
X(s) = x0 2 2
+ · 2 (4.17)
s + ωn ωn s + ωn2
Next, using the inverse Laplace transform x(t) = L [X(s)]−1 ] and the two following Laplace
transforms (#5 and #6):
s
f (t) is cos(ωt) when F(s) is 2 (4.18)
s + ω2
ω
f (t) is sin(ωt) when F(s) is 2 (4.19)
s + ω2
Therefore, we can obtain the solution for the system output X(s) as:
v0
ωn
The same procedure can be used to calculate the under-damped and forced responses. How-
ever, when calculating these responses the algebraic solution for X(s), s often contains quotients
of polynomials. These Polynomial ratios may not be found in simple Laplace tables and must
be solved using the method of partial fractions. An example of this procedure can be found in
Appendix B of Inman.
4.2.2 Impulse Response Function

Shock loads on mechanical systems represent a very common source of vibration. These short-
duration forces are also called an impulse. An impulse excitation is defined as a force that is
applied for a very short, or infinitesimal, length of time. An impulse is a nonperiodic force that
is represented by the symbol δ . The response of a system to an impulse load is the same as the
system’s free response provided that the correct initial conditions are applied. This is illustrated
in the following where the applied force F(t) is impulsive (i.e., large magnitude over a very short
time).
90
Figure 4.4: An impulse function with the impulse at t = 0.
The impulse response function can be solved for analytically, however, we will solve it using
the transfer function approach. Here we will consider the under-damped spring-mass system. First,
assume that the system is at rest (no initial conditions). Next, we write the EOM as:
mẍ + cẋ + kx = δ (t) (4.21)
Taking the Laplace transform of both sides of the equation yields
m s2 X(s) − sx(0) − ẋ(0) + c sX(s) − x(0) + kX(s) = 1

(4.22)
note that the L [δ ] = 1 per #1 in the transform table. However, if we assume zero initial conditions
(a system at rest when the impulse happens), the equation simplifies too.
ms2 X(s) + csX(s) + kX(s) = 1 (4.23)
or
(ms2 + cs + k)X(s) = 1 (4.24)
Solving this equation for X(s):
1 1
X(s) = · 2 (4.25)
m s + 2ζ ωn s + ωn2
Again, the mass is extracted to develop a formulation that can be found in the Laplace tables.
Setting the constraint that ζ < 1 and consulting #10 in the table for Laplace transforms results in:
1 −ζ ωnt
x(t) = e sin(ωd t) (4.26)
mωd
where this is the general solution for a damped system subjected to an impulse loading function.
For the undamped case, a solution can be obtained by setting ζ = 0. This results in the following
form for the undamped case:
1
x(t) = sin(ωnt) (4.27)
mωn
91
Below is a typical response for both an undamped and underdamped 1-DOF system subject to an
impulse response at t = 0 seconds.
Figure 4.5: Temporal responses from underdamped and undamped 1-DOF systems to an impulse
response function.
4.2.3 Unit Step function

Now consider a unit step function, denoted with a capital Greek Phi Φ:
Figure 4.6: A Step function.
A step function is a common loading situation and can represent the dropping of a load into a
truck, a car going over a curve, or a motor starting up.
The Laplace transform of the function, for a unit step function Φ, is:
e−∞ e−0 1
Z ∞
L [Φ(t)] = e−st dt = − + =
0 s s s
92
This also lines up with Laplace Transform #3 from the Laplace table. This would be expected as
Φ is used to represent the unit step function (i.e. a step function with a displacement of 1). As we
consider linear systems in this class, we can scale the magnitude of the response by the magnitude
of the impulse after the transform is performed.
4.2.4 Undamped Spring-mass System

For a spring-mass system subjected to a unit step, assuming both initial conditions are zero, the
solution can be obtained using the transform method. First, the EOM is
mẍ(t) + kx(t) = Φ(t) (4.28)
Taking the Laplace transform of both sides and assuming zero initial conditions yields:
1
ms2 X(s) + kX(s) = (4.29)
s
Next, this equation is solved for X(s) as:
1
X(s) = (4.30)
s(ms2 + k)
This can be rearranged as:
1 1
X(s) = · 2 (4.31)
m s(s + ωn2 )
where m1 will pass through the Laplace function. Therefore, taking the inverse Laplace transform
using #9 of the provided Laplace transforms yields:
1 1 1
x(t) = · 2 1 − cos(ωnt) = 1 − cos(ωnt) (4.32)
m ωn k
4.2.5 Under Damped Spring-mass System

For a spring-mass-damper system subjected to a unit step, assuming both initial conditions are
zero, the solution can be obtained using the transform method. First, the EOM is:
mẍ(t) + cẋ(t) + kx(t) = Φ(t) (4.33)
Converting to the standard form results in:
1
ẍ(t) + 2ζ ωn ẋ(t) + ωn2 x(t) = · Φ(t) (4.34)
m
taking the Laplace transform of both sides and assuming zero initial conditions yields:
1 1
s2 X(s) + 2ζ ωn sX(s) + ωn2 X(s) = · (4.35)
m s
Next, this equation is solved for X(s) as:
1 1 1
X(s) = · · (4.36)
s2 + 2ζ ωn s + ωn2 m s
93
ωn2
multiplying the right-hand-side of this equation by ωn2
results in:
1 ωn2
X(s) = · (4.37)
mωn2 s(s2 + 2ζ ωn s + ωn2 )
1
Again, the mω 2 will pass through the Laplace function. Therefore, taking the inverse Laplace
n
transform using #11 on the Laplace transform sheet yields:
1 ωn −ζ ωnt
x(t) = · 1 − e sin(ωd t + φ ) , where φ = cos−1 (ζ ), where ζ < 1 (4.38)
mωn2 ωd
After obtaining equations for the undamped and under-damped cases, the responses for the unit
step, solved with the transform method, can be plotted as:
Figure 4.7: Temporal responses from underdamped and undamped 1-DOF systems subjected to an
impulse response function.
Note that the system will settle out around F0 /k where F0 Φ is a scaling factor for the step
loading.
Example 4.1 A load of dirt m is dumped into the back of a dump truck. Assuming the dirt not
move, the bed of the truck can be modeled as a spring-mass-damper system (of values k, m,
and c, respectively). Next, the load of dirt is modeled as a force F(t) = mg that is applied to the
spring-mass-damper system, as illustrated in the following figure. The simplification of the
system allows for the analysis of the truck bed’s vibrations as a simple spring-mass-damper
system. First, assume that the truck’s damper is broken, how does the maximum dynamic
displacement compare to the static displacement? What would happen to the maximum dis-
placement if the damper was repaired on the truck?
94
Figure 4.8: Dump truck being loaded with dirt showing (a) dirt going into the truck beda ; and
(b) the single-degree-of-freedom vibration model.
Solution:
First, setting the load applied to the truck as 1 unit, it can be seen that this is a unit step
loading condition and a broken damper represents an undamped case. To obtain a rough idea
about the nature of static and dynamic displacement, the undamped displacement is
1 mg
1 − cos(ωnt) =
x(t) = 1 − cos(ωnt) (4.39)
k k
This equation has a maximum amplitude when the cos(ωt ) = −1, resulting in:
mg
x(t) = 1 − (−1) (4.40)
k
This can be rearranged for the maximum displacement value xmax as:
mg
xmax = 2 (4.41)
k
Note that the dynamic displacement of the truck bed is twice that of the static displace-
ment. Therefore, if the truck manufacturer designed the truck to only take the static load of
rosters (i.e., if the dirt were placed gently into the truck bed), the frame of the truck would
be damaged when the dirt were loaded into the truck dynamically. From this, it can be un-
derstood that it is important to consider the dynamic responses of a system during the design
phase.
a SGT Marvin Lynchard, A dump truck is filled with dirt, by members of the 459th Civil Engineering Flight,
for use in repairing a damaged runway during Exercise Prime Beef ’82, Public Domain, via picryl.com
Vibration Case Study 4.1 A smokestack or chimney stack is used to exhaust combustion
gases into the outside air. The design of large stacks poses considerable challenges from a
structural dynamics perspective. As high winds pass over the tower creating a combination
of oscillating wind currents and complex vortex shedding that load the tower with a variety
of wind-induced frequencies. This are called vortex-induced vibrations a . This wide band-
95
Open Vibrations 4.3 System Response to Arbitrary Inputs
width of excitation results in a stack that is loaded near its resonant frequency. To mitigate
this, stack designers have designed stacks with changing diameters to ensure that different
parts of the stack have different resonant frequencies. Also, wind bands in the forms of
protruding bricks or helical strakes are added to the stacks to prevent vortex shedding which
reduces the loading on the tower.
Sadly, engineers understanding of vortex shedding and structural dynamics lagged be-
hind the development of these structures; leading to multiple wind-induced collapses of
smokestacks during the industrial revolution. The use of the transfer function approach gives
the practitioner the ability to easily model the complex response of smokestack excited with
a wide bandwidth of excitation.
Figure 4.9: Methods used to reduce vortex-induced vibrations in smokestacks, showing: (a)
helical steel strakes on a chimney stack b , and; (b) Tapered chimney with wind bands at a
Weaving Factory in the UK c .
a Wang, Lei, and Xing-yan Fan. “Failure cases of high chimneys: A review.” Engineering failure analysis
105 (2019): 1107-1117.
b tromBer, CC BY-SA 3.0 <https://creativecommons.org/licenses/by-sa/3.0>, via Wikimedia Commons
c P Flannagan / Large Chimney Stack of the disused Weaving Factory, Donaghcloney. / CC BY-SA 2.0
4.3 System Response to Arbitrary Inputs

The time-domain response of a system to an arbitrary input force in time can be calculated using a
series of impulses as shown in figure 4.10. This method allows the practitioner to easily calculate
the response of an arbitrary input to a system using a single expression executed in a “for loop”.
This type of analysis is often more efficient in terms of programming than more direct methods
such as the transfer functions shown in this text.
96
Figure 4.10: Generalized response showing that any signal can be represented as a series of impulse
signals.
To solve for a generalized response to arbitrary inputs, consider that the forces applied in a loop
can be expressed as an open-ended summation such that
δ x(t) = ∑ F(τ)g(t − τ)∆τ (4.42)
By letting δt → 0 the summation can be transferred into the continuous integration
Z t
δ x(t) = F(τ)g(t − τ)dτ (4.43)
0
Knowing the solution to an impulse load, as expressed in equation 4.26, and inserting this into the
prior equation at g(t) we can show that for a system of impulses, the total response is
Z t
1
F(t)e−ζ ωn (t−τ) sin ωd (t − τ) dτ

x(t) = (4.44)
mωd 0
Again, this represents the total system response (without initial conditions) for an arbitrary exci-
tation F(t). These equations are called the Duhamel integral. While in many cases, the form of
the function F(t) allows for explicit integration of equations 4.43 and 4.44. However, numerical
evaluation is always possible and many times easier given the simplicity of coding.
Example 4.2 In testing, a hammer is used to excite a 1-DOF system with an impact (i.e. im-
pulse), however, the hammer impacts the system twice by ascendant (a double hit). The first
impact has a force of 0.2 N, while the second has a force of 0.1 N and happens 0.1 seconds
after the first impact. Plot the response for the double impact. The system has the parameters
m = 1 kg, c = 0.5 kg/s, k = 4 N/m.
Solution:
First, we can define the forcing function as:
F(t) = 0.2δ (t) + 0.1δ (t − τ) (4.45)

where τ is the offset between the first and second impacts. Next, considering that the unit
impulse has a magnitude of 1 we can obtain solutions for the first impact by first writing its
EOM:
97
mẍ(t) + cẋ(t) + kx(t) = 0.2δ (t) (4.46)

Taking the Laplace transform of both sides of the equation yields
m s2 X(s) − sx(0) − ẋ(0) + c sX(s) − x(0) + kX(s) = 0.2

(4.47)
However, assuming zero initial conditions, the equation simplifies to.
(ms2 + cs + k)X(s) = 0.2 (4.48)
Solving this equation for X(s):
0.2 1
X(s) = · 2 (4.49)
m s + 2ζ ωn s + ωn2
Again, consulting #10 in the table for Laplace transforms results in:
0.2 −ζ ωnt
x1 (t) = e sin(ωd t) (4.50)
mωd
where this is the general solution for a damped system subjected to an impulse loading func-
tion. The second impact can now be solved for using the same method. However, now the
time (t) must be offset by (τ) to allow the impact to still be located at t = 0 in terms of the
second impact. This results in:
0.2 −ζ ωnt
x1 (t) = e sin(ωd t) (4.51)
mωd
0.1 −ζ ωn (t−τ)
x2 (t) = e sin ωd (t − τ) (4.52)
mωd
Next, using the knowledge that the systems are linear and that the Laplace transform of a linear
combination of two transforms is the same as the linear transformation of these functions we
can build the piecewise function:
(
0.2 −ζ ωn t
mωd e sin(ωd t) if t < τ
x(t) = 0.2 −ζ ω t
n sin(ω t) + 0.1 e−ζ ωn (t−τ) sin ω (t − τ)

mωd e d mωd d if τ ≤ t
For the mass, damping, and stiffness values given above can be plotted as:
98
Example 4.3 Consider the base exciton as shown below subjected to an arbitrary base exci-
tation. Derive an equation (Duhamel integral) for its displacement (z), when the displacement
is expressed at the relative displacement of the mass such that z = x − y.
Figure 4.11: A base excited 1-DOF spring-mass-damper system.
Defining the EOM for the system results in
mẍ + cẋ + kx = cẏ + ky (4.53)
Using z = x − y, this can be simplified to
mz̈ + cż + kz = −mÿ (4.54)
Given that we can replace −mÿ with F, this is the same equation as
mz̈ + cż + kz = F (4.55)
meaning that the solutions for an arbitrary force-excited problem can transfer to a base-excited
problem if we consider the relative displacement of the mass. Therefore, we can write the
equation for the relative displacement of the mass as
Z t
1
ÿ(t)e−ζ ωn (t−τ) sin ωd (t − τ) dτ

z(t) = − (4.56)
mωd 0
99
Open Vibrations 4.4 Transfer Function for Response to Random Inputs
4.4 Transfer Function for Response to Random Inputs

Consider the following system
Figure 4.12: Generic block diagram of a system H(s) subjected to an input F(s) and its corre-
sponding output X(s) where the (s) denotes that the considered system is in the s-plane.
where F(s) is the input, H(s) is the system, and X(s) is the output from the system. This for-
mulation is called the transfer-function approach and is commonly used for the formulation and
solution of dynamic problems in the control literature. It can also be used for solving various
forced-vibration problems including those from complex or stochastic inputs.
4.4.1 Defining the Transfer Function H(s)

Again, consider the generic system represented in figure 4.12. For this system representation,
F(s) is the Laplace of the transform of the driving force, and H(s) is the Laplace transform of the
response of the system h(t).
We need to define the transfer function H(s) for a generic system. To do this let us show the
reasoning behind the transfer function. Here we will show that the output of any system (x(t)) can
be related to the input of the system ( f (t)) through a series of polynomial coefficients (a and b).
Consider the general nth -order linear, time-invariant differential equation that governs the behavior
of the dynamic system.
d n x(t) d n−1 x(t) d m f (t) d m−1 f (t)

an + a n−1 + ... + a0 x(t) = bm + bm−1 + ... + b0 f (t) (4.57)
dt n dt n−1 dt m dt m−1
where x(t) is the output and f (t) is the input. Note that this is similar to the formulation we have
had before for the EOM. Taking the Laplace transform of both sides of the above equation yields
an sn X(s) + an−1 sn−1 X(s) + ... + a0 X(s) + initial condition for x(t) = (4.58)
bm sm F(s) + bm−1 sm−1 F(s) + ... + b0 F(s) + initial condition for f (t)
It can be seen that this equation is a purely algebraic expression. If we assume the initial conditions
to be zero, the equation reduces to the following:
(an sn + an−1 sn−1 + ... + a0 )X(s) = (bm sm + bm−1 sm−1 + ... + b0 )F(s) (4.59)
if werearrange equation 4.59 to solve for the relationship between the Laplace variables X(s) and
F(s) and the algebraic expressions we get:
X(s) bm sm + bm−1 sm−1 + ... + b0

= (4.60)
F(s) an sn + an−1 sn−1 + ... + a0
100
this shows that the ratio of the input algebraic expressions over the output algebraic expressions
is equal to the ratio of the output Laplace variable over the input Laplace variable. This shows
that we can relate the Laplace variables to the algebraic expressions. Therefore, we can define the
transfer function H(s) as:
X(s)
H(s) = (4.61)
F(s)
In a more formal term, the transfer function is defined as: “The ratio of the Laplace transforms of
the output or response function to the Laplace transform of the input or forcing function assuming
zero initial conditions”.
Equation 4.61 can be rearranged to show that the output of the system X(s), can be obtained if
we know the input F(s) and the transfer function H(s):
X(s) = H(s)F(s) (4.62)
4.4.2 Transfer Function Method (Steady-state Solution)

Considering the forced system:
Figure 4.13: A spring-dashpot-mass model of a 1-DOF system with external excitation.
that can be expressed as the equation of motion

mẍ(t) + cẋ(t) + kx(t) = F0 cos(ωt) (4.63)
Here F0 cos(ωt), is used at the input but any input will develop the same transfer function as the
transfer function is bounded to the system and not the input. From the #6 in the table for Laplace
Transforms, we know that
s
L [cos(ωt)] = 2 (4.64)
s + ω2
Therefore,
F0 s
F(s) = 2 (4.65)
s + ω2
Ignoring the initial conditions, and therefore considering only the particular solution, and taking
the Laplace transform of the EOM equation yields:
F0 s
(ms2 + cs + k)X(s) = (4.66)
s + ω2
2
where X(s) denotes the Laplace transform of the unknown function x(t) and s is the complex
transform variable. Rearranging the above equation for X(s) yields:
F0 s
X(s) = (4.67)
(ms + cs + k)(s2 + ω 2 )
2
101
Now that we have F(s) and X(s) we can obtain H(s) as
X(s) F0 s s2 + ω 2 1
H(s) = = 2 2 2
· = 2 (4.68)
F(s) (ms + cs + k)(s + ω ) F0 s ms + cs + k
or
1
H(s) = (4.69)
ms2 + cs + k
This ratio is termed the transfer function of a system and is an important tool in vibration analysis.
Sometimes, how the system responds to an input with certain frequency components is impor-
tant in understanding the system in general, therefore, we want to solve for the frequency response
function of the system. The frequency response function is denoted as H( jω) where the complex
number s is replaced by the frequency component of the system while considering the imaginary
portion in the complex plane (i.e., s = jω). Therefore, the frequency response function of the
system becomes:
1 1
H( jω) = = (4.70)
m( jω)2 + c jω + k −mω 2 + c jω + k
rearranging into a standard form yields:
1
H( jω) = (4.71)
k − mω 2 + cω j
recall that j2 = −1. This is the frequency response function of the system. Therefore, it can be seen
that the frequency response function of the system is the transfer function of the system evaluated
along the imaginary axis s = jω. However, this expression contains imaginary values (that help
to account for the phase in the system) and therefore can be challenging to work with. As the
amplitude |H( jω)| of the response (the real portion of the equation) is useful to the practitioner, it
is prudent to consider the special case of amplitude response while neglecting the phase response.
Consider that:
H( jω) = R + I j (4.72)
so p
|H( jω)| = R2 + I2 (4.73)
multiplying H( jω) by 1 that is represented by its unit complex conjugate yields:
k − mω 2 − cω j

1
H( jω) = (4.74)
k − mω 2 + cω j k − mω 2 − cω j
k − mω 2

−cω
= j (4.75)
(k − mω 2 )2 (cω)2 (k − mω 2 )2 (cω)2
k−mω 2 −cω
therefore, R = (k−mω 2 )2 (cω)2
and I = (k−mω 2 )2 (cω)2
. Now, calculating the amplitude of H( jω) we
102
get:
H(ω) = |H( jω)| (4.76)

p
= R2 + I2 (4.77)
s
(k − mω)2 + (−cω)2
= 2 (4.78)
(k − mω 2 )2 + (cω)2 )
s
1
= (4.79)
(k − mω )2 + c2 ω 2
2
1
=p (4.80)
(k − mω 2 )2 + c2 ω 2
where H(ω) represents only the amplitude of the frequency response function and therefore drops
the j term from the expression.
Review 4.2 To recap, for a single DOF damped spring-mass system the transfer function is:
1
H(s) = (4.81)
ms2 + cs + k
And the frequency response function is:
1
H( jω) = (4.82)
k − mω 2 + cω j
While the amplitude of the frequency response is:
1
H(ω) = |H( jω)| = p (4.83)
(k − mω 2 )2 + c2 ω 2
Example 4.4 Considering the forced system:
Set the forcing function to be F0 sin(ωt) and calculate the transfer function.
103
Solution: The equation of motion for the system is:
mẍ + cẋ + kx = F0 sin(ωt) (4.84)
From the #6 in the table for Laplace Transforms, we know that:

ω
L [sin(ωt)] = (4.85)
s2 + ω 2
Therefore,
F0 ω
F(s) = (4.86)
s + ω2
2
Ignoring the initial conditions and taking the Laplace transform of the EOM equation yields:
F0 ω
(ms2 + cs + k)X(s) = (4.87)
s + ω2
2
Solving algebraically for the X(s) yields:

F0 ω
X(s) = (4.88)
(ms + cs + k)(s2 + ω 2 )
2
Now that we have F(s) and X(s) we can obtain H(s) as
X(s) F0 ω s2 + ω 2 1
H(s) = = · = (4.89)
F(s) (ms2 + cs + k)(s2 + ω 2 ) F0 ω ms2 + cs + k
or
1
H(s) = (4.90)
ms2 + cs + k
This is identical to the solution obtained using F0 cos(ωt) as would be expected because the
transfer function is related to the system and not to the input.
Review 4.3 The frequency domain is a mathematical representation of a signal or data in

terms of its frequency components, as opposed to its temporal or time-based representation.
The frequency domain provides a different perspective on the signal by decomposing it into
its constituent sinusoidal signals at discrete frequencies and their respective magnitudes. A
3D rerensation of this process is shown in figure 4.15. The transformation between the time
domain and the frequency domain is typically achieved using mathematical techniques such
as the Fourier Transform or the Fast Fourier Transform (FFT).
104
Figure 4.15: 3D visualization of time and frequency domains where a temporal signal is
decomposed into constituent sinusoidal signals.
4.4.3 Response to Random Inputs

The transfer and frequency response functions can be very useful for determining the system’s
response to random inputs. Up to this point, we have solved for deterministic input.
• Deterministic-For a known time t, the value of the input force F(t) is precisely known.
• Random For a known time t, the value of the input force F(t) is known only statistically.
To expand, a random signal is a signal with no obvious pattern. For these types of signals, it is
not possible to focus on the details of the input signal, as is done with a deterministic signal, rather
the signal is classified and manipulated in terms of its statistical properties.
Randomness in vibration analysis can be thought of as the result of a series of results obtained
from testing a system’s repeatability for various inputs under varying conditions. In these cases,
one record or time history is not enough to describe the system. Rather, an ensemble of various
tests are used to describe how the system will respond to the various inputs.
First, let us consider two inputs, a deterministic input (typical sin wave), and a random input
(white noise). These inputs are shown in figure 4.16.
105
Figure 4.16: Two arbitrary inputs: (a) sinusoidal; and (b) uniform random noise.
One of the first factors to consider is the mean of the random signal x(t), defined as:
Z T
1
E[x] = x̄ = lim x(t)dt (4.91)
T →∞ T 0
where T is the length in time of the data collected. However, for random signals, we often want to
consider signals with an average mean of zero (i.e. x̄(t) = 0). For signals not centered around zero,
we can obtain a zero-centered signal if the signal is stationary and we subtract the mean value x̄
from the signal x(t). This can be written as:
x′ (t) = x(t) − x̄ (4.92)
where the x′ (t) is now centered around zero. As mentioned before, it is important to consider
whether or not the input signals are stationary. A signal is stationary if its statistical properties
(usually expressed by its mean) do not change with time. Here, it can be seen that for our inputs
considered the signals are stationary if a long enough time period is considered.
Another important variable is the variance (or mean-square value) of the random variable x(t)
defined as:
1 T
Z
2
E[(x − x̄) ] = lim (x(t) − x̄)2 dt (4.93)
T →∞ T 0
and provides a measure of the magnitude of the fluctuations in the signal x(t). If the signal has an
expected value of zero, or E[x] = 0, this simplifies to.
Z T
2 1
E[x ] = x2 = lim x2 (t)dt (4.94)
T →∞ T 0
This expression leads to the calculation of the root-mean-square (RMS) of the signal:
p
xrms = x2 (4.95)
106
Considering a nonstationary signal, an important measure of interest is how fast the value of
the variables changes. This is important to understand as it provides context for how long a signal
must be sampled before a meaningful representation of the signal can be calculated in a statistical
sense. One way to quantify how fast the values of signal change is the autocorrelation function:
Z T
1
Rxx (τ) = lim x(t)x(t + τ)dt (4.96)
T →∞ T 0
The subscript xx denotes that this is a measure of the response for the variable xx, τ is the time
difference between the values at which the signal x(t) is sampled. The autocorrelation for the two
inputs considered above is expressed in figure 4.17.
Figure 4.17: Responses from the autocorrelation function for the inputs shown in figure 4.16 show-
ing: (a) a sinusoidal; and (b) uniform random noise.
Note that the value of τ selected in the autocorrelation function greatly affects its response for the
sinusoidal input. This is because the values for the sinusoidal are highly correlated. To expand, the
value at any time t is greatly affected by the values immediately before and after it. This is not the
case for the random input where the signal is not correlated and therefore there is little difference
in changing the value of τ on the response of the autocorrelation function.
Next, if we take the Fourier transform of the autocorrelation function we obtain the power
spectral density (PSD) defined as:
1
Z ∞
Sxx (ω) = Rxx (τ)e− jωτ dτ (4.97)
2π −∞
where the integral of Rxx (τ) changes the real number τ into the frequency-domain value ω. The
frequency spectrum is denoted with S and the subscript of the considered variable (e.g., Sxx (ω)).
The frequency spectrum for the two input cases considered are plotted in figure 4.18.
107
Figure 4.18: Power spectral density plots for the inputs shown in figure 4.16 showing: (a) a sinu-
soidal; and (b) uniform random noise.
where the flat frequency response for the random input denotes that the random input is white
noise input. This flat frequency response in the frequency domain can be denoted S0 , such that
S f f (ω) = S0 or Sxx (ω) = S0 , depending on whether the frequency spectrum of the input ( f f ) or
output (xx) is being considered. While a true white noise input would be perfectly flat, white noise
is really just a theoretical concept as all real-world data will have some variation in the frequency
domain as diagrammed in figure 4.18(b).
Recall that Sxx is the spectrum of the response of the system. For the one-DOF system consid-
ered here, we can express the arbitrary input as a series of impulse inputs as discussed in section
4.3. This knowledge, along with the frequency response function can be used to relate the spectrum
of the input S f f (ω) to the output through the transfer function as:
" Z #
1 ∞
Sxx (ω) = |H( jω)|2 R f f (τ)e− jωτ dτ (4.98)
2π −∞
This can also be expressed in symbolic form as:
Sxx (ω) = |H( jω)|2 S f f (ω) (4.99)
where R f f denotes the autocorrelation function of F(t) and S f f denotes the PSD of the forcing
function F(t). The notation |H( jω)|2 is the square of the magnitude of the complex frequency
response. A more detailed derivation can be found in Raoa , Inmanb , or Newlandc , but here it is
more important to study the results rather than the derivations.
a Singiresu, S. Rao. “Mechanical vibrations”. Boston, MA: Addison Wesley, 1995.

b Inman, Daniel J., and Ramesh Chandra Singh. “Engineering vibrations”. Vol. 3. Englewood Cliffs, NJ: Prentice
Hall, 1994.
c Newland, David E. “Random vibrations, spectral & wavelet analysis.” Longman Scientific & Technical (1993).
108
Example 4.5 Consider the following system
Calculate the PSD of the response x(t) given that the PSD of the applied force S f f (ω) is
white noise.
Solution: From the system we know that the EOM is
mẍ(t) + cẋ(t) + kx(t) = F(t) (4.100)
The frequency response function for this system is

1
H( jω) = (4.101)
k − mω 2 + cω j
while the amplitude of the response is:
1
H(ω) = |H( jω)| = p (4.102)
(k − mω 2 )2 + c2 ω 2
Applying the equation that relates S f f (ω) to Sxx (ω) we get:

2
1
Sxx (ω) = |H( jω)|2 S f f (ω) = S f f (ω) (4.103)
k − mω 2 + cω j
White noise means the forcing function S f f (ω) is constant across the frequency spectrum,
therefore, S f f (ω) = S0 . Additionally as:
2
1 1
|H( jω)|2 = = (4.104)
k − mω 2 + cω j (k − mω 2 )2 + c2 ω 2
where the absolute value is the amplitude of the system. Therefore, we obtain:
1 S0
Sxx (ω) = |H( jω)|2 S0 = S0 = (4.105)
(k − mω 2 )2 + c2 ω 2 (k − mω )2 + c2 ω 2
2
Using various values for the elements in the system, the PSD for the system considered looks
like:
109
Figure 4.20: Response for considered 1-DOF systems subjected to a white noise input.
Another useful quantity to consider is the expected output, in terms of its mean and variance,
for a given input. Working within the constraint that the system will oscillate about zero, E[x] = 0,
the mean-square value can be directly related to the PSD function as:
Z ∞
2
E[x ] = x2 = |H( jω)|2 S f f (ω)dω (4.106)
−∞
For a constant input S0 , as diagrammed in figure 4.18(b), the mean-square value can be expressed
as: Z ∞
E[x2 ] = x2 = S0 |H( jω)|2 dω (4.107)
−∞
After inspecting Rthe above equation, it becomes clear that to obtain the square of the expected
value, a solution for −∞ ∞
|H( jω)|2 dω must be obtained. For cases where S f f (ω) = S0 and as such
S f f (ω) can be pulled out of the integral, these integrals have beenRsolved [Random Vibrations,
Spectral & Wavelet Analysis, Newland (1993)]. For example, given −∞ ∞
|H( jω)|2 dω:
2
B0 πB20
Z ∞
dω = (4.108)
−∞ A0 + jωA1 A0 A1
and
2
B0 + jωB1 π(A0 B21 + A2 B20 )
Z ∞
dω = (4.109)
−∞ A0 + jωA1 − ω 2 A2 A0 A1 A2
When combined with equation 4.107, these integrals allow for the easy calculation of the expected
values.
Example 4.6 For the system below, calculate the mean-square response of the system given
that the spectrum of the input force F(t) is a perfect theoretical white noise.
110
Solution: Again, as the forcing function S f f (ω) is constant across the frequency spectrum
S f f (ω) = S0 the mean-square response can be calculated as:
Z ∞
2
E[x ] = x2 = S0 |H( jω)|2 dω (4.110)
−∞
Using the already tabulated response:

2
B0 + jωB1 π(A0 B21 + A2 B20 )
Z ∞
dω = (4.111)
−∞ A0 + jωA1 − ω 2 A2 A0 A1 A2
and the frequency response function for the system as derived in equation 4.71:
1
H( jω) = (4.112)
k − mω 2 + cω j
when B0 = 1, B1 = 0, A0 = k, A1 = c, and A2 = m. Therefore, using the tabulated expression
we can show that:
πm S0 π
E[x2 ] = S0 = (4.113)
kcm kc
111
Table of Laplace Transforms for Vibrations
This is a partial list of important Laplace transforms for vibrations and assumes
zero initial conditions, 0 < t, and ζ < 1.
f (t) L [ f (t)] = F(s) f (t) L [ f (t)] = F(s)
δ (t) 1 (1) 1 1
3
ωt − sin(ωt) (17)
ω s2 (s2 + ω 2 )
δ (t − t0 ) e−st0 (2)
1
3
sin(ωt) − ωt cos(ωt) . . .
1 2ω
1 (3) 1
s (18)
(s2 + ω 2 )2
1
eat (4) t s
s−a sin(ωt) (19)
2ω (s2 + ω 2 )2
ω
sin(ωt) (5) 2ωs
s2 + ω 2 t sin(ωt) (20)
(s2 + ω 2 )2
s
cos(ωt) (6)
s2 + ω 2 s2 − ω 2
t cos(ωt) (21)
ω (s2 + ω 2 )2
sinh(ωt) (7)
s − ω2
2
ω
eat sin(ωt) (22)
s (s − a)2 + ω 2
cosh(ωt) (8)
s − ω2
2
s−a
eat cos(ωt) (23)
1 1 (s − a)2 + ω 2
1 − cos(ωt) (9)
ω2 s(s2 + ω 2 )
ω
eat sinh(ωt) (24)
1 −ζ ωt 1 (s − a)2 − ω 2
e sin(ωd t) (10)
ωd s + 2ζ ωs + ω 2
2
s−a
eat cosh(ωt) (25)
ω (s − a)2 − ω 2
1 − e−ζ ωt sin(ωd t + φ ), φ = cos−1 (ζ ) . . .
ωd
ω2 1 1
(11) sin(ω2t) − sin(ω1t) . . .
s(s2 + 2ζ ωs + ω 2 ) ω2 ω1
ω12 − ω22
(26)
t n−1 1 (s2 + ω12 )(s2 + ω22 )
, n = 1, 2, . . . (12)
(n − 1)! sn
s(ω12 − ω22 )
cos(ω2t) − cos(ω1t) (27)
n! (s2 + ω12 )(s2 + ω22 )
t n , n = 1, 2, . . . (13)
sn+1
eat f (t) F(s − a) (28)
n!
t n eωt , n = 1, 2, . . . (14)
(s − ω)n+1 f (t − a)Φ(t − a) e−as F(s) (29)
1 1
(1 − e−ωt ) (15) e−as
ω s(s + ω) Φ(t − a) (30)
s
1 −ωt 1
(e + ωt − 1) (16) f ′ (t) sF(s) − f (0) (31)
ω2 s2 (s + ω)
Open Vibrations
5 Multiple Degree-of-freedom Systems

Until now we have only considered and modeled systems that can require one coordinate system
to describe their motion. In this chapter, we will develop the mathematical tools required to model
multiple degree-of-freedom systems that require multiple independent coordinates to describe their
motion. As before, the equations that describe the motion of rigid bodies in space are developed
from Newton’s second law of motion. However, unlike before, there exists an independent equation
for each body in motion. These equations are therefore coupled by the system and are often
expressed in matrix notation such that the mass, damping, and stiffness matrices are easily defined.
Review 5.1 Linear Algebra

Linear algebra allows for the efficient solving of these coupled equations. In this text,
matrices are expressed as bold capital letters (X), vectors are denoted with an arrow (⃗x), and
scalars/variables with italic letters (x). However, given the range of notation needed, it is not
always possible to strictly follow this formulation.
The dot product allows us to multiply matrices and is defined as:

a b e ae + b f
= (5.1)
c d f ce + d f
Another arrangement of the same principle, in a format more related to vibrations, is:

a1 + a2 b e (a1 + a2 )e + b f
= (5.2)
c d f ce + d f
The transpose of a matrix is an operator which flips a matrix over its diagonal. For a
matrix A, the transpose AT can be written as:
 
a b
T a c e
A = c d  → A = (5.3)
b d f
e f
A matrix is symmetric if A = AT . Therefore, the symmetric matrix must be square and

can be written as:
   
a b c a d g
A = d e f  = AT = b e h , where b = d, c = g, f = h (5.4)
g h i c f i
The determinant of a matrix is a scalar value that is a function of the entries of a square
matrix. The determinant characterizes the matrix and its linear map. The determinant is
often writted as det(A), det A, or |A|. For a 2 × 2 matrix this is defined as:

a b
det(A) = ad − bc, when A = (5.5)
c d
The inverse of a square matrix is such that AA−1 = A−1 A = I where I is the identity
matrix:

1 0
I= (5.6)
0 1
Open Vibrations 5.1 General Discussion on Mode Shapes
and the inverse of a 2 × 2 matrix is defined as:

−1 1 d −b a b
A = , when A = (5.7)
det(A) −c a c d
A matrix that does not have an inverse is called a singular matrix.
5.1 General Discussion on Mode Shapes

Studying and characterizing the natural frequencies of a system allows for the detailed investigation
of the system response. Modern vibration analysis relies heavily on the concepts of mode shapes
for various engineering tasks. Practical applications of the study of mode shapes (often called
experimental modal analysis) include
• Correlation Finite Element Analysis with structures
• Structural Dynamic Modification
• Reduction of Finite Element Analysis models
• Forced Response Prediction
• Active Vibration Control
Vibration Case Study 5.1 In automotive engineering, the requirements for safe and com-
fortable vehicles necessitate the need for a thorough understanding of the vehicle’s dynamic
properties and how any design changes affect its dynamics. Experimental modal analysis is
an important troubleshooting and model-updating tool in the study of vehicle noise and vi-
bration harshness (NVH). Oftentimes, experimental modal analysis is performed on a “body
in white” or a sub-frame structure to develop a better understanding of the dynamics of the
structure. Overall, experimental modal analysis is an important tool used in improving a
vehicle’s NVH performance.
114
Open Vibrations 5.1 General Discussion on Mode Shapes
Figure 5.1: Experimental modal analysis of an automotive (Jaguar) body in white, typically
done to reduce vehicle noise and vibration harshness a .
a Cjp24, CC BY-SA 3.0 <https://creativecommons.org/licenses/by-sa/3.0>, via Wikimedia Commons
Mode shapes are not the displacement of a system, rather they describe the configurations
into which a structure will naturally displace at a given frequency. For example, consider the 4-
DOF system shown in figure 5.2 that represents a pole (i.e. cantilever beam). Assuming that the
system experiences a linear response and using the mode-superposition method we can see that the
displaced shape ⃗x is a function of all of the mode shapes ui and their corresponding participation
factors qi . Note that the mode shapes associated with the lower frequencies tend to provide the
greatest contribution to structural response. As the frequencies that excite the modes increase, the
mode shapes contribute less, are predicted less reliably, and are harder to measure. Therefore, the
analysis of the system is often truncated after the first few modes and rarely exceeds the 10th mode.
Figure 5.2 shows a structure with N degrees of freedom that therefore had N corresponding
mode shapes. Each mode shape is independent and normalized such that the maximum displace-
ments are the same. The summation of the mode shapes multiplied by their corresponding partici-
pation factors (qi ) yields the deflection of the structure.
115
Open Vibrations 5.2 Modeling Undamped Two Degree of Freedom Systems
Figure 5.2: Deflection of a vertical cantilever,⃗x, is a function of the considered mode shapes ui and
their corresponding participation factors qi .
5.2 Modeling Undamped Two Degree of Freedom Systems

Consider the undamped 2-DOF systems presented in figure 5.3. This system with a single mass
capable of moving in two directions. To expand, figure 5.3(a) reports a mass that can move horizon-
tally or vertically in space. However, this mass does not rotate during its movements. Moreover,
figure 5.3(b) presents a system that rotates about the spring and displaces vertically. These are
examples of 2 DOF systems because each system has two independent coordinate systems that
express the movement of the mass.
Figure 5.3: Examples of single mass 2-DOF systems that: (a) displaces in the vertical and hori-
zontal directions, and; (b) rotates about the spring and displaces in the vertical direction.
Another example of a 2-DOF system with two masses, each with their own independent coordi-
nate system, is presented in figure 5.4. The two coordinates that describe the system’s movements
are x1 and x2 .
116
Figure 5.4: 2-DOF system with two masses and two independent coordinate systems x1 and x2 .
5.2.1 Solution for the Two-Degree-of-Freedom System

Before we derive a model for undamped 2-DOF systems, let us first consider the solution to the
system shown in figure 5.4. The solution consists of two equations, one for each mass. This
solution will be derived in section 5.2.2 and is expressed by the coupled equations:
x1 (t) = A1 sin(ω1t + φ1 )u11 + A2 sin(ω2t + φ2 )u12 (5.8)
x2 (t) = A1 sin(ω1t + φ1 )u21 + A2 sin(ω2t + φ2 )u22 , ω1 or ω2 ̸= 0

These two equations can be written as a single equation in matrix form as:
⃗x(t) = A1 sin(ω1t + φ1 )⃗u1 + A2 sin(ω2t + φ2 )⃗u2 , ω1 or ω2 ̸= 0 (5.9)
Where the bold text denotes vectors. Therefore, the vectors ⃗u1 and ⃗u2 are the mathematical expres-
sions that “couple” or tie the equations together. Expanding these vectors shows:

x1 (t) u11 u
⃗x(t) = , ⃗u1 = , ⃗u2 = 12 , (5.10)
x2 (t) u21 u22
The four key components of the solution expressed in equation 5.9 are:
1. ω1 and ω2 are the natural frequencies of the system. They are not the frequencies of the
masses. The solution states that each of the masses oscillates at the two frequencies ω1 and
ω2 . Moreover, consider the special case where the initial conditions are selected to force
A2 = 0, in this case, each mass would only oscillate at only one frequency, ω1 .
2. A1 and A2 are the constants of integration and determine the amplitude of the system.
3. φ1 and φ2 represent the phase shift of the system
4. ⃗u1 and ⃗u2 are the first and second mode shapes of the system and couple the system together.
5.2.2 Deriving the Solution for the Two-Degree-of-Freedom System

To derive this solution for the system under consideration an FBD for figure 5.4 can be constructed
for the forces acting on each mass. First we have to make the assumption that x1 < x2 , this allows
us to say that m2 pulls on m1 and results in:
117
Figure 5.5: Free body diagram for the 2-DOF system presented in figure 5.4.
Applying Newton’s second law and summing the forces on each mass in the horizontal direction
yields:
m1 ẍ1 = −k1 x1 + k2 (x2 − x1 ) (5.11)

m2 ẍ2 = −k2 (x2 − x1 )
These equations can be rearranged in terms of x1 and x2 as:
m1 ẍ1 + (k1 + k2 )x1 − k2 x2 = 0 (5.12)

m2 ẍ2 − k2 x1 + k2 x2 = 0
where these are two coupled second-order differential equations that each require two initial con-
ditions to solve. These initial conditions can be obtained from the displacement and velocity terms
as:
x1 (0) = x10 (5.13)

ẋ1 (0) = ẋ10 = v10
x2 (0) = x20
ẋ2 (0) = ẋ20 = v20
As before, these initial conditions will be the constants of integration used to solve the two second-
order differential equations. This solution will provide the free response of each mass in the system.
There is a multitude of ways to solve these two coupled second-order differential equations, how-
ever, here we will just consider a matrix notation solution. This matrix notation solution is used as
this formulation is readably solved using computers and is expandable to more than 2 DOF.
To initiate the solution, let us first develop the matrix formulation of the two coupled ODEs:

m1 0 x¨1 k1 + k2 −k2 x1 0
+ = (5.14)
0 m2 x¨2 −k2 k2 x2 0
This equation can also be expressed as the vector equation:
M⃗x¨ + K⃗x = 0 (5.15)
and is known as the EOM in vector form. In this formulation, the mass matrix (M) is defined as:

m1 0
M= (5.16)
0 m2
118
while the stiffness matrix (K) is:

k1 + k2 −k2
K= (5.17)
−k2 k2
along with the displacement, velocity, and acceleration matrices:

x1 ẋ1 ẍ
⃗x = ˙
, ⃗x = , ⃗x = 1
¨ (5.18)
x2 ẋ2 ẍ2
Beyond these equations, we can write the initial conditions as:

x1 (0) ˙ ẋ1 (0)
⃗x0 = , ⃗x0 = (5.19)
x2 (0) ẋ2 (0)
This simple connection between vibration analysis and matrix analysis allows computers to be
used to solve large and complicated vibration problems quickly.
Recall that the 1-DOF version of the equation of motion was solved by calculating the values
of the constants in an assumed harmonic solution. The same approach is applied here in order to
solve for the displacement of the two-DOF system. This time, the solution is assumed in the form:
⃗x(t) = ⃗ue jωt (5.20)

where ⃗u is a vector of constants to be demerited and can be written as:

u
⃗u = 1 (5.21)
u2
√
From before, ω is also a constant to be determined. Again, j = −1. In the same manner as
before, e jωt represents harmonic motion as e jωt = cos(ωt) + j sin(ωt). Taking the derivatives of
⃗x(t) = ⃗ue jωt yields:
˙ = jω⃗ue jωt
⃗x(t) (5.22)
¨ = −ω 2⃗ue jωt
⃗x(t) (5.23)
Substituting this into the EOM in vector form (M⃗x¨ + K⃗x = 0) yields:
−ω 2 M⃗ue jωt + K⃗ue jωt = 0 (5.24)
or
(−ω 2 M + K)⃗ue jωt = 0 (5.25)
As e jωt ̸= 0 for any value of t and not allowing ⃗u to be zero it can be demerited that (−ω 2 M + K)
must satisfy the vector equation. Therefore,
(−ω 2 M + K)⃗u = 0, ⃗u ̸= 0 (5.26)
This forms a homogeneous set of algebraic equations. To be useful, these equations have a
nonzero solution for the system must exist. For this to be true, the inverse of the coefficient matrix
(−ω 2 M + K) must not exist. To expand, assume that the inverse of (−ω 2 M + K) does exist, by
119
multiplying both sides of the equation by (−ω 2 M + K)−1 yields ⃗u = 0. This is a trivial solution
(it is not useful) as no motion in the system is implied. Therefore, the logical connection can be
drawn between the solution of the equation and the inverse of the coefficient matrix (−ω 2 M + K).
Applying the singularity condition to the coefficient matrix of equation (−ω 2 M +K)⃗u = 0, ⃗u ̸=
0 results a nonzero solution of ⃗u. However, for this to exist the following must be true:
det(−ω 2 M + K) = 0 (5.27)
Solving this expression results in one algebraic equation with one unknown (ω). Expanding the
above equation to consider the values for the matrices M and K results in:
−ω 2 m1 + k1 + k2

−k2
det =0 (5.28)
−k2 −ω 2 m2 + k2
Using the definition of the determinant yields that the unknown quantity ω 2 must satisfy:
m1 m2 ω 4 − (m1 k2 + m2 k1 + m2 k2 )ω 2 + k1 k2 = 0 (5.29)
This expression is called the characteristic equation for the system and is used to determine the
constants ω1,2 , in the assumed form of the solution given by the assumed solution ⃗x(t) = ⃗ue jωt ,
once the values of the physical parameters m1 , m2 , k1 , and k2 are known. Note that ω1,2 are not in
the characteristic equation, therefore, solving for ω1,2 will be done by factoring the equation above
to obtain two solutions ω1 and ω2 . The characteristic equation is in the form of the quadratic
formula if you set x = ω 2 , as:
ax2 + bx + c = 0 (5.30)
After finding the value of ω1,2 using the characteristic equation, the values in ⃗u can be found
using equation (−ω 2 M + K)⃗u = 0, ⃗u ̸= 0 for each value of ω 2 . That is, for both ω1 and ω2 there
is a vector ⃗u that satisfies the equation. These solutions can be written as:
(−ω12 M + K)⃗u1 = 0 (5.31)
and
(−ω22 M + K)⃗u2 = 0 (5.32)
The direction of the vectors ⃗u1 and ⃗u2 can be obtained by solving the above expressions, however,
the information regarding the magnitude of is not contained in this expression. To verify this,
assume that ⃗u1 satisfies the equation, therefore, the vector a⃗u1 also satisfies the equation where a
is any nonzero number. Hence the vectors satisfying the above are of arbitrary magnitude.
The values obtained for ⃗u1 and ⃗u2 can now be combined with the assumed solution:
⃗x(t) = ⃗ue jωt (5.33)
to form a set of solutions:
⃗x(t) = ⃗u1 e− jω1t , ⃗u1 e jω1t , ⃗u2 e− jω2t , ⃗u2 e jω2t (5.34)
Since the equation to be solved is linear, the solution is the sum of these solutions. This results in:
⃗x(t) = (ae jω1t + be− jω1t )⃗u1 + (ce jω2t + de− jω2t )⃗u2 (5.35)
120
where a, b, c, and d are the arbitrary constants of integration to be determined by the initial con-
ditions. Applying Euler’s formulas for the sin functions (where ω1 or ω2 ̸= 0) reorganizes this
equation as:
⃗x(t) = A1 sin(ω1t + φ1 )⃗u1 + A2 sin(ω2t + φ2 )⃗u2 , ω1 or ω2 ̸= 0 (5.36)
Another way to write this equation is in the form:

x1 (t) A1 sin(ω1t + φ1 )
= ⃗u1 ⃗u2 , ω1 or ω2 ̸= 0 (5.37)
x2 (t) A2 sin(ω2t + φ2 )
Where the values for A1 and A2 can be obtained by setting applying the boundary conditions and
taking the derivatives of the equations as done in the 1-DOF problems.
The final form of the equation provides physical insight into the solution of the system. It
states that each mass in the system oscillates at both of the natural frequencies of the system (ω1
and ω2 ). Furthermore, the importance of the initial conditions can be understood. Assume that
initial conditions are chosen that result in A2 = 0, this cancels out the second natural frequency
such that each mass oscillates at only one frequency, ω1 . Moreover, the positions of the masses
can be determined by the values of the vector ⃗u1 at any given time. For this reason, ⃗u1 is termed
the first mode shape of the system. Likewise, if the opposite initial conditions are chosen such that
A1 = 0, then both system coordinates (e.g., masses in the systems we have studied) will oscillate
at ω2 and again, the positions can be obtained from the vector ⃗u2 . Where ⃗u2 is termed the second
mode shape. The interactions between mode shapes and natural frequencies are very important
and form the basis of several areas in the field of vibrations.
Example 5.1 Considering the following system:
Figure 5.6: 2-DOF system with two masses and two independent confidante systems x1 and
x2 .
Calculate response for the system if m1 =9 kg, m2 =1 kg, k1 = 24 N/m, and k2 = 3 N/m with
the initial conditions x10 = 1 mm, v10 = 0 mm/s, x20 = 0 mm, and v20 = 0 mm/s.
Solution:
We have already obtained a characteristic equation for this system. This is shown in
Equation 5.29 and is given as:
m1 m2 ω 4 − (m1 k2 + m2 k1 + m2 k2 )ω 2 + k1 k2 = 0 (5.38)
121
Substituting our values into this obtains:
9 · 1ω 4 − (9 · 3 + 1 · 24 + 1 · 3)ω 2 + 24 · 3 = 0 (5.39)
or
ω 4 − 6ω 2 + 8 = 0 (5.40)
This can then be factored into:
(ω 2 − 2)(ω 2 − 4) = 0 (5.41)
This results in solutions of ω12 = 2 and ω22 = 4. Leading to:
√
ω1 = ± 2 rad/sec, ω2 = ±2 rad/sec (5.42)
Next, we need to obtain solutions for ⃗u1 and ⃗u2 . Having

√ solved for ω1 and ω2 we can obtain
T
this. First, knowing ⃗u1 = [u11 u21 ] and using ω1 = 2 and the following equation:
(−ω12 M + K)⃗u1 = 0 (5.43)
yields simplified to

9 0 24 + 3 −3 u11 0
−2 + = (5.44)
0 1 −3 3 u21 0
simplified to
27 − 9 · 2 −3 u11 0
= (5.45)
−3 3 − 2 u21 0
or
9 −3 u11 0
= (5.46)
−3 1 u21 0
Taking the dot product of the matrix equation yields:
9u11 − 3u21 = 0, and − 3u11 + u21 = 0 (5.47)
Both of these equations yield the same equation, that is:

u11 1
= (5.48)
u21 3
As mentioned before, only the ratio of the elements is determined here. To show this is true it
is easily seen that:
1 1
u11 = u21 → au11 = au21 (5.49)
3 3
To obtain a numerical value, we arbitrarily assign a value to one of the elements. Here, let
u21 = 1 so let u11 = 1/3. Therefore,
1

⃗u1 = 3 (5.50)
1
122
The same processes can be used for obtaining ⃗u2 using ω2 = 2, this results in:

−9 −3 u12 0
= (5.51)
−3 −1 u22 0
−9u12 − 3u22 = 0, and − 3u12 − u22 = 0 (5.52)
Both of these equations yield the same equation, that is:

u12 1
=− (5.53)
u22 3
Again, assuming u22 = 1 this can be rearranged into ⃗u2 as:
1
−3
⃗u2 = (5.54)
1
Where ⃗u1 and ⃗u2 represent only the directions and shape of the mode shapes and not the
magnitude of the mode shapes. Now that we have the mode shapes, we can solve for the
initial conditions A1 and A2 . To do this, let us use the following formulation of the solution:

x1 (t) A1 sin(ω1t + φ1 )
= ⃗u1 ⃗u2 , ω1 or ω2 ̸= 0 (5.55)
x2 (t) A2 sin(ω2t + φ2 )
Adding our values for the problem at t = 0 this becomes:

1
− 31 A1 sin(φ1 )

1 3
= (5.56)
0 1 1 A2 sin(φ2 )
and after applying the dot product:

1 1

1 A1 sin(φ 1 ) − A 2 sin(φ2 )
= 3 3 (5.57)
0 A1 sin(φ1 ) + A2 sin(φ2 )
Next we can differentiate the equation for x(t) to obtain the velocity solution. Adding our
values for the problem at t = 0 obtains:
" √2 2
#
ẋ1 (0) v 0 A cos(φ1 ) − 3 A2 cos(φ2 )
= 10 = = √3 1 (5.58)
ẋ2 (0) v20 0 2A1 cos(φ1 ) + 2A2 cos(φ2 )
Now that we have 4 equations for 4 unknowns we can use these equations to solve for A1 , A2 ,
φ1 , and φ2 . The 4 equations are:
3 = A1 sin(φ1 ) − A2 sin(φ2 ) (5.59)
123
0 = A1 sin(φ1 ) + A2 sin(φ2 ) (5.60)

√
0 = 2A1 cos(φ1 ) − 2A2 cos(φ2 ) (5.61)
√
0 = 2A1 cos(φ1 ) + 2A2 cos(φ2 ) (5.62)
Setting these last two equations equal to each other yields:
√ √
0 = 2A1 cos(φ1 ) + 2A2 cos(φ2 ) = 2A1 cos(φ1 ) − 2A2 cos(φ2 ) (5.63)
or:
0 = −4A2 cos(φ2 ) (5.64)
√
For this equation to be true, φ2 = π2 . Therefore, applying this to 0 = 2A1 cos(φ1 )+2A2 cos(φ2 )
results in: √
0 = 2A1 cos(φ1 ) (5.65)
where again, for this equation to be true, φ1 = π2 . Now the first two equations become:
3 = A1 − A2 (5.66)
0 = A1 + A2 (5.67)
Where this shows us that A1 = 23 and A2 = − 32 . Therefore, now that we have the initial
conditions we can find a solution for the temporal response of each mass. Using the equations
from before:
And applying our obtained values
√

3 π 1 3 π 1
x1 (t) = sin( 2t + ) + − sin(2t + ) − (5.70)
2 2 3 2 2 3
√

3 π 3 π
x2 (t) = sin( 2t + ) + − sin(2t + ) (5.71)
2 2 2 2
results in:
√

1 π π
x1 (t) = sin( 2t + ) + sin(2t + ) (5.72)
2 2 2
√

3 π π
x2 (t) = sin( 2t + ) − sin(2t + ) (5.73)
2 2 2
These results can be plotted as:
124
Figure 5.7: Temporal response for each of the rigid bodies in the 2-DOF system.
Example 5.2 Mode shapes can be better understood through a graphical representation. To
do this, consider the 2-DOF system presented in figure 5.8(a). Assuming that x1 < x2 the FBD
for the system is expressed in figure 5.8(b).
Figure 5.8: (a) 2-DOF system with two masses arranged in a vertical configuration; and (b)
FBD of system.
For simplicity, all masses, and spring stiffness are considered equal and that m = 1 and k = 1.
Solution:
125
From the previous investigations in this text, we know that the forces caused by gravity
will cancel out. Therefore, the EOM for the system can be written as:
mẍ1 = −kx1 + k(x2 − x1 ) (5.74)

mẍ2 = −k(x2 − x1 ) − kx2
These equations can be written in matrix notation as:

m 0 x¨1 2k −k x1 0
+ = (5.75)
0 m x¨2 −k 2k x2 0
Substituting the values of the matrices M and K into this expression det(−ω 2 M + K) = 0
yields:
−ω 2 m + 2k

−k
det =0 (5.76)
−k −ω 2 m + 2k
The determinant yields that the unknown quantity, ω 2 , must satisfy:
m2 ω 4 − 4kmω 2 + 3k2 = 0 (5.77)
using the quadratic formula we obtain
3k √
r r
k
ω1 = ± = 1 rad/sec, ω2 = ± = 3 rad/sec (5.78)
m m
Now, we need to obtain solutions for ⃗u1 and ⃗u2 . Knowing (−ω12 M + K)⃗u1 = 0 yields:

1 −1 u11 0
= (5.79)
−1 1 u21 0
u11 − u21 = 0, and − u11 + u21 = 0 (5.80)
Setting u11 = 1 results in u21 = 1 . The same processes can be performed for ⃗u2 to show that
if we set u12 = 1, u22 = −1. Therefore, the mode shapes can be expressed as:

u11 u12 1 1
⃗u1 ⃗u2 = = (5.81)
u21 u22 1 −1
The displacement of the masses as a function of time and the general mode shape plots are
graphically represented in figure 5.9. In the 2-DOF system considered here, the second mode
shape has a spot at the center of the middle spring that does not move (i.e. has zero displace-
ment). This point is called a node. Nodes correspond to points in the mode shape where the
126
Open Vibrations 5.3 Explicit method for Solving 2-DOF Systems
displacement is always zero. Furthermore, the displacement of the node points remain zero at
all times, as diagrammed in the top-right of figure 5.9.
Figure 5.9: Modes of vibration for the system shown in figure 5.8 showing the: (a) first mode;
and (b) second mode.
5.3 Explicit method for Solving 2-DOF Systems
Figure 5.10: Forced 2-DOF damped system showing: (a) system, and; (b) FBD.
127
As before, we can use explicit methods for solving multiple degree of freedom problems. Cramer’s
rule is an explicit formula for the solution of a system of linear equations with as many equations
as unknowns. Cramer’s rule is valid whenever the system has a unique solution and can be used as
a more generalized approach to solving for the temporal solution to a 2-DOF. Consider the 2-DOF
systems shown in figure 5.10, where x2 displaces more than x1 . The two coupled equations of
motion are expressed as:
m1 ẍ1 + (c1 + c2 )ẋ1 − c2 ẋ2 + (k1 + k2 )x1 − k2 x2 = 0 (5.82)
m2 ẍ2 + (c2 + c3 )ẋ2 − c2 ẋ1 + (k2 + k3 )x2 − k2 x1 = 0
As before, taking the Laplace of the EOM (while ignoring the initial conditions) changes the
equation from the temporal domain to the complex s-plane. This yields:
m1 s2 X1 (s) + (c1 + c2 )sX1 (s) − c2 sX2 (s) + (k1 + k2 )X1 (s) − k2 X2 (s) = F1 (s) (5.83)
m2 s2 X2 (s) + (c2 + c3 )sX2 (s) − c2 sX1 (s) + (k2 + k3 )X2 (s) − k2 X1 (s) = F2 (s)
these equations can be rearranged in terms of X1 and X2 as follows:
[m1 s2 + (c1 + c2 )s + (k1 + k2 )]X1 (s) − [c2 s + k2 ]X2 (s) = F1 (s) (5.84)
[m2 s2 + (c2 + c3 )s + (k2 + k3 )]X2 (s) − [c2 s + k2 ]X1 (s) = F2 (s)
These equations show two linear equations in terms of X1 and X2 that can be solved for using
Cramer’s rule, resulting in the expression:
D1 (s)
X1 (s) = (5.85)
D(s)
D2 (s)
X2 (s) =
D(s)
where:
F1 (s) −(c2 s + k2 )
D1 = 2 (5.86)
F2 (s) m2 s + (c2 + c3 )s + (k2 + k3 )
= [m2 s2 X2 (s) + (c2 + c3 )s + (k2 + k3 )]F1 (s) + (c2 s + k2 )F2 (s)
m1 s2 + (c1 + c2 )s + (k1 + k2 ) F1 (s)

D2 = (5.87)
−(c2 s + k2 ) F2 (s)
= [m1 s2 + (c1 + c2 )s + (k1 + k2 )]F2 (s) + (c2 s + k2 )F1 (s)
m1 s2 + (c1 + c2 )s + (k1 + k2 ) m2 s2 + (c2 + c3 )s + (k2 + k3 )

D= (5.88)
−(c2 s + k2 ) −(c2 s + k2 )
= m1 m2 s4 + [m2 (c1 + c3 ) + m1 (c2 + c3 )]s3
+[m2 (k1 + k2 ) + m1 (k2 + k3 ) + c1 c2 + c2 c3 + c3 c1 ]s2
+[(k1 + k2 )(c2 + c3 ) + c1 k2 + c1 k3 − c2 k2 + c2 k3 ]s
+(k1 k2 + k2 k3 + k3 k1 )
128
The denominator, D(s) is a 4th polynomial in s and is the characteristic polynomial of the
system. The system is considered a 4th order system because the characteristic polynomial of the
system is of order 4.
Vibration Case Study 5.2 Multi-span concrete bridges like the Trigno V bridge (figure 5.11)
over the Trigno river in Italy have repeating segments that make up the bridge decks. The
structural components of the segmented bridge decks are separate components sitting on
bearing pads and piers where the only connecting material between decks is the overlay that
is added to provide a contentious road surface. This configuration forms what is known as a
partially-connected bridge decka .
Figure 5.11: The Trigno V bridge over the Trigno river that carries SS650 north of Trivento
Italy is made up of seven repeating concrete bridge deckb .
This system can be modeled as a multi-degree of freedom problem. However, the chal-
lenge is that with so many nearly identical bridge components, the natural frequencies of
each bridge deck will be close, but not identical. This results in a clustering of natural
frequencies as shown in figure 5.12 where the frequencies of the 1st and 2nd modes of the
various bridge deck components are clustered in groups and therefore hard to distinguish.
Moreover, obtaining the characteristic structural dynamics of any particular bridge deck
section would be difficult as the decks are coupled through the pavement overlay making it
challenging to isolate the dynamic measurements of just one bridge section.
129
Open Vibrations 5.4 Eigenvalue-based Solution for Natural Frequencies and Mode Shapes
Figure 5.12: Measured acceleration signal for the bridge showing the estimated 1st and 2nd
frequencies for the obtained using the method outlined in Tomassini et al.a .
a Tomassini E., GarcÃa-MacÃas E., Reynders E., Ubertini F., Modal analysis for damage identification of
partially continuous multi-span bridges, Journal of Physics: Conference Series, Eurodyn 2023: XII Interna-
tional Conference on Structural Dynamics (2023).
b Imagery 2003 Google and Maxar Technologies used in accordance with their general guidelines on sharing
(2003) and likely falls under free use in the U.S.
5.4 Eigenvalue-based Solution for Natural Frequencies and Mode Shapes

The process of calculating the mode shapes presented in section 5.2 is long and tedious. Therefore,
methods that can be easily deployed on computers are of great interest to the practitioner. An
eigenvalue-based solution that takes advantage of the symmetry in the M and K matrices and can
be easily implemented on a computer is discussed in this section.
Review 5.2 Eigenvalues and Eigenvectors

In linear algebra, eigenvalues (λ ) and eigenvectors (⃗v) are concepts that appear promi-
nently in the analysis of linear transformations. By definition, if⃗v is a vector (in vector space
V over a field F) and T is a linear transformation into itself, then ⃗v is an eigenvector of T if
T (⃗v) is a scalar multiple of ⃗v:
T (⃗v) = λ⃗v (5.89)
where λ is a scalar in the field F, known as the eigenvalue associated with the eigenvector⃗v.
If the linear transformation is expressed in the form of an n × n matrix A, then the eigenvalue
equation for a linear transformation above can be rewritten as the matrix multiplication
Av = λ v (5.90)
where v is a n × 1 matrix of the eigenvectors. For the matrix A, eigenvalues and eigenvectors
130
can be used to decompose the matrix.
Figure 5.13: Matrix A acts by stretching the vector v, not changing its direction, so v is an
eigenvector of A.
The generalized eigenvalue problem is an important formulation for the study of vibra-
tions and is written as
Av = ⃗λ Bv (5.91)
where A and B are real matrices. As written, this expression maps a general space A into
B using ⃗λ and v. In the study of vibrations, the general generalized eigenvalue problem is
used to link mass (M) and stiffness (k) matrices such that
Kv = ⃗λ Mv (5.92)
5.4.1 Deriving the Eigenvalue-based Solution

To derive an eigenvalue-based solution for calculating the natural frequencies and mode shapes
in a computationally efficient way, we need to merge our mass and stiffness into one expression;
termed the mass normalized stiffness K e matrix that we mathematically define later. First, let us
consider that the vast majority of mass (M) and stiffness (K) matrices are symmetric and positive
definite due to the physical meaning of these matrices. Therefore, M can be factored into two
terms using the Cholesky decomposition:
M = LLT (5.93)
Review 5.3 Cholesky Decomposition The Cholesky decomposition of a real positive-definite

matrix A is a decomposition of the form:
A = LLT (5.94)
131
where L is the lower triangular matrix of A. A matrix is positive definite if the scalar xT Ax
is positive for any non-zero vector x comprised of real numbers:
xT Ax > 0 (5.95)
For the unique case diagonal mass matrices (all the mass values lie along the diagonal of the
matrix) the Cholesky decomposition (L) is defined as:
√
1/2 m1 0
L=M = √ (5.96)
0 m2
While a special case that is not always true, it is a commonly encountered mass matrix formulation
due to the nature of mass matrices. Moreover, the example considered within this test all consists
of a diagonal mass matrix. For the special case diagonal mass matrices, equation 5.96 factors into:
M = M 1/2 M 1/2 (5.97)
Moreover, the inverse of the diagonal matrix (M 1/2 ) is denoted as M −1/2 and defined as:
" 1 #
√ 0
−1 −1/2 m1
L =M = (5.98)
0 √1
m2
Now, let us consider the previously derived EOM for an undamped 2-DOF system:
M⃗x¨ + K⃗x = 0 (5.99)
This expression can be transformed into a symmetric eigenvalue problem, allowing us to leverage
the strengths of symmetric eigenvalue mathematics and computer solvers. To solve the perform this
transform, we set ⃗x = M −1/2⃗q and multiply the equation by M −1/2 such that the EOM becomes:
M −1/2 MM −1/2⃗q¨ + M −1/2 KM −1/2⃗q = 0 (5.100)
As M −1/2 MM −1/2 is equal to the identity matrix I and defining M −1/2 KM −1/2 as the mass nor-
malized stiffness K
e yields the simplified expression:
I⃗q¨ + K⃗
eq = 0 (5.101)
where Ke = M −1/2 KM −1/2 is equivalent to the expression k/m from the 1-DOF system as the are
both mass-normalized stiffness values.
As before, a solution is found by assuming a solution, taking the derivatives of the solution,
and substituting it into the EOM. Following these steps and assuming a solution of:
⃗q = ve jωt (5.102)
where v in an n × n matrix for a system with n degrees of freedom. Adding this assumed solution
to the EOM results in the form:
−vω 2 e jωt + Kve
e jωt = 0 (5.103)
132
driving out the nonzero scaler e jωt and rearranging the above expression results in:
e = ω 2v
Kv (5.104)
Knowing that v ̸= 0, as a matrix of zeros would mean no motion is present in the system, this
equation can be expressed in a typical eigenvalue formulation:
e = λv
Kv (5.105)
where v is a column matrix made up of the eigenvectors (v = [v1 , v1 , · · · , vn ]) and λ is a square
matrix with eigenvalues on the diagonal. As K e is symmetric, this is a symmetric eigenvalue prob-
lem.
An important attribute of eigenvectors to note is that the eigenvectors only encode information
about the direction of the transformation while information on the magnitude is captured by the
eigenvalue. Therefore, different values within an eigenvector may be used to represent the same
direction. A challenge for the entry-level practitioner is that different software systems may re-
turn different eigenvectors for the same problem. For example, MATLABa returns eigenvectors
such that the 2-norm of each is 1. However, when solved symbolicallyb non-normalized eigen-
vectors are returned. Various other engineering-focused applications may return normalized or
non-normalized eigenvectorsc . Therefore, it is helpful for practitioners to normalize computed
eigenvectors to unit norm eigenvectors to allow for comparison between different computational
tools.
Review 5.4 Vector Norms

The Euclidean norm of a vector (also termed as 2-norm, Euclidean length or the vector
magnitude) is defined as: s
n √
||v|| = ∑ (v2i ) = vT v (5.106)
i=1
If ||v|| = 1 it is a “unit norm”. If ||v|| is not a a unit norm vector, in can be converted to one
in by applying a scalar α such that such that αv = 1. In general, a nonzero vector v of any
length can be normalized to vnormalized using the following expression:
1
vnormalized = √ v (5.107)
vT v
Example 5.3 Normalize the vector ⃗v1 = [1/3 1]T .
Solution: p √
First, let’s check the Euclidean norm of⃗v1 , this is 12 + 1/32 = 1.11 = 1.05; therefore,
the unit vector is not unit norm.
To normalize the vector⃗v1 , a scalar (α) is calculated to make αv = 1. Therefore, following
a MATLAB 2023a eig function.

b MATLAB 2023a Symbolic Math Toolbox.
c LAPACK.
133
the definition of an orthogonal vector:
(α⃗v1 )T (α⃗v1 ) = 1 (5.108)
or: 1
α[1/3 1]α 3 = α 2 (1/9 + 1) = 1 (5.109)
1
√
Therefore, α = 3/ 10. Resulting in a normalized unit vector of
" 1 #
√
⃗v1−normalized = α⃗v1 = 10 (5.110)
√3
10
q √ √
as (1/ 10)2 + (3/ 10)2 = 1
√
Eigenvalues from the EOM are equal to ω 2 . Or more importantly, ωi = λi . Moreover, we
can relate the eigenvectors to the modes shapes by a factor of the mass matrix:
⃗u1 = M −1/2⃗v1 (5.111)
The important thing to remember is that the natural frequencies are the square root of the eigen-
values and the mode shapes are related to the eigenvectors through the mass matrix. Expanding on
equation 5.111, one can go from the mode shapes to the eigenvector through:
⃗v1 = M 1/2⃗u1 (5.112)
therefore, it can be seen that the eigenvectors and mode shapes are related through the mass nor-
malization process.
Example 5.4 Consider the system presented in example 5.1 and repeated below where m1 =9
kg, m2 =1 kg, k1 = 24 N/m, and k2 = 3 N/m with the initial conditions x10 = 1 mm, v10 = 0
mm/s, x20 = 0 mm, and v20 = 0 mm/s. Calculate the natural frequencies and the mode shapes
using the eigenvalue solution.
x2 .
Solution:
Writing the mass and stiffness matrix of the system as:
134

9 0
M= (5.113)
0 1
and
27 −3
K= (5.114)
−3 3
we can compute K
e using the following expression:
e = M −1/2 KM −1/2
K (5.115)
where KM −1/2 is computed first to maintain symmetry. This results in:
27 −3 13 0

−1/2 9 −3
KM = = (5.116)
−3 3 0 1 −1 3
and: 1
e = M −1/2 KM −1/2 = 3 0 9 −3 3 −1
K = (5.117)
0 1 −1 3 −1 3
Now a solution must be obtained for the eigenvalue problem:
e = λv
Kv (5.118)
While this can be obtained using computers for such a simple case it is more appropriate to
solve this expression by had. Therefore, the above expression can be rewritten as:
e − λ I)v = 0
(K (5.119)
However, as v ̸= 0 the matrix must be singular, the determinant of the (K

e − λ I) matrix must
equal zero. Or:
3 − λ −1
det =0 (5.120)
−1 3 − λ
This can be expanded to the characteristic equation:
λ 2 − 6λ + 8 = 0 (5.121)
with the roots (eigenvalues):

λ1 = 2 and λ2 = 4 (5.122)
√
Therefore, ω1 = 2 and ω2 = 2. These are the same values computed in example 5.1. The
eigenvectors for λ1 are computed as:
e − λ1 I)v = 0
(K (5.123)
or:
3 − 2 −1 v11 0
= (5.124)
−1 3 − 2 v21 0
135
This results in two dependent scalar equations:
v11 − v21 = 0 and − v11 + v21 = 0 (5.125)
That show us that v11 = v21 or ⃗v1 = [1 1]T . However, ⃗v1 so let’s apply a scalar α to normalize
it to a unit vector. Therefore, using (α⃗v1 )T (α⃗v1 ) = 1 we obtain:

1
α[1 1]α = α 2 (2) = 1 (5.126)
1
√
or α = 1/ 2. This allows us to normalize the vector knowing α⃗v1 = 1, resulting in a normal-
ized vector of:
1 1 0.71
α⃗v1 = √ = (5.127)
2 1 0.71
A similar process is followed for λ2 = 4 that leads to the normalized vector

1 −1 −0.71
α⃗v2 = √ = (5.128)
2 1 0.71
Lastly, the normalized eigenvectors can be converted to mode shapes using u = M −1/2 v. Re-
sulting in:
1
0 0.71 0.24
⃗u1 = 3 = (5.129)
0 1 0.71 0.71
and:
1
0 −1
3 −0.24
⃗u2 = = (5.130)
0 1 1 0.71
While these mode shapes are correct, it is common practice to report them normalized with a
maximum value of 1, therefore, ⃗u1 = [1/3 1] and ⃗u2 = [−1/3 1]. While unit-normalization
of the eigenvectors was not required in this example to obtain the right solution, it is good
practice. Note that these are the same mode shape vectors as computed in example 5.1.
Example 5.5 Eigenvalue Approach Solved using MATLAB
x2 .
136
Open Vibrations 5.5 Transfer-function Method
Using the Eigenvalue approach and MATLAB, determine the natural frequencies and
mode shapes of the system shown in figure 5.15, where m1 =9 kg, m2 =1 kg, k1 = 24 N/m,
and k2 = 3 N/m with the initial conditions x10 = 1 mm, v10 = 0 mm/s, x20 = 0 mm, and
v20 = 0 mm/s with the EOM expressed as:

m1 0 x¨1 k1 + k2 −k2 x1 0
+ = (5.131)
0 m2 x¨2 −k2 k2 x2 0
Solution:
Using the eigenvalue method, the following MATLAB code will solve for the natural
frequencies and mode shapes:
Listing 3: MATLAB code to find the frequencies and mode shapes of a 2-DOF system.
% define the M and K matrix
M = [9 0; 0 1]
K = [24+3 -3; -3 3]
% build the M inverse square - root and mass normalized stiffness matrix
M_inv_sqr = sqrt ( inv ( M ) )
K_mass_norm = M_inv_sqr * K * M_inv_sqr
% Using , K_mass_norm * v = lambda * v

[v , lambda ] = eig ( K_mass_norm )
% Solve for natural frequencies

omega_1 = sqrt ( lambda (1 ,1) )
omega_2 = sqrt ( lambda (2 ,2) )
% solve for the mode shapes

u_1 = M_inv_sqr * v (: ,1)
u_2 = M_inv_sqr * v (: ,2)
where ω1 =1.41 rad/sec and ω2 =2 rad/sec while ⃗u1 = [0.333 1] and ⃗u2 = [−0.333 1].
5.5 Transfer-function Method

As in 1-DOF systems, transfer functions can be used to solve for the temporal response of 2-DOF
systems under a variety of inputs. Again, the transfer function of a differential equation is defined
as the ratio of the Laplace transform of the output (system response) to the Laplace transform of
the input (forcing function). Moreover, the procedure for using the Laplace transform to solve the
equations of motion is the same and follows three steps:
1. Take the Laplace transform of both sides of the EOM while treating the time derivatives
symbolically.
2. Solve for X(s) in the obtained equation.
3. Apply the inverse transform x(t) = L [X(s)]−1
137
Open Vibrations 5.5 Transfer-function Method
Example 5.6 2-DOF System Subjected to Impulse

Two masses are connected through a spring, as shown in figure 5.16.
Figure 5.16: 2-DOF system subjected to an impulse showing: (a) system, and (b) FBD.
Solution:
Assuming that x1 displaces more than x2 , the equations of motion are:
m1 ẍ1 + k(x1 − x2 ) = F0 δ (t) (5.132)

m2 ẍ2 + k(x2 − x1 ) = 0
Taking the Laplace of both equations (step 1) yields:
(m1 s2 + k)X1 (s) − kX2 (s) = F0 (5.133)

−kX1 (s) + (m2 s2 + k)X2 (s) = 0
solving these two equations for X1 and X2 (step 2) results in:
F0 (m2 s2 + k)
X1 (s) = 2 (5.134)
s [m1 m2 s2 + k(m1 + m2 )]
F0 k
X2 (s) = 2
s [m1 m2 s2 + k(m1 + m2 )]
Using partial fractions, or a symbolic toolbox in MATLAB or Python, these expressions can
138
Open Vibrations 5.6 Multiple Degrees of Freedom
be rewritten as:

F0 1 m2 ω
X1 (s) = + (5.135)
m1 + m2 s2 ωm1 s2 + ω 2

F0 1 1 ω
X2 (s) = +
m1 + m2 s2 ω s2 + ω 2
where:
2 1 1
ω =k + (5.136)
m1 m2
Taking the inverse transform of the expressions for X1 (s) and X2 (s) (step 3) results in expres-
sions in the time domain and yields:

F0 m2
x1 (t) = t+ sin(ωt) (5.137)
m1 + m2 ωm1

F0 1
x2 (t) = t + sin(ωt)
m1 + m2 ω
Considering a system where F0 = 10 N, m1 = 1000 kg, m2 = 1000 kg, and k = 1500 N/m the
temporal response is annotated in figure 5.17.
Figure 5.17: Temporal response for the considered 2-DOF system subjected to an impact load.
5.6 Multiple Degrees of Freedom

This chapter introduces methodologies for the solving of systems with more than 2-DOF. As shown
in Chapter 5, 2-DOF systems can be solved analytically using 2 EOM coupled through their mode
shapes. However, these methods become tedious when extended to systems with damping or even
beyond the 2-DOF system.
139
Example 5.7 Multiple Mode Shapes
Figure 5.18: A Beechcraft Baron in flighta along with the Free-Free 3-DOF model simplified
as a mass-spring model.
Modeling the vibrations of a twin-engine airplane as a three-degree-of-freedom system can

be done as shown in figure 5.18 where the fuselage is a center mass, and the engines are point
masses suspended by cantilevers from the center mass. The stiffness of the wing corresponds
to the modulus of the wing E and its moment of inertia I. Assuming that m1 = m3 = 1m,
m2 = 3m, and k = 3EIl , the EOM can be written as:
       
1 0 0 ẍ1 3 −3 0 x1 0
EI 
m 0 3 0 ẍ2  + −3 6 −3 x2  = 0 (5.138)
l
0 0 1 ẍ3 0 −3 3 x3 0
calculate the natural frequencies and mode shapes of the system and plot the mode shapes in
relation to the considered Beechcraft Baron.
Solution using the mass normalized stiffness matrix K̃:

Solving for the modes shapes using the mass normalized stiffness matrix K̃ requires solv-
ing for M −1/2 and K̃ such that:
 
1 0 0
M −1/2 = 0 0.577 0 (5.139)
0 0 1
 
3 −1.732 0
3EI 
K̃ = M −1/2 KM −1/2 = −1.732 2 −1.732 (5.140)
l
0 −1.732 3
Then, the eigenvalue problem, formulated as:
K̃v = λ v (5.141)
140
is solved for the eigenvalues and normalized eigenvectors using a computer, resulting in:
λ1 = 0, λ2 = 1.73, λ3 = 2.23 (5.142)

     
0.447 −0.707 0.548
⃗v1 = 0.775 , ⃗v2 =  0.0  , ⃗v3 = −0.632 (5.143)
0.447 0.707 0.548
When the eigenvalue problem √ is solved using the mass normalized stiffness matrix K̃ the
natural frequencies are ωi = λi while the mode shapes are derived from the eigenvectors as
⃗u = M −1/2⃗v. This results in:
ω1 = 0 rad/sec, ω2 = 1.414 rad/sec, ω3 = 1.826 rad/sec (5.144)

     
0.447 −0.707 0.548
⃗u1 = 0.447 , ⃗u2 =  0.0  , ⃗u3 = −0.365 (5.145)
0.447 0.707 0.548
Next, normalizing the mode shapes by the max of the vector results in:
     
1 −1 1
⃗u1 = 1 , ⃗u2 =  0  , ⃗u3 = −0.667 (5.146)
1 1 1
these mode shapes can than be plotted as:
Figure 5.19: The unit vector normalized displacement of the mode shapes solved for using
the mass normalized stiffness matrix K̃.
141
Open Vibrations 5.7 Modal Analysis
Solution using the generalized eigenvalue approach:

The mode shapes can also be solved for using the generalized eigenvalue approach where
the eigenvalue problem is written as:
Kv = ⃗λ Mv (5.147)
solving for the eigenvalues and eigenvectors yields:
λ1 = 0, λ2 = 1.73, λ3 = 2.23 (5.148)

     
−0.577 0.707 0.639
v1 = −0.577 , v2 =  0.0  , v3 = −0.426 (5.149)
−0.577 −0.707 0.639
Note that the eigenvalues are the same as those solved for using the normalized stiffness
matrix approach while the eigenvectors appear to be different (mode 2). software tools and
computing languages do not all follow the same standards in terms of returning eigenvectors
as the information stored in the eigenvectors is just the direction of the transform. However,
mode 2 reported here is still correct as only the shape of the eigenvalue matters.
Figure 5.20: The unit vector normalized displacement of the mode shapes solved for using
the generalized eigenvalue approach.
a “A Beechcraft Baron 58 in flight” by San Diego Air & Space Museum Archives, Public Domain.
5.7 Modal Analysis

Modal analysis is the study of a system’s dynamic properties and is done in the frequency domain.
Consider a system with n degrees of motion, modal analysis allows for the uncoupling of the
EOM into n single-degree-of-freedom system (represented as 2nd -order DOF systems) where the
142
displacements of the masses are expressed as the linear summations of the normal modes of the
system. If every mode shape is considered, the solution is equivalent to the solution obtained from
the original nth -degree-of-freedom system.
Consider the generic multidegree-of-freedom system under external forces, expressed as:
M⃗x¨ + K⃗x = ⃗F (5.150)

where damping is not considered and the vector ⃗F is a set of deterministic inputs. To expand
this equation by modal analysis, the eigenvalue problem must first be solved. The generalized
eigenvalue problem is written at:
λ Mv = Kv (5.151)
For the nth -degree-of-freedom, the generalized eigenvalue problem can be simplified to:
ωi2 M⃗vi = K⃗vi (5.152)
Considering that the total displacement of the system, expressed as ⃗x(t) , is the summation of
the displacement of each of the noncontributing modes; assuming a linear system, the temporal
response of the system can be written as:
⃗x(t) = q1 (t)⃗v1 + q2 (t)⃗v2 + q3 (t)⃗v3 + · · · + qn (t)⃗vn (5.153)
where the time-dependent generalized scalars q1 (t), q1 (t), · · · , q1 (t) are the modal participation
coefficients (also called principal coordinates). Defining the modal matrix P as:
P = [⃗v1 ⃗v2 ⃗v3 · · · ⃗vn ] (5.154)
where ⃗v1 = M −1/2 u1 and are the orthonormal eigenvectors of K̃ (i.e. the mass normalized eigen-
vectors) and not the eigenvectors of the original system formation shown in equation 5.151. Note
that the modal matrix is made of of the eigenvectors of K̃ and not the mode shapes of the system;
for context see review 5.5.
Review 5.5 Modal Matrix A modal matrix is a mathematical concept taken from linear
algebra and not specific to vibrations or structural dynamics. This is why the modal matrix
does not contain the modes of the system but rather eigenvectors.
From linear algebra, the modal matrix B for the matrix A is a matrix of size n × n consist-
ing of the eigenvectors of A as columns in B. It is used in the definition of matrix similarity
such that
C = B−1 AB (5.155)
where C is a n × n diagonal matrix with the eigenvalues of A on the main diagonal (zeros
elsewhere). D is the spectral matrix of A. The eigenvalues must appear in the diagonal (top-
left to bottom-right in the same order as their corresponding eigenvectors are arranged in B
(column-wise left to right).
143
The linear combination of the normal modes (equation 5.153) can be more concisely written
as:
⃗x(t) = P⃗q(t) (5.156)
where ⃗q(t) = [q1 q2 q3 · · · qn ]T . Next, the relationship that relates the physical space to the modal
space for the acceleration component is written as:
¨ = P⃗q(t)
⃗x(t) ¨ (5.157)
combining these two terms results in the EOM that can be written as:
¨ + KP⃗q(t) = ⃗F
MP⃗q(t) (5.158)
To convert the EOM into the standard form, first the PT is multiplied through the equation as:
¨ + PT KP⃗q(t) = PT ⃗F
PT MP⃗q(t) (5.159)
If the modes are normalized, the following is true:
PT MP = I (5.160)
where I is the identity matrix and

 
↖ 0 0
PT KP =  0 ω 2 0  (5.161)
0 0 ↘
⃗
Next we define Q(t) ⃗
as vector of generalized forces in the modal space such that Q(t) = PT ⃗F.
This results in a EOM in the modal space expressed as:
 
↖ 0 0
¨ +  0 ω 2 0 ⃗q(t) = Q(t)
⃗q(t) ⃗ (5.162)
0 0 ↘
For a system with n degrees of freedom, this equation can be broken down into:
q̈i (t) + ωi2 qi (t) = Qi (t), i = 1, 2, · · · , n (5.163)
This expression is the same ODE that we have solved multiple times in this text. Therefore, we
know the solution to be:
q̇i0
qi (t) = qi0 cos(ωit) + sin(ωit) (5.164)
ωi
Lastly, to solve for a solution in the modal space, the initial conditions that were given in the phys-
ical space must be converted to the modal space. This can be done by generalizing the velocities
in terms of the modal matrix:
⃗q(0) = PT M⃗x(0) (5.165)
˙
⃗q(0) = PT M⃗v(0) (5.166)
144
Example 5.8 Free vibration response

Solve for the free vibration response of the 2-DOF presented in figure 5.21 using modal
analysis. Show the temporal response for the entire system for its first 20 seconds using the
full modal reconstruction and the reconstruction truncated to just include the first mode. Also,
plot the variations in the modal participation coefficients through time. Apply the parameters,
f1 = 0 N, f2 = 0 N, m1 = 10 kg, m2 = 1 kg, k1 = 30 N/m, k2 = 5 N/m, k3 = 1 N/m, x1 (0) = 1
mm, x2 (0) = 0 mm, v1 (0) = 0 mm/s, and v2 (0) = 0 mm/s.
Figure 5.21: Forced 2-DOF damped system showing: (a) system, and (b) FBD.
Solution:
The equations of motion that couple the system are:
m1 ẍ1 + (k1 + k2 )x1 − k2 x2 = 0 (5.167)

m2 ẍ2 + (k2 + k3 )x2 − k2 x1 = 0
In matrix form, these become:

m1 0 ẍ1 k1 + k2 −k2 x1 0
+ = (5.168)
0 m2 ẍ2 −k2 k2 + k3 x2 0

1 0
⃗x(0) = , ⃗v(0) =
0 0
The natural frequencies and mode shapes can then be obtained by solving the eigenvalue
problem. Setting up the generalized eigenvalue problem:
Kv = λ Mv (5.169)
145
and solving yields:

0.55
λ1 = 2.73, ⃗v1 = (5.170)
0.84

−0.15
λ2 = 6.76, ⃗v2 =
0.99
this is then related to the natural frequency and mode shapes as:

p 0.55
ω1 = λ1 = 1.65 rad/s, ⃗v1 =⃗v1 α1 = α (5.171)
0.84 1

p −0.15
ω2 = λ2 = 2.60 rad/s, ⃗v2 =⃗v2 α2 = α2 (5.172)
0.99
recall that the eigenvalues only contain information about the direction of the linear transform,
and therefore, their magnitudes are arbitrary. Therefore, they must be scaled proportionally
to each other. For this reason, the scalars α1 and α2 are added. By orthogonalizing the modal
vectors with respect to the mass matrix, the values of α1 and α2 are found as:
1 =⃗vT1 M⃗v1 (5.173)

2
10 0 0.55
1 = α1 0.55 0.84 (5.174)
0 1 0.84
and:
1 =⃗vT2 M⃗v2 (5.175)

2
10 0 −0.15
1 = α2 −0.15 0.99 (5.176)
0 1 0.99
therefore, α1 = 0.52 and α2 = 0.91. Now, applying the proper scaling values to the modal
vector, the modal matrix becomes:

0.284 −0.14
P = [⃗v1 ⃗v2 ] = (5.177)
0.43 0.900
Next, check that the normal modes in the modal matrix (P) are normalized, per equation 5.160.
This yields,
T 1 −2.775e − 16
P MP = ≈I (5.178)
−2.775e − 16 1
which is close enough to I. Considering that ⃗x(t) = P⃗q(t), the EOM for the system can be
expressed as: 2
¨ ω1 0 ⃗ 0
⃗q(t) + 2 ⃗q(t) = Q = (5.179)
0 ω2 0
rewriting this in scalar form for each modal coefficient yields:
q̈i (t) + ωi2 qi (t) = 0, i = 1, 2 (5.180)
146
where the solution for this ODE is:

q̇i0
qi (t) = qi0 cos(ωit) + sin(ωit) (5.181)
ωi
where qi0 and q̇i0 are the initial conditions in modal space. Therefore, the given initial condi-
tions must be transferred into modal space as:

T 0.284 0.43 10 0 1 2.85
⃗q(0) = P M⃗x(0) = = (5.182)
−0.14 0.900 0 1 0 −1.378

˙ T 0.28 0.43 10 0 0 0
⃗q(0) = P M⃗v(0) = = (5.183)
−0.14 0.90 0 1 0 0
therefore,
q1 (t) = 2.85 · cos(1.65t) (5.184)

q2 (t) = −1.34 · cos(2.6t)
converting back into the time domain is done knowing ⃗x(t) = P⃗q(t), therefore,

0.28 −0.14 2.85 · cos(1.65t)
⃗x(t) = P⃗q(t) = (5.185)
0.43 0.90 −1.38 · cos(2.6t)
This is further simplified into:
x1 (t) = 0.81 · cos(1.65t) + 0.19 · cos(2.6t) (5.186)

x2 (t) = 1.24 · cos(1.65t) − 1.24 · cos(2.6t)
These results are plotted in figure 5.22.
147
Figure 5.22: Temporal response for the 2-DOF reconstructed using just all the modal coordi-
nates.
Next, the truncated response can be computed by only considering the first mode response
for the system (i.e. ⃗x(t) = q1 (t)⃗v1 ). This is obtained as:

0.28 −0.14
⃗x(t) = P⃗qtruncated (t) = 2.85 · cos(1.65t) (5.187)
0.43 0.90
This is further simplified into:
x1 (t) = 0.81 · cos(1.65t) (5.188)

x2 (t) = 1.24 · cos(1.65t)
These results are plotted in figure 5.23. Note that this only considers the response of the
system that is a function of the first mode. Note that this captures some of the “general” idea
of the system while missing out on the finer points that the 2nd mode contributes.
148
Figure 5.23: Truncated temporal response for the 2-DOF reconstructed using just the first
modal coordinates.
Lastly, the participation of the two modes can be plotted from the time series response of
equation 5.185.
Figure 5.24: Modal participation coefficients.
5.8 Numerical Methods

Numerical methods can be used to solve the response of multi-degree of freedom system subjected
to forced vibrations. While not the most computationally efficient method, the EOM is an ODE
that can be solved directly while considering the initial directions to obtain the response of the
system.
149
Example 5.9 Directly Solving the ODE of the EOM for a 2-DOF system.
Consider the system presented in figure 5.25(a) where m1 =2 kg, m2 =1 kg, k1 = 20 N/m,
k2 = 10 N/m, c1 = 0.5 kg/s, and c2 = 1 kg/s; initially at rest. m1 is subjected to the ramp and
hold load shown in figure 5.25(b). Using MATLAB, solve the EOM for the temporal response
2-DOF system using a numerical ODE solver.
x2 .
Solution:
Assuming x1 < x2 , the matrix form of the system is

m1 0 x¨1 c1 + c2 −c2 ẋ1 k1 + k2 −k2 x1 R(t)
+ + = (5.189)
0 m2 x¨2 −c2 c2 ẋ2 −k2 k2 x2 0
where R(t) is the piecewise ramp function shown in figure 5.25(b). This expression can be
re-arranged to:
⃗x¨ = M −1 (Ft −C ·⃗x˙ − K ·⃗x) (5.190)
which is the format required by MATLAB’s ode45 solver. Thereafter, the code in listings 4
and 5 can be used to develop the results shown in figure 5.26.
Listing 4: MATLAB code for solving the EOM of the two-degree-of-freedom system.
% Time span for simulation
tspan = [0 , 10]; % Start time and end time
% Initial conditions [ x1 , x1 ' , x2 , x2 ']

in it ia l_c on di tio ns = [0 , 0 , 0 , 0];
% Use ode45 to solve the system of ODEs

[t , y ] = ode45 (@ equations_of_motion , tspan , i ni ti al_ co nd iti on s ) ;
% Extract displacements of masses
150
x1 = y (: , 1) ;
x2 = y (: , 3) ;
Listing 5: Functions for Matlab code.

% Equations of motion for the system
function [ dydt ] = e q ua t io n s _o f _m o t io n (t , y )
% Setup the system parameters

m1 =2; m2 =1; k1 =20; k2 =10; c1 =0.5; c2 =1;
% Build the Mass , Damping , and Stiffnes matrices

M = [ m1 , 0; 0 , m2 ];
C = [ c1 + c2 , - c2 ; -c2 , c2 ];
K = [ k1 + k2 , - k2 ; -k2 , k2 ];
% Unpack the state variables

x = y (1:2) ;
x_dot = y (3:4) ;
% Get the force excitation vector at time t

F_t = f o r c e _ e x c i t a t i o n _ v e c t o r ( t ) ;
% Equations of motion
x_dotdot = inv ( M ) * ( F_t - C * x_dot - K * x ) ;
% Pack the derivatives into the output vector dydt

dydt = [ x_dot ; x_dotdot ];
end
% Define the force excitation vector F ( t )

function F_t = f o r c e _ e x c i t a t i o n _ v e c t o r ( t )
if t <1 % Ramp load from 0 to 1 second

f1_t = t ;
else % constant load after 1 second
f1_t =1;
end
% Force vector , with no load on f2

F_t = [ f1_t ; 0];
end
151
Figure 5.26: Displacement response of the 2-DOF system.
152
Open Vibrations
Part II
Applied Topics
The tuned mass damper inside Taipei 101. The largest damper ball in the world, weights
660-metric-tons and consists of 41 circular steel plates that are 125 mm (4.92 in) thick. a
a Armand du Plessis, CC BY 3.0 <https://creativecommons.org/licenses/by/3.0>, via Wikimedia Commons
153
Open Vibrations
6 Vibration Control
Throughout this text, we have studied various aspects related to analyzing and modeling vibrating
systems. Therefore, it becomes prudent to look at methods for reducing or eliminating unwanted
vibrations. However, before vibrations in a system can be effectively reduced they must be better
understood in terms of their effects on the system under study. For this reason, this chapter first
introduces the vibration Nomograph, which is then followed by vibration isolation, absorption, and
active suppression.
6.1 Vibration Nomograph

There exist various methods and standards for measuring and describing acceptable levels of vibra-
tions in systems, these include ISO/AWI 2631 for the evaluation of human exposure to whole-body
vibrations and ISO 4866 for the measurement and effects of vibrations on structures. A common
way to present the acceptable limit of vibration is in a vibration nomograph. A vibration nomo-
graph is a simplified way to express the acceptable limits on a system while considering the dis-
placement, velocity, acceleration, and frequency of a system. A typical nomograph with various
limits is presented in figure 6.1.
A vibration nomograph is a logarithmic plot that allows us to easily express the relationships
between displacement, velocity, acceleration, and frequency of a system. The vibration nomo-
graph presented in figure 6.1 considers an undamped 1-DOF system with constant amplitude (A)
experiencing harmonic motion that can be modeled as:
x(t) = A sin(ωt) (6.1)
Therefore, the velocity and acceleration terms can be found by taking the derivatives of the dis-
placement expression to yield:
ẋ(t) = Aω cos(ωt) (6.2)
and:
ẍ(t) = −Aω 2 sin(ωt) (6.3)
These equations are converted from a circular frequency in rad/sec to a linear frequency ( f ) in Hz,
such that ω = 2π f . Therefore, equations 6.1-6.3 become:
x(t) = A sin(ωt) (6.4)
v(t) = ẋ(t) = 2π f A cos(ωt) (6.5)

2 2
a(t) = ẍ(t) = −4π f A sin(ωt) (6.6)
Thereafter, the maximum values for velocity vmax and acceleration amax are related to amplitude
through:
vmax = 2π f A (6.7)
amax = −4π 2 f 2 A = −2π f vmax (6.8)
by taking the natural log of both sides of equation 6.7 we obtain:
ln vmax = ln(2π f ) + ln A (6.9)
154
Open Vibrations 6.1 Vibration Nomograph
Figure 6.1: Vibration nomograph showing the acceptable limits of vibration for various applica-
tions.
155
Open Vibrations 6.1 Vibration Nomograph
doing the same for equation 6.8 leads to:
ln amax = − ln(2π f ) − ln vmax (6.10)
It can be seen that both of these expressions are linear.

The nomograph sets the x−axis as frequency in Hz and the y−axis as velocity in mm/s. Equa-
tion 6.9 tells us that For a constant amplitude of displacement (A), ln vmax is linearly proportional
to ln(2π f ), at a rate of 2π. As the x−axis in a nomograph is frequency, measured in Hz and
thereby accounting for the 2π, ln(2π f ) is a straight line with a positive slope of 1 with respect to
the frequency axis (i.e. x−axis). Therefore, a line on the nomograph that represents a constant
displacement is at a 45◦ angle from the x-axis.
For a constant value of velocity, (vmax ), equation 6.10 shows that acceleration (ln amax ) is lin-
early proportional to − ln(2π f ), at a rate of 2π. Again, as the x−axis in a nomograph is frequency,
measured in Hz, acceleration is represented by a straight line that varies with − ln(2π f ), therefore,
ln amax is a straight line with the slope of -1. This is also represented by a line of constant accelera-
tion set at a -45◦ angle from the x-axis. These equations are expressed in the vibration nomograph
plot of figure 6.1 where each point on the plot represents a specific sinusoidal (harmonic) vibration
for a 1-DOF system.
Vibration Case Study 6.1 Vibration mitigation can be achieved through better system-level
design. For example, early rubber tires consisted of uniform tread patterns, but as cars
became faster the noise, vibration, and harshness (NVH) generated by the tire would be
constrained around a single frequency; thereby amplifying NVH felt by the driver and pas-
sengers. This challenge led to the development of tires with irregular tread patterns to spread
the energy created at the road/tire interface out over a wider bandwidth of excitation; thereby
reducing NVH felt and heard by the passengers.
156
Open Vibrations 6.2 Vibration Isolation
Figure 6.2: Illustration for the 1935 US patent which proposed the use of irregular tread
patterns to control the pitch of road noisea .
a US patent number US2006197A, inventors Elliott S Ewart and Arthur W. Bull; Public Domain
6.2 Vibration Isolation

To mitigate vibrations in a system the ideal approach would be to limit the source of vibrations.
However, this is not always applicable to the system you are considering. Therefore, isolating the
157
Open Vibrations 6.3 Vibration Absorption
system from vibrations is the next best step. One approach to this is to design systems around limit-
ing the force and displacement transmissibility discussed prior; where both force and displacement
transmissibility are considered isolation problems.
One way to do this is to track the transmissibility ratio which is denoted as T.R. and defines the
ratio of the magnitude of the transmitted (FT ) to applied force (F0 ).
s
FT 1 + (2ζ r)2
T.R. = = (6.11)
F0 (1 − r2 )2 + (2ζ r)2
6.3 Vibration Absorption

Vibration absorbers, also termed dynamic vibration absorbers, are a class of mechanical devices
that seek to reduce unwanted vibrations in a system. In contrast to a traditional dash-pot style
damper, these systems seek to “redirect” the vibrations from the system to another mass connected
to the system. In this way, the main system is protected from the bandwidth of vibrations that
the vibration absorbers are tuned for. As the vibration absorbers must be tuned for the system,
it is generally limited to devices that operate at a fixed frequency like industrial equipment or
cables suspended in the air and subjected to wind loading. Figure 6.3 presents a Stockbridge and a
Dogbone damper designed to remove certain bandwidths of excitation from wind-excited cables.
Figure 6.3: Vibration absorbers deployed on wind excited cables showing: (a) a Stockbridge
damper on a high-power transmission linea , and; (b) a dogbone damper on a suspender cable
of a suspension bridgeb .
6.3.1 Vibration Absorption for Undamped Systems

Vibration absorbers are most often designed to shift the resonance frequency of the first mode of the
system away from the expected excitation frequency. This is done by adding an additional degree-
of-freedom in the form of a mass (the vibration absorber) connected to the system with a spring
a “Stockbridgedampers installed on high voltage power lines” by Badics CC BY-SA 3.0
b “Dogbone dampers on the road-support cables of the Severn Bridge” by Bassaar CC BY-SA 4.0
158
to alter the natural frequency of the combined system away from the original excitation frequency.
Dashpots may also be added in parallel to the spring element if additional energy dissipation is
needed beyond that provided by the original system.
Figure 6.4: A vibration absorber (m2 ) for mitigating unwanted dynamics in a device (m1 ).
The tuning of a 2-DOF system can be done by setting the displacement of the mass to be
controlled to zero and solving for the mass and stiffness of the vibration absorber. Consider the
system presented in figure 6.4, here m1 and k1 are the mass and stiffness of the system while m2
and k2 are the mass and stiffness of the vibration absorber. A good assumption to make when
designing a vibration absorber is that the mass of the vibration absorber should be between 1% and
5% of the mass of the system to be damped. Therefore, for this case let m1 = 20 kg, m2 = 1 kg, and
k1 = 20 kN. Assuming a sinusoidal input where F0 = 1 kN, the equations of motion are:
m1 x¨1 + k1 x1 + k2 (x1 − x2 ) = F0 sin(ωt) (6.12)
m2 x¨2 + k2 (x2 − x1 ) = 0 (6.13)

Assuming the temporal solution is of a harmonic form, the following is true:
xi (t) = Xi sin(ωt), i = 1, 2 (6.14)
using the transfer function approach and assuming no initial conditions, the following steady-state
solution can be obtained for m1 and m2 :
(k2 − m2 ω 2 )F0
X1 = (6.15)
(k1 + k2 − m1 ω 2 )(k2 − m2 ω 2 ) − k22
k2 F0
X2 = (6.16)
(k1 + k2 − m1 ω 2 )(k2 − m2 ω 2 ) − k22
p
Next, the natural frequency of m1 (ω1 ) can be solved for as ω1 = k1 /m1 . In order to eliminate
movement for m1 at a given driving frequency ω, the numerator of equation 6.15 should be set to
zero. Note that setting F0 to zero is a trivial solution and provides no benefit to the system in terms
of vibration control. Therefore:
k2 = m2 ω 2 (6.17)
159
note that this will force

pthe frequency of the tuned vibration absorber to match that of the system,
therefore ω1 = ω2 = k2 /m2 . Next, normalizing the input force F0 by the stiffness of the main
system k1 yields:
F0
δst = (6.18)
k1
using this term, equations 6.15 and 6.16 can be rearranged as:
2
X1 1 − ωω2
=h ih i (6.19)
ω 2 ω 2
δst 1 + kk21 − − kk21

ω1 1− ω2
X2 1
=h 2
ih i (6.20)
ω 2
δst 1 + kk21 − − kk21
ω

ω1 1− ω2
Figure 6.5 reports the normalized displacement of the system over a frequency range for the system
with and without a vibration absorber. Note that at ω = 1 the original system is in resonance
while the system with the vibration absorber has no displacement. However, no system is without
compromise. From equation 6.20 it can be seen that at ω = ω1 = ω2 the second mass needs a
displacement equal to:
k1 F0
X2 = − δst = − (6.21)
k2 k2
or 1 m using the given parameters. Therefore, the mass and stiffness values of the vibration ab-
sorber should be selected based on the allowable travel of the vibration absorber (i.e. X2 ), among
other factors. Moreover, from this equation it can be seen the force exerted by the second mass
operates in the direction opposite the original force (−F0 − k2 X2 ), thereby canceling it. Lastly, note
that the addition of the vibration absorber creates two resonate frequencies of the system, termed
Ω1 and Ω2 . These resonate frequencies represent the roots of the system and care should be taken
to limit the time the system spends at these frequencies (i.e. on startup). The locations of these
roots can be solved analytically by setting the denominators of equation 6.19 to zero.
Figure 6.5: Frequency response of the undamped system with and without the vibration absorber.
160
6.3.2 Vibration Absorption for Damped Systems

As shown in section 6.3.1 and figure 6.5 in particular, vibration absorption for undamped systems
only shift the resonance response from part of the spectrum to another. However, it is communally
desired to limit the resonance response of the system while also absorbing vibration energy at
one specific frequency. While not derived here, the frequency response of a damped vibration
absorber like that shown in figure 6.4 (but with the addition of a damper) can be expressed as the
dimensionless amplitude of the response of the primary mass:
v
2 + r2 − β 2 2
u
X1 u (2ζ r)
=t 2 2 2 (6.22)
δst 2ζ r r2 − 1 + µr2 + µr2 β 2 − (r2 − 1)(r2 − β 2 )
This expression requires four design variables [µ, β , r, ζ ] to be set by the practitioner. First, ω1
is the natural frequency of the primary mass onto which the vibration absorber is attached, and is
defined as: p
ω1 = k1 /m1 (6.23)
p
which leads to a similar expression for the frequency of the vibration absorber ω2 = k2 /m2 .
Next, we define the ratio of natural frequencies β = ω2 /ω1 and the ratio of the masses µ = m2 /m1 .
However, to function as a vibration absorber, it is often desired to set β = 1. Next, we build an
expression for the “mixed damping ratio”:
c
ζ= (6.24)
2m2 ω1
where c is the damping value of the added damper. Again, we define r = ω/ω1 to create a variable
of the driving frequency to the frequency of the system. Figure 6.6 shows how the selection of the
four design variables [µ, β , r, ζ ] results in different spreads of the response of the primary system
mass. As the mixed damping ratio increases, the response of the system converges on that of a
system without a vibration absorber as the system mass and vibration absorber are more tightly
coupled.
Figure 6.6: Frequency response of a 1-DOF system with various damped vibration absorbers.
161
Open Vibrations 6.4 Active Vibration Suppression
Vibration Case Study 6.2 Tuned mass dampers for vibration mitigation in tall structures
Figure 6.7: Illustration of Taipei 101’s main tuned mass damper. a
A tuned mass dampers (TMD), also known as a harmonic absorber or seismic damper,
are devices that are designed into a structure to mitigate structural vibration. The mass
is typically a block of steel or concrete and is mounted on suspended cables to create a
pendulum and damped in relation to the structure. By tuning the oscillating frequency of the
damping system to be near the same natural frequency of the structure, energy is transferred
to the mass and extracted through the dampers. Thereby reducing vibration which prevents
discomfort or damage. While discussed here in the context of tall buildings, tuned mass
dampers are also frequently found in automobile components and power transmission lines.
a Someformofhuman, CC BY-SA 4.0 <https://creativecommons.org/licenses/by-sa/4.0>, via Wikimedia
Commons
6.4 Active Vibration Suppression

Vibrations in systems can be mitigated through a number of active systems; typically it’s easiest
to consider this as an actuator that adds energy to the system at the correct time to cancel out
vibrations.
6.4.1 Metrics for Vibration Control

There are various performance indicators that one can use to judge the performance of an active
control scheme. They depend on the system order (1st , 2nd ) and the excitation experienced by
the system. For simplicity in this introductory text, we will define four performance indicators
162
subjected to a step response, each shown in figure 6.4.1. The performance indicators are:
• peak time (tp ) is the time to the first peak.
• peak value (xp ) is the maximum value experiences by the system
• settling time (ts ) is the time it takes the system to get within an error (±ε%) of the steady-
state displacement (xss ) and stay there.
• max percentage overshoot (Mp ) is defined as M p = (x p /xss − 1) · 100.
6.4.2 Position-Derivative (PD) Control

Active vibration control adds energy to the system in order to mitigate the vibrations in the system.
As depicted in figure 6.8(a), an active vibration control system requires a sensor to acquire data
from the system, control hardware, and algorithms to process this data, and an actuator to exert
physical control on the system. These systems together are called a feedback loop, as a movement
in the mass results in a controlled force ( fu ) being exerted on the system. This control force is
diagrammed in figure 6.8(b).
163
Figure 6.8: Active vibration control system showing: (a) the system with a feedback loop that
takes a signal from the sensor, converts it to a control signal, and drives the actuator; and (b) the
free body diagram.
Adding the control force to the EOM for the 1-DOF system presented in figure 6.8 results in:
mẍ + cẋ + kx = F(t) = f + fu (6.25)
A common method for providing control for vibration suppression is called position and derivative
control or PD-control. A PD-controller is a state-variable feedback controller as it uses velocity and
displacement obtained from the measured acceleration, assuming that the acceleration is properly
integrated. PD-control measures the position and velocity of the mass and uses these to compute
the control force needed to mitigate the vibration to an acceptable level. A simple way to code a
PD-controller is to provide a control force proportional to the displacement velocity (derivative of
displacement) of the mass such that:
fu = −g1 x − g2 ẋ (6.26)
where g1 and g2 are the proportional gains of the systems. The control gains can be constants
determined by the designer or variables updated through time by an algorithm. Here we will
consider the gains to be constant, therefore, the EOM for the closed-loop system in figure 6.8
becomes:
mẍ + (c + g2 )ẋ + (k + g1 )x = F(t) = f (6.27)
This formulation lets g1 act as additional stiffness while g2 acts as additional damping. This closed-
loop EOM can be used to solve for the effective natural frequency of the system, given by:
r
k + g1
ωn = (6.28)
m
and the effective damping ratio of the system
c + g2
ζ= p (6.29)
2 m(k + g1 )
164
Review 6.1 Continuous control systems have been widely used for centuries. For example,
consider that the centrifugal governor which uses spinning weights was used by Christiaan
Huygens in the 1600s in the Netherlands to regulate the gap between millstones in windmills
or by James Watt who famously linked a stem regulator to a centrifugal governor to control
steam turbines.
Arguably, the Russian American engineer Nicolas Minorsky was the first to develop the
theoretical analysis for the three-term control we now call PID. This was done in 1922 while
he was researching and designing automatic ship steering for the US Navy. He based his
work on watching how a ship’s helmsman responds to wave loading on a ship, with a delayed
input to the helm that not only considered the current ship course but also past errors and
the desired rate of change for the ship. For a helmsman, the goal is stability, not absolute
control, which simplifies how one thinks about the challenge of control.
Figure 6.9: Historical perspective of PID control showing: (a) Portrait of Nicolas Minorsky
a and (b) the battleship USS New Mexico (BB-40) of the United States Navy which was the
first to implement PID control in its steering b .

a Peter Minorsky, grandson of Nicolas Minorsky, CC BY-SA 1.0
b U.S. Navy, Public domain, via Wikimedia Commons
6.4.3 Proportional-Integral-Derivative (PID) Control

Proportional-Integral-Derivative (PID) Control is a three-term controller that employs feedback
that is widely used in continuous control systems, including for the control of structural systems. A
PID controller seeks to minimize the measured error value e(t) between a desired setpoint (SP) and
a measured process variable (PV) by applying corrections based on the proportional (P), integral
(I), and derivative (D) terms (denoted P, I, and D respectively), from which it gets its name.
165
Figure 6.10: Generalized PID controller for a system with feedback, where r(t) is the desired
setpoint (SP) and y(t) is the measured process value (PV).
The overall control equation is defined as

Z t
de(t)
u(t) = Kp e(t) + Ki e(τ) dτ + Kd (6.30)
0 dt
where Kp , Ki , and Kd are non-negative coefficients for the proportional, integral, and derivative
terms, respectively. The PID controller is diagrammed in figure 6.10 for a system with feedback
control, such as that shown in figure 6.8. Moreover, in the Laplace-derived s domain, the transfer
function of the PID controller is defined as
Ki
L [s] = Kp + + Kd s (6.31)
s
where s is the complex frequency. A temporal response for the 1-DOF shown in figure 6.8 when
controlled with a PID controller is reported in figure 6.11.
Figure 6.11: System response for a 1-DOF system controlled with a PID.
166
Open Vibrations
7 Experimental Vibrations
Experimental testing requires the practitioner to understand the basics of testing hardware and
digital signal processing. An understanding of how to acquire and process vibration data is key to
being able to apply one’s knowledge of vibrations to real-world systems.
Vibration Case Study 7.1 On August 14, 2018, the Ponte Morandi viaduct in Genoa Italy
collapsed, killing 43 and displacing hundreds of people from their homes. The Morandi
viaduct was a cable-stayed bridge with uniquely few stays, typically only two per span. The
Stays were a hybrid of steel cables overlaid with concrete. The concrete overlay made the
direct inspection of the stays impossible.
While the exact cause may never be known, is suspected that one of the stay cables
within the concrete failed due to corrosion and poor maintenance causing a bridge with very
little redundancy in its design to faila .
In 2017, researchers from the Polytechnic University of Milan instrumented and stud-
ied the vibration characteristics of the bridge and noted that the modal frequencies of the
stays on pillar 9 (the one that collapsed) were more than 10% different than other stays on
the bridge. While it’s always hard to draw conclusions from one test, comparing modal
frequencies between two similar structures can be useful for tracking damage.
Figure 7.1: The Ponte Morandi bridge, showing the bridge: a) before the collapseb , and; b)
after collapsec .
a Rymsza, Janusz. “Causes of the Collapse of the Polcevera Viaduct in Genoa, Italy.” Applied Sciences 11,
no. 17 (2021): 8098. https://doi.org/10.3390/app11178098.

b Davide Papalini, CC BY-SA 3.0 <https://creativecommons.org/licenses/by-sa/3.0>, via Wikimedia Com-
mons
c Michele Ferraris, CC BY-SA 4.0 <https://creativecommons.org/licenses/by-sa/4.0>, via Wikimedia
Commons
167
Open Vibrations 7.1 Sensing and Data Acquisition
Figure 7.2: Key components for performing experimental modal analysis.
The measurement of vibrating systems requires specialized hardware. While a variety of ven-
dors sell vibration measurement systems in a number of form factors, the general hardware re-
quirements remain constant. The basic hardware requirements are: Exciter - A system to provide
a measurable input to the system, Transducers - Sensors used for converting the mechanical move-
ments of the structure to signals, and data acquisition - Hardware for digitizing the signal generated
by the transducers. Figure 7.2 shows some of the key systems required for vibration testing and
their interactions.
7.1 Sensing and Data Acquisition

The output of a vibrating system is measured through a combination of sensors and data acquisition
systems.
7.1.1 Accelerometers
Accelerometers are by far the most common type of sensor used for measuring vibrations. Various
types of Accelerometers exist, including Micro-electromechanical systems (MEMS) based systems
that are commonly found in cell phones, piezo-resistive-based systems used for high acceleration
loading (greater than 10,000 gn ), or piezo-electric sensors commonly deployed in industrial set-
tings. In terms of dedicated vibration testing, piezo-electric sensors are the most common sensor
system.
Piezo-electric sensors use a piezo-electric material to convert small movements into a small
electrical charge (measured in Coulomb) in and out of the piezo-electric material. On its own,
the signal encoded by this charge is hard to measure and susceptible to electromagnetic noise
if run over medium to long wires. Therefore, amplifiers are added to the sensors to assist in
transferring this signal back to the data acquisition; thereby creating Integrated Electronics Piezo-
Electric (IEPE) sensors. Figure 7.3(a) shows the cross section of a common IEPE sensor. Through
168
tuning the piezo-electric material and packaging, IEPE sensors can be made to measure a variety
of applications (figure 7.3(b)). Table 1 reports specifications for five different IEPE sensors that
are used to measure a range of applications from the structural motion of buildings to packages
subjected to high-shock loading (e.g., missals, plane crashes).
Figure 7.3: Integrated Electronics Piezo-Electric (IEPE) accelerometers, showing: (a) the cross-
section of a typical IEPE) accelerometer with key components annotated, and; (b) selection of
IEPE accelerometers for various applications.
Table 1: Specifications for various IEPE accelerometers.

specifications accelerometers
model number PCB 393B31 PCB 393B04 PCB 352C67 PCB 352A21 PCB 352A92
Sensitivity(± 10 %) 10.0 V/g 1000 mV/g 100 mV/g 10 mV/g 0.25 mV/g
Measurement Range ± 0.5 g pk ± 5 g pk ± 50 g pk ± 500 g pk ± 20 kg pk
Frequency Range(± 5 %) 0.1 to 200 Hz 0.06 to 450 Hz 0.5 to 10 kHz 1.0 to 10 kHz 1.2 to 10 kHz
Resonant Frequency >700 Hz >2.5 kHz >35 kHz >50 kHz >100 kHz
Non-Linearity ≤1% ≤1% ≤1% ≤1%
Transverse Sensitivity ≤5% ≤5% ≤5% ≤5 %
169
Figure 7.4: Graph of a generic frequency roll-off where the cutoff frequency is at -3.01 dBa .
The frequency range of a sensor reports the frequency of vibration the sensor is designed to
acquire. The upper limit is defined as the cutoff frequency of the sensor which is typically de-
fined as the frequency at which the sensor experiences -3.0 dB (relative unit) of signal power loss.
However, at -3.0 dB the power of the measured signal is about half as strong as that of the ideal
signal. Figure 7.4 shows a generic frequency roll-off chart. For the practitioner, it is important to
note that the signal starts to die off well before the cutoff frequency. Therefore, it is important to
be cognizant of a sensor’s frequency response when trying to obtain measurements of signals near
the upper end of an accelerometer.
7.1.2 Data Acquisition

Acquiring data from sensors was once a challenging task that has been greatly simplified by the
development of modern hardware and software systems that are far too numerous to list here. One
important item to note here is that IEPE sensors require IEPE data acquisition system, as shown
in figure 7.5. As IEPE sensors generate very small signals, they require amplification before they
can be accurately measured. Therefore, an IEPE data acquisition system consists of a front-end
amplifier that is specifically designed to amplify these small signals, a signal conditioning circuit
that filters and shapes the signals, and an analog-to-digital converter (ADC) that converts the analog
signals to digital signals to be analyzed by a computer.
a PDerivative work: KrishnavedalaOriginal: Omegatron, CC BY-SA 3.0
170
Open Vibrations 7.2 Controlled Force Excitation
Figure 7.5: Integrated Electronics Piezo-Electric (IEPE)-based measurement system showing the:
(a) simplified circuit schematica ; and (b) IEPE data acquisition systems in various form factors.
7.2 Controlled Force Excitation

To develop an accurate understanding of vibrating systems, it is important to understand the energy
input. Moreover, it is important to test systems under the vibration inputs they may encounter
during transportation, operation, or storage as this can help identify potential weaknesses or design
flaws in the product.
7.2.1 Modal Hammers
Figure 7.6: Model hammer, showing: (a) instrumented hammer with interchangeable tips, and; (b)
the temporal and frequency response from various tips.
A modal hammer is an instrumented hammer used to impart a measured impact force into the
structure. The frequency range of the resulting vibrations is determined by the duration of the
a “IEPE sensor connected to the input of an instrument” by JanBurg CC BY-SA 4.0
171
Open Vibrations 7.3 Digital Signal Processing
impact. A shorter impact duration leads to higher frequencies being excited. To achieve varying
frequency bandwidths with the same amount of impact energy, special hammer tips of different
stiffnesses can be utilized. Figure 7.6 shows a model hammer with interchangeable tips and the
responses generated by the tips. In general, there are three factors to be considered when selecting
the proper modal hammer.
Factor 1: Frequency bandwidth. Softer hammer tips result in longer pulse durations and nar-
rower frequency bandwidths, while harder tips lead to shorter pulse durations and broader fre-
quency bandwidths. However, when using a hard tip, the power spectral density of the excitation
may be insufficient to excite vibration modes in the system. In such cases, increasing the impact
force by swinging the hammer harder or adding a head extender may be attempted, but there is
a risk of overloading the IEPE force transducer. An alternative solution is to switch to a hammer
model with a larger measurement range or use a softer tip to concentrate the impact energy at
lower frequencies. Moreover, the duration of the pulse may also be affected by the stiffness of the
specimen being impacted.
Factor 2: Energy of the impact. The diverse shapes, masses, and material properties (e.g., stiff-
ness or damping) of objects being tested necessitate a range of force pulses with varying parameters
to achieve optimal excitation. Compact objects typically have higher resonance frequencies and
require less energy to be excited than larger objects. As a result, a short-duration force pulse can
be generated using small or medium-sized hammers. In contrast, larger structures require higher-
energy impacts, which are typically concentrated in a low-frequency bandwidth. Modal hammers
are available in various masses with measurement ranges ranging from 100 N to 20 kN, allowing
practitioners to deliver force pulses with different energies without requiring large swings. As large
swings make it difficult to control the force and angle of the hammer tip’s impact on the structure.
Factor 3: Tests repeatability. Hammer impacts performed by the practitioner during testing may
vary in terms of impact energy, the frequency bandwidth of excited vibrations, and the angle of
impact. Therefore, it is common practice to average multiple results obtained during testing to
develop high-quality and consistent data.
7.2.2 Electrodynamic Shakers

An electrodynamic shaker can generate a wide range of frequencies and amplitudes that simulate
different vibration environments to simulate real-world vibration environments. This is in contrast
to a modal hammer that only subjects the measured item to an impulse force. An electrodynamic
shaker is shown in figure 7.2. It consists of a strong electromagnet that generates a magnetic field
and a moving coil. When an alternating electrical current is passed through the coil, it generates
a magnetic field that interacts with the magnet, causing the instrumented system to vibrate. The
vibration produced by the electrodynamic shaker can be controlled by adjusting the frequency,
amplitude, and waveform of the electrical signal applied to the coil through an amplifier.
7.3 Digital Signal Processing

An analog signal is a continuous-time signal that can take any value within a certain range, while
a digital signal is a discrete-time signal that takes on only a finite number of values at discrete
time intervals. Digitization in signal processing is the process of converting an analog signal into
a digital form.
172
7.3.1 Sampling and Quantization

The process of digitization involves two main steps: sampling and quantization and is visualized in
figure 7.7 for a sinusoidal and a more complex signal. In the sampling step, the continuous analog
signal is measured at regular time intervals, known as the sampling rate; measured in samples-per-
second (S/s). The resulting discrete-time signal is a sequence of samples that represent the value
of the analog signal at each sampling instant. In the quantization step, each sample is converted
from its continuous value to a digital value that can be represented using a finite number of bits.
The accuracy of the digitized signal depends on the number of bits used for quantization; the more
bits used, the more accurately the signal can be represented.
Figure 7.7: Digitization of two continuous time-series signals sampled at 5 S/s.
Review 7.1 Harry Nyquist (February 7, 1889 - April 4, 1976) was a Swedish physicist and
electronic engineer. His parents emigrated to the U.S. in 1907. He attended the University
of North Dakota starting in 1912 where he obtained a B.S. in 1914 and an M.S. in 1915, both
in electrical engineering (entry to M.S. was 3 years!). Thereafter, he went to Yale University
where he received a Ph.D. in physics in 1917.
173
Figure 7.8: Picture of Harry Nyquist from the American Institute of Physics.a
a Fair use, via Wikimedia Commons
7.3.2 Aliasing
In signal processing, aliasing is an effect that causes different signals to become indistinguishable
from each other, as shown in figure 7.9. In this way, the signals become an alias of one another
when sampled. Aliasing accounts for the development of distortion or artifact in a reconstructed
signal when compared to the original continuous signal. Aliasing occurs when a continuous-time
signal is sampled at a rate that is too low, resulting in a higher-frequency component in the signal
being incorrectly represented as a lower-frequency component due to undersampling.
174
Figure 7.9: Aliasing of a 3 Hz signal that is sampled at 5 S/s where the 3 Hz signal folds back on
itself to create a 2 Hz signal.
The Nyquist-Shannon sampling theorem is a theorem in the field of signal processing that
defines the sample rate that permits a discrete sequence of samples to sample a continuous-time
signal of finite bandwidth. It states that a signal must be sampled at a rate at least twice its highest
frequency component to be accurately the frequency domain of the signal; this is known as the
Nyquist limit. Otherwise, the higher-frequency components of the signal will “fold” back into the
lower-frequency range, resulting in a distorted representation of the signal. Moreover, the signal
must be sampled at twice its highest frequency component with one additional sample to accurately
reproduce the temporal domain of the signal. For example, suppose a sine wave with a frequency
of 3 Hz is sampled at a rate of 5 S/s; as diagrammed in figure 7.9. According to the Nyquist-
Shannon sampling theorem, the signal should be sampled at a rate of at least 6 S/s plus one sample
to rebuild the signal in the temporal domain. Because the sampling rate is lower than the Nyquist
rate, the higher-frequency component of the signal (3 Hz) will be aliased to a lower frequency
(2 Hz), resulting in the distorted representation of the signal shown by the dashed orange line in
figure 7.9.
Rebuilding discretely sampled continuous signals requires much more than just sampling at the
Nyquist limit of 2× the desired frequency content of the signal plus one additional data point. This
is because the Nyquist limit only applies to rebuilding perfect sinusoidal signals and real-world
signals are complex. A good rule of thumb is that a signal must be sampled at least 10 times per
cycle to accurately rebuilt the temporal response of the signal.
175
7.3.3 Time-Frequency Analysis
Figure 7.10: Spectrogram of a 0-10 Hz chirp signal.
The frequency components of a signal can change over time, requiring time-frequency techniques
to analyze. Of these, a spectrogram such as that shown in figure 7.10 is a visual representation of
the spectrum of frequencies of a signal over time. The spectrogram is created by dividing the signal
into short time windows and computing the Fourier transform of each window. By applying the
Fourier transform to each time window of the signal, the spectrogram displays the variation of the
frequency content of the signal over time. Spectrograms can be used for a variety of purposes, such
as identifying and analyzing patterns in the frequency content of a signal, detecting and visualizing
changes in the frequency content over time, and identifying specific frequency components that
may be associated with modes of the vibrating system.
176
Open Vibrations
8 Structural Dynamics
The dynamic response of civil infrastructures, including buildings, bridges, and towers, can be
studied by applying fundamental vibrations concepts studied in the previous chapters.
8.1 Single-story frame

Let’s start by considering the single-story frame shown in Fig. 8.1 (a) in free vibration (no external
load is applied to the structure). The frame has height H and bay width L. As shown in Fig. 8.1,
the frame consists of two columns with a modulus of elasticity E and moment of inertia (second
moment of the cross-sectional area) I. The columns are fixed at the base. The frame in Fig. 8.1 (a)
can be modeled as a single-degree of freedom (SDOF) system under the following assumptions:
• Shear building: flexible columns (EI ̸= 0), beam infinitely rigid (EIb = ∞), axial deforma-
tions of beams and columns negligible (EA = 0);
• Lumped mass system: floor-mass concentrated at the floor level.
Fig. 8.1 (b) illustrates a SDOF with mass m and stiffness k that can be used to model the
dynamic behavior of the single-story frame considering no damping (ζ = 0).
Figure 8.1: (a) Single story frame; (b) undamped single degree of freedom system.
The response of a SDOF system can be written in general notation as:

v0
ωn
where ωn is the natural frequency of the frame, x0 and v0 are the initial conditions. In order to
find ωn , we need to calculate the stiffness of the system. The mass is usually given.
The stiffness of the system can be found by applying Hooke’s law: F = kx. To find k, let’s
imagine applying an arbitrary lateral force F to the frame and analyzing a single column. At the
top, the column will be subjected to a force F and to a moment M0 , as schematically shown in Fig.
8.2 (a). Applying the equilibrium equations to the column, it can be found that M0 = FH 2 .
177
Open Vibrations 8.1 Single-story frame
Figure 8.2: Single column subjected to: (a) force and moment; (b) force only; (c) moment only.
Since the system is linear, we can calculate the effects of F and M0 separately and then summed
them together (superposition principle). The maximum deflection due to F occurs at the top of the
column, as shown in Fig. 8.2 (b), and it is equal to:
FH 3
xmax,F = (8.2)
3EI
while the maximum deflection caused by M0 (Fig. 8.2 (c)) is:
M0 H 2
xmax,M0 = (8.3)
2EI
The displacements in Eq. (8.2) and (8.3) were found using engineering tables. The total dis-
placement x at the top of the column is obtained from the sum of the two displacements:
FH 3 M0 H 2
x= − (8.4)
3EI 2EI
where the xmax,M0 is negative in sign because the displacement caused by M0 goes in opposite
direction to xmax,F . Replacing M0 = FH2 in Eq. (8.4):
FH 3 FH 3 FH 3
x= − = (8.5)
3EI 4EI 12EI
Applying Hooke’s law:
FH 3
F = kc x = kc (8.6)
12EI
where kc is the stiffness of the column. Therefore:
12EI
kc = (8.7)
H3
Since the frame has two columns, the total stiffness of the SDOF system will be:
12EI
k= ∑ kc = ∑ (8.8)
columns 2 H3
where k is also called lateral stiffness. Note that the lateral stiffness of the frame is independent
on the length of the bay L and it depends only on the properties of the columns (E, I, and H). It is
possible at this point to calculate the natural frequency of the frame:
178
Open Vibrations 8.1 Single-story frame
s
∑2 12EI
r
k H3
ωn = = (8.9)
m m
If the columns have same properties, Eq. (8.9) becomes:
r r
k 24EI
ωn = = (8.10)
m H 3m
Finally, the response of the system to initial conditions x0 and v0 can be obtained:
v0
ωn
Example 8.1 Let’s consider the single-story frame shown in Fig. 8.1 with mass m = 0.15 kip
s2 /ft, L = 12 ft, EI= 1800 kip f t 2 . a) Determine the EOM and the natural period of the frame;
b) assume that the moment of inertia of the right column is 2I. Will the EOM change?
Solution a) :
The frame can be modeled as a single degree of freedom in free vibration. Therefore, the
EOM is:
mẍ + kx = 0 (8.12)
The lateral stiffness of the system is:
12EI 24EI
k= ∑ kc = ∑ = (8.13)
columns 2 H3 H3
Thus, the natural frequency and period are:

r r r
k 24EI 1800 rad
ωn = = 3
= 3
= 12.91 (8.14)
m mH 0.15 · 12 s
2π
Tn = = 0.48s (8.15)
ωn
Solution b) :
The EOM won’t change, but the lateral stiffness of the system will be:
12EI 24EI 36EI
k= ∑ kc = + 3 = (8.16)
columns H3 H H3
The same principle can be applied to a single-story frame with damping ratio ζ ̸= 0. In this
case, the displacement of the frame will be given by:
179
Open Vibrations 8.2 Duhamel’s Integral

(v0 + x0 )ωn
(−ζ ωn t)
x(t) = e cos(ωd t) + x0 sin(ωd t) (8.17)
ωd
where ωd is the damped natural frequency of the system:
q
ωd = ωn 1 − ζ 2 (8.18)
8.2 Duhamel’s Integral

In Chapter 4, the frequency response method was used to solve the EOM of an SDOF system
subjected to an arbitrary force. Here, an alternative method widely employed in structural dynam-
ics to find the solution of the EOM is presented. This method exploits a specif integral, named
Duhamel’s integral.
Let’s consider an underdamped SDOF system subjected to an arbitrary force F(t). The EOM
is:
mẍ + cẋ + kx = F(t) (8.19)

Let’s assume that the system is at rest: x(0) = 0 and ẋ = 0. The assumption underlying
Duhamel’s integral method is that a generic force F(t) can be expressed as a sequence of impulses
of very small duration and the response of the system as the sum of the response to individual unit
impulses.
An impulsive force can be defined as a very large force applied in a very short time interval.
Fig. 8.3 (a) shows an impulsive force F(t) = ε1 applied at time t = τ. Assuming to apply an
impulsive force to a generic mass m and applying Netwon’s second law:
mẍ = F(t) (8.20)

and integrating both sides between two generic time instants t1 and t2 yields:
Z t2
F(t)dt = m(ẋ1 − ẋ2 ) (8.21)
t1
where the left-hand side of the equation represents the magnitude of the force and the right-
hand side the change in momentum.
In the limit case in which ε tends to
R t2
0, F(t) tends to 1 and the impulsive force is called unit
impulse. In the case of a unit impulse, t1 F(t)dt = 1 and t1 tends to t2 . Therefore, the velocity of
the mass can be found as:
1
ẋ(τ) = (8.22)
m
A similar concept applies to an SDOF system. Since the impulse is applied in a very short time
interval, the spring and the damper do not have the time to react. When we apply a unit impulse to
an underdamped SDOF, the system will start vibrating with velocity ẋ(τ) given by Eq. (8.22) and
displacement x(τ) = 0. The response of the system is given by the following equation:
1 −ζ ωn (t−τ)
x(t) = h(t − τ) = e sin(ωd (t − τ)) (8.23)
mωd
180
Open Vibrations 8.2 Duhamel’s Integral
where τ is the time instant at which the impulse is applied. Note that the Dirac delta function
δ (t − τ) mathematically defines a unit impulse centered at t = τ.
Figure 8.3: (a) Impulsive force; (b) arbitrary force decomposed in a series of impulses.
Let’s now consider a force F(t) varying arbitrarily with time. As shown in Fig. 8.3 (b), F(t)
can be represented as a sequence of infinitesimaly short impulses. The response of a linear system
to F(t) can be therefore expressed as the response to a series impulses, following:
Z t
x(t) = p(τ)h(t − τ)dτ (8.24)
0
where h(t −τ) is the response to a unit impulse and p(τ) is the magnitude of the actual impulse.
For the case of an underdamped SDOF system, Eq. (8.24) can be re-written as:
1 t
Z
x(t) = p(τ)e−ζ ωn (t−τ) sin(ωd (t − τ))dτ (8.25)
mωd 0
Eq. (8.25) represents the Duhamel’s integral.
Similarly, the response of an undamped SDOF system to an arbitrary force can be expressed
through the Duhamel’s integral as:
Z t
1
x(t) = p(τ)sin(ωn (t − τ))dτ (8.26)
mωn 0
If F(t) is characterized by a simple function, Duhamel’s integral can be evaluated in closed
form. If the equation of F(t) is complicated, Duhamel’s integral can be solved with numerical
methods.
Note that Eq. (8.25) and (8.26) apply when the initial conditions are zero (the system is at rest).
If the initial conditions are different than zero, we need to add the free vibration response of the
system to Eq. (8.25) and (8.26), respectively.
Example 8.2 Let’s consider an undamped SDOF system subjected to a step function force
with constant amplitude F0 , as schematically represented in Fig. 8.4. Assume that the system
is at rest (initial conditions: x(0) = ẋ(0) = 0) and compute the system response x(t).
181
Open Vibrations 8.3 Two-story frame
Figure 8.4: (a) Step function force; (b) undamped SDOF system.
Solution:
The system is undamped, therefore we can use the Duhamel’s integral in Eq. (8.26) to find
x(t):
Z t
1
x(t) = F0 sin(ωn (t − τ))dτ (8.27)
mωn 0
Considering that F0 is constant:
F0 cos(ωn (t − τ)) t

F0
x(t) = = [1 − cos(ωnt)] (8.28)
mωn ωn 0 mωn2
Reminding that ωn2 = k/m, x(t) becomes:
F0
x(t) = [1 − cos(ωnt)] (8.29)
k
where Fk0 is the displacement that the system would undergo if the force F0 was applied
statically.
In the case of underdamped SDOF system, the response becomes:
" !#
F0 ζ
x(t) = 1 − e−ζ ωnt cos(ωd t) + p sin(ωd t) (8.30)
k 1−ζ2
8.3 Two-story frame

The concepts discussed in Sec. 1 can be extended to the 2-story frame represented in Fig. 8.5. In
fact, a 2-story frame can be modeled as a 2-DOF system under the following assumptions:
• shear building: flexible columns (EI ̸= 0), beam infinitely rigid (EIb = ∞), axial deformations
of beams and columns negligible (EA = 0);
• lumped mass system: floor-mass concentrated at the floor level.
Under such assumptions and free vibrations, we expect that the building moves following the
182
deformed shape reported in Fig. 8.5 (dotted line). Let’s call the degrees of freedom of the frame
x1 (t) and x2 (t).
Figure 8.5: 2-story frame with lumped masses.
The forces acting on the 2-DOF system are reported in Fig. 8.6. It follows that the equation of
motion of the two masses are:
m1 ẍ1 + k1 x1 + k2 (x2 − x1 ) + c1 ẋ1 + c2 (ẋ2 − ẋ1 ) = 0 (8.31)

m2 ẍ2 − k2 (x2 − x1 ) − c2 (ẋ2 − ẋ1 ) = 0
Figure 8.6: (a) 2-DOF system used to model the 2-story frame; (b) free body diagram of the two
masses.
In matrix notation, these two equations become:
183

m1 0 x¨1 k1 + k2 −k2 x1 c1 + c2 −c2 ẋ1 0
+ + = (8.32)
0 m2 x¨2 −k2 k2 x2 −c2 c2 ẋ2 0
where we can define the mass matrix M as:

m1 0
M= (8.33)
0 m2
the stiffness matrix K as:

k1 + k2 −k2
K= (8.34)
−k2 k2
and the damping matrix C as:

c1 + c2 −c2
C= (8.35)
−c2 c2
While mass and damping of a frame are usually given, the stiffness values k1 and k2 need to be
calculated as a function of the columns properties (EI) and geometry (h). As demonstrated in Sec.
1, the stiffness of a column with clamped ends can be determined as:
12EI
kc = (8.36)
h3
The lateral stiffness of each floor can be computed as the sum of the stiffness of the columns at
that floor:
12EI
k= ∑ kc = ∑ (8.37)
columns 2 h3
Therefore, for the frame in Fig. 8.5, the stiffness values are:
24EI
k1 = k2 = (8.38)
h3
The solution of the EOM in Eq.(8.32) was derived in Chapter 5 and can be summarized as:

x1 (t) A1 sin(ω1t + φ1 )
= u1 u2 , ω1 or ω2 ̸= 0 (8.39)
x2 (t) A2 sin(ω2t + φ2 )
where u1 and u2 are eigenvectors (or mode shapes), ω1 and ω2 are the natural frequency of vibra-
tion, φ1 , φ2 , A1 , and A2 are constants that can be found based on the initial conditions (see Chapter
5 for more details).
184
Example 8.3 Consider the frame in Fig. 8.7. Determine the natural frequency of vibration
and mode shapes of the system.
Figure 8.7: Example of a 2-story frame with floors with different dynamic properties.
Solution:
Assumption: the frame can be modeled as a shear building with mass lumped at the floor
levels. The lateral stiffness at the first floor is:
12(2EI) 48EI
k1 = 2 = 3 (8.40)
h3 h
The lateral stiffness at the second floor is:
12(EI) 24EI
k2 = 2 = 3 (8.41)
h3 h
Therefore, the stiffness matrix can be written as:
48EI
+ 24EI − 24EI

h3 h3 h3
24EI 3 −1
K= = 3 (8.42)
− 24EI
h3
24EI
h3
h −1 1
The EOM of the system is:

2m 0 x¨1 24EI 3 −1 x1 p
+ 3 = 1 (8.43)
0 m x¨2 h −1 1 x2 p2
In order to determine the natural frequency of vibration and the mode shapes of the system,
we need to solve the characteristic equation:
185
det(−ω 2 M + K) = 0 (8.44)
leading to:
2m2 ω 4 + 5kmω 2 + 2k = 0 (8.45)

This equation has two solutions:
k
ω12 = (8.46)
2m
2k
ω22 = (8.47)
m
Therefore, the two natural frequencies of vibration of the system are:
r
k
ω1 = (8.48)
r2m
2k
ω2 = (8.49)
m
where k = 24EI
h3
. The mode shapes of the frame can be found by solving the following equation:
(−ω12 M + K)u1 = 0 (8.50)

Replacing mass and stiffness matrix the equation becomes:

k 2m 0 3 −1 u11 0
− +k = (8.51)
2m 0 m −1 1 u21 0
simplified to

2k −k u11 0
= (8.52)
−k 2k u21 0
leading to two equations:
k
2ku11 − ku21 = 0, and − ku11 + u21 = 0 (8.53)
2
It follows that:
1
2u11 = u21 , and u11 = u21 (8.54)
2
To obtain a numerical value, we arbitrarily assign a value to one of the elements. Here, let
u21 = 1 so let u11 = 1/2. Therefore,
1

u1 = 2 (8.55)
1
186
Similarly, u2 can be obtained by solving the following equation:
(−ω22 M + K)u2 = 0 (8.56)

leading to:

−1
u2 = (8.57)
1
Fig. 8.8 represents the two-mode shapes of the building.
Figure 8.8: Mode shapes of the 2-story frame.
The temporal response of the system is given by:

x1 (t) A1 sin(ω1t + φ1 )
= u1 u2 , ω1 or ω2 ̸= 0 (8.58)
x2 (t) A2 sin(ω2t + φ2 )
where u = [u1 , u2 ] is the time invariant part of the equation.
187

Open Vibrations

Uploaded by

Copyright:

Available Formats

Open Vibrations

Uploaded by

Document Information

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

Open Vibrations

Uploaded by

Copyright:

Available Formats

Open Vibrations

Austin R.J. Downey

December 24, 2023

3.6 Numerical Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83

4 Transfer Function Approach 86

5 Multiple Degree-of-freedom Systems 113

II Applied Topics 153

7 Experimental Vibrations 167

8 Structural Dynamics 177

1 Basic Concepts in Vibrations

Review 1.1 Newton’s three laws of motion:

1. In an inertial frame of reference, an object either remains at rest or continues to move

1.1 Single Degree-of-Freedom Systems

3. door swinging on axis

4. throttle (gas pedal)

1.1.1 Spring-Mass Model

“All models are wrong, but some are useful ”

Figure 1.3: A single-degree-of-freedom (1-DOF) spring-mass model showing: (a) annotated

1.1.2 Linear Springs

Figure 1.4: Force-displacement plot for a linear spring.

1.2 Equivalent Stiffness

It follows that the change in length ∆l can be expressed as:

1.2.1 Equivalent Stiffness of Structural Systems

This expression can be expanded into the form:

Figure 1.9: Compressed cantilever rod.

1.2.2 Springs in Series and Parallel

Figure 1.11: Two springs k1 and k2 combined in series.

where the total displacement is

Figure 1.12: Two springs k1 and k2 combined in parallel.

The displacement in both springs is the same, so the total displacement is

The forces in the direction of spring elongation sum to zero, therefore:

this simplifies to:

Example 1.2 Calculate the equivalent stiffness of the following system:

Figure 1.13: Equivalent stiffness for springs in series and parallel.

1.3 Equation of Motion for an Oscillating System

2. Apply Newton’s second law to both FBDs ( equilibrium and displaced).

aẍ + bẋ + cx = 0 (1.24)

Example 1.3 Considering the system:

equilibrium state displaced state

The equation for the equilibrium state is:

and in the displaced state:

mẍ = −kx (1.28)

equilibrium state displaced state

The equation for the equilibrium state is:

and in the displaced state:

2.1 Mathematical Modeling of Free Vibration

Figure 2.3: 1-DOF spring-mass system.

can be modeled expressed with the following EOM

mẍ(t) + kx(t) = 0 (2.1)

Figure 2.4: Summary of the temporal response for a 1-DOF system.

2.1.1 Solve for the Natural Frequency (ωn ) of the System

• System oscillates at different speed → use a parameter to adjust ωn in rad/s.

• Systems have different amplitudes → use a parameter to adjust A in meters.

• System has different starting points → use a parameter to adjust φ in rad.

Using these four constituent components, an equation can be proposed:

x(t) = Asin(ωnt + φ ) (2.2)

Take the derivative to get velocity:

ẋ(t) = Aωn cos(ωnt + φ ) (2.3)

Take the derivative again to get acceleration:

ẍ(t) = −Aωn2 sin(ωnt + φ ) (2.4)