2019

An Introduction to Shock and Vibration Response Spectra

Revision A

BY TOM IRVINE
1 CONTENTS
1 CONTENTS ............................................................................................................................................. 1
3 TOM IRVINE BACKGROUND .................................................................................................................. 8
4 VARIABLES ............................................................................................................................................. 9
5 INTRODUCTION ..................................................................................................................................... 1
5.1 NATURAL FREQUENCY PHOTO GALLERY....................................................................................... 1
5.2 SOME HISTORY .............................................................................................................................. 3
5.3 MECHANICAL FAILURE .................................................................................................................. 5
5.3.1 Failure & Damage Photo Gallery ........................................................................................... 5
5.3.2 Historical Failures .................................................................................................................. 6
5.3.3 Aerospace and Other Industries ........................................................................................... 7
5.3.4 Failure Modes ....................................................................................................................... 7
5.4 MACHINE HEALTH MONITORING.................................................................................................. 8
5.5 VIBRATION ENERGY HARVESTING ................................................................................................ 9
5.6 HUMAN VIBRATION ...................................................................................................................... 9
5.7 VIBRATION EXCITATION SOURCES .............................................................................................. 10
5.7.1 Source Characterization ...................................................................................................... 10
5.7.2 Initial Displacement or Velocity .......................................................................................... 10
5.7.3 Applied Force ...................................................................................................................... 10
5.7.4 Base Excitation .................................................................................................................... 10
5.7.5 Self-Excited Vibration .......................................................................................................... 10
5.7.6 Source Photo Gallery, Base Excitation ................................................................................ 11
5.7.7 Source Photo Gallery, Applied Force .................................................................................. 15
5.7.8 Source Photo Gallery, Self-Excited Vibration ...................................................................... 18
6 SINGLE-DEGREE-OF-FREEDOM SYSTEMS & BASIC CONCEPTS ............................................................ 19
6.1 FREE VIBRATION MODEL............................................................................................................. 19
6.2 CONSERVATION OF ENERGY METHOD ....................................................................................... 20
6.3 RAYLEIGH METHOD ..................................................................................................................... 21
6.4 NATURAL FREQUENCY EXAMPLE ................................................................................................ 22
6.5 OCTAVE RULE .............................................................................................................................. 22
6.6 SOME MATHEMATICS ................................................................................................................. 23
6.6.1 Complex Variables............................................................................................................... 23
6.6.2 Euler’s Equation .................................................................................................................. 23

1
7 TWO-DEGREE-OF-FREEDOM SYSTEMS ............................................................................................... 24
7.1.1 Grounded Two-DOF System, Equation of Motion .............................................................. 24
7.1.2 Grounded Two-DOF System, Eigenvalues & Vectors .......................................................... 25
7.1.3 Grounded Two-DOF System Example, Free Vibration ........................................................ 26
7.1.4 Grounded Two-DOF System, Participation Factors & Effective Modal Mass ..................... 28
7.1.5 Two-DOF System Example, Response to Initial Displacement ........................................... 29
7.1.6 Semi-definite, Two-DOF System ......................................................................................... 31
8 CONTINUOUS SYSTEMS WITH LONGITUDINAL VIBRATION ................................................................ 33
8.1 ROD LONGITUDINAL NATURAL FREQUENCY .............................................................................. 33
8.1.1 One-Dimensional Longitudinal Vibration Equation of Motion ........................................... 33
8.1.2 Rod Longitudinal Vibration Energy Formulas ..................................................................... 36
8.1.3 Beam Longitudinal Vibration Example ................................................................................ 36
8.2 COIL SPRING SURGE .................................................................................................................... 37
8.3 RING FREQUENCY........................................................................................................................ 40
9 BERNOULLI-EULER BEAMS .................................................................................................................. 42
9.1.1 Beam Equations .................................................................................................................. 42
9.1.2 Common Boundary Conditions ........................................................................................... 44
9.1.3 Beam Bending Fundamental Frequencies .......................................................................... 45
9.1.4 Beam Bending Participation Factors & Effective Modal Mass............................................ 47
9.1.5 Bending Wave Speed & Wavelength .................................................................................. 48
9.1.6 Beam Bending Energy Formulas ......................................................................................... 49
9.1.7 Beam Example, Wind Chime ............................................................................................... 50
10 THIN PLATE VIBRATION ................................................................................................................... 55
10.1.1 Rectangular Plate, Bending Vibration ................................................................................. 55
10.1.2 Circular Plate, Bending Vibration ........................................................................................ 58
10.1.3 Honeycomb Sandwich Plate ............................................................................................... 62
11 DAMPING ........................................................................................................................................ 66
11.1 DAMPING METRICS ..................................................................................................................... 66
11.2 FOUR DAMPING CATEGORIES ..................................................................................................... 67
11.3 DAMPING MECHANISMS ............................................................................................................ 69
11.4 HUNTSVILLE, ALABAMA PEDESTRIAN BRIDGE DAMPING ........................................................... 70
11.5 PEGASUS LAUNCH VEHICLE DAMPING ....................................................................................... 72
11.5.1 Drop Transient .................................................................................................................... 72

2
 11.5.2 Drop Transient Damped Sine Curve-fit ............................................................................... 73
11.5.3 Drop Transient Logarithmic Decrement Method ............................................................... 74
11.6 TRANSAMERICA BUILDING DAMPING ........................................................................................ 75
12 SINUSOIDAL VIBRATION.................................................................................................................. 77
12.1 DESCRIPTIVE STATISTICS ............................................................................................................. 77
12.2 SINE TESTING .............................................................................................................................. 80
12.3 SAMPLE SINE CALCULATIONS ..................................................................................................... 80
12.4 BEAT FREQUENCY EFFECT ........................................................................................................... 81
13 DYNAMIC MODELS FOR EXTERNAL EXCITATION ............................................................................ 83
13.1 TWO MODELS ............................................................................................................................. 83
13.2 APPLIED FORCE ........................................................................................................................... 83
13.2.1 SDOF Response Equations .................................................................................................. 83
13.2.2 Input & Response Example, Applied Force ......................................................................... 85
13.2.3 Six Transfer Functions for Applied Force and Response ..................................................... 87
13.2.4 Helicopter Ground Resonance ............................................................................................ 88
13.2.5 Two-Degree-of-Freedom System Subjected to Applied Forces .......................................... 90
13.2.6 Two-Degree-of-Freedom System, Frequency Domain Example ......................................... 96
13.3 BASE EXCITATION ........................................................................................................................ 97
13.3.1 SDOF Response Equations .................................................................................................. 97
13.3.2 Input & Response Example, Base Excitation ..................................................................... 100
13.3.3 Alternate Base Excitation Method using Seismic Mass .................................................... 102
13.3.4 Two-DOF Response Equations, Time Domain................................................................... 105
13.3.5 Two-DOF Response Equations, Time Domain Example .................................................... 107
13.3.6 Two-DOF Response Equations, Frequency Domain .......................................................... 109
13.3.7 Two-DOF Response Equations, Frequency Example......................................................... 111
14 SINE SWEEP TESTING .................................................................................................................... 112
14.1 SINE SWEEP INTRODUCTION .................................................................................................... 112
14.2 SINE SWEEP TIME DOMAIN EQUATIONS .................................................................................. 113
14.2.1 Linear Sweep Rate............................................................................................................. 113
14.2.2 Logarithmic Sweep Rate.................................................................................................... 113
14.2.3 Sweep Type Comparison ................................................................................................... 114
14.3 SINE SWEEP SPECIFICATION EXAMPLE ..................................................................................... 114
14.4 SINE SWEEP INPUT & RESPONSE EXAMPLE, SDOF SYSTEM...................................................... 116

3
 14.5 SINE SWEEP INPUT & RESPONSE EXAMPLE, CONTINUOUS BEAM ........................................... 119
14.6 RANDOM EXCITATION SOURCES .............................................................................................. 122
14.7 STATISTICAL PARAMETERS........................................................................................................ 123
14.8 RANDOM VIBRATION HISTOGRAMS ......................................................................................... 125
14.9 WHITE & PINK NOISE ................................................................................................................ 129
14.10 WHITE NOISE INPUT & RESPONSE EXAMPLE ........................................................................ 130
14.11 CREST FACTOR STUDY ........................................................................................................... 133
14.12 NON-GAUSSIAN RANDOM VIBRATION ................................................................................. 134
14.13 HALT/HASS TESTING WITH NON-GAUSSIAN RANDOM VIBRATION ..................................... 135
15 SAMPLE RATE & ALIASING ............................................................................................................ 141
15.1 DATA COLLECTION ACCURACY .................................................................................................. 141
15.2 FIRST REQUIREMENT ................................................................................................................ 141
15.3 SECOND REQUIREMENT ............................................................................................................ 142
15.4 ANTI-ALIASING FILTER ............................................................................................................... 143
15.5 ALIASING EXAMPLE ................................................................................................................... 145
16 FOURIER TRANSFORMS................................................................................................................. 148
16.1 FOURIER TRANSFORM EQUATIONS .......................................................................................... 148
16.2 TEXTBOOK FOURIER TRANSFORM ............................................................................................ 148
16.2.1 Discrete Fourier Transform ............................................................................................... 150
16.2.2 Useful Discrete Fourier Transform .................................................................................... 151
16.3 LEAKAGE & HANNING WINDOW .............................................................................................. 153
16.3.1 Leakage Error .................................................................................................................... 153
16.3.2 Hanning Window............................................................................................................... 157
16.4 FAST FOURIER TRANSFORM...................................................................................................... 160
16.4.1 FFT Characteristics ............................................................................................................ 160
16.4.2 FFT Example ...................................................................................................................... 160
16.4.3 FFT Algorithm Background ................................................................................................ 161
16.5 ELECTRICAL TRANSFORMER HUM ............................................................................................ 163
16.6 APACHE HELICOPTER FLYOVER DATA ....................................................................................... 166
16.7 TURBOPROP ENGINE................................................................................................................. 169
16.8 WATERFALL FFT......................................................................................................................... 171
16.8.1 Introduction ...................................................................................................................... 171
16.8.2 Solid Rocket Motor Example ............................................................................................. 171

4
 16.8.3 Wind Chime Acoustics....................................................................................................... 174
16.9 WHITE NOISE FOURIER TRANSFORM ........................................................................................ 174
17 POWER SPECTRAL DENSITY........................................................................................................... 177
17.1 DECIBELS & SLOPES ................................................................................................................... 177
17.1.1 Differences ........................................................................................................................ 177
17.1.2 Slopes ................................................................................................................................ 177
17.2 PSD SPECIFICATIONS ................................................................................................................. 178
17.3 PSD INTEGRATION IN THE FREQUENCY DOMAIN ..................................................................... 180
17.4 PSD CALCULATION METHODS .................................................................................................. 182
17.5 PSD CALCULATION EXAMPLES ................................................................................................. 185
17.5.1 White Noise PSD ............................................................................................................... 185
17.5.2 Taurus Auto PSD................................................................................................................ 186
17.6 PSD TIME HISTORY SYNTHESIS .................................................................................................. 189
17.7 MILES EQUATION FOR BASE EXCITATION ................................................................................. 198
17.8 GENERAL METHOD FOR BASE EXCITATION .............................................................................. 198
18 VIBRATION RESPONSE SPECTRUM................................................................................................ 201
18.1 VRS BASICS ................................................................................................................................ 201
18.2 VRS COMPARISON OF TWO PSDS ............................................................................................. 203
18.3 PSD ENVELOPE DERIVATION VIA VRS ....................................................................................... 205
19 SHOCK TESTING & ANALYSIS ......................................................................................................... 207
19.1 SHOCK ISOLATOR PHOTO GALLERY .......................................................................................... 207
19.2 SIMPLE DROP SHOCK ................................................................................................................ 208
19.3 CLASSICAL SHOCK...................................................................................................................... 210
19.4 HALF-SINE SHOCK EXAMPLE ..................................................................................................... 214
19.5 RESPONSE TO ARBITRARY EXCITATION .................................................................................... 221
19.5.1 Shock Response Spectrum Test Specification Objective .................................................. 222
19.6 SEISMIC SHOCK ......................................................................................................................... 224
19.6.1 Seismic Waveforms ........................................................................................................... 224
19.6.2 Seismic Response Spectrum Method................................................................................ 226
19.6.3 El Centro Earthquake ........................................................................................................ 227
19.6.4 Golden Gate Bridge ........................................................................................................... 231
19.6.5 Vandenberg AFB................................................................................................................ 232
19.6.6 Seismic Testing .................................................................................................................. 234

5
 19.7 PYROTECHNIC SHOCK ............................................................................................................... 235
19.7.1 Flight Events ...................................................................................................................... 235
19.7.2 Pyrotechnic Device Photo Gallery ..................................................................................... 236
19.8 PYROTECHNIC SHOCK DATA...................................................................................................... 239
19.8.1 Initial SRS Slopes ............................................................................................................... 239
19.8.2 Re-entry Vehicle Separation, Flight Data .......................................................................... 240
19.9 WATER IMPACT SHOCK ............................................................................................................. 242
19.10 SHOCK RESPONSE SPECTRUM SYNTHESIS ............................................................................ 245
19.10.1 Synthesis Objectives ..................................................................................................... 245
19.10.2 Wavelet Synthesis ......................................................................................................... 246
19.10.3 Damped-Sine Synthesis ................................................................................................. 251
19.10.4 Modal Transient Finite Element Analysis for Uniform Base Excitation ........................ 256
19.10.5 Shock Fields ................................................................................................................... 257
19.10.6 Shock Mitigation ........................................................................................................... 259
19.10.7 Pyrotechnic Shock Testing Methods ............................................................................. 261
19.10.8 Shock Failure Modes ..................................................................................................... 263
20 PSEUDO VELOCITY & THE STRESS-VELOCITY RELATIONSHIP ........................................................ 264
20.1 PSEUDO VELOCITY..................................................................................................................... 264
20.2 STRESS-VELOCITY RELATIONSHIP .............................................................................................. 265
20.3 BEAM STRESS ............................................................................................................................ 267
20.4 SHOCK VELOCITY THRESHOLD RULE-OF-THUMB...................................................................... 268
20.5 VELOCITY THRESHOLD FROM YIELD STRESS ............................................................................. 268
20.6 MORSE CHART........................................................................................................................... 269
21 DIGITAL FILTERING ........................................................................................................................ 271
21.1 FILTERING INTRODUCTION ....................................................................................................... 271
21.2 HIGH-PASS & LOW-PASS FILTERS .............................................................................................. 271
21.3 BUTTERWORTH FILTER CHARACTERISTICS ............................................................................... 271
21.4 PHASE CORRECTION.................................................................................................................. 276
21.5 PSD CALCULATION VIA DIGITAL FILTER..................................................................................... 277
21.6 SEISMIC DATA FILTERING .......................................................................................................... 280
21.7 INTEGRATION AND DIFFERENTIATION OF ACCELEROMETER SIGNALS .................................... 286
21.8 BAND LIMITED WHITE NOISE .................................................................................................... 290
21.9 LAUNCH VEHICLE COUPLED-LOADS ANALYSIS FILTERING ........................................................ 291

6
22 THE GREAT SPACECRAFT BASE INPUT VIBRATION TEST DEBATE ................................................. 292
22.1 DEBATE INTRODUCTION ........................................................................................................... 292
22.2 Arguments against Spacecraft Vibration Testing .............................................................. 293
22.3 Arguments For Spacecraft Vibration Testing ..................................................................... 294
22.4 SPACECRAFT TESTING EXPERIENCES ......................................................................................... 294
22.5 SINE SWEEP CONTROL .............................................................................................................. 295
23 TACOMA NARROWS BRIDGE FAILURE .......................................................................................... 296
23.1 BRIDGE INTRODUCTION ............................................................................................................ 296
23.2 BRIDGE COLLAPSE ..................................................................................................................... 297
23.3 FAILURE THEORIES .................................................................................................................... 298
23.3.1 Candidates......................................................................................................................... 298
23.3.2 Random Turbulence .......................................................................................................... 298
23.3.3 Vortex Shedding ................................................................................................................ 298
23.3.4 Aerodynamic Instability .................................................................................................... 299
23.4 BRIDGE’S STROUHAL FREQUENCY ............................................................................................ 302
23.5 REPLACEMENT BRIDGE ............................................................................................................. 302
24 BIBLIOGRAPHY .............................................................................................................................. 304

7
3 TOM IRVINE BACKGROUND

Tom Irvine earned B.S. and M.S. degrees in engineering science from Arizona State University, in 1985 and
1987, respectively. He wandered into the specialty of vibration serendipitously and found that vibration’s
quirky traits were well-suited to his unconventional personality.

He worked at Orbital Sciences Corporation in Chandler, Arizona, for 11 years as a full-time employee, and
for 9 years as a consultant. He derived the shock and vibration levels and performed hands-on testing for
the first Pegasus launch vehicle, which successfully flew in 1990. He has analyzed flight accelerometer
data from Pegasus and numerous launch vehicles, as well as data from ground separation tests. He has
posted hundreds of tutorial and software programs at his Vibrationdata website, covering shock,
vibration, acoustics, structural dynamics, and signal processing.

He is currently a scientist at Dynamic Concepts, Inc., in Huntsville, Alabama, where he works on the NASA
Space Launch System (SLS) vehicle. He is also an industry representative for the NASA Engineering &
Safety Center (NESC) Loads and Dynamics Technical Disciplines Team.

Tom can be contacted via Email: [email protected]

Readers are welcome to visit his blog: https://vibrationdata.wordpress.com/

Tom either performed or verified the calculations in this e-book using his Vibrationdata Matlab GUI
package which is freely available at the above blog link. He hopes this book will help bridge the gap
between shock and vibration analysis as taught in universities and as practiced in industry.

Also note that this is expected to be the first of a series. Future works will include fatigue, finite element
analysis, vibroacoustics, etc.

8
4 VARIABLES

1. Amplitude m mass
A
2. Cross-sectional area
Generalized displacement
a Acceleration N
coordinate in frequency domain
Bandwidth of an ideal rectangular 1. Number of octaves
B n
filter 2. Number of points
1. Damping coefficient matrix 3. Frequency index
C
2. Crest factor PE Potential energy
1. Damping coefficient
c p Probability density function
2. Speed of sound in material
CL Longitudinal wave speed 1. Amplification factor
Q
2. Eigenvector matrix
cp Bending wave phase speed 1. Eigenvector
cg Bending wave group speed q 2. Generalized displacement
coordinate in time domain
D Plate stiffness factor S Strouhal number
d Diameter s Laplace variable
E Elastic modulus T Period or duration
F Applied force in frequency domain t Time
1. Frequency
f U Free stream air velocity
2. Applied force in time domain
fn Natural Frequency (Hz) u Longitudinal displacement
fr Ring frequency V Shear force
fs Strouhal frequency v Velocity

G, g Gravity acceleration, 9.81 m/sec2 W Wavelet or damped sine amplitude

H Transfer function w Hanning window function

1. Impulse response function X, x Displacement amplitude
h
2. Plate thickness
1. Amplitude
1. Area moment of inertia Y
I 2. Beam transverse displacement
2. Identity matrix
KE y Base displacement amplitude
Kinetic energy
1. Stiffness Z Plate transverse displacement
k
2. Wavenumber
z Relative displacement
L Length
1. Mass matrix
M
2. Bending moment

9
 Beam bending wavenumber
1. Logarithmic decrement0

2. Dirac Delta Function
 Loss factor

 Kurtosis

 Wavelength

 1. Poisson ratio
2. Mean
 Viscous damping ratio
 Mass density (mass/volume)
L Mass density (mass/length)
1. Standard deviation

2. Stress
Time delay in the convolution

integral
 Phase angle

 Eigenvalue matrix

 Angular frequency (rad/sec)

d Damped angular natural frequency
n Angular natural frequency
 “If you want to find the secrets of the universe, think in terms of
energy, frequency and vibration.”

- Nikola Tesla

5 INTRODUCTION

5.1 NATURAL FREQUENCY PHOTO GALLERY

Figure 5.1. Mother Earth

The Earth experiences seismic vibration. The fundamental natural frequency of the Earth is 309.286 micro
Hertz, equivalent to a period of approximately 54 minutes [1].

1
Figure 5.2. IRVE 2 Launch & Free-Free Beam Analogy

A rocket vehicle behaves as a free-free beam during flight. The vehicle’s body bending modes can be
excited by wind gusts, aerodynamic buffeting, thrust offset, maneuvers, etc. The image shows the IRVE 2
launch from Wallops Island. The vehicle is a Black Brant 9, which has a Terrier Mark 70 first stage. Flight
accelerometer data showed that the fundamental bending frequency began at about 8 Hz and then swept
up to 13 Hz as propellant mass was expelled. Autopilot guidance and control algorithms need to account
for the body bending mode to maintain stability.

2
Figure 5.3. Tuning Fork, First Mode, A4 Note, Fundamental Frequency 440 Hz

The mode shape from a finite element model is shown exaggerated. The two tines undergo an in-plane
bending mode, 180 degrees out-of-phase with one another. The stem also participates in this mode, but
its displacement is so relatively small that it is not apparent in the figure.

5.2 SOME HISTORY

The Ancient Greek philosopher Pythagoras (570–495 BC) studied the vibration of stringed instruments
and developed a theory of harmony. He first identified that the pitch of a musical note is in inverse
proportion to the length of the string that produces it, and that intervals between harmonious sound
frequencies form simple numerical ratios.

The ancient Greeks believed that al objects were composed of some combination of the basic elements:
earth, water, air, and fire. Aristotle (384–322 BC) attempted an explanation of earthquakes based on
natural phenomena. He postulated that winds within the earth whipped up the occasional shaking of the
earth's surface. He noted that earthquakes sometimes caused the water to burst forth in what would
now be called a tsunami.

Galileo Galilei (1564-1642) performed numerous experiments with oscillating pendulums during the
Renaissance. He discovered that the period of swing of a pendulum was independent of its amplitude.

3
He used his own pulse as a time measurement because there were no watches at that time. Christiaan
Huygens (1629-1695) successfully built a pendulum clock based upon Galileo’s work.

Robert Hooke (1635-1703) developed the law of linear spring stiffness. This law has been generalized to
the elasticity principle that strain in a body is proportional to the applied stress.

Sir Isaac Newton (1643-1727) derived his laws of motion, which are the foundation of mechanical
vibration analysis. He published his findings in a book known by its abbreviated title Principia.

Jacob Bernoulli, Daniel Bernoulli, and Leonhard Euler derived the equation of motion for beam vibration
circa 1750.

Lord Rayleigh, John William Strutt, published Theory of Sound in two volumes during 1877-1878. Volume
I covered harmonic vibrations, systems with one degree of freedom, vibrating systems in general,
transverse vibrations of strings, longitudinal and torsional vibrations of bars, vibrations of membranes and
plates, curved shells and plates, and electrical vibrations. Volume II covered aerial vibrations, vibrations
in tubes, reflection and refraction - of plane waves, general equations, theory of resonators, Laplace’s
functions and acoustics, spherical sheets of air, vibration of solid bodies, and facts and theories of
audition.

The modern field of mechanical vibration analysis has been built upon the foundation of these authors’
works and has further developed as the result of failures and disasters, as well as the need to design,
analyze and test structures and component to withstand dynamic environments. Vibration analysis is also
important in other areas including semiconductor manufacturing, human exposure, energy harvesting,
machine health monitoring, etc.

4
5.3 MECHANICAL FAILURE

5.3.1 Failure & Damage Photo Gallery

Figure 5.4. Turbine Blade Fatigue Failure

A Qantas Boeing 747-400 aircraft was flying from Sydney to Singapore on May 9, 2011 when the crew
noticed abnormalities from the aircrafts No. 4 engine during a climb. The indications included an increase
in both the exhaust gas temperature and vibration levels. The plane continued to Singapore for a safe
landing and disembarkation of the passengers and crew. The Australian Transport Safety Bureau (ATSB)
determined that the cause was a broken turbine blade.

5
Figure 5.5. Washington Monument Crack

A one-inch wide, four-foot long crack formed in the Washington Monument, near the top of the 555-
foot obelisk, due to the Mineral, Virginia earthquake on August 23, 2011.

5.3.2 Historical Failures

Railcar axles were failing under repeated “low level” cyclic stress, in the mid nineteenth century. These
stresses puzzled engineers because the levels were much lower than the material yield stress. The
Versailles rail accident occurred on May 8, 1842 after the leading locomotive broke an axle, and the
carriages behind piled into it and caught fire causing many deaths. This prompted German scientist
August Wöhler to develop the S-N curves used in fatigue analysis. The term “fatigue” was chosen to
describe the “tired metal” in the axles.

The collapse of the Tacoma Narrows Bridge “Galloping Gertie” was captured on film on November 7, 1940.
This failure is often referred to as the classic “resonant vibration” failure. But it was more properly a “self-
excited” or “flutter” response, as discussed in Section 23.

6
5.3.3 Aerospace and Other Industries

Aerospace has myriad examples of potential vibration problems. Helicopters may undergo “ground
resonance” prior to takeoff. Launch vehicles may have pogo oscillations in liquid engine propulsion
systems. Solid rocket boosters may have thrust oscillations. Both high performance aircraft and launch
vehicles must withstand random vibration due to turbulent boundary layers and shock waves as they
accelerated through the transonic velocity and encounter the maximum dynamic pressure condition.

Ships, automobiles, machine tools, buildings, nuclear reactors, and other mechanical systems and
structures all have their own vibration concerns and failure modes.

5.3.4 Failure Modes

There are many types of potential failure modes including yielding, buckling, ultimate stress, fatigue,
fretting, fastener loosening, relay chatter, and loss of sway space. Engineers must understand these
hazards so that the components, systems, and structures may be designed and tested accordingly.

Resonant excitation and is often a factor in these failures. Engineers thus have a responsibility to identify
equipment natural frequencies through analysis and testing. The natural frequency is the frequency at
which the system would oscillate if it were given an initial displacement and then allowed to vibrate freely.
Resonance occurs when the excitation frequency is at or near the system’s natural frequency. Damping
values, mode shapes, effective modal mass and other dynamic parameters are also needed for this
analysis.

7
5.4 MACHINE HEALTH MONITORING

Figure 5.6. Vibration Monitoring for Predictive Maintenance

Machines, pumps, and other equipment with reciprocating or rotating parts may experience vibration due
to blade passing frequencies, gear mesh frequencies, magnetostriction motor hum, shaft misalignment,
rotating imbalance, etc. In addition, fluid handling machines, like fans and pumps will experience
broadband turbulence.

Some vibration is normal. But excessive vibration may cause accelerated wear and premature failure.
High vibration levels could also indicate a bearing defect. Equipment manufacturers should provide
acceptable limits in terms of amplitude and frequency. The amplitude specification may be in terms of
displacement, velocity or acceleration.

Vibration monitoring is thus needed to reduce maintenance costs, extend life, and improve safety.
Machine vibration can be measured with accelerometers which are permanently mounted on the
machine and monitored continuously with a wireless network. Or a technician may use a handheld
device with an accelerometer that can be temporarily mounted to the machine via a magnetic base, as a
periodic check.

8
5.5 VIBRATION ENERGY HARVESTING

Figure 5.7. Cantilever Piezoelectric Energy Harvesting Kit

Piezoelectric transducers or electromagnetic induction devices can be used to convert ambient vibration
energy to electrical current to charge batteries or to power wireless sensor networks. As an example,
consider a cooling pump in a factory. The pump’s vibration could be used to power its own health
monitoring sensor.

Piezoelectric crystals are asymmetric. The asymmetry in the unit cell of the material sets up the
mechanism whereby deforming the crystal leads to a small potential difference. The piezoelectric
harvesters are typically cantilever beam structures. A mass may be added to the free end of the beam to
tune the device to the source’s dominant vibration frequency. Note that steady vibration with a dominant
frequency is best suited for harvesting. Outputs of 10 to 50 milliwatts are possible, depending on the
vibration amplitude and frequency.

Another potential benefit of an energy harvesting device is that it removes energy from a system, thus
providing some damping which can extend the system’s life.

5.6 HUMAN VIBRATION

There are also biomedical concerns for the case where humans are exposed to vibration. Motion sickness
and seasickness may result from vibration exposure in the 0.1 Hz to 0.5 Hz domain. This may occur
because the brain has difficulty processing apparently conflicting data from the eyes, ear canals, and other
somatosensory organs.

Furthermore, operators of farm equipment, busses, and trains may suffer spinal damage due to long-term
exposure. The human spinal column natural frequency is 10 to 12 Hz. Each organ and bodily part has its

9
own natural frequency. Whole body vibration is addressed in ISO 2631 [2]. Hand-arm vibration is another
concern for operators of power tools.

5.7 VIBRATION EXCITATION SOURCES

5.7.1 Source Characterization

Engineers must also determine the characteristics of the shock or vibration excitation so that they can
properly analyze and test the affected components or structures. Design modifications can then be made
to avoid dynamic coupling between the excitation frequencies and the structure’s natural frequencies.
The excitation sources can be grouped into four types.

5.7.2 Initial Displacement or Velocity

A common example is a guitar string which is given an initial displacement and then suddenly released.
Another is the Pegasus drop transient which is discussed in Section 11.5. An initial velocity example is a
drop shock test machine.

5.7.3 Applied Force

Wind, acoustic pressure, and ocean currents impacting structures are all cases of applied force. So is
footfall excitation on a pedestrian bridge. Shafts with rotating imbalance also impart an applied force
within a machine. Another example is the force from a hammer striking an object. Impulse hammers are
used to purposely excite a structure in modal testing to identify the structure’s natural frequencies,
damping ratios and mode shapes. The force may also be applied by a small shaker attached to the
structure via a stinger rod.

5.7.4 Base Excitation

Seismic excitation is the classic example of base excitation. Another is an automobile traveling over a
speed bump or down a washboard road. Aircraft and launch vehicle avionics components are subjected
to base excitation from their mounting structures during flight. These components are tested on large
shaker tables to verify that they can withstand flight vibration prior to installation into the vehicle.

5.7.5 Self-Excited Vibration

Self-excited vibration is a special case where the alternating force that sustains the motion is created or
controlled by the motion itself. This source type is also referred to as negative damping. Airfoil and bridge

10
flutter are two examples in this category. See the Tacoma Narrows Bridge failure in Section 23. The pogo
oscillation in a launch vehicle with liquid engines is another case.

5.7.6 Source Photo Gallery, Base Excitation

Figure 5.8. Washboard Road, Base Excitation

Ripples can form on gravel and dirt roads with dry, granular road materials. This pattern creates an
uncomfortable ride for the occupants of traversing vehicles and hazardous driving conditions for vehicles
that travel too fast to maintain traction and control. The resulting vibration may also damage suspension
components.

11
Figure 5.9. Loma Prieta Earthquake 1989, Base Excitation

The earthquake caused the Cypress Viaduct to collapse, resulting in 42 deaths. The Viaduct was a raised
freeway which was part of the Nimitz freeway in Oakland, which is Interstate 880. The Viaduct had two
traffic decks. Resonant vibration caused 50 of the 124 spans of the Viaduct to collapse. The reinforced
concrete frames of those spans were mounted on weak soil. As a result, the natural frequency of those
spans coincided with the excitation frequency of the earthquake ground motion. The Viaduct structure
thus amplified the ground motion. The spans suffered increasing vertical motion. Cracks formed in the
support frames. Finally, the upper roadway collapsed, slamming down on the lower road. The remaining
spans which were mounted on firm soil withstood the earthquake.

12
Figure 5.10. Shaker Table Testing, Lateral Base Excitation

A small satellite is mounted to a slip table which is driven by a large electromagnetic shaker. The purpose
of the test is to verify that the satellite can withstand the flight vibration that will be imparted by the
launch vehicle.

13
Figure 5.11. Shaker Table Testing, Vertical Base Excitation

An equipment rack is mounted via wire rope isolators to an expander head which is mounted in turn on
an electromagnetic shaker.

14
5.7.7 Source Photo Gallery, Applied Force

(Image courtesy of the Modal Shop)

A small shaker applies a force

excitation to an automobile fender via
a stinger rod. The applied force and
resulting acceleration at the input
point are measured by an impedance
head transducer. Response
accelerometers may be mounted at
various locations on the vehicle.

Figure 5.12. Automobile Modal Test, Applied Force

Figure 5.13. Thumper Truck, Applied Force

A thumper truck is a vehicle-mounted ground impact system which can be used to provide a seismic
source to perform both reflection and refraction seismic surveys for oil, natural gas and mineral
exploration. A heavy weight is raised by a hoist at the back of the truck and dropped about three meters,
to impact the ground. The resulting ground waveforms are measured with geophones. Some thumpers
use a technology called "Accelerated Weight Drop" (AWD), where high pressure gas is used to accelerate
a heavy weight Hammer (5,000 kg) to hit a base plate coupled to the ground.

15
Figure 5.14. Space Shuttle Orbiter Ku-Band Antenna Dithering, Applied Force to Orbiter

The Ku-band antenna is the disk below the lower Atlantis decal. Dither is a vibration employed in some
mechanical systems to avoid stiction and to ensure smooth motion. Stiction is short for static friction.
The antenna was dithered via a command signal at a frequency of 17 Hz to maintain its ability to smoothly
search for Tracking and Data Relay Satellite System (TDRSS) satellites. There was a concern that the
antenna’s resulting vibration would interfere with sensitive microgravity experiments such as crystal
growth.

16
Figure 5.15. Space Shuttle & Ares I Solid Rocket Boosters

Solid rocket boosters have elongated internal combustion cavities which act as organ pipes. Vortex
shedding within the hot exhaust gases drives standing pressure waves inside these cavities. These
pressure waves have a fundamental thrust oscillation frequency with integer harmonics. The cavities can
be modeled as closed-closed pipes because the nozzle throat diameter is relatively small compared to the
cavity diameter. The Space Shuttle boosters’ thrust oscillation frequency was 15 Hz. The Ares I vehicle
was projected to have a 12 Hz oscillation, but this vehicle was cancelled prior to what would have been
its first flight. Note that the speed of sound in the internal hot exhaust gas is about 3500 ft/sec.

17
5.7.8 Source Photo Gallery, Self-Excited Vibration

(Image: https://space.skyrocket.de/doc_lau_det/titan-2-
glv.htm )

Figure 5.16. Titan II Gemini Vehicle, Pogo Oscillation

Rocket vehicles with liquid engines, such as the Titan II, may experience combustion instability, which
causes excessive vibration forces. This is a potential source of self-excited vibration, whereby the elastic
vehicle structure and the propulsion system form a closed-loop feedback system.

There are several types of combustion instability vibration effects. The most common effect is “Pogo,”
which is similar to Pogo stick motion. In this case, a low frequency oscillation in the combustion chamber,
or propellant feed system, excites the longitudinal vibration mode of the entire rocket vehicle, or some
other structural mode. This may create a cyclical energy exchange between the vibration mode and the
propulsion system oscillation. The problem may also be initiated when a wind gust or some other
perturbation excites the vibration mode. This vibration in turn causes an oscillation in the propulsion
system, which further excites the vibration.

Placing accumulators in the fuel and oxidant lines to damp out the pressure fluctuations solves this Pogo
problem. The accumulator contains a volume of gas that acts like a soft spring to reduce the propellant
frequency to well below that of critical structural frequencies. The accumulator volume must be carefully
selected to meet this goal.

18
6 SINGLE-DEGREE-OF-FREEDOM SYSTEMS & BASIC CONCEPTS

6.1 FREE VIBRATION MODEL

x
x
m m

k c kx

Figure 6.1. Single-degree-of-Freedom System, Spring-Mass Model

The mass value is m. The viscous damping coefficient is c. The spring stiffness from Hooke’s law is k. The
displacement is x. The velocity is x . The equation of motion is derived using Newton’s law.

F = mx (6.1)

m x = − kx − cx (6.2)

mx + c x + kx = 0 (6.3)

The resulting equation of motion is a second-order, ordinary differential equation, linear, homogenous
with constant coefficients. The angular natural frequency in radians/sec is

n = k m (6.4)

The natural frequency in cycles/second or Hz is

f n =  n 2 (6.5)

The natural frequency is the frequency at which the system would oscillate if it were given an initial
displacement and then allowed to vibrate freely. The period T is the inverse of the natural frequency.

T = 1 / fn (6.6)

Furthermore, the damping coefficient divided by mass can be represented as

19
 c / m = 2  n (6.7)

The corresponding amplification factor is

Q = 1 / ( 2 ) (6.8)

The damped natural frequency is

(6.9)
 d =  n 1 − 2

The free vibration solution to equation (6.3) can found using a Laplace transform and applying the initial
displacement and velocity per Reference [3]. The resulting displacement is


  x ( 0 ) + n x ( 0 )  

x(t) = exp ( −n t )   x ( 0 )  cos ( d t ) +   sin ( d t )  (6.10)

  d  


The resulting time history is a damped sinusoid, similar to the Pegasus drop transient flight accelerometer
data shown later in this document in Figure 11.7.

6.2 CONSERVATION OF ENERGY METHOD

The conservation of energy method is based on the time derivative of the total energy of an undamped
system

d
( K E + PE ) = 0 (6.11)
dt

Apply this method to the SDOF system in Figure 6.1 to derive the equation of motion.

d 1 2 1 2
 mx + kx  = 0 (6.12)
dt  2 2 

mxx +kxx = 0 (6.13)

Dividing through by the velocity term yields the equation of motion for the undamped SDOF system.

m x+kx=0 (6.14)

20
6.3 RAYLEIGH METHOD

The Rayleigh method can be used to determine the fundamental frequency of an undamped system by
setting the maximum kinetic energy K E to the maximum potential energy PE .

 K E  max =  PE  max (6.15)

The kinetic energy for the SDOF system in Figure 6.1 is

1
KE = m x2 (6.16)
2

The potential energy for the SDOF system is

1
PE = k x2 (6.17)
2

Assume a displacement of

x(t) = A sin ( n t ) (6.18)

The corresponding velocity is

x ( t ) = n A cos ( n t ) (6.19)

The energy terms become

1
KE = m n 2 A 2 cos 2 ( n t ) (6.20)
2

1
PE = k A 2 sin 2 ( n t ) (6.21)
2

Substitute the two energy terms into equation (6.15).

1 1
m n 2 A 2 = k A 2 (6.22)
2 2

Algebraic simplification yields the expected natural frequency

21
 n = k m (6.23)

6.4 NATURAL FREQUENCY EXAMPLE

An avionics component is modeled as a solid mass per Figure 6.1. The goal is to mount the component
via elastomeric isolator bushings so that the natural frequency is 30 Hz. The isolators will filter out high
frequency vibration energy from the base excitation. Calculate the required isolator stiffness using
equations (6.4) and (6.5).

k = n 2 m = ( 2  f n ) m
2
(6.24)

2  1 lbf sec 2 /in 

  cycles   lbf
k =  2   30  ( 5 lbm )   = 460

(6.25)
  sec    386 lbm  in

The stiffness k in equation (6.25) is the total isolator stiffness. Now assume that the component will be
mounted via four isolators in parallel. The individual isolator stiffness would then be 115 lbf/in. The next
step would be to determine whether isolators with that stiffness value are commercially available.
Otherwise, some adjustment or compromise would be needed.

6.5 OCTAVE RULE

A one octave increase in frequency means that the higher frequency is twice the lower frequency.
Conversely, the lower frequency is one-half the higher one. The number of octaves n between any two
frequencies f 1 and f 2 is

log ( f 2 f 1 )
n= (6.26)
log (2)

Note that a piano keyboard has steps of one-twelfth octave counting both the black and white keys.

The octave rule-of-thumb in mechanical vibration analysis states that there should be at least a one octave
separation between two frequencies to mitigate dynamic coupling. For example, consider an SDOF
system subjected to a harmonic force or base excitation. The system’s natural frequency should be tuned
to less than one-half or greater than twice the excitation frequency.

Isolation or low tuning would be the case where the system’s natural frequency was at most one-half the
excitation frequency. High tuning occurs when the natural frequency is at least twice the excitation

22
frequency. The system is said to be hard-mounted if the natural frequency significantly greater than the
excitation frequency.

As another example, consider two SDOF systems that are to be joined together. The respective natural
frequencies should be separated by on octave prior to mounting one with other. The Pegasus launch
vehicle in Section 11.5 has a natural frequency of about 10 Hz. The payload’s own natural frequency
should be high-tuned to 20 Hz or more to reduce dynamic coupling effects. Note that tuning below 5 Hz
is a poor choice due to the possibility of high relative displacement of the payload within the vehicle
fairing, as well as possible interference with the autopilot control algorithm stability. Further information
about this frequency requirement can be found in the Pegasus Payload User Guide.

6.6 SOME MATHEMATICS

6.6.1 Complex Variables

Complex variables with real and imaginary components are used extensively in signal analysis and
structural dynamics. The imaginary component is expressed as a scale factor applied to j = − 1 .

Let C be a complex variable with real component A and imaginary component B .

C = A + jB (6.27)

The variable C is a vector in the complex, Euclidean plane. The magnitude is the norm.

C = A 2 + B2 (6.28)

The phase angle φ is

 = arctan ( B / A ) (6.29)

6.6.2 Euler’s Equation

Euler’s equation is a complex exponential function used for Fourier transforms and for structural modal
response.

23
 exp (  j  ) = cos (  )  j sin (  ) (6.30)

Complex exponentials can simplify trigonometry, because they are easier to manipulate than their
sinusoidal components.

7 TWO-DEGREE-OF-FREEDOM SYSTEMS

7.1.1 Grounded Two-DOF System, Equation of Motion

m1

Figure 7.1. Two-Degree-of-Freedom System, Spring-Mass Model

24
Newton’s law can be applied to the system in Figure 7.1 to derive the equations of motion. The steps are
omitted for brevity. The resulting pair of ordinary differential equations can be represented in matrix
form and are coupled via the damping and stiffness matrices.

 m1 0   x1  c1 + c2 −c2   x1   k1 + k 2 −k 2   x1  0

 0 m   x  +  −c +
c2   x 2   −k 2
=
k 2   x 2  0
(7.1)
 2 2  2

A shorthand form is

Mx + Cx + K x = 0 (7.2)

The matrices and vectors are

m 0  c + c −c2  k + k 2 −k 2   x1 
M= 1  , C= 1 2  , K= 1  , x=  (7.3)
 0 m2   −c 2 c2   −k 2 k2  x 2 

Consider the undamped, homogeneous form of equation (7.2).

M x+Kx = 0 (7.4)

7.1.2 Grounded Two-DOF System, Eigenvalues & Vectors

Seek a harmonic solution of the form

x = q exp ( jt ) (7.5)

x = j q exp ( jt ) (7.6)

x = −2 q exp ( jt ) (7.7)

The q vector is the generalized coordinate vector. The undamped, homogeneous equation is transformed
through substitution and algebraic manipulation into the generalized eigenvalue problem.

25
  K − 2M  q = 0 (7.8)

The eigenvalues 2 can be found by setting the determinant equal to zero.

 
det K − 2 M = 0
(7.9)

There is an eigenvalue for each degree-of-freedom. Each angular natural frequency is then calculated
from the square root of the respective eigenvalue. The corresponding eigenvectors represent orthogonal
mode shapes. The eigenvector for each mode is found via

 K −  i2 M  q i = 0 , i = 1, 2
(7.10)

An eigenvector matrix Q can be formed. The eigenvectors are inserted in column format.

q q12  (7.11)
Q =  q1 | q2  =  11 
q 21 q 22 

The coefficient q i j represents the modal displacement of mass i for mode j .

Each eigenvector can be multiplied by an arbitrary scale factor. A mass-normalized eigenvector matrix
Q̂ can be calculated such that the following orthogonality relations are obtained.

ˆT MQ
ˆ =I (7.12)
Q

ˆT KQ
ˆ = (7.13)
Q

The superscript T represents matrix transpose. The identity matrix is I , a diagonal matrix of ones. The
diagonal matrix of eigenvalues  .

7.1.3 Grounded Two-DOF System Example, Free Vibration

Now define a modal coordinate (t ) in terms of the normalized eigenvector matrix such that the
displacement vector is

26
 ˆ ˆ
ˆ  =  q11 q12  
x=Q qˆ  (7.14)
 21 qˆ 22 

Substitute equation (7.14) into equation (7.2).

ˆ  + CQ
MQ ˆ  + KQ
ˆ = 0 (7.15)

Premultiply by the transpose of the normalized eigenvector matrix.

ˆ T MQ
ˆ +Q
ˆ T CQ
ˆ +Q
ˆ T KQ
ˆ =Q
ˆ T0 (7.16)
Q

The orthogonality relationships yield

ˆ T CQ
I +Q ˆ  +  = 0 (7.17)

Furthermore, the following simplifying assumption is made for the damping matrix.

ˆ  =  2 1 1
ˆ T CQ 0  (7.18)
Q  0 2 2 2 


The equation of motion can now be written as

1 0   1   2 1 1 0   1  12 0   1  0  (7.19)
+
0 1     0 +   = 
   2  2 2 2  2   0 22  2  0 


The two equations are now uncoupled in terms of the modal coordinates. The modal displacement for
free vibration response to initial conditions is found via Laplace transforms as shown in Reference [4].
The modal displacement for dof i is

 
( )
(
 i (t) = exp − i  i t  i (0)cos d,i t +)
1
(
 i  i  i (0) +  i (0)  sin d,i t
  ) (7.20)

 d,i 


27
The damped angular natural frequency is

(7.21)
 d, i =  i 1 −  i 2

The physical displacement can then be found by substituting equation (7.20) into (7.14). The physical
displacements can be expressed in terms of the eigenvector coefficients and modal coordinates as

x 1(t) = qˆ11 1(t) + qˆ12  2(t) (7.22)

x 2(t) = qˆ 21 1(t) + qˆ 22  2(t) (7.23)

The modal displacement initial conditions are required for the complete solution. These can be found
from the following transformations.

( 0) = Q
ˆ T M x (0) (7.24)

( 0) = Q
ˆ T M x (0) (7.25)

7.1.4 Grounded Two-DOF System, Participation Factors & Effective Modal Mass

The participation factors and effective modal masses are indicators of how excitable the modes are given
uniform base excitation to the given system. Some modes are more readily excited than others, and some
cannot be excited at all for this excitation type.

The system’s generalized mass matrix m̂ is given by

(7.26)
m̂ = QT M Q

Again, the generalized mass will be the identity matrix if the eigenvectors are mass normalized. Let r be
the influence vector which represents the displacements of the masses resulting from static application
of a unit ground displacement. Define a coefficient vector L as

(7.27)
L = QT M r

28
The modal participation factor matrix  i for mode i is

i = Li m
ˆ ii (7.28)

The effective modal mass m eff , i for mode i is

(7.29)
m eff , i = L i 2 ˆ ii
m

7.1.5 Two-DOF System Example, Response to Initial Displacement

The system in Figure 7.1 has the following parameters.

Table 7.1. Grounded Two-DOF System, Parameters

Variable Value
m1 2 lbm
m2 1 lbm
k1 150 lbf/in
k2 100 lbf/in
x1 0.01 in
x2 0.005 in

Assume 5% damping for each mode. The initial velocity is zero for each mass.

Apply the mass and stiffness values from Table 7.1 into equation (7.1) to form the following matrices.

 2 / 386 0   250 −100

M=  , K=  (7.30)
 0 1/ 386  −100 100 

The 1/386 factor is needed to convert lbm to lbf sec^2/in.

The eigenvalues and vectors are found by inserting these matrices into equations (7.9) and (7.10). The
results are shown in Table 7.2.

29
Table 7.2. Grounded Two-DOF System Modal Parameters

Mode Angular Natural Natural Participation Modal Mass

Frequency (rad/sec) Frequency (Hz) Factor Ratio
1 125.25 20.0 0.085 0.935
2 266.72 42.5 -0.022 0.065

The mass-normalized mode shapes in matrix format are

 8.93 −10.64 
Q̂ =   (7.31)
15.05 12.63 

The response time histories for the initial conditions are shown in Figure 7.2.

Figure 7.2. Grounded, Two-DOF System, Response to Initial Displacement

30
7.1.6 Semi-definite, Two-DOF System

7.1.6.1 Semi-definite, Two-DOF System Equation of Motion

x1

x1 x2
m1 k ( x1 - x2 )
k x2

m1 m2
k ( x1 - x2 ) m2

Figure 7.3. Ungrounded, Two-DOF System

Consider the system in Figure 7.3. This system could represent a simple, two-node finite element model
of a rod’s longitudinal vibration. A characteristic of this system is that the fundamental mode is a rigid-
body mode at zero frequency. Both masses move in unison for the rigid-body mode. Newton’s law can
be applied to the system to derive the equations of motion. The steps are omitted for brevity. The
resulting pair of ordinary differential equations can be represented in matrix form.

 m1 0   x1   k −k   x1  0
 0 m   x  +  −k k   x  = 0  (7.32)
 2  2    2  

The eigenvalues and vectors are found using the method shown in Section 7.1.2. Each angular natural
frequency is then calculated from the square root of the respective eigenvalue.

1= 0 (7.33)

k ( m1+ m 2 )
 2= (7.34)
m1m 2

The mass-normalized eigenvectors, representing the mode shapes, are

31
 1 1
q1 = 1 (7.35)
m1 + m2 

1  m2 
q2 =   (7.36)
m1m 2( m1 + m2 )  −m1 

Again, the two masses move in unison for the first mode. They move 180 degrees-out-phase for the
second mode. The semi-definite system is revisited in Section 13.3.3.

7.1.6.2 Semi-definite, Two-DOF System Example

Consider the system in Figure 7.3 with the values in Table 7.3.

Table 7.3. Semi-definite, Two-DOF System, Parameters

Variable Value
m1 2 lbm
m2 1 lbm
k 2000 lbf/in

The mass and stiffness matrices are

 2 / 386 0   2000 −2000

M=  , K=  (7.37)
 0 1/ 386  −2000 2000 

The 1/386 factor is needed to convert lbm to lbf sec^2/in.

The eigenvalues and vectors are found by inserting these matrices into equations (7.9) and (7.10). The
angular natural frequencies are

1= 0 (7.38)

 2 = 1076 rad / sec (171.3 Hz) (7.39)

The mass-normalized eigenvectors in column format are

32
 11.343 −8.021 (7.40)
Q =  q1 | q2  = 
 11.343 16.042 

The total mass sum is 3 lbm (0.00777 lbf sec^2/in) which is also the effective modal mass of the rigid-body
modes. The second mode has zero effective modal mass.

8 CONTINUOUS SYSTEMS WITH LONGITUDINAL VIBRATION

8.1 ROD LONGITUDINAL NATURAL FREQUENCY

8.1.1 One-Dimensional Longitudinal Vibration Equation of Motion

u (x, t)

Figure 8.1. Slender Rod

The longitudinal displacement u ( x, t ) in a rod in undamped free vibration is governed by the second
order, partial differential equation

  u   2u
EA(x)  = L (x) 2
x 
(8.1)
x  t

The term EA(x) is the product of the elastic modulus and cross-sectional area. The equation for a
uniform rod is

2 u 2
   u
=  (8.2)
x 2  E  t 2

Note that

33
  = L A (8.3)

Equation (8.2) is a common formula for one-dimensional wave propagation. Similar equations govern the
propagation of sound in a pipe and the torsional vibration in a shaft. Note that a wave is the phenomenon
in which physical energy propagates through space relative to a medium. Wave propagation is discussed
throughout this document, as shown for seismic waves in Section 19.6.1 as for example.

The method of separation of variables can be applied to equation (8.2) as

u(x, t) = U(x) exp ( j t ) (8.4)


u(x, t) = j  U(x)exp ( j t )
t (8.5)

2
u(x, t) = −2 U(x)exp ( j t ) (8.6)
t 2

Substitute equations (8.4) and (8.6) into (8.2). The resulting spatial equation after some simplification is

d2 
U(x) + 2   U(x) = 0 (8.7)
dx 2
E

Longitudinal waves are non-dispersive. The group and phase wave speeds are equal. The wave speed is
calculated from the elastic modulus and mass density.

CL = E (8.8)

The speed of sound is related to the frequency and wavelength by

CL = f  (8.9)

The spatial equation can be expressed as

d2 2
U(x ) + U( x) = 0 (8.10)
dx 2 CL 2

34
The spatial solution is

 x   x  (8.11)
U(x) = a1 sin   + a 2 cos  
 CL   CL 

The a i values in equation (8.4) are coefficients that depend on the boundary conditions.

The boundary condition for a fixed left end is

U ( 0) = 0 (8.12)

The boundary condition for a free left end is

d (8.13)
U(x) =0
dx x =0

The same boundary conditions could be applied at the right end. The longitudinal natural frequencies are
given in Table 8.1.

Table 8.1. Longitudinal Beam or Rod Natural Frequencies

Configuration Natural Frequencies Mode Shapes

 ( 2n − 1)  x 
Fixed-Free fn =
( 2n − 1 ) CL , n = 1, 2,3
U(x) = sin  
4 L  2L 

fn =
n CL
, n = 1, 2,3  n x 
Free-Free U(x) = cos  
2 L  L 

fn =
n CL
, n = 1, 2,3  n x 
Fixed-Fixed U(x) = sin  
2 L  L 

The mode shapes are normalized to a value of one. Note that the Free-Free beam also has a rigid-body
mode at zero frequency. Further solution details are given in Reference [5]. Also, longitudinal vibration
is revisited in Section 20.2.

The previous method of separation of variables is a modal solution approach. As an alternative, a wave
solution can be applied to equation (8.2) as

35
 u(x, t) = A sin ( kx −  t −  ) (8.14)

The wavenumber k is related to the angular frequency, longitudinal wave speed and wavelength as

k =  c = 2  (8.15)

8.1.2 Rod Longitudinal Vibration Energy Formulas

Natural frequency and modes shapes can also be derived from a rod’s energy terms via the Rayleigh-Ritz
or finite element method. The following formulas are given for reference.

The total strain or potential energy P E of a uniform rod is

2
1 L  du 
2 0
PE = EA   dx (8.16)
 dx 

The total kinetric energy K E of a uniform rod is

2
1 L  du 
2 0
KE = EA   dx (8.17)
 dx 

8.1.3 Beam Longitudinal Vibration Example

Pendulum
Hammer
Test Item
Beam

Figure 8.2. Longitudinal Beam for Shock Testing

The suspended beam can be idealized as a free-free beam for longitudinal vibration. The pendulum
hammer is raised to some initial angular displacement and then released. The hammer strikes the beam’s

36
end plate on the left end, delivering a force impulse. The longitudinal modes of the beam amplify the
input energy and deliver a base input shock pulse to the test article mounted near the right end.

Bai and Thatcher described this shock test method in Reference [6]. This method is still used today,
although the pendulum hammer is often replaced by a pneumatic gun firing a projectile into the beam’s
end plate.

The goal was to choose the beam length so that its fundamental frequency would match the “knee
frequency” of a shock response spectrum specification at 2000 Hz. (Further information on this type of
specification is given in Section 19.) The length for an aluminum beam to meet this goal is calculated as

1 CL 200, 000 in / sec

L= = = 50 in  4 ft (8.18)
2 f n 2 ( 2000 cycles / sec )

Note that the longitudinal wave speed in aluminum is about 200,000 in/sec.

8.2 COIL SPRING SURGE

Figure 8.3. Automotive Valve Springs & Lifters

37
Spring surge is another form of longitudinal vibration, which arises because springs have both stiffness
and inertia properties. This is in contrast to the usual simplifying assumption in vibration analysis that
springs are massless. Springs used in high-speed machinery must have natural frequencies well in
excessive of the frequency of motion that they control. Otherwise, the spring itself may resonate,
resulting in loss of engine performance or even catastrophic failure. This can occur in automotive value
springs when the engine camshaft RPM speed is increased well above the normal operating speed. Race
car engines are at particular risk for this problem. Spring surge is also a potential concern when coil
springs are used as isolation mounts for equipment.

The following spring natural frequencies are taken from Reference [7]. Note the similarities with the
longitudinal rod formulas in Table 8.1.

Figure 8.4. Fixed-Free Spring (hand drawn)

The fundamental frequency of the fixed-free spring is

 2n − 1  k
f n =  , n = 1, 2,3 (8.19)
 4  L L

Its corresponding normalized mode shape is

 ( 2n − 1)  x 
U(x) = sin   (8.20)
 2L 

The wave speed in the spring is

38
 CL = k L L (8.21)

The wave speed equation applies for other spring boundary conditions as well.

Figure 8.5. Fixed-Fixed Spring

The fundamental frequency of the fixed-fixed spring is

n k
fn= , n = 1, 2,3 (8.22)
2 L L

The corresponding normalized mode shape is

 n x 
U(x) = sin   (8.23)
 L 

Note that the Fixed-Fixed and Free-Free cases have the same natural frequency equation.

39
 L

Figure 8.6. Spring with End Mass

The simplified fundamental frequency formula for a spring with and end mass M is

1 k
fn= (8.24)
2 M + 0.33 L L

A more rigorous formula is

 k  M 
fn= n , cot (  n ) =    n , n = 1, 2,3 (8.25)
2 L L  L L 

The corresponding mode shape is

sin (  n x / L ) , n = 1, 2,3,...
(8.26)

8.3 RING FREQUENCY

40
Figure 8.7 Ring Mode from Finite Element Model

Consider a thin ring with a rectangular cross section and with completely free boundary conditions. The
ring frequency corresponds to the mode in which all points move radially outward together and then
radially inward together. This is the first extension mode. It is analogous to a longitudinal mode in a rod.
The ring frequency is the frequency at which the longitudinal wavelength in the skin material is equal to
the vehicle circumference.

CL
fr = (8.27)
d

The ring frequency is an important concept for launch vehicles. The front end of a launch vehicle is
typically composed of cylindrical shell segments which house avionics components and the payload.
These shells tend to have notable vibration responses at their respective ring frequencies due to external
acoustic environments and stage separation shock events. Furthermore, there is a rule-of-thumb in
statistical energy analysis (SEA) that a cylinder tends to behave as a flat plate above its ring frequencies.
This is a simplifying assumption. Also note that the ring frequency is an idealized concept for a cylindrical
shell. In practice, cylindrical shells tend to have a high modal density near the ring frequency.

Consider an aluminum launch vehicle cylindrical shell with a diameter of 72 inches. Again, the longitudinal
wave speed in aluminum is approximately 200,000 in/sec. The ring frequency is

f r = ( 200,000 in/sec ) (  60 in ) = 1060 Hz (8.28)

41
9 BERNOULLI-EULER BEAMS
9.1.1 Beam Equations

y(x, t)
EI, L

Figure 9.1. Cantilever or Fixed-Fixed Beam

y(x, t)
EI, L

Figure 9.2. Simply-Supported or Pinned-Pinned Beam

The governing equation for beam bending free vibration is a fourth order, partial differential equation.

2 y  2 y  2 y
−  EI(x)  =  L (9.1)
 x 2   x 2  t2

The term EI(x) is the stiffness which is the product of the elastic modulus and area moment of inertia.
The equation for a uniform beam is

4 y 2 y
−EI = L (9.2)
x4 t2

The method of separation of variables can be applied as

42
 y(x, t) = Y(x)exp ( j t ) (9.3)

y(x, t) = j Y(x)exp ( j t ) (9.4)
t

2
y(x, t) = −2 Y(x)exp ( j t ) (9.5)
t 2

Substitution into equation (9.2) leads to spatial ordinary differential equations

d4 2  L 
Y(x) −    Y(x) = 0 (9.6)
dx 4  EI 

Define the wavenumber β as

  (9.7)
4 = 2  L 
 EI 

The differential equation can be expressed as

d4
4
Y(x) − 4 Y(x) = 0 (9.8)
dx

The spatial solution has the form

Y(x) = a1 sinh ( x ) + a 2 cosh ( x ) + a 3 sin ( x ) + a 4 cos ( x ) (9.9)

Substitute the spatial solution into the differential equation.

4 a1 sinh ( x ) + a 2 cosh (x ) + a 3 sin (x ) + a 4 cos (x )


− 2  a1 sinh ( x ) + a 2 cosh (x ) + a 3 sin (x ) + a 4 cos (x ) = 0
 EI  (9.10)

Equation (9.10) is satisfied by the wavenumber relationship in equation (9.7), which gives credibility to
the assume spatial solution. The a i values in equations (9.9) and (9.10) are coefficients that depend on
the boundary conditions discussed in 9.1.2. The spatial solution gives eigen function modes shapes. It
also gives eigenvalue roots from which the natural frequencies are calculated.

43
The eigenvalues are represented as n L . The angular natural frequencies n are

n = n2 EI L (9.11)

9.1.2 Common Boundary Conditions

9.1.2.1 Cantilever Beam

Recall the cantilever or fixed-free beam in Figure 9.1. The displacement and slope at the fixed end are
both zero.

d
Y(0) = 0 , Y(0) = 0 (9.12)
dx

The moment and shear force at the free end are both zero.

 d2  d  d2 
M =  EI(x) 2 Y(x)  =0 , V=  EI(x) Y(x)  =0
 dx  dx  dx 2  (9.13)
x =L x =L

Application of these boundary conditions to the spatial displacement in equation (9.9) and its derivatives
yields the following equation for finding the a i coefficients.

 sin ( L ) + sinh ( L ) cos (L ) + cosh (L )   a1  0 

    
   =  
cos ( L ) + cosh ( L ) − sin (L ) + sinh (L )  a 2  0 
(9.14)

By inspection, equation (9.14) can only be satisfied if a1 = 0 and a 2 = 0 . Set the determinant of the
matrix to zero in order to obtain a nontrivial solution. This yields the following transcendental equation
for finding the roots by the Newton-Raphson method. There are multiple roots which satisfy equation. A
subscript is thus added.

cos ( n L ) cosh ( n L ) = −1 (9.15)

The fundamental root for the cantilever beam is n L = 1.87510 . The natural frequencies can be found
by substituting each root in to equation (9.11).

44
9.1.2.2 Pinned-Pinned Beam

Recall the pinned-pinned beam in Figure 9.2. The boundary conditions at the left end of the pinned-pinned
beam are

 d2 
Y(0) = 0 , M =  EI(x) 2 Y(x)  =0 (9.16)
 dx  x =0

The boundary conditions at the right end are similar.

Application of these boundary conditions to the spatial displacement in equation (9.9) and its derivatives
yields the following equation for finding the a i coefficients.

sinh ( L ) sin ( L )   a1  0 
    
   =  
sinh ( L ) − sin ( L )  a 2  0 
(9.17)

By inspection, equation (9.14) can only be satisfied if a1 = 0 and a 2 = 0 . Set the determinant of the
matrix to zero in order to obtain a nontrivial solution. This yields the following transcendental equation.
There are multiple roots which satisfy equation. A subscript is thus added.

sin (  n L ) sinh ( n L ) = 0 (9.18)

The equation is satisfied if

n L = n, n = 1,2,3,.... (9.19)

The natural frequencies can be found by substituting each root in to equation (9.11).

9.1.3 Beam Bending Fundamental Frequencies

The solutions to the beam equation for various boundary conditions are derived in Reference [8]. A
summary of beam formulas is given in

45
Table 9.1. Note that the free-free beam also has a rigid-body mode with zero frequency.

46
Table 9.1 Beam Bending Fundamental Frequencies

Configuration Fundamental Frequency

1  3.5156 
EI L
Cantilever (Fixed-Free) 2  L2 

1 3EI
Cantilever with end mass 2 ( 0.2235 m beam + mend ) L3

2
1 
Pinned-Pinned   EI L
2  L 

1  22.373 
EI L
2  L2 
Free-Free & Fixed-Fixed

1 15.418 
EI L
2  L2 
Fixed-Pinned

The pinned-pinned beam has integer harmonics as follows.

2
1 n  EI
fn = , n = 1,2,3... (9.20)
2  L  

The beam bending frequencies for other configurations have non-integer harmonics. The free-free beam
higher modal frequencies are shown in

47
Table 9.2. The free-free formulas also apply to the fixed-fixed beam. The fixed-free beam higher modal
frequencies are shown in Table 9.3.

48
Table 9.2 Free-Free & Fixed-Fixed Beam Bending Frequencies, Elastic Modes

Mode Natural Frequency (Hz)

1  22.373  EI
f1 =
2  L2  

f2 = 2.757 f 1

f3 = 5.404 f 1

f4 = 8.402 f 1

Table 9.3 Fixed-Free Beam Bending Frequencies, Elastic Modes

Mode Natural Frequency (Hz)

1  3.5156 
f1 = EI L
2  L2 
f2 = 6.268 f 1

f3 = 17.456 f 1

f4 = 34.386 f 1

9.1.4 Beam Bending Participation Factors & Effective Modal Mass

Again, the mode shapes Yn (x) and their corresponding mode shapes are found by applying the boundary
conditions to the displacement shape in equation (9.9). The participation factors are then

L
 n =  L Yn (x) dx
(9.21)
0

49
The effective modal mass is
2
 L  Y (x)dx 
 0 L n 
m eff , n =L (9.22)
 L  Yn (x) dx
2
0

9.1.5 Bending Wave Speed & Wavelength

Figure 9.3. Traveling Bending Wave Packet Dispersion, Snapshots at Three Different Times

Both the red and blue circles oscillate in the vertical axis in Figure 9.3. The red circle travels at the phase
velocity. The blue circle propagates twice as fast at the group velocity and is always at or near the positive
or negative peak in the wave packet. The phase speed is the more important of the two speeds for
vibroacoustics analyses where the bending waves are excited by the external sound field, or vice versa.

50
Bending waves are dispersive as a result of the governing fourth-order partial differential equation. The
wave speed varies with the angular natural frequency. The phase speed for a given frequency is

cp =   EI L 
1/4 (9.23)

The group speed is twice the phase speed.

cg = 2 c p (9.24)

Demonstrating bending wave propagation in terms of these two speeds is best done with an animation,
but the still images in Figure 9.3 may be useful. Furthermore, the phase speed is related to the frequency
and bending wavelength by

cp = f  (9.25)

9.1.6 Beam Bending Energy Formulas

Natural frequency and modes shapes can also be derived from a beam’s energy terms via the Rayleigh-
Ritz or finite element method. The following formulas are given for reference.

The total strain or potential energy P E of a uniform beam is

2
1 L  d2 y 
PE =
2 0 EI  2  dx
 dx 
(9.26)
 

The total kinetric energy K E of a uniform beam is

1 2 L
n    y  dx
2
KE = (9.27)
2 0

51
9.1.7 Beam Example, Wind Chime

Figure 9.4. Author’s Wind Chimes

The wind chimes in the figure are made from anodized aluminum. The tallest one is 34.7 inch long, 1.25
inch diameter, and 0.094 inch wall thickness. The chimes are like free-free beams undergoing bending
vibration when struck. The vibration then generates sound waves at discrete frequencies matching the
bending frequencies. The fundamental frequency of the tallest chime is calculated via the following steps.

  4 1 
d o − di 4  = (1.25 in ) − (1.16 in )  = 0.0321in 4
4 4
I= 
64   64    (9.28)

 lbm   1 lbf sec2 / in      2 lbf sec2

  ( ) ( ) 
2
L =  0.1 3   1.25in − 1.16in = 4.59e-05
in   386 lbm   4  
(9.29)
 in 2

1  22.373  
fn =   1.0e + 07
2  ( 34.7in )  
2
lbf 
in 
2  (
0.0321in 4 ) 4.59e-05
lbf sec 2
in 2
= 247 Hz (9.30)
 

52
The tallest wind chime was struck separately from the others. The resulting sound was recorded and
analyzed. The sound tone represented bending modes. The theoretical value from equation (9.30) agrees
closely with the measured frequency 244 Hz. A comparison of the first four bending modes and their
respective theoretical values from

53
Table 9.2 is given in Table 9.4. There are two orthogonal planes where transverse vibration occurs. Each
mode actually represents a pair of modes. The pair should ideally have the same frequency but may be
slightly different due to manufacturing and material imperfections, as well as two holes in each for the
connecting cord.

Table 9.4. Tallest Chime Elastic Mode Natural Frequencies

Mode Calculated (Hz) Measured (Hz) Musical Note

1 238 244 B

2 657 664 E

3 1288 1273 D#

4 2127 2050 C

54
Figure 9.5. Wind Chime Elastic Mode Shapes

The Y-axis represents transverse displacement which has been normalized to have a maximum absolute
value of one. The points where each curve crosses the zero baseline are nodal points.

55
Figure 9.6. Tallest Wind Chime Time History

The time history shows some beat frequency effect because each mode is really a pair of orthogonal
modes with closely-spaced frequencies.

56
Figure 9.10. Tallest Wind Chime Fourier Transform

The highest response is at the fundamental bending mode at 244 Hz. The fourth bending modes at 2050
Hz is barely excited. Note that Fourier transforms are covered in Section 16.

57
10 THIN PLATE VIBRATION

10.1.1 Rectangular Plate, Bending Vibration

10.1.1.1 Rectangular Plate Equations

y a

0 x

Figure 10.1. Rectangular Plate

The fourth-order governing equation of motion for the out-of-plane displacement is

 4z  4z  4z   2z
D  4 + 2 2 2 + 4  + h 2 = 0
 x x y y  t
(10.1)


The plate stiffness factor is

Eh 3
D=
(
12 1 −  2 ) (10.2)

The kinetic energy in the plate is

58
 h  2 b a 2
2 0 0
KE = Z dx dy (10.3)

The potential or strain energy in the plate is

 2  2  2  2  2  2   2  
2
D b a
PE =     Z  +   Z  + 2   Z   Z  + 2 (1 −  )   Z   dxdy
2 0 0  x 2   y 2   x 2  y 2   x y   (10.4)
         

Assume a harmonic displacement response for the partial differential equation (10.1).

z(x, y, t) = Z(x, y)exp( jt) (10.5)

Substitution of the assumed response leads to the spatial equation

 4Z 4Z 4Z 

D  4 + 2 2 2 + 4  − h2 Z = 0
 x x y y 
(10.6)


The modal solution equation (10.6) for most boundary conditions requires an approximation method such
as the Rayleigh-Ritz or finite element method. Further information on plate vibration solutions is given in
Reference [9].

10.1.1.2 Rectangular Plate, Simply-Supported along All Edges

The case of plate with simply-supported boundary conditions on all four edges can be solved for exactly.
The boundary conditions for this plate are

Z(x, y) = 0, M x ( x, y) = 0 for x = 0, a (10.7)

Z(x, y) = 0, M y (x, y) = 0 for y = 0, b
(10.8)

The displacement and bending moment along each edge is equal to zero.

59
Assume the following displacement function which satisfies the boundary conditions, where m and n are
indices. The mass-normalized modes shapes for the plate are

2  m x   n y 
Zmn = sin   sin  
a b h  a   b  (10.9)

The angular natural frequencies are

D   m   n  
2 2
 mn =   +  
h  a   b   (10.10)


10.1.1.3 Rectangular Plate Example

Consider a rectangular circuit board made from G10 material. The modulus of elasticity is 2.7e+06 lbf/in^2
with Poisson ratio of 0.12. The dimensions of the circuit board are 4 in x 2 in x 0.063 in. The board is
simply-supported about its perimeter. Assume that the electronic components do not add any stiffness.

The board has a uniform mass distribution. The total mass is 0.113 lbm. This includes the G10 board and
the electronic components.

The mass per area is

 h = 0.0141 lbm / in 2 = 3.65e − 05 lbf sec2 / in 3 (10.11)

The circuit board stiffness factor is

D=
( 2.7e + 06 lbf / in 2 ) ( 0.063 in )
3
= 57.1 lbf in
12 (1 − 0.12 2 )
(10.12)

The fundamental angular natural frequency is

57.1 lbf in    2   2 
11 =   +   = 3875 rad / sec
3.65e − 05 lbf sec 2 / in 3   4 in  2 in   (10.13)
  

60
The fundamental frequency is equivalent to 613.7 Hz. The mode shape is shown in Figure 10.2.

Figure 10.2. Circuit Board, Fundamental Bending Mode Shape

10.1.2 Circular Plate, Bending Vibration

10.1.2.1 Circular Plate Equations


X
0

Figure 10.3. Circular Plate

The governing fourth-order spatial equation of motion for the out-of-plane displacement is

61
 4 Z ( r,  ) − 4 Z ( r,  ) = 0 (10.14)

The biharmonic operator is

 2 1  1 2   2 1  1 2 
 4 =  2 2 =  2 + + 2 2  2 + + 2 2 
 r r r r   r r r r   (10.15)
  

The wave number is

1/4
 2 h 
=  (10.16)
 D 

The kinetic energy in the plate is

h 2 2 R 2
2 0 0
KE = Z r dr d (10.17)

The potential or strain energy in the plate is

 2
D 2 R   2 Z 1 Z 1  2 Z   2 Z  1 Z 1  2 Z 
( )
2 0 0   r 2 r  r r 2 2 
PE =  + +  − 2 1 −  
2  r r
+ 2 2
  r  r  

2
   1 Z  
+ 2 (1 −  )     r dr d
  r  r    


(10.18)
The solution of equation (10.14) proceeds by assuming separation of the spatial variables.

Z(r, ) = R(r) () (10.19)

Further information on this solution method is given in Reference [9].

62
10.1.2.2 Circular Plate, Simply-Supported

The boundary conditions for the simply-supported plate are

Z(a, ) = 0 (10.20)

Mr = 0 at r = a (10.21)

The mass-normalized, fundamental mode shape is

Z00 =
1.945
m
 J0 ( 2.2215 r / a ) − 0.03686 I0 ( 2.2215 r / a )  (10.22)

Note that I 0 and J 0 in equation (10.22) are Bessel functions. This fundamental mode shape is sometimes
referred to as an “oil can” mode. The corresponding fundamental frequency is

4.979 D
00 = , for  = 0.33 (10.23)
a 2 h

A characteristic equation for variable mode numbers and Possion ratios is given in Reference [9]. The
roots of this equation determine the numerical coefficient on the righthand side of equation (10.23).

10.1.2.3 Circular Plate Example, Simply-Supported

Assume a 48 inch diameter, 0.5 inch thick, aluminum circular plate, with a simply-supported
circumference. Calculate the fundamental frequency and mode shape.

The plate stiffness factor is

D=
(1.0e + 07 lbf / in 2 ) ( 0.5 in )
3
= 1.169e+05 lbf in
12 (1 − 0.332 )
(10.24)

The mass per area is

h = ( 0.000259 lbf sec2 /in 4 ) ( 0.5 in ) =1.295e-04 lbf sec2 /in3 (10.25)

63
The angular natural frequency is

4.979 1.169e+05 lbf in

00 = = 259.7 rad / sec
( 24 in ) 2
1.295e-04 lbf sec2 /in 3 (10.26)

The angular frequency is equivalent to 41.33 Hz.

Figure 10.4. Circular Plate Fundamental Mode Shape

64
10.1.3 Honeycomb Sandwich Plate

10.1.3.1 Honeycomb Sandwich Plate Characteristics

Figure 10.5. Sample Honeycomb Sandwich Plates

t1

Top Skin

Core

Bottom Skin

t2

Figure 10.6. Honeycomb Sandwich Geometry

Honeycomb sandwich structures are designed to have a high stiffness-to-mass ratio. The stiff, strong face
sheets carry the bending loads, while the core resists shear loads. The face sheets are typically made from
aluminum or carbon fiber with epoxy resin. The honeycomb core material is usually aluminum for
aerospace applications. Other core materials include Nomex aramid or Kevlar para-aramid fiber sheets
saturated with a phenolic resin. In addition, closed cell foams such as Rohacell are substituted for
honeycomb in some sandwich panel designs.

65
Table 10.1. Honeycomb Sandwich Plate Normal Modes

Domain Characteristic

Low Frequencies Bending of the entire structure as if were a thick plate

Mid Frequencies Transverse shear strain in the honeycomb core governs the behavior

High Frequencies The structural skins act in bending as if disconnected

The modal response of honeycomb sandwich plates has the characteristic shown in Table 10.1. The lower
natural frequencies may be calculated using the plate formulas in Sections 10.1.1 and 10.1.2 by properly
accounting for the material and geometric properties.

Assume

1. The skin elastic modulus is much greater than the core modulus.
2. Each skin has the same material.
3. Each skin is "thin" relative to the core.

The honeycomb sandwich plate bending stiffness from Reference [10] is

 E   t t h2  1
D= f 2  1 2  , h = h c +  t1 + t 2  (10.27)
1 −    t1 + t 2  2

The elastic modulus of each skin or face sheet is E f .

66
10.1.3.2 Honeycomb Sandwich Plate Example, Part I

A suborbital rocket vehicle has a circular bulkhead made from honeycomb. The material is aluminum for
both core and skins. The bulkhead properties are shown in Table 10.2. Calculate the fundamental
frequency of the bare bulkhead.

Table 10.2 Bulkhead Parameters

Parameter Value
Boundary Condition Simply Supported
Diameter 40 inch
Core Thickness 1.0 inch
Thickness of Each Skin 0.063 inch
Total Thickness 1.125 inch
Skin Elasticity 10.0e+06 lbf/in^2
Core Elasticity Negligible
Skin Density 0.10 lbm/in^3
Core Density 0.01 lbm/in^3
Poisson's Ratio 0.3

The plate stiffness is calculated in two steps.

1
h = 1.0 in +  0.063 + 0.063 in = 1.063 in (10.28)
2

10e + 06 lbf / in 2   (0.063in)(0.063in) (1.063in )2 

D=    = 3.99e + 05 lbf in (10.29)
 1 − 0.332   0.063in + 0.063in 
 

The total skin mass is


( )
mskin = 2   ( 40in ) ( 0.063 in ) 0.1 lbm / in 3 = 15.83 lbm
4
2
(10.30)

The core mass is

67
 
( )
mcore =   ( 40in ) (1in ) 0.01 lbm / in 3 = 12.57 lbm
4
2
(10.31)

The total bulkhead mass is

m total = 15.83 lbm + 12.57 lbm = 28.4 lbm (10.32)

The total mass per area is

28.4 lbm
A = = 0.0226 lbm / in 2 = 5.854e − 05 lbf sec2 / in 3
(  4 )( 40in ) 2 (10.33)

The angular fundamental frequency of the bare bulkhead for μ=0.33 is

4.979 3.99e + 05 lbf in

00 = = 1028 rad / sec
( 20 in ) 2
5.854e − 05 lbf sec2 /in 3 (10.34)

The angular frequency is equivalent to 164 Hz. The mode shape is qualitatively similar to that in Figure
10.4.

10.1.3.3 Honeycomb Sandwich Plate Example, Part II

Now assume that 100 lbm of avionics components are added to the bulkhead and that this nonstructural
mass is spread uniformly across the bulkhead. Neglect the geometry and stiffness of the avionics. Repeat
the natural frequency calculation.

The total mass per area is

A =
(100 + 28.4 ) lbm = 0.1022 lbm / in 2 = 2.647e − 04 lbf sec2 / in3
(  4 )( 40in )2 (10.35)

The angular fundamental frequency of the mass-loaded bulkhead for μ=0.33 is

68
 4.979 3.99e + 05 lbf in
00 = = 483.3 rad / sec
( 20 in ) 2
2.647e − 04 lbf sec2 /in 3 (10.36)

The equivalent frequency is 76.9 Hz.

11 DAMPING

11.1 DAMPING METRICS

The following metrics are demonstrated by example throughout this document.

The viscous damping ratio is defined in terms of the damping coefficients as

 = c cr (11.1)

The critical damping coefficient is

cr = 2 k m (11.2)

The amplification factor is

Q = 1 / ( 2 ) (11.3)

The amplification factor can be calculated from measured frequency response function data via the half-
power method.

Q = f f (11.4)

where f is the peak center frequency and f is the difference between the two -3 dB points on either
side of the peak.

The loss factor is

 = 2 (11.5)

69
The logarithmic decrement is

2 (11.6)
=
1 − 2

11.2 FOUR DAMPING CATEGORIES

The four damping types are summarized in Table 11.1. The description is for the free vibration of a single-
degree-of-freedom system due to initial displacement or velocity.

Table 11.1. Damping Types

Type Value Response Description

Undamped =0 Simple harmonic motion, sinusoidal response

Underdamped 0   1 Damped sine response with exponential decay

Critically damped  =1 Border between the overdamped and underdamped cases

Overdamped  1 Sum of two decaying exponentials with no oscillation

The first three of these types are shown by examples in Figure 11.1. The underdamped type emphasized
throughout this document.

70
Figure 11.1. SDOF Response to Initial Displacement for Three Damping Cases

71
11.3 DAMPING MECHANISMS

Damping occurs as vibration energy is convert to heat, sound or some other loss mechanism. Damping is
needed to limit the structural resonant response. Common sources are:

1. viscous effects
2. Coulomb damping, dry friction
3. aerodynamic drag
4. acoustic radiation
5. air pumping at joints
6. boundary damping

The dominant source for assembled structures is usually joint friction. This damping mechanism may be
nonlinear due to joint microslip effects. The damping value tends to increase at higher excitation levels.

Damping may also decrease as the natural frequency increases, such that the amplification factor
increases with natural frequency. The equation per Reference [11] for approximating Q for an electronic
system subjected to a sine base input is
0.76
 fn 
Q = A  (11.7)
 ( G )0.6 
 in 

where

A = 1.0 for beam-type structures

= 0.5 for plug-in PCBs or perimeter supported PCBs
= 0.25 for small electronic chassis or electronic boxes

Gin Sine Base Input (G)

Beam structures: several electronic components with some interconnecting wires or cables.

PCB: printed circuit board well-populated with an assortment of electronic components.

Small electronic chassis: 8-30 inches in its longest dimension, with a bolted cover to provide access to
various types of electronic components such as PCBs, harnesses, cables, and connectors.

72
11.4 HUNTSVILLE, ALABAMA PEDESTRIAN BRIDGE DAMPING

Figure 11.2. University Drive Pedestrian Bridge

Figure 11.3. Slam Stick X, Triaxial Accelerometer & Data Logger, Shown on Book Shelf

The bridge in Figure 11.2 is near the University of Alabama in Huntsville. The author walked to the center
of the bridge and mounted the triaxial accelerometer in Figure 11.3 on the deck floor using double-sided
tape. He was the only person on the bridge. The ambient vibration response of the bridge was negligibly
low in each of the axes, while he remained standing. The wind was very light on this day. The author then
jumped up-and-down on the bridge to excite its vibration modes. The response in the vertical axis was
significant. The response in each of the lateral axes remained negligibly low.

73
Figure 11.4. Signal Identification via Damped Sine Curve-fit

The jumps were performed near the 80 second mark. The accelerometer data was band-passed filtered
from 1 to 20 Hz. The response after 82 seconds represents the free vibration decay. A damped sine curve-
fit synthesis was performed on the vertical acceleration time history as shown in Figure 11.4, using trial-
and-error with convergence. The bridge’s fundamental frequency is 2.2 Hz with 0.16% damping. The
damping ratio is very low. But note that the fundamental frequency and damping may be nonlinear. A
pedestrian’s vertical forces correspond to each footfall, and typically occur at 2.0 Hz. This is very close to
the bridge’s 2.2 Hz natural frequency.

Also note that there is a potential for pedestrians to synchronizing their steps with the bridge motion and
with one another. This behavior is instinctive rather than deliberate. Pedestrians find that walking in
synchronization with the motion of a bridge is more comfortable, even if the oscillation amplitude is
initially very small. This cadence makes their interaction with the movement of the bridge more
predictable and helps them maintain their balance. But the synchronization also causes the pedestrians’
gait to reinforce the bridge’s oscillation in a resonance-like manner. These sorts of problems occurred
after the opening of the London Millennium Bridge in 2000. Both passive and tuned mass dampers were
added to the bridge for vibration control. Retrofitting the University Drive Bridge with dampers is
unnecessary due to its low pedestrian traffic volume. But this could be an interesting project for the
nearby engineering students.

74
11.5 PEGASUS LAUNCH VEHICLE DAMPING
11.5.1 Drop Transient

Figure 11.5. Pegasus Launch Vehicle, Drop & Stage 1 Burn

Figure 11.6. Pegasus Fundamental Bending Mode Shape, Exaggerated

A modified L-1011 aircraft carries the Pegasus vehicle up to an altitude of nearly 40,000 ft and a speed of
Mach 0.8, as shown in Figure 11.5. Pegasus is suspended underneath the aircraft by hooks, where it
develops an initial displacement due to gravity. The strain energy is suddenly released at the onset of the
drop transient, causing Pegasus to oscillate nearly as a free-free beam. This is a significant “coupled loads”
event for the payload which is enclosed in the fairing at the front end of the vehicle. The Pegasus first
stage then ignites. The payload is eventually delivered into a low earth orbit. A certain Pegasus/payload
configuration was analyzed via a finite element modal analysis. The resulting wire mesh model of the
fundamental bending mode is shown in Figure 11.6. Note that this was a different Pegasus configuration
than that represented by the flight data in Figure 11.7.

75
11.5.2 Drop Transient Damped Sine Curve-fit

Figure 11.7. Pegasus Drop Transient, Flight Accelerometer Data, Free Vibration Response

The data was measured in the transverse axis at the payload interface. The response is nearly a textbook
quality damped sinusoid, with an exponential decay. The natural frequency and damping ratio are
identified via a synthesized damped sine, curve-fit.

76
11.5.3 Drop Transient Logarithmic Decrement Method

Figure 11.8. Pegasus Drop Transient, Flight Accelerometer Data, Logarithmic Decrement

The logarithmic decrement for a starting peak and a peak n cycles later is

1  xo 
= ln   (11.8)
n  xn 

The logarithmic decrement for the two peaks shown in Figure 11.8 is

1  0.7473  (11.9)
 = ln   = 0.0773
8  0.4026 

The logarithmic decrement value is equivalent to 1.2% damping in agreement with the previous damped
sine curve-fit. The logarithmic decrement method has been included for historical reasons. The damped
sine curve-fit method is more robust.

77
11.6 TRANSAMERICA BUILDING DAMPING

Figure 11.9. Transamerica Pyramid

The Transamerica Pyramid is built from a steel frame, with a truss system at the base. The height is 850
ft (260 m). Reference [12] gives natural frequency and damping as obtained in the 1989 Loma Prieta
earthquake and from ambient vibration. The ambient vibration was presumably due to wind, low level
micro-tremors, mechanical equipment, outside street traffic, etc.

78
Table 11.2. Transamerica Pyramid, Modal Parameters

Direction Loma Prieta Earthquake Ambient Vibration

fn (Hz) Damping fn (Hz) Damping

North-South 0.28 4.9% 0.34 0.8%

East-West 0.28 2.2% 0.32 1.4%

The results show non-linear behavior with an increase in damping during the severe earthquake relative
to the benign ambient vibration.

79
12 SINUSOIDAL VIBRATION

12.1 DESCRIPTIVE STATISTICS

Figure 12.1. Sine Function, One-Cycle

A sine function is characterized by its frequency and amplitude. The frequency is the inverse of the period,
where the period is the duration of one cycle. The amplitude may be represented by any of the three
parameters shown on the right-hand-side.

The zero-to-peak and RMS amplitudes for sinusoidal oscillations are related by

A Peak = 2 A RMS (12.1)

80
Note that

A RMS2 = A mean 2 + A12 (12.2)

where A1 is the standard deviation

The RMS value is thus equation to the standard deviation for the case of zero mean. A pure sine function
always has zero mean.

Figure 12.2. Sine Function & Histogram

A sine function and its corresponding histogram are shown in Figure 12.2. The histogram divides the
amplitude points into bins and then counts the points in each bin. The number of points in each bin are
then plotted as a bar graph. Some of the histograms in this document are turned 90 degrees clockwise
and shown adjacent to the corresponding time history for educational purposes.

A sine function has a bathtub histogram. The amplitude tends to dwell at either the peak or valley with a
fairly rapid transition in between. Histograms are important for fatigue analysis.

Now let x(t) be a sinusoidal displacement function. The displacement equation is

81
 x(t) = Xsin ( t ) (12.3)

where

X is the zero-to-peak amplitude

The corresponding velocity v(t) is

v(t) =  X cos ( t ) (12.4)

The acceleration a(t) is

a(t) = −2 Xsin ( t ) (12.5)

Table 12.1. Peak Values Referenced to Peak Displacement

Parameter Value

displacement X
velocity X
acceleration 2 X

Table 12.2. Peak Values Referenced to Peak Acceleration

Parameter Value
acceleration A
velocity A/
displacement A / 2

Let A be the peak acceleration.

A = 2 X (12.6)

Conversely,

82
 X = A / 2 (12.7)

12.2 SINE TESTING

Table 2 has enormous consequences for vibration testing. Commercial shaker tables, such as those in
Figure 5.10 and Figure 5.11, are typically limited to a few inches peak-to-peak displacement unless they
are large-scale tables specifically designed for earthquake simulation. Vibration tests are usually specified
in terms of acceleration. The test specification starting frequency is typically 10 to 20 Hz to limit the
displacement.

12.3 SAMPLE SINE CALCULATIONS

Astronaut Michael Collins wrote [13]:

The first stage of the Titan II vibrated longitudinally, so that someone riding on it would be
bounced up and down as if on a pogo stick. The vibration was at a relatively high frequency,
about 11 cycles per second, with an amplitude of plus or minus 5 Gs in the worst case.

The equivalent angular frequency for 11 Hz is 69.1 rad/sec. Calculate the corresponding displacement X
using equation (12.7).

X = (5 G )(386 in/sec2 / G) / ( 69.1 rad/sec)2 = 0.40 inch zero-to-peak (12.8)

M.J. Griffin [14] cited a study where 5% of the passengers on a ship experienced vomiting while the ship
underwent a 0.2 m/sec^2 acceleration amplitude during a 3 hour. Assume a 0.2 Hz frequency (1.26
rad/sec^2). Calculate the displacement.

X = (0.2 m/sec^2)(100 cm/m) / ( 1.26 rad/sec)2 = 12.7 cm zero-to-peak (12.9)

83
12.4 BEAT FREQUENCY EFFECT

Figure 12.3. Channel Beam Free Vibration Response at Free End to Force Impulse

Consider the case of two closely-space vibration frequencies generating sound. The envelope of the two
signals oscillates at a frequency which is the differences of the two frequencies. This beat frequency is
detected by the human hearing and brain system, but it is a psychoacoustic effect which does not affect
the mechanical hardware.

A channel beam was mounted in a cantilever configuration. The free end was struck with an impulse
hammer. The resulting acceleration time history at the free end is shown in Figure 12.3. A synthesized
damped sine curve-fit yielded natural frequencies at 583 and 691 Hz. The beat frequency is the difference
which is 108 Hz.

84
Figure 12.4. Channel Beam Free Vibration Response Fourier Transform

A Fourier transform of the beat frequency data is shown in Figure 12.4. Peaks occurs 583 and 691 Hz but
not at the 108 Hz beat frequency. Fourier transforms are covered in Section 16.

85
13 DYNAMIC MODELS FOR EXTERNAL EXCITATION

13.1 TWO MODELS

There are two basic SDOF spring-mass models used for external excitation in shock and vibration analysis.
One is the applied force model, and the other is the base excitation model. An introduction to these types
and images were given in Section 5.7. The purpose of this section is to demonstrate solution methods for
each of these type in both the time and frequency domains.

13.2 APPLIED FORCE

13.2.1 SDOF Response Equations

The applied force model and its free-body diagram are

f f
x
x
m m

k c kx

Figure 13.1. Single-degree-of-freedom Model, Applied Force

The equation of motion for applied force F is derived using Newton’s law.

F = mx (13.1)

m x = f (t) − kx − cx (13.2)

86
 m x + cx + kx = f (t) (13.3)

1 (13.4)
x + 2  n x +  n 2 x = f (t)
m

The resulting equation of motion is non-homogenous. It can be solved via Laplace transforms if the force
is deterministic such as a sine function. A convolution integral is needed if the force varies arbitrary with
time. The resulting displacement for an applied force is

0 f () exp −n (t − ) sin d (t − )d

1 t
x(t) = (13.5)
md

The convolution integral in equation (13.5) can be solved in terms of a ramp invariant digital recursive
relationship [15].

Alternate methods for numerically solving equation (13.4) include the Runge-Kutta and Newmark-beta
methods. The Runge-Kutta method may unstable for stiff systems, however.

The displacement convolution integral has an embedded impulse response function.

h d (t) =
1 
md 
(  )
exp −n t  sin ( d t )  (13.6)

The Laplace transform of the impulse response function is the steady-state frequency response function
(FRF).

1 1 
Hd (s) =  2  (13.7)
m  s + 2 n s + n 2 

This FRF can be represented in the frequency domain by setting s = j  , where j = −1 .

 
1 1 
H d () = (13.8)
 n ( )
m   2 − 2 + j2  
n 


The FRF in equation (13.8) is the Fourier transform of the impulse response function in equation (13.7).
Fourier transforms are covered in Section 16. An equivalent FRF is

87
  
1 n 2  (13.9)
H d () =
(
 n )
k   2 − 2 + j2  
n 


13.2.2 Input & Response Example, Applied Force

The single-degree-of-freedom system in Figure 13.1 has a mass of 20 kg, stiffness of 78960 N/m, natural
frequency of 10 Hz and an amplification factor Q=10 for this example. It is excited into resonance by a
20 N applied force sinusoidal excitation at 10 Hz. The exact response can be calculated via a Laplace
transform solution. The result is shown in Figure 13.2.

Figure 13.2. Sine Applied Force & Response

The peak response is 2.53 mm. Note that the following ratio is approximately equal to the Q value.

kx
=
( 78960 N/m )( 0.00253 m )  10 (13.10)
F 20 N

88
Figure 13.3. Transfer Function for Steady-State Force Excitation

The response calculation could be repeated for a number of frequency and amplification factor cases to
determine the maximum steady-state response for each case. The transfer function represents the
displacement divided by the applied force as normalized by the stiffness such that the magnitude is
unitless. The transfer function magnitude is

kx 1
= ,  = f / fn
( )
F 2 2
(13.11)
+ ( 2 )
2
1− 

The phase angle by which the displacement lags the force is

 (
 = arctan  2  1 − 2 
 ) (13.12)

The applied force transfer function in Figure 13.3 can roughly be divided into three zones. The function is
controlled by stiffness at the low frequency ratio end. Its magnitude converges to the simple static
Hooke’s law as the frequency ratio approaches zero. The function is controlled by damping for frequency
ratios at or near the resonance condition where the ratio is equal to 1. The function is controlled by mass
at higher ratios, above √2.

89
13.2.3 Six Transfer Functions for Applied Force and Response

Recall equation (13.8). The following equation represents the ratio of displacement to force in terms of
their respective Fourier transforms.

X ( )
= Hd () (13.13)
F ( )

Similar frequency response ratios can be formed for the velocity and acceleration responses relative to
force. Inverses of these functions can also be taken. The names of these transfer functions are given in
Table 13.1 and Table 13.2.

Table 13.1. Transfer Function Nomenclature

Displacement / Force Velocity / Force Acceleration / Force

Admittance,
Accelerance,
Compliance, Mobility
Inertance
Receptance

Table 13.2. Reciprocal Transfer Function Nomenclature

Force / Displacement Force / Velocity Acceleration / Force

Apparent Mass,
Dynamic Stiffness Mechanical Impedance
Dynamic Mass

Furthermore, the receptance frequency response function can be modeled in terms of a block diagram.
Similar diagrams could be drawn for the other transfer functions.

Input Force Transfer Function Displacement Response

Figure 13.4. FRF Diagram with Input & Output

90
13.2.4 Helicopter Ground Resonance

f
x

k c

Figure 13.5. Helicopter Preparing for Takeoff and SDOF System Analogy

A sketch of a generic helicopter is shown in Figure 13.5. The helicopter may be modeled as a single-
degree-of-freedom system excited by an applied force.

A new helicopter design must undergo thorough testing to ensure the reliability of the design with respect
to vibration. It may encounter severe vibration while it is on the ground, preparing for takeoff. A similar
problem may occur immediately after landing.

Many helicopters have fully-articulated rotor systems. The advancing blades flap upward and the
retreating blades flap downward in order to maintain uniform lift during forward flight. A helicopter with
a fully-articulated rotor has some measure of an unbalanced rotational force, which causes an oscillating
force. This force may excite one of the helicopter’s natural frequencies, as the rotor accelerates to its full
operating speed prior to takeoff. This condition is called resonant excitation.

John O. Emmerson wrote in Reference [16]:

In the mid-1950s, the emerging helicopter introduced a frightening problem called ground
resonance. A perfectly sound production helicopter firmly supported on its landing gear with
the rotor turning at or near flight rpm would suddenly initiate a “tramping” on the gear and
then, almost immediately, completely disintegrate. The Sikorsky company, which was leading
the pack in developing the helicopter, quickly gathered the most data on these ground-based
disasters and thus shouldered the lion’s share of corrective understanding and action. The
company established that the lead/lag pivots on the rotor blades, the pylon elastic stiffness, and
the elastic components of the landing gear, including the tire and shock struts, all working
together as a single structural entity, had an unfortunate resonant frequency near the rotor
operating rpm. The engineers were able to show that a fearless pilot, entering the resonant

91
 phase, could save the day be pulling up the collective pitch control and unloading the landing
gear (i.e., getting airborne). Lowering the resonant frequency and adding damping offered a
permanent solution.

Today, all of us in the helicopter development business, with the apprehension borne of
substantial ignorance, approach initial ground whirl testing of each new helicopter with
tremulous respect for the demon of ground resonances.

Modern helicopters avoid ground resonance by using dampers on the blades and on the landing gear.
The landing gear may have shock absorbed struts.

Figure 13.6. TH-55 Osage, Military Version of the Hughes 269A

Figure 13.7. TH-55 Osage, Destroyed by Ground Resonance

An intact TH-55 helicopter is shown in Figure 13.6. A similar TH-55 helicopter was destroyed by ground
resonance in March 1967, at Ft. Wolters, TX, as shown in Figure 13.7.

92
13.2.5 Two-Degree-of-Freedom System Subjected to Applied Forces

13.2.5.1 Two-Degree-of-Freedom System Equations of Motion, Time Domain

Figure 13.8. Two-DOF System with External Forces

Recall the grounded, two-DOF system from Section 7.1.1. It is shown again in Figure 13.8 with added
forces. The coupled equations of motion are

 m1 0   x1  c1 + c2 −c2   x1   k1 + k 2 −k 2   x1   f 1 
 0 m   x  +  −c +
c2   x 2   −k 2
= 
k 2   x 2  f 2 
(13.14)
 2 2  2

Note that each of the coefficient matrices is symmetric. This is typical for linear structural dynamics
systems. A shorthand form of equation (13.4) is

93
 Mx + Cx + K x = f (13.15)

The matrices and vectors are

m 0  c + c −c 2  k + k 2 −k 2  f 1   x1 
M= 1  , C= 1 2  , K= 1  , f = , x=  (13.16)
 0 m2   −c 2 c2   −k 2 k2  f 2  x 2 

The natural frequencies and mode shapes are found from the undamped, homogeneous equation as was
shown in Section 7.1.1. Now define a modal coordinate (t ) in terms of the normalized eigenvector matrix
Q̂ such that the displacement vector is

ˆ ˆ
ˆ  ( t ) =  q11 q12   ( t )
x(t) = Q qˆ  (13.17)
 21 qˆ 22 

Repeat the steps in Section 7.1.1 with the added force term.

ˆ T MQ
Q ˆ +Q
ˆ T CQ
ˆ +Q
ˆ T KQ
ˆ =Q
ˆ Tf (13.18)

The orthogonality relationships yield

ˆ T CQ
I +Q ˆ  +  = Q
ˆ TF (13.19)

The equation of motion can now be written with the modal damping simplification as

1 0   1   2 1 1 0   1  12 0   1   qˆ11 qˆ 21   f 1 
+
0 1     0 +   = 
   2  2 2 2   2   0 22   2   qˆ12 qˆ 22   f 2 
 (13.20)

The two equations are now uncoupled in terms of the modal coordinates. Now consider the case where
the forcing functions are deterministic time domain functions. The modal displacement for the response
to the harmonic forces can be found via Laplace transforms. Alas, the resulting modal displacement
equations even for sinusoidal forces are too unwieldy for inclusion in this book but are given in Reference
[17]. The physical displacements are then found via equation (13.17). The solution for the case of

94
arbitrary forces can be calculated using the ramp invariant digital recursive filtering relationship in
Reference [15].

13.2.5.2 Two-Degree-of-Freedom System, Time Domain Example

The system in Figure 13.8 has the parameters previously shown in Table 7.1, with 5% damping for each
mode. Its natural frequencies are again 20.0 and 42.5 Hz. The mass-normalized mode shapes in matrix
format are

 8.93 −10.64 
Q̂ =   (13.21)
15.05 12.63 

It is subjected to a single sinusoidal force applied to the top mass, with zero initial conditions. The
amplitude is 1 lbf. The forcing frequency is 20.0 Hz, which is the same as the fundamental frequency. The
results are shown in Figure 13.9. Now consider the steady-state amplitude ratios between the mass
responses. The responses are very nearly in proportion to the eigenvector coefficients in the first column
of equation (13.21). The response of mass 2 is thus nearly (15.05/8.93) times greater than that of mass
1. The measured ratio is slightly inexact became the 20.0 Hz excitation also excites the 42.5 Hz mode, but
the second mode’s participation can be considered as negligibly low due to the greater than one-octave
frequency separation.

The analysis is repeated with the forcing frequency at 42.5 Hz, as shown in Figure 13.10. The steady-state
response of mass 2 is nearly (12.63/10.64) times than that of mass 1, per the second column of equation
(13.21). The negative sign in the second column indicates that the responses are 180 degrees out-of-
phase with one another.

95
Figure 13.9. Displacement Responses to 20.0 Hz Excitation

Figure 13.10. Displacement Responses to 42.5 Hz Excitation

96
13.2.5.3 Two-Degree-of-Freedom System Equations of Motion, Frequency Domain

Equation (13.20) can be used in the frequency domain for the case of steady-state harmonic forces.
Define force and modal displacement coordinates as follows

f ( t ) = Fexp ( jt ) (13.22)

 ( t ) = N exp ( jt ) (13.23)

Substitute these harmonic terms into Equation (13.20). This process is akin to taking a Fourier transform
of the equation. The result for modal displacement i after simplification is

(
 i )
  2 − 2 + j2    N = qˆ F + qˆ F
i i  i 1i 1 2i 2 (13.24)

qˆ 1i F1 + qˆ 2i F 2 (13.25)
Ni =
( )
  2 − 2 + j2   
 i i i 

The modal displacement equation can be generalized to a system with n degrees-of-freedom as

 
n  q̂ p i Fp 
Ni =   
( )
(13.26)
   2 − 2 + j2    
p =1  i i i 
 

The physical displacements are then found from

ˆ ˆ
ˆ N ( ) =  q11 q12  N ( )
X ( ) = Q qˆ  (13.27)
 21 qˆ 22 

The physical displacement at degree-of-freedom k is

97
 n
Xk =  qˆ ki N i (13.28)
i =1

The principle of linear superposition holds in this case. The displacement can be calculated for an
individual force. Then the total displacement is the sum of the displacements for the individual forces.

Assume that the only force acting on the system is at degree-of-freedom p . The modal displacement at
degree-of-freedom i is

q̂ p i Fp
Ni = (13.29)
( )
  2 − 2 + j2   
 i i i 

The physical displacement at degree-of-freedom k for the single force case is

 
n  q̂ ki q̂ p i Fp 
Xk =   
( )
(13.30)
   2 − 2 + j2    
i =1  i i i 
 

The receptance function for the single force applied becomes

 
Xk n  q̂ ki q̂ p i 
Hd, kp () = =   (13.31)
Fp  2
( ) 
i =1   i −  + j2  i  i  

2


The numerator in the series shows the relationship between the receptance function and the mass-
normalized eigenvector terms.

98
13.2.6 Two-Degree-of-Freedom System, Frequency Domain Example

Figure 13.11. Two-DOF System Receptance FRF

The system in Figure 13.8 has the parameters previously shown in Table 7.1, with 5% damping for each
mode. Its natural frequencies are again 20.0 and 42.5 Hz. Its receptance frequency response function
from equation (13.31) is given in Figure 13.11 for four permutations. The first index in the legend is the
response location. The second is the applied force location, where each location is a degree-of-freedom.
The phase angles could also be plotted but are omitted for brevity.

99
13.3 BASE EXCITATION

13.3.1 SDOF Response Equations

This case is also referred to as support motion or seismic excitation. The base excitation model and its
free-body diagram are

m m

k c

k(y-x)

Figure 13.12. Single-degree-of-freedom Model, Base Excitation

The base displacement is y, and base velocity is y . The equation of motion for base excitation is derived
using Newton’s law.

F = mx (13.32)

m x = −k ( y − x ) − c ( y − x ) (13.33)

m x + cx + kx = cy + ky (13.34)

The challenge with equation (13.34) is that the base displacement and base velocity are usually unknown,
but the base acceleration is known. The workaround is to define a relative displacement z as

z=x−y (13.35)

Through substitution and algebraic manipulation, the following ordinary differential equation is derived
in terms of the relative displacement.

100
 m z + cz + kz = − m y (13.36)

The resulting equation of motion is non-homogenous. The m y term on the right-hand side can be
regarded as a virtual force. The equation for base excitation can be further simplified as

z + 2 n z + n 2 z = − y (13.37)

The absolute acceleration is then calculated as

x =z+y (13.38)

Equation (13.37) can be solved via Laplace transforms if the base acceleration is deterministic such as a
sine function. A convolution integral is needed if the excitation varies arbitrarily with time. The resulting
relative displacement equation for acceleration base excitation is

0 y() exp −n (t − ) sin d (t − )d

1 t
z(t) = − (13.39)
d

Equation (13.37) can also be solved in the frequency domain. Time and frequency domain analyses should
yield the same results per Parseval's theorem.

The relative displacement convolution integral has an embedded impulse response function.

h z (t)=
1 
d 
(  )
exp −n t  sin ( d t )  (13.40)

The Laplace transform of the impulse response function is the steady-state transmissibility function for
relative displacement.

 1 
Hz (s)=  2  (13.41)
 s + 2n s + n 
2

This transmissibility function can be represented in the frequency domain by setting s = j  , where
j= −1 .

101
  
H z () =  1 
( )
  2 − 2 + j2   (13.42)
 n n 


The relative displacement transmissibility function in equation (13.42) is the Fourier transform of the
impulse response function in equation (13.40). Fourier transforms are covered in Section 16.

The absolute acceleration response convolution integral is

0 y() exp −n (t − ) ( d ) 

1 t
x(t) = 2
− 2n 2 sin d (t − ) + 2n d cos d (t − ) d
d

(13.43)

The convolution integral in equation (13.43) is computationally inefficient, however. An alternative is to

use the Smallwood ramp invariant digital recursive filtering relationship [18]. The recursive filtering
algorithm is fast and is the numerical engine used in almost all shock response spectrum software. It is
also accurate assuming that the data has a sufficiently high sample rate and is free from aliasing. The
algorithm is discussed further in Section 19.5.

The embedded impulse response function in equation (13.43) is

h a (t) =
1 
d 
( ) ( )
exp −n t   d 2 − 2n 2 sin ( d t ) + 2n d cos ( d t )

(13.44)

The Laplace transform of the impulse response function is the steady-state transmissibility function for
absolute acceleration.

 2 n s + n 2 
Ha (s) =  2  (13.45)
 s + 2 n s + n 
2

The transmissibility function in the frequency domain is

 
 n 2 + j2n  
H a ( ) =
( )
(13.46)
  2 − 2 + j2  
 n n 


102
13.3.2 Input & Response Example, Base Excitation

The single-degree-of-freedom system in Figure 13.12 has a natural frequency of 10 Hz and an amplification
factor Q=10 for this example. It is excited into resonance by a 1 G sinusoidal base excitation at 10 Hz. The
exact response can be calculated via a Laplace transform solution. The result is shown in Figure 13.13.

Figure 13.13. Sine Base Input & Response

103
Figure 13.14. Transmissibility Function for Steady-State Sine Base Input

The response calculation could be repeated for a number of frequency and amplification factor cases to
determine the maximum steady-state response for each case. The transmissibility results for three
amplification factor cases are shown in Figure 13.14. The frequency ratio is the excitation frequency f
divided by the natural frequency fn. The magnitude is the steady response divided by the input
amplitude. The transmissibility tends to unity gain when the excitation frequency is much less than the
natural frequency. Resonant amplification occurs when the frequency ratio is at or very near one.
Isolation occurs when the excitation frequency is greater than √2 times the natural frequency.

Note that the Q value is approximately equal to the peak transmissibility at resonance in Figure 13.14.
Also, the resonant time domain response in Figure 13.13 reached a steady-state peak of Q times the base
input. This suggests that sine vibration testing could be a way for measuring the amplification factor and
its corresponding damping ratio by noting the peak response either in the frequency or time domain. But
this simple relationship should be considered only for the special case of a single-degree-of-freedom

104
system. Multi-degree-of-freedom systems and continuous systems have responses which depend on
damping as well as on mode shape parameters.

Consider a cantilever beam subjected to base excitation exciting its fundamental bending mode. An
accelerometer is mounted on the beam’s free end, and another at its midpoint. The peak response at
free end will be higher than that at the midpoint due to the respective location-dependent modal gain
factors. The half-power bandwidth method is a much more reliable damping identification method for
measured response data from this beam and from general structures. This method would yield the same
damping value from either of the two accelerometer data results. See the examples in Sections 14.4 and
14.5.

The transmissibility function plotted in Figure 13.14 was introduced in Section 13.3. The magnitude H (  )
can be represented as

1 + ( 2 )
2
H () = ,  = f / fn
( )
(13.47)
2 2
+ ( 2 )
2
1− 

The phase angle by which the response lags the input is

 
23
 = arctan  
( )
(13.48)
 1 +  2 4 2 − 1 
 

13.3.3 Alternate Base Excitation Method using Seismic Mass

Figure 13.15. Seismic Mass for Base Excitation, Semi-definite System

Recall the semi-definite system in Section 7.1.6 and the two-DOF system subjected to an applied force in
Section 13.2.5. A hybrid of these systems can be configured for base excitation by setting the bottom

105
mass to a value vastly larger than the top mass by a factor of, say, 1e+04. Then an external force is applied
to the base to produce the desired acceleration. This approach is unnecessary for SDOF system base
excitation but becomes very useful for complex, multi-degree-of-freedom systems. The equation of
motion is

 m1 0   x1   c −c   x1   k −k   x1  f   m1 y 
 0 m   x  +  −c c   x  +  − k k   x  =  0  =  0  (13.49)
 2  2   2   2    

The hybrid equation is solved using the methods in Section 13.2.5. There is a rigid-body mode at zero
frequency. The second modal frequency is

k ( m1+ m 2 )
 2= (13.50)
m1m 2

The second modal frequency becomes

2  k m 2 , f or m 1  m 2 (13.51)

The second frequency is essentially equal to the top mass and the spring as a single-degree-of-freedom
system.

Repeat the example in Section 13.3.2. The natural frequency was 10 Hz in this example with Q=10. Use
of the seismic mass method requires explicit knowledge of the mass and stiffness values, as shown for
this example in Table 13.3. The top mass is excited into resonance by a 1 G sinusoidal base excitation at
10 Hz via an equivalent 1e+05 lbf force applied to the bottom, seismic mass. The acceleration response
is shown in Figure 13.16. The result is the same as that in Figure 13.13 with the response time histories
in both plots reaching 10 G steady-state.

Table 13.3. Seismic Mass, Hybrid System, Parameters

Variable Value
m1 1.0e+05 lbm
m2 10 lbm
k1 102.2 lbf/in

106
Figure 13.16. Sine Base Input & Response, Seismic Mass Method

Note that in some cases it may be necessary to ground the seismic mass with a very soft spring for
numerical stability purposes, such that the fundamental frequency is just slightly above the zero frequency
but is still much less than the frequency in equation (13.51).

107
13.3.4 Two-DOF Response Equations, Time Domain

Figure 13.17. Two-DOF System Subjected to Base Excitation

Consider the system in Figure 13.17, which is sometimes referred to as two-stage isolation per Reference
[19]. The equation of motion can be derived using Newton’s law.

 m1 0   x1  c1 + c2 −c2   x1   k1 + k 2 −k 2   x1  c1 y + k1 y 

 0 m   x  +  −c +
c2   x 2   −k 2
=
k 2   x 2    (13.52)
 2 2  2 0 

The familiar challenge with equation (13.52) is that the base displacement and base velocity are usually
unknown, but the base acceleration is known. The workaround is to define a relative displacement z
such that

108
 zi = xi − y (13.53)

Substitution of the relative displacements into (13.52) yields the following equation after algebraic
manipulation.

 m1 0   z1  c1 + c2 −c2   z1   k1 + k 2 −k 2   z1   −m1y 

 0 m   z  +  −c +
c2   z 2   −k 2
=
k 2   z 2   −m2 y 
(13.54)
 2 2  2

The natural frequencies and modes shapes are calculated using the method from Section 7.1.2. The modal
decoupling process is essentially the same as that in Section 13.2.5, with force vector changed to

 −m y 
F= 1  (13.55)
 −m2 y 

The relative displacement is represented in terms of the mass-normalized eigenvector matrix and the
modal displacement.

ˆ ˆ
ˆ  =  q11 q12   ( t )
z =Q qˆ  (13.56)
 21 qˆ 22 

The uncoupled equations of motion are

1 0   1   2 1 1 0   1   2 0   1   qˆ11 qˆ 21   − m1y 
0 1    +  0   +  1   = 
   2  2 2 2   2  0 2  
2  2   qˆ12 qˆ 22   − m 2 y 


(13.57)

The modal displacement can be calculated via Laplace transforms for a deterministic base acceleration.
Modal accelerations can be calculated by twice differentiating the modal displacements.

A ramp invariant digital filtering relationship can be used for a base excitation which varies arbitrarily with
time. Such equations are available for direct calculation of both the modal displacement and acceleration.
The filtering equations are taken from Reference [15] given the virtual force and displacement
characteristics of equation (13.57), rather than the SDOF formulas for base excitation given in References
[18], [20].

The physical relative displacement can then be calculated by applying the modal displacements to
equation (13.56). The physical relative acceleration can then be determined from the modal acceleration

109
via the same eigenvector matrix transformation in equation (13.56). The physical absolute acceleration
is then found by adding the base acceleration to the physical relative acceleration.

13.3.5 Two-DOF Response Equations, Time Domain Example

The system in Figure 13.17 has the parameters previously shown in Table 7.1 with 5% damping for each
mode. Its natural frequencies are again 20.0 and 42.5 Hz. It is subjected to sinusoidal base excitation,
with zero initial conditions. The base acceleration amplitude is 1 G with frequency at 20.0 Hz, which is
the same as the fundamental frequency. The results are shown in

Figure 13.18.

The analysis is repeated with the excitation frequency at 42.5 Hz, as shown in

Figure 13.19. The responses are 180 degrees out-of-phase with one another.

110
Figure 13.18. Acceleration Responses to 20.0 Hz Excitation

111
Figure 13.19. Acceleration Responses to 42.5 Hz Excitation

13.3.6 Two-DOF Response Equations, Frequency Domain

Equation (13.57) can be used in the frequency domain for the case of steady-state harmonic base
excitation. Define base acceleration and modal displacement coordinates as follows

y ( t ) = YA exp ( jt ) (13.58)

 ( t ) = N exp ( jt ) (13.59)

112
Substitute these harmonic terms into Equation (13.57). The result for modal displacement i after
simplification is

(
 i )
  2 − 2 + j2    N = − qˆ m + qˆ m Y
i i  i 1i 1 2i 2 A ( ) (13.60)

Ni =
(
− qˆ 1i m1 + qˆ 2i m 2 YA ) (13.61)

( i 2
)
− 2 + j2  i  i 

The physical relative displacements are then found from

ˆ ˆ
ˆ N ( ) =  q11 q12  N ( )
Z ( ) = Q qˆ  (13.62)
 21 qˆ 22 

The physical relative displacement at degree-of-freedom k is

2 2  q̂ ki qˆ 1i m1 + qˆ 2i m 2 
( )
Zk =  qˆ ki N i = − YA  
)

(
(13.63)
i =1   i −  + j2  i  i  
2 2
i =1
 

The transmissibility function relating physical relative displacement to base excitation is

Zk 2  (
 q̂ ki qˆ 1i m1 + qˆ 2i m 2 ) 

H rd ( ) = = − 
)

(
(13.64)
i =1   i −  + j2  i  i 
YA 2 2 
 

The transmissibility function relating physical relative acceleration to base excitation is

113
 2  q̂ ki qˆ 1i m1 + qˆ 2i m 2
( ) 

H ra ( ) = 2  
)

(
(13.65)
i =1   i −  + j2  i  i 
2 2 
 

The transmissibility function relating physical absolute acceleration to base excitation is

2  q̂ ki qˆ 1i m1 + qˆ 2i m2
( ) 

Haa ( ) = 1 + 2  2 2
)

(
(13.66)
i =1   i −  + j2  i  i  
 

13.3.7 Two-DOF Response Equations, Frequency Example

Figure 13.20. Two-DOF System Transmissibility

114
The system in Figure 13.17 has the parameters previously shown in Table 7.1, with 5% damping for each
mode. Its natural frequencies are again 20.0 and 42.5 Hz. The transmissibility function for each mass
relative to the base is given in

Figure 13.20 as calculated per equation (13.66). The damping could be determined by apply the half-
power bandwidth method to the curves if this were measured data.

14 SINE SWEEP TESTING

14.1 SINE SWEEP INTRODUCTION

Figure 14.1. Sample Sine Sweep Time History with Very Brief Duration for Educational Purposes Only

Structures and components may be mounted on shaker tables, similar to those in Figure 5.10 and Figure
5.11, and subjected to sine sweep base excitation. The essence of a sine sweep test is that the base
excitation input consists of a single frequency at any given time. The frequency itself, however, is varied
with time. The sine sweep test may begin at a low frequency and then sweep to a high frequency, or vice-

115
versa. Some specifications require several cycles, where one cycle is defined as from low to high
frequency and then from high back to low frequency.

The specification might require either a linear or a logarithmic sweep rate. The sweep will spend greater
time at the lower frequency end if the sweep is logarithmic. The example in Figure 14.1 had a logarithmic
sweep rate and a constant amplitude. Nevertheless, the specification might require that the amplitude
vary with frequency.

Some potential test objectives are

1. Identify natural frequencies and amplification factors or damping ratios

2. Perform sine sweep before and after random vibration test to determine if any parts
loosened, etc.
3. Check for linearity of stiffness and damping by varying the input amplitude
4. Workmanship screen for defective parts and solder joints
5. Represent an actual environment such as a rocket motor oscillation
6. NASA/GSFC typically uses sine sweep vibration for spacecraft testing

14.2 SINE SWEEP TIME DOMAIN EQUATIONS

14.2.1 Linear Sweep Rate

The normalized amplitude function for a linear sweep is

  t 
Y(t) = sin  ( f 2 − f1 )   + f1   (14.1)
  T 

The test duration is T. The starting and ending frequencies are f1 and f 2 , respectively.

14.2.2 Logarithmic Sweep Rate

The normalized amplitude function for a logarithmic sweep is

( )(
  f −1 + 2Rt
 )  
Y(t) = sin 2 
1
(14.2)
 
 R ln ( 2 ) 

   

116
The sweep rate in octaves per time is R . Recall the octave formula Section 6.5. A typical sweep rate is
1 octave/minute for sine sweep tests intended for natural frequency identification.

14.2.3 Sweep Type Comparison

Figure 14.2. Sweep Rate Comparison Example

A comparison of the frequency time variation for two sweep rate types is shown in Figure 14.2. The
sweeps were performed over the frequency domain from 20 to 2000 Hz over a 400 second duration. This
is equivalent to a rate of 1 octave/minute for the logarithmic sweep.

14.3 SINE SWEEP SPECIFICATION EXAMPLE

A vendor has a product that must withstand sinusoidal vibration with an amplitude of 12 G, over the
frequency domain from 10 Hz to 2000 Hz, with a logarithmic sweep rate.

The shaker table has a displacement limit of 1.0 inch peak-to-peak, or 0.5 inch zero-to-peak. Recall that
the displacement limit is a constraint at low frequencies. How should the test be specified? The answer
is to use a specification with two amplitude segments as shown in Figure 14.3. The first segment is a

117
constant displacement ramp, 0.5 inch zero-to-peak. The second segment is a constant acceleration
plateau at 12 G.

Figure 14.3. Sample Sine Sweep Specification

The cross-over frequency is calculated from the square root of the ratio of the peak acceleration to the
peak displacement. The respective peaks should be in zero-to-peak format.

1 X peak
f cross = (14.3)
2 X peak

118
The cross-over frequency for the specification in Figure 11.2 is

 in / sec 2 
(12 G )  386 
1  G  (14.4)
fcross = = 15.3 Hz
2 0.5in

This is the frequency at which the 12 G acceleration corresponds to a 0.5 in zero-to-peak displacement. It
is the knee frequency in Figure 14.3.

14.4 SINE SWEEP INPUT & RESPONSE EXAMPLE, SDOF SYSTEM

The SDOF system in Figure 13.12 has a natural frequency of 40 Hz and amplification factor Q=10. It is
subjected to a base input of 1 G from 10 to 80 Hz, with a logarithmic rate of 1 octave/min. The input and
response time histories are shown in Figure 14.4. The peak response is 10 times the peak input
acceleration. This corresponds to Q=10, but this is the special of a single-degree-of-freedom excited into
resonance. A better method for damping identification is the half-power bandwidth method, which is
particularly needed for multi-degree-of-freedom and continuous systems due to modal gain factors.

A transmissibility function can be made relating the response to the input in the frequency domain. This
is done by tracking the frequency in terms of the zero-crossing rate. The result is shown in Figure 14.5
along with the half-power bandwidth damping calculation.

Each of the data points in the Figure 14.5 are approximations constrained by the resolution and plot
cursor limitations. The points at 37.85 and 41.91 Hz are each half-power points, 3 dB less than the
resonant peak. The amplification factor is the frequency difference between the half-power points divided
by the center frequency.

f 40.05 Hz
Q= = = 9.86  10 (14.5)
f ( 41.91 − 37.85) Hz

Recall that Q = 10 corresponds to 5% damping.

119
Figure 14.4. Sine Sweep Input & Response, Total Duration & Close-up Views

120
Figure 14.5. SDOF Transmissibility Function from Input & Response Time Histories

121
14.5 SINE SWEEP INPUT & RESPONSE EXAMPLE, CONTINUOUS BEAM

y(x, t)

w(t)

Figure 14.6. Cantilever Beam, Base Excitation

Table 14.1. Cantilever Beam Parameters

Material Aluminum
Length 24 inch
Cross-section Circular, 1 inch diameter
Area Moment of Inertia 0.04909 in^4
Elastic Modulus 1.0e+07 lbf/in^2
Mass/Length 0.0785 lbm/in (2.035e-04 lbf sec^2/in^2)
Amplification Factor Q=10 for all modes

Consider the cantilever beam in Figure 14.6 with the properties in Table 14.1. The fundamental frequency
of the cantilever beam using the formula in Table 9.3.

1  3.5156 
f1 =  
(1.0e + 07lbf / in 2 )( 0.04909in 4 ) = 47.7 Hz (14.6)
2  ( 24 in )2  2.035e − 04 lbf sec / in
2 2
 

The first four natural frequencies and their associated parameters are shown in Table 14.2.

122
Table 14.2. Cantilever Beam Modal Parameters

Mode Natural Participation Modal Mass

Frequency (Hz) Factor Ratio
1 47.7 2.99E-03 0.646
2 299 9.20E-04 0.198
3 837 3.16E-04 0.068
4 1641 1.62E-04 0.035

Figure 14.7. Cantilever Beam Response Time History at Free End

The cantilever beam is subjected to a 1 G sine sweep from 20 to 2000 Hz at a logarithmic rate of 1
octave/minute. The acceleration response is shown in Figure 14.7, as calculated using Reference [21].
The corresponding transmissibility function is shown in
Figure 14.8.

123
Figure 14.8. Cantilever Beam Transmissibility at Free End

Again, each of the modes had an amplification factor of 10, but each of the four peaks in
Figure 14.8 has an amplitude other than 10. The peaks depend on both the amplification factor and the
mode shape parameters. Nevertheless, the half-power bandwidth method can be used to extract the
proper amplification values as show for the first and second modes in the following equations.

The amplification factor for the first bending mode is

f 47.86 Hz
Q= =  10 (14.7)
f ( 50.15 − 45.36 ) Hz

124
The amplification factor for the second bending mode is

f 299.9 Hz
Q= =  10
f ( 314.1 − 284 ) Hz (14.8)

The half-power bandwidth calculation could be repeated for the third and fourth modes.

14.6 RANDOM EXCITATION SOURCES

Figure 14.9. Nonstationary Flight Accelerometer Data

The liftoff acoustic event occurs during the first two seconds. The vibration then tapers down
to a lower lever during the subsonic portion of flight. The envelope then increases as the vehicle
accelerates through the transonic velocity and encounters its maximum dynamic pressure
condition in the 10 to 20 second segment. The level then gradually decreases throughout the

125
 remainder of the first stage burn. Some transient events occur beginning at 60 seconds due to
the attitude control system and other sources.

Consider the following vibration sources:

1. Turbulent airflow passing over an aircraft wing

2. Oncoming turbulent wind against a building
3. Rocket vehicle liftoff acoustics
4. Earthquake excitation of a building

One common characteristic of these examples is that the motion varies randomly with time. Thus, the
amplitude cannot be expressed in terms of a "deterministic" mathematical function.

Dave Steinberg wrote [11]:

The most obvious characteristic of random vibration is that it is nonperiodic. A knowledge of the
past history of random motion is adequate to predict the probability of occurrence of various
acceleration and displacement magnitudes, but it is not sufficient to predict the precise magnitude
at a specific instant.

Random vibration analysis is thus a matter of probability and statistics. The parameters include: mean,
standard deviation, histogram, power spectral density, etc. Stationary random vibration is where the key
statistical parameters remain constant over time. The goal of most shaker table vibration tests is to apply
a stationary base input to the test unit. This would be the case where the specification was a power
spectral density. But real-world vibration tends to be nonstationary as show in Figure 14.9.

14.7 STATISTICAL PARAMETERS

The mean value μ of the time history is calculated as

1 n
=  Yi
n i=1
(14.9)

where Yi is the instantaneous amplitude and n is the total number of points.

The variance 2 of the time history is

1 n
2 =   Yi − 2 (14.10)
n i=1

126
The root-mean-square value YRMS of the time history is

1 n 2
YRMS =  Yi
n i=1
(14.11)

A normal distribution is symmetric about its mean value, but some time histories may have skewed
histograms. The skewness  is calculated as

n
  Yi − 
3

i =1
(14.12)
 =
n 3

Kurtosis is a measure of the "tailedness" of the histogram. The kurtosis  is calculated as

n
  Yi − 
4

i =1
(14.13)
 =
n4

Sinusoidal time histories have  = 1.5. Broadband random, stationary, Gaussian time histories have
 = 3.0. Time histories with higher kurtosis values usually have higher sigma peaks than those expected
from a Gaussian distribution.

The crest factor is the ratio of the absolute peak divided by the RMS.

A peak
C= (14.14)
A RMS

The denominator on the right-hand side may be replaced by the standard deviation for the case of zero
mean. The crest factor for pure sine vibration is always √2. The crest factor for random vibration is
typically assumed to be 3, but higher crest factors can occur even for closed-loop stationary broadband
vibration test performed on shaker tables [22].

Now consider a single-degree-of-freedom system subjected to a stationary Gaussian broadband random

applied force or base excitation. The expected crest factor Ce per Reference [23] is

127
 0.5772
Ce = 2 ln ( fn T ) + (14.15)
2 ln ( fn T )

Note that equation (14.15) can be extended to a multi-degree-of-freedom response by replace the natural
frequency with the positive-slope zero-crossing frequency.

14.8 RANDOM VIBRATION HISTOGRAMS

Again, the histogram divides the time history amplitude points into bins and then counts the points in
each bin. The fundamental histogram for random vibration is the Gaussian or Normal distribution, also
known as a bell curve. This function represents the distribution of the instantaneous points for a
stationary random vibration environment, which may either be broadband or narrowband. White noise
is a special case of broadband vibration which is present later in this document. The Gaussian histogram
can be normalized as probability density function such that the total area under the curve is equal to one.
These functions are useful in vibration analysis for estimating peak responses and fatigue damage.

Figure 14.10. Gaussian Probability Density Function

Again, the mean value is μ and the standard deviation is σ. The mean value in random
vibration is usually zero, either because the vibrating object has zero net displacement
or because the accelerometer signal has been AC coupled.

128
The Gaussian distribution equation is

1  − ( x −  )2 
(14.16)
p(x) = exp  
 
 2
2
2 2


The probability that the variable x will be between any two given limits is obtained by integrating the area
under the curve, which requires the gamma function. Probability examples are shown in Table 14.3 and
Table 14.4.

Table 14.3. Gaussian Distribution, Probability of Occurrence Between Limits

Probability
Statement Percent
Ratio
- < x < + 0.6827 68.27%
-2 < x < +2 0.9545 95.45%
-3 < x < +3 0.9973 99.73%

Table 14.4. Gaussian Distribution, Probability of Occurrence Outside Limits

Probability
Statement Percent
Ratio
|x|> 0.3173 31.73%
| x | > 2 0.0455 4.55%
| x | > 3 0.0027 0.27%

129
Figure 14.11. Rayleigh Distribution

The Rayleigh distribution is a continuous probability distribution for positive-valued random variable with
the equation

A  −A 2  (14.17)
p(A) = exp  2  , A  0
2  2 

The X-axis represents a scale factor times the standard deviation, as denoted by A which is used to denote
the absolute value of x.

The area under the curve from zero to a given A value represents probability of occurrence. The total
area under the curve is one. The probability of exceeding the A value is one minus the area as shown in
Table 14.5. The probability that a given peak will be within zero and 1σ is 39.35%. The probability of
exceeding 1σ is 60.65%.

An SDOF response time history can be reduced to a subset of points representing local peaks. A peak
occurs at a point where the initial slope has the opposite polarity of the ending slope. A peak can be
negative or positive regardless of whether it is a local minimum or maximum. The peak distribution will
tend to be a Rayleigh distribution for random excitation, whether by applied force or base input.
Characterizing the distribution is important for fatigue analysis.

130
 Table 14.5. Rayleigh Distribution, Probability of Occurrence Outside Limits

 Prob [ A >  ]
0.5 88.25 %
1.0 60.65 %
1.5 32.47 %
2.0 13.53 %
2.5 4.39 %
3.0 1.11 %
3.5 0.22 %
4.0 0.034 %

Thus, 1.11 % of the peaks will be above 3σ for a signal whose peaks follow the Rayleigh distribution.

131
14.9 WHITE & PINK NOISE

Figure 14.12. White Light Passing through Prism

White light passed through a prism produces a spectrum of colors. White noise is a type of analogous
random vibration.

White noise and pink noise are two special cases of random vibration.

White noise is a random signal which has a constant power spectrum for a constant frequency bandwidth.
It is thus analogous to white light, which is composed of a continuous spectrum of colors. Static noise
over a non-operating TV or radio station channel tends to be white noise. In addition, Commercial white
noise generators are designed to produce soothing random noise which masks household noise as a sleep
aid.

Pink noise is a random signal which has a constant power spectrum for each octave band. This noise is
called pink because the low frequency or “red” end of the spectrum is emphasized. Pink noise is used in
acoustics to measure the frequency response of an audio system in a particular room. It can thus be used
to calibrate an analog graphic equalizer. Waterfalls and oceans waves may generate pink-like noise.

132
14.10 WHITE NOISE INPUT & RESPONSE EXAMPLE

Figure 14.13. White Noise Base Input Time History and SDOF Response with Histograms

Stationary white noise base input excitation is applied to the SDOF model in Figure 13.12. The system has
a natural frequency of 600 Hz and an amplification factor of Q equal to 10, or 5% damping. The input time
history is shown in Figure 14.13, with a Gaussian distribution and a standard deviation of 10 G. The
response is calculated using the Smallwood ramp invariant digital recursive filtering relationship [18]. The
response likewise has a Gaussian distribution.

133
Figure 14.14. SDOF Response to White Noise Base Input, Close-up View

Table 14.6. SDOF Response to White Noise Statistics

Parameter Value Parameter Value

Mean 0G Kurtosis 3.0

Std Dev 21.8 G Crest Factor 4.28

RMS 21.8 G Maximum 93.3 G

Skewness 0 Minimum -90.6 G

A brief segment from the input and response are given in Figure 14.14. The white noise input is broadband
random. The response tends to vibrate at 600 Hz natural frequency, but the amplitude envelope varies
randomly with time. The response is narrowband random. The response statistics for the entire duration
are given in Table 14.6.

The response mean is zero. The RMS and standard deviation values are thus equal per equation (17) and
have a value of 21.8 G. The highest input absolute peak is 93.3 G. The ratio of the peak relative to the

134
standard deviation is called the “crest factor,” which is 4.28 in this example. This value can also be
represented as 4.28σ. A common assumption is that the crest factor is 3σ, which is three times the
standard deviation. But 4σ and 5σ peaks are common in stationary random vibration.

Figure 14.15. SDOF Response to White Noise, Absolute Peak Values with Histogram

The histogram of the absolute local peaks has a Rayleigh distribution.

The expected crest factor for the sample problem is

0.5772
C= 2 ln ( 600 Hz )( 60 sec ) + = 4.71 (14.18)
2 ln ( 600 Hz )( 60 sec ) 

The crest factor for sample problem was 4.28, which is 10% lower than the expected value.

135
14.11 CREST FACTOR STUDY

As a numerical experiment, the white noise synthesis and response calculation steps for the previous
example were repeated 20,000 times. A crest factor was calculated for each trial. A histogram of the
crest factor is shown in Figure 14.16. The histogram has a shape somewhat similar to a Rayleigh
distribution, but further investigation is needed to properly characterize it. The crest factors ranged from
4.14 to 6.62, with a mean of 4.76. This experimental mean value was 1.1% higher than the theoretical
4.71 value. Again, stationary random vibration time history with Gaussian distributions typically have
peaks ranging from 4σ to 5σ, but higher crest values may also occur. This is true both for broadband
inputs and narrowband responses. Designing structures and components for 3σ vibration loads
is thus insufficient, even though this is a common rule-of-thumb.

Figure 14.16. Histogram of Response Peaks for 20,000 Trials

136
14.12 NON-GAUSSIAN RANDOM VIBRATION

Stationary Gaussian random time histories can be generated in closed-loop shaker table tests. But time
histories measured in the field usually depart from the Gaussian ideal. A theoretical example is shown in
Figure 14.17. Nevertheless, the response of an SDOF system to a non-Gaussian random input will have
some tendency to be Gaussian due to the “central limit theorem,” as shown in Figure 14.17.

The central limit theorem in probability theory establishes that, in some situations, when independent
random variables are added, their properly normalized sum tends toward a normal distribution even if
the original variables themselves are not normally distributed.

Figure 14.17. Non-Gaussian Broadband Random Time History & Response

The base input has standard deviation=1.0, kurtosis=9.3 and crest factor=11.1. The response statistics
are shown in Table 14.7.

137
Table 14.7. SDOF Response to Non-Gaussian Time History

Parameter Value Parameter Value

Mean 0G Kurtosis 3.5
Std Dev 2.12 G Crest Factor 6.9
RMS 2.12 G Maximum 12.7 G
Skewness 0 Minimum -14.6 G

14.13 HALT/HASS TESTING WITH NON-GAUSSIAN RANDOM VIBRATION

Figure 14.18. HALT/HASS Testing of Circuit Boards, Combined Thermal Cycling and Vibration

Highly Accelerated Life Testing (HALT) is a process for ruggedization of preproduction products. Highly
Accelerated Stress Screening (HASS) is the production screen for products once they have been
characterized in HALT, typically performed on 100% of units per lot. The random vibration environment

138
is nonstationary, non-Gaussian vibration with high kurtosis. The overall GRMS of the vibration is
controlled, but the spectral content is uncontrolled. The thermal cycling and vibration may be applied
separately or simultaneously depending on the test phase. Note that HALT/HASS is the offspring of the
classical Environmental Stress Screening (ESS) process.

Figure 14.19. HALT/HASS Chamber Vibration Impactors

The underside of the chamber is connected to hammer impactors which create a repetitive shock time
history.

139
Figure 14.20. HALT Time History Simulation with Histogram

A time history was synthesized to represent a typical HALT chamber vibration environment, as shown in
Figure 14.20. The kurtosis is 8.9 and the crest factor is 9.3. The time history appears stationary, but the
close-up view in Figure 14.21 reveals that it is composed of a series of transient pulses.

140
Figure 14.21. HALT Time History Simulation, Close-up View

The vibration is a series of repetitive shock pulses due to the hammer strikes underneath the table.

141
Figure 14.22. Electronics Life Cycle Bathtub Curve

A successfully screening process for mass-produced electronics will precipitate latent defects to eliminate
the “infant mortality” failures shown in Figure 14.22. The defects could be bad solder joints, loose
connectors or screws, bent clamps, flawed pieces parts, etc. An example of a bad solder joint is shown in
Figure 14.23.

142
Figure 14.23. Cracked Solder Joints in Piece Part with “J Leads”

The purpose of subjecting electronic components to HALT/HASS is to uncover flaws such as bad solder
joints. The electronic components can be powered and monitored before, during, and after the tests to
verify that they are functional.

143
15 SAMPLE RATE & ALIASING

15.1 DATA COLLECTION ACCURACY

The wagon-wheel effect, or stroboscopic effect, is an optical illusion in which a spoked wheel appears to
rotate at a different speed than its true rotation. The wheel can appear to rotate more slowly than the
true rotation, appear stationary, or even appear to rotate in the opposite direction from the true rotation.
The wagon-wheel effect is seen in film or television depictions of stagecoaches or wagons in Western
movies, but the effect may also be seen in videos of helicopter rotors and aircraft propellers.
This optical illusion is a form of aliasing. It can occur in a video recording where the frame rate is too slow
relative to the rotational or oscillation frequency of the moving object. A similar problem can occur when
accelerometer data is digitized unless precautions are taken.
Engineers collect accelerometer data in a variety of settings both in the field and in lab tests. The
accelerometers measure the data in analog form. The accelerometer may have an integral mechanical
lowpass filter. Furthermore, the signal conditioning unit may have an analog lowpass filter. Lowpass
filtering of the analog signal is necessary to prevent aliasing. Eventually, the data is passed through an
analog-to-digital converter.

The proper lowpass frequency and sampling rate must be selected to ensure that the digitized data is
accurate. There are two rules-of-thumb which must each be met.

15.2 FIRST REQUIREMENT

The first requirement is that the sampling rate must be greater than the maximum analysis frequency.
Industry guidelines are summarized in Table 15.1, where

(minimum sampling rate) > ( N )( maximum analysis frequency ) (15.1)

Table 15.1. Sample Rate First Requirement

Analysis Type N
Frequency Domain 2
Time Domain 10

The frequency domain requirement is based on the need for at least two time-domain coordinates per
cycle to resolve a sine wave. The frequency domain analysis thus extends up to the Nyquist frequency
which is one-half the sample rate. Note that some conservative references specify an N of 2.5 for
frequency domain calculation.

144
A sampling rate of 100 kHz is thus required for a shock response spectrum (SRS) analysis extending to 10
kHz. The SRS function is covered in Section 19. Recall that the SRS is calculated in the time domain.

The IES Handbook for Dynamic Data Acquisition and Analysis [24] gives the following guidelines:

Unlike other spectral quantities evolving from the discrete Fourier transform computations, the
SRS is essentially a time domain quantity. Hence, the digital sampling rate given by Rs=1/Δt
introduces errors beyond those associated with aliasing about the Nyquist frequency. Thus, Rs
must be high enough to accurately describe the response of the SRS oscillators. To minimize
potential error, it is recommended that the SRS computations be performed with a sampling
rate of Rs > 10 fh, where fh is the highest natural frequency of the SRS computation.

15.3 SECOND REQUIREMENT

Shannon’s sampling theorem states that a sampled time signal must not contain components at
frequencies above the Nyquist frequency.

Again, the Nyquist frequency is equal to one-half the sampling rate. Shannon's theorem applies to
frequency domain analysis.

Thus, the second requirement is that the sampling rate must be greater than the maximum frequency
present in the source energy at the measurement location. This requirement is independent of the
maximum analysis frequency.

The guidelines are summarized in Table 15.2, where

(minimum sampling rate) > ( M )( maximum frequency in source energy ) (15.2)

Table 15.2. Sample Rate Second Requirement

Analysis Type M
Frequency Domain 2
Time Domain 10

Note the similarity between Table 15.1 and Table 15.2.

145
15.4 ANTI-ALIASING FILTER

Accelerometer Signal Conditioner Analog-to-Digital

& Lowpass Filter Converter

Figure 15.1. Data Acquisition Flow Chart

The challenge is that the maximum frequency in the source energy maybe unknown, which is the usual
case for field data. The solution is to use an analog lowpass filter prior to digitization. This filter is also
referred to as an anti-aliasing filter. The filter cut-off frequency is typically set at, or slightly above, the
maximum analysis frequency.

The IES Handbook for Dynamic Data Acquisition and Analysis gives the following guidelines.
Let
fc be the cutoff frequency
fN be the Nyquist frequency

• A lowpass anti-aliasing filter with a cutoff rate of at least 60 dB/octave should be used for the
analog-to-digital conversion of all dynamic data.
• With a 60 dB/octave cutoff rate, the half-power point cutoff frequency of the filter should be set
at fc < 0.6 fN.
• If the anti-aliasing filter has a more rapid cutoff rate, a higher cutoff frequency can be used, but
the bound fc < 0.8 fN should never be exceeded.

146
Figure 15.2. Pyrotechnic Shock Test of an Avionics Component

Mild detonating cord is source shock and is coiled on the underside of the plate. Accelerometers are
mounted on the plate to record the shock level. The component is powered and monitored during the
test. The source shock may have energy components exceeding 100 kHz. Use of an analog anti-aliasing
filter is critical. Otherwise, high frequency energy will be erroneously folded down to lower frequencies
in the shock response spectrum.

147
15.5 ALIASING EXAMPLE

Figure 15.3. A 20 Hz Sine Function Sampled at 200 Hz

There are 10 points per period. The Nyquist frequency is 100 Hz.

Figure 15.4. An 1800 Hz Sine Function Samples at 2000 Hz

148
The red markers show that there are 1.11 points per 1800 Hz cycle. The Nyquist frequency is 1000 Hz.
Aliasing occurs!

Figure 15.5. Two Signals Digitized at 2000 Hz

The two signals both appear as 200 Hz although only one truly is.

149
Figure 15.6. Fourier Magnitude of Each Signal Sampled at 200 Hz

The 180 Hz signal is folded about the Nyquist frequency which is 100 Hz. The resulting energy is deposited
at 20 Hz. Fourier transforms are covered in Section 16.

As an example, consider a pyrotechnic shock test with a time domain maximum analysis frequency of 10
kHz. The sample rate should be at least 100 kHz. The Nyquist frequency is 50 Hz. Assume a lowpass filter
with a 60 dB/octave cutoff rate. The maximum cutoff frequency should thus be 30 kHz. These are values
to satisfy basic requirements. A better approach would be to set the sample rate at 250 kHz or more if
the measurement is going to be in the “near field” relative to the source. This higher sample rate can be
readily achieved in a ground test, but not in flight where the telemetry bandwidth is limited.

As another example, consider a vibration test with a maximum analysis frequency of 2 kHz in the time
domain. The sample rate should be at least 20 kHz. The Nyquist frequency is 10 kHz. The maximum
cutoff frequency should thus be 6 kHz.

150
16 FOURIER TRANSFORMS

16.1 FOURIER TRANSFORM EQUATIONS

16.2 TEXTBOOK FOURIER TRANSFORM

The Fourier transform is a method for identifying the spectral content of a time history. This transform
has already been used previously in this document, including for the channel beam vibration response in
Figure 12.3 and Figure 12.4. This section will give background equations starting with the rarefied
mathematician approach and finishing with the practical engineering version.

The textbook Fourier transform X(f ) for a continuous time series x(t) is defined as


X(f ) =  x(t)exp ( − j2 f t ) dt (16.1)
−

Note that j = − 1 . And both the time and frequency domains are considered to extend from minus to
plus infinity in their respective spheres.

There is also an inverse Fourier transform defined as


x(t) =  X(f )exp ( + j2 f t ) df (16.2)
−

As an aside, the corresponding Euler’s equation is

exp (  j2 f t ) = cos ( 2 f t )  jsin ( 2 f t ) (16.3)

There are several impractical quirks begat from this transform and inverse pair. One is that a signal as
simple as a pure sine function will have a spectral component at a negative of its true frequency, as shown
in Figure 16.1. How is an engineer to design and test for a negative frequency? In addition, measured
time histories have a distinct start and stop time with a finite duration. Next, integration is analogous to

151
multiplication. Assume that x(t) represents an acceleration time history. The resulting Fourier transform
from equation (16.1) would return a transform with amplitude dimension [acceleration-time]. What
would this be? Velocity?

The last peculiarity is that measured time histories are strictly real. The resulting Fourier transform of a
measured signal from equation (16.1) would be complex, with real and imaginary components. Now take
this Fourier transform and insert it into the inverse formula in (16.2). The resulting time history will be
complex. Fortunately, this concern is easily resolved because the imaginary components will all be zero.

Imaginary X(f)

f=0

Figure 16.1. Mathematicians’ Fourier Transform of a Sine Function (hand drawn)

The time history behind the imagery Fourier transform equation (16.1) is

(
x(t) = A sin 2fˆ t ) (16.4)

The real Fourier transform is omitted because it is zero for a pure sine function with zero phase angle.
The symbol δ represents the Dirac delta function, a mathematical oddity unsuitable for engineering design
and test purposes.

Note the following traits for any given signal:

1. The real Fourier transform is symmetric about the f=0 line.

2. The imaginary Fourier transform is anti-symmetric about the f=0 line.

152
16.2.1 Discrete Fourier Transform

The discrete Fourier transform rectifies some, but not all, impracticalities of the textbook Fourier
transform. Engineers collect data using data acquisition systems which return digitized outputs for
measured analog inputs. The integral in equation (16.1) thus needs to be exchanged for a series, which is
easily done. The resulting discrete Fourier transform is

1 N− 1   2   (16.5)
Xk =   x n exp − j nk   , f or k = 0,1,..., N − 1
N n =0   N 

Its corresponding inverse is

N− 1
  2  (16.6)
xn =   Xk exp + j N nk  , f or n = 0,1,..., N − 1
n =0

Eureka! Equation (16.5) does away with negative frequencies. The starting frequency will now be zero.
Also, the Fourier transform amplitude dimension will now be the same that of the time history. The time
domain dimension of acceleration will thus carry directly over to the Fourier transform. And both the
time and frequency domain functions may now be defined over finite domains. Alas, this discrete pair
will need further tweaking before the equations are ready for practical use, but progress has been made.

Recall the beloved sine function from equation (16.4). Set its amplitude at 1 G and its frequency at 1 Hz.
Digitize it with a sample rate of 32 samples per second with a 512 second duration. The frequency step
in the Fourier transform will be (1/512) Hz. The imaginary discrete Fourier transform is shown in Figure
16.2.

153
Figure 16.2. Discrete Fourier Transform of a Sine Function

The sine function has a frequency of 1 Hz. Again, the Fourier transform of a pure sine function is imaginary
and antisymmetric. The antisymmetric reflection occurs about 16 Hz which is the Nyquist frequency, one-
half the sample rate. Does the discrete Fourier transform imply that that a hypothetical component
subjected to this environment needs to be designed and tested for 31 Hz as well as 1 Hz? How would the
peak at 31 Hz change if the sample rate were increased to 64 samples per second? Is this some
manifestation of Heisenberg uncertainty principle? And why is the absolute amplitude at 1 Hz equal to
0.5 but the sine function amplitude is 1?

The discrete Fourier transform in equation (16.5) will not be abandoned. It reappears in certain
calculations such as going back-and-forth between impulse response functions and transfer functions,
convolution applications, inverse Fourier transforms, etc. These intermediate concepts will be covered
in a future document. But a more practical transform is needed presently.

16.2.2 Useful Discrete Fourier Transform

Sidestepping the quandaries raised in the previous section, let us modify the discrete Fourier transform
into a practical tool. Henceforth, the Fourier transform will only be shown up to the Nyquist frequency.
And the Fourier transform will be represented in terms of magnitude rather than real and imaginary
components. The phase may also be shown with magnitude if needed. And the amplitudes at all nonzero

154
frequencies will be doubled. A more elegant of way of expressing this is that the energy components
above the Nyquist frequency will be folded about the Nyquist frequency to lower frequencies. The energy
at 31 Hz in Figure 16.2 would thus be deposited at 1 Hz. Note that the spectral magnitude at zero
frequency is a special case since it remains the simple mean value of the signal.

The full amplitude Fourier transform up to the Nyquist frequency is

 1 N −1


 x n , f or k = 0
N n =0

Xk = 
 N −1
 2  2  N
N  x n exp  − j N n k  , f or k = 1,

,
2
−1
 n =0 (16.7)

Figure 16.3. Full Amplitude, Discrete Fourier Transform of a Sine Function, Up to Nyquist Frequency

The beloved sine function’s Fourier magnitude is plotted in

155
Figure 16.3 using equation (16.6). A clean spectral peak of 1 G appears at 1 Hz. The phase angle is omitted
because it would only be valid at the sinusoid’s own frequency. It would be noise at other frequencies.
There is a need for phase angle in other situations, such as modal testing, comparing the lag from input
to response, etc.

16.3 LEAKAGE & HANNING WINDOW

16.3.1 Leakage Error

There are a number of error sources associated with the Fourier transform, including "leakage,” which is
a smearing of energy throughout the frequency domain. Leakage results when both of the following
conditions are present:

1. The signal is taken over a finite duration

2. The signal is "non-periodic" in the time record

Both these conditions are usually present in engineering data, so leakage is common. For example, leakage
occurs if a Fourier transform is calculated for a non-integral number of sine function cycles.

156
Figure 16.4. Sine Function and Memory Buffer Example 1

Consider that a data acquisition system is used to monitor a sine function over a 10 second duration. The
sine function has an amplitude of 1 G and a frequency of 1 Hz, as shown in the top subplot. The sample
rate is 32 samples per second. Now assume that the data acquisition system has a limited memory buffer
and measures only three cycles as shown in the middle subplot. Note that the time history amplitude is
zero at the start and end of the record, with an integer number of cycles. In essence, the Fourier transform
will correctly assume that the original signal is a series of three-cycle segments as shown in the bottom
subplot, such that no leakage occurs as shown in Figure 16.6.

157
Figure 16.5. Sine Function and Memory Buffer Example 2

Repeat the previous example but with a tighter memory buffer that only measures 2 1/2 cycles as shown
in the middle subplot. Leakage will occur due to the non-integer number of cycles. The Fourier transform
will assume that the original signal is a distorted series of these segments as shown in the bottom subplot,
as shown in Figure 16.6.

158
Figure 16.6. Fourier Magnitude Comparison, Both Examples

The top subplot shows a clean spectral line at 1 Hz for the case of three whole cycles. The bottom subplot
shows that the energy is smeared across a number of bands as a result of leakage error from the non-
integer number of cycles.

On a related note, the top subplot has a spectral line available at exactly 1 Hz. But the bottom subplot
lacks a line at 1 Hz and instead has lines at 0.8 and 1.2 Hz. Recall that the frequency step is the inverse of
the duration. The three-cycle spectral increment is thus 0.333 Hz, whereas the two and half cycle step is
0.4 Hz.

159
16.3.2 Hanning Window

Figure 16.7. Hanning Window

The leakage error can be reduced by subjecting the time history to an appropriate window. Two common
types of windows are the rectangular window and the Hanning window.

The rectangular, or flat, window leaves the time history data unmodified and is equivalent to no window
at all. A rectangular window is appropriate for transient data or nonstationary data, including the type
of shock time histories presented in Section 19. It does not reduce leakage.

One of the most common windows is the Hanning window, or the cosine squared window. It works best
with stationary data. This window tapers the time history data so that the amplitude envelope decreases
to zero at both the beginning and end of the time segment. The Hanning window w(t) can be defined the
is plotted in Figure 16.7.

 t 
w(t) = 1 − cos 2   , 0  t  T (16.8)
T

A normalization factor of the square root of (8/3) may be applied to the Hanned data to compensate for
lost energy, with the goal of maintaining the original overall RMS level.

160
Figure 16.8. Sample Sine Function, with Non-Integer Number of Cycles

A sine function before and after Hanning window application is show in Figure 16.8. The non-integer
number of cycles will cause a leakage error which is mitigated by the Hanning window although not
eliminated. The Fourier transforms of each are shown in Figure 16.9.

161
Figure 16.9. Fourier Spectra, Sine Function Before and After Hanning Window

Ideally, the Fourier transform would have a single, discrete line at 1 Hz with an amplitude of 1 G. A
challenge is that there are spectral lines available at 0.945 Hz and 1.05 Hz but none at 1.0 Hz. Both the
rectangular and Hanning Fourier transforms have some leakage error as a result. The rectangular window
produces more leakage error than the Hanning window. Thus, the Hanning window is recommended for
stationary data.

162
16.4 FAST FOURIER TRANSFORM

16.4.1 FFT Characteristics

Table 16.1. Suitable Record Lengths for the FFT Calculation

2 256 32,768
4 512 65,536
8 1024 131,072
16 2048 262,144
32 4196 524,288
64 8192 1,048,576
128 16,384 2,097,152

The discrete Fourier transform requires a tremendous amount of calculations for digital signals with large
numbers of points. A time history with M coordinates would require M2 complex multiplication steps.
The discrete Fourier transform can be carried out by a Fast Fourier transform (FFT) method, however. The
method is based on a time series with a number of points equal to 2N, where N is an integer. The FFT
requires M log 2 M complex multiplication steps, where M = 2N . Some length examples are shown in
Table 16.1.

16.4.2 FFT Example

Now consider a time history with 1,000,000 points. A regular Fourier transform would require 10 12
complex multiplication steps. On the other hand, an FFT would only require approximately 2(10 7) steps.
Thus, the FFT achieves the calculation in 1/50,000th of the time.

The above example is not quite correct. Again, the FFT is based on a time series with 2N coordinates.
Note that 219 = 524,288 and 220 = 1,048,576. But a time history with 1,000,000 points falls between
these two cases.

There are some options for dealing with a time history that is not an integer power of 2. One option is to
truncate the time history. This should be acceptable if the data is stationary. In the above example, the
time history would thus be truncated to 524,288 points. The second option is to pad the time history with
trailing zeroes to bring its length to an integer power of 2. A problem with this option is that it artificially

163
reduces the amplitude of the Fourier transform spectral lines. A third option is to carefully perform linear
interpolation on the data so that its length becomes 2N. Or the problem can be circumvented altogether
by configuring the data acquisition system to collect exactly 2N points. This was a feature of some vintage
systems.

16.4.3 FFT Algorithm Background

The Fourier transform requires multiplying the time points by complex weighting functions. These
functions have a certain symmetry for the case of time series with 2N points. The calculation can thus be
simplified by taking advantage of this symmetry using concepts such as the unit circle and the butterfly
algorithm. A few introductory steps are given here. Define a weighting factor W n .

 2 n 
W n = exp  − j  (16.9)
 N 

The discrete Fourier transform becomes

1 N−1
X(k) = 
N n =0
 
x(n)W nk , for k = 0,1,..., N-1
(16.10)

The matrix representation is

 X(0) 
W0 W0 W0 W 0(N −1)   x(0) 
 
 X(1)   0 W1 W2 W1(N −1)   x(1) 
  1 W  (16.11)
 X(2)  =  W 0 W2 W4 W 2(N −1)   x(2) 
  N  
     
 X(N − 1)  W0  − 
W (N −1)(N −1)  
x(N 1)
 W1(N −1) W 2(N −1)

Note that the weighting matrix in equation (16.11) is symmetric. Also note

W0 = 1 (16.12)

The Fourier transform equation in matrix form simplifies to

164
  X(0)  1 1 1 1   x(0) 
 X(1)   
W1 W2 W1(N −1)   x(1) 
 1
1
  
 X(2)  = 1  (16.13)
W2 W4 W 2(N −1)   x(2) 
  N  
    
 X(N − 1)   
1 W1(N −1) (N −1)(N −1)   − 
W 2(N −1) W 
x(N 1)

The weights can be arranged geometrically as show in Figure 16.10.

j W6

W5 W7

W4 W0

-1 0 1

W3 W1

-j W2

Figure 16.10. Unit Circle for N=8, Polar Coordinates with Real and Imaginary Axes (hand drawn)

The units circle shows some emerging symmetries. The real part of W1 is equal to that of W 7 . The
imaginary part of W1 is equal to that of W 3 . Additional symmetries can be identified by noting that
some values are the negative of others. Further development of the FFT algorithm is beyond the scope
of this document but is given in References [25] and [26].

165
16.5 ELECTRICAL TRANSFORMER HUM

Figure 16.11. Transformer Utility Box

Magnetic fields from alternating current induce vibration in certain materials, such as iron, steel, and
ferrite alloys. This effect is called magnetostriction. Transformer hum is a common example of
magnetostriction. This hum is a potential nuisance, and it is an energy loss mechanism.

Figure 16.12. Step-up Transformer

166
A basic transformer has two sets of windings, the primary and the secondary. Certain transformers may
have secondary windings or taps. The primary connects to the power source. The secondary connects to
the load usually at a lower voltage. The windings are wrapped around a core, which may be iron or ferrite,
depending on the design requirements. An iron core may be laminated with layers of iron and non-
conducting material stacked together. A laminated iron core is used to reduce the creations of eddy
currents in the iron core that would dissipate the energy being transferred from the primary coil to the
secondary coil in the form of heat.

An alternating current is input to the primary winding. The current produces a magnetic flux, or magnetic
lines of force in the core. The resulting magnetic field induces an alternating current in the secondary
winding. The number of turns on each winding determines the output voltage from the transformer. The
output voltage from the secondary is proportional to the ratio of the turns on the windings.

Figure 16.13. Magnetostriction Cycle

The magnetic field stresses the core material. The stress causes the core volume to expand regardless of
the polarity. The core contracts to its equilibrium volume as the magnetic flux returns to zero. The
material expands and contracts twice per each magnetic cycle, as shown in Figure 16.25. The core thus
vibrates at 120 Hz in response to the changing magnetic field from 60 Hz AC power. The deformation
amplitude may be so small that it is imperceptible to the unaided human eye. Nevertheless, it may be
large enough to generate a hum. The hum frequency is 120 Hz, with integer harmonics. A hypothetical
strain gage mounted on the core would output a rectified version of the curve in Figure 16.25.

167
Figure 16.14. Transformer Hum Fourier Spectra

The Fourier magnitude from a sound output of the transformer in Figure 16.11 is shown in Figure 16.14,
with its 120 Hz peak and integer harmonics.

168
16.6 APACHE HELICOPTER FLYOVER DATA

Figure 16.15. Apache AH-64 Helicopter apache_helicopter.emf

The Apache is a twin-engine army attack helicopter, as shown in Figure 16.15. It entered service with the
US Army in 1984.

Helicopter noise consists of a complex mixture of sounds. The repeating impulse noise from the rotor
blades is the dominant source of distinct spectral tones. The fundamental tone occurs at a characteristic
frequency depending on the number of blades and their rotation speed. The blade-passing frequency
(BPF) is the shaft rotation frequency times the number of rotor blades. Noise occurs at the BPF as well as
at integer harmonics. Both the main rotor and the tail rotor generate this type of noise.

Rotor-vortex interaction noise is another source. This source is also referred to as blade-vortex interaction
(BVI) noise. This source is significant in descent or level flight at low and medium velocities. Vortices form
in the wake of the blades. A given rotor blade can run into the tip vortex shed by a preceding blade. This
causes a distinctive type of noise: the annoying "blade-slap" of helicopters with slow-turning rotors, or a
sharp fluttering noise for rotorcraft with fast-turning rotors.

Furthermore, the velocity at the tips of the advancing blades in high-speed forward flight may approach
Mach 1. Shock waves begin to form at this condition, leading to shock induced flow separation. The shock
waves and turbulent flow generate further noise. This source is referred to as highspeed impulsive (HSI)
noise.

In addition, most helicopters are powered by turbine engines. Engine and gearbox noise consists of
distinct spectral peaks as well as broadband random noise.

169
Figure 16.16. Apache Flyover Acoustic Time History

A segment of an acoustic time history from an Apache flyover is shown in Figure 16.16. The signal is a
summation of sine tones, mostly from the main and tail rotor blade passing frequencies.

Figure 16.17. Apache Flyover Acoustic Fourier Spectra with Main Rotor BPFs Identified

170
The measured blade passing frequency is 21 Hz with integer harmonics. The Apache’s main rotor has four
blades. The apparent main hub frequency is thus 5.25 Hz. The actually main hub frequency is  4.84 Hz
per published data in Reference [27]. The apparent frequency is higher due to a Doppler shift with the
Apache flying toward the microphone.

Figure 16.18. Apache Flyover Acoustic Fourier Spectra with Tail Rotor BPFs Identified

The Apache tail rotor has four blades. The blades, however, are not oriented 90° (perpendicular) from
each other as in most helicopters. Specifically, one set in front of the other at a 55° angle. The
supplementary angle is 125°. This unusual arrangement is required because the two sets of blades use a
"Delta-Hinge" which allows the blades to simultaneously flap and feather. The four blades appear to
behave as two for the tail rotor blade pass frequency.

The purpose of this design is to render the blades more efficient and to reduce noise. The main rotor hub
and the tail rotor hub are driven by the same source. Thus, there is a fixed ratio between the spin rates of
each rotor. The tail rotor spin rate is about 4.86 times higher than the main rotor rate per Reference [27].
Again, the measured, Doppler-shifted, main rotor spin rate is 5.25 Hz, or 315 rpm. The tail rotor spin rate
is 25.5 Hz, or 1530 rpm, per the fixed ratio. The unique configuration of the blades yields a 2X blade passing
frequency at 51 Hz. The blades behave as two pairs rather than as four individual blades. Tail rotor integer
harmonics of this frequency occur at 102, 153, and 204 Hz. Recall that 153 Hz was the dominant peak in
the Fourier transform in Figure 16.17, due to the tail rotor’s interaction with the main rotor wake.

171
16.7 TURBOPROP ENGINE

Figure 16.19. Bombardier Q400 Aircraft

Figure 16.20. Turboprop Engine Noise Fourier Spectra

The author once flew from Spokane to Seattle, Washington on an Alaska/Horizon Air Bombardier Q400
aircraft, similar to the one in Figure 16.19. He made an audio recording using his Android phone. The

172
Fourier magnitude is shown in Figure 16.20. This model aircraft has two Pratt & Whitney Canada PW150A
turboprop engines. The PW150A engine/propeller rotation rate during takeoff and climb is 1020 RPM but
is throttled back at cruise altitude to 850 RPM, or 14.17 Hz. There are six blades on each engine, so the
blade passing frequency is 85 Hz. This frequency and its integer harmonics are clearly visible in the spectral
analysis.

173
16.8 WATERFALL FFT

16.8.1 Introduction

A waterfall FFT is a 3D plot with axes of time, frequency, and amplitude. It is useful for visualizing how
spectral peak amplitudes and frequencies change with time. A spectrogram is a related 2D plot that
shows this same data but with amplitude represented in terms of contour colors.

16.8.2 Solid Rocket Motor Example

Figure 16.21. Solid Rocket Motor

Solid rocket motors may have standing pressure oscillations which form in the combustion chamber
cavities, caused by vortex-shedding and other effects. This condition is called “Resonant Burn” or “Thrust
Oscillation.” The sinusoidal oscillation frequency may sweep downward as the cavity volume increases
due to the conversion of propellant to exhaust gas. A typical upper-stage motor susceptible to this
problem in shown in Figure 16.21. Flight data from an accelerometer mounted on a bulkhead adjacent
to this motor during its resonance burn is shown in Figure 16.22. The sinusoidal character is apparent in
the close-up view.

174
Figure 16.22. Solid Motor Resonance Burn, Accelerometer Data, Full and Close-up Views

175
Figure 16.23. Solid Motor Resonance Burn, Accelerometer Data, Waterfall FFT

The resonant burn spectral peaks begin at 520 Hz and then sweep downward to 450 Hz. The
corresponding amplitude varies with both frequency and time probably due to structural resonance
effects as measured by the bulkhead accelerometer. The causes of peaks at 280 and 350 Hz are unknown.

176
16.8.3 Wind Chime Acoustics

Figure 16.24. Wind Chime, Unscaled Sound Pressure, Waterfall FFT

Recall the wind chime example from Section 9.1.7. The bending frequencies for the tallest chime were
238, 657, 1288 and 2127 Hz. A waterfall of the chime’s sounding is shown in Figure 16.24. The
fundamental tone is very persistent with low damping. The second mode makes a very small contribution.
The third mode initially has the highest generated sound, but the tone decays quickly. The fourth mode’s
spectral peaks are barely visible in the plot.

16.9 WHITE NOISE FOURIER TRANSFORM

The Fourier transform is an excellent tool for resolving both the frequency and amplitude of a pure
sinusoidal signal, aside from leakage concerns. Alas, it is a poor choice for random vibration. This is the
motivation for representing random vibration in terms of power spectral density as presented in Section
17. The problem with the Fourier transform approach for random vibration relates to the dependence of
the amplitude on the frequency resolution as shown in the following example. Consider a stationary
white noise signal. An engineer may choose the entire signal or a representative segment for the Fourier
transform calculation. A longer segment duration gives a finer frequency resolution. There would thus
be more points in the Fourier transform magnitude to represent the signal’s energy. The corresponding
amplitude of each point would be, on average, lower than if a shorter segment were used with a wider
frequency resolution and fewer points.

177
Figure 16.25. White Noise Time History

A stationary white noise signal is shown in Figure 16.25. The standard deviation of the entire time history
is 1 G. The standard deviation of the segment from 0 to 1 second is also 1 G, approximately. Separate
Fourier transforms are taken from the entire duration and for the 0 to 1 second segment, as shown in
Figure 16.26.

178
Figure 16.26. White Noise Fourier Transform

Each curve represents a standard deviation of 1 G, equivalent to 1 GRMS given zero mean. The overall
GRMS value can be calculated by taking the “square root of the sum of the squares” from each curve. The
blue curve has fewer points, by a factor of one-tenth, to carry the energy than the red curves. The blue
curve’s average amplitude is thus higher than the denser red curve.

179
17 POWER SPECTRAL DENSITY

17.1 DECIBELS & SLOPES

17.1.1 Differences

Power spectral density functions are sometimes specified in terms of decibels. The dB difference
between two levels A & B in units of G2 or G2/Hz is

dB = 10 log (A / B) (17.1)

A 6 dB increase would quadruple a G2/Hz level. A 6 dB decrease would reduce the G2/Hz level by a ratio
of one-fourth.

The dB difference between two levels C & D in units of G or GRMS is

dB = 20 log (C / D) (17.2)

A 6 dB increase would double a G level. A 6 dB decrease would reduce the G level by a ratio of one-half.

17.1.2 Slopes

Consider a frequency domain plot with a straight-line segment in log-log format with coordinates ( f1, y1 )
and ( f 2 , y 2 ) .

The slope N between the coordinate pair is

log  y2 / y1 
N= (17.3)
log  f 2 / f1 

The line equation for this pair is

(
y(f ) = y1 f N f1N ) (17.4)

180
The dB/octave slope is for units of G2 or G2/Hz is

dB / oct = 10 N log 2 (17.5)

The dB/octave slope is for units of G or GRMS is

dB / oct = 20 N log 2 (17.6)

17.2 PSD SPECIFICATIONS

0.1
Overall Level = 6.0 grms
2
0.04 g / Hz

+3 dB / octave -3 dB / octave
PSD ( g / Hz )
2

0.01

0.001
20 80 350 2000

FREQUENCY (Hz)

Figure 17.1. Navmat P-9492 Acceleration PSD

The power spectral density in Figure 17.1 is taken from Reference [28]. The lowercase g convention is
carried over from this source, but uppercase G is used elsewhere in the present document This level is
intended as a base input for a shaker table test to screen components for parts and workmanship defects.
It is representative of similar levels found in various military and NASA standards. The Y-axis unit is more
properly grms2/Hz. The specification consists of three straight line segments in log-log format. The overall
grms value is the “square root of the area under the curve.” The integration approach is to take the grms2
area under each of the three segments. Then sum the areas and take the square root. This requires
special integration formulas due to the log-log format which depend on the slope from equation (17.3).

181
The area a i for a segment is

  
  yi   1  f N +1 − f N +1  , for N  −1
  f N   N + 1   i +1 i 
  i 

ai = 

   f 
  yi f i  ln  i +1   , for N = −1
    f 
   i 

(17.7)

The overall level for m total segments is then

m
a RMS =  ai (17.8)
i =1

Performing these steps for the Navmat P-9492 specification in Figure 17.1 requires some preliminary
work, because the specification did not explicitly give the first and last amplitudes. The first coordinate
amplitude is 0.01 g2/Hz by inspection of the graph. Identifying last coordinate amplitude requires a two-
step calculation. The slope is

−3dB
N=  −1.0 (17.9)
10log 2

The final coordinate amplitude is

 0.04 g 2 / Hz 
y=  ( 2000 Hz )−1 = 0.007 G 2 / Hz (17.10)
 ( 350 Hz )−1 
 

The final amplitude is whimsically referred to as “Bond.” The overall grms calculation is summarized in
Table 17.1.

182
Table 17.1. Navmat P-9492 Overall Level Calculation

Frequency
Slope N grms2
Domain
20 to 80 1 1.5
80 to 350 0 10.8
350 to 2000 -1 24.4
Total 36.7

The overall level is the square root of 36.7 grms2 which is 6.06 grms.

17.3 PSD INTEGRATION IN THE FREQUENCY DOMAIN

An acceleration PSD has a corresponding velocity PSD and a displacement PSD.

Let

APSD = Acceleration PSD

VPSD = Velocity PSD
DPSD = Displacement PSD

The integration formulas are

VPSD(ω) = APSD(ω) / ω2 (17.11)

DPSD(ω) = APSD(ω) / ω4 (17.12)

The resulting velocity and displacement PSDs for the Navmat specification are shown in Figure 17.2 and
Figure 17.3, respectively. These plots are useful to determine whether a given acceleration PSD test can
be performed on a certain shaker table given the table’s velocity and displacement limits.

183
Figure 17.2. Navmat P-9492 Velocity PSD

Figure 17.3. Navmat P-9492 Displacement PSD

184
17.4 PSD CALCULATION METHODS

Power spectral density functions of measured data may be calculated via three methods:

1. Measuring the RMS value of the amplitude in successive frequency bands, where the signal
in each band has been bandpass filtered
2. Taking the Fourier transform of the autocorrelation function
3. Taking the limit of the Fourier transform X(f) times its complex conjugate divided by its
period T as the period approaches infinity.

The filtering option is most instructive method and is covered in Section 21. The autocorrelation option
is the Wierner-Khintchine approach and is omitted from this document. A modified form of the third
option is covered in this section and is the preferred method for computational efficiency.

The one-sided power spectral density function X PSD (f ) is calculated from the discrete Fourier transform
X(f ) as

lim  1 X (f ) X * (f ) 
X PSD (f ) =   (17.13)
f → 0 2 f 

The Fourier transform in equation extends from zero frequency to the Nyquist frequency, which is one-
half the sample rate. Each Fourier transform has peak dimension. The one-half factor is needed to convert
the amplitude from peak2/Hz to rms2/Hz.

The frequency step is finite in practice and is the inverse of the total measured duration.

f = 1 / T (17.14)

This frequency step is the smallest sine wave frequency which can be resolved, where this sine wave has
one period equal to the total duration. The frequency step is linked to the number of statistical-degrees-
of-freedom. The reliability of the power spectral density data is proportional to the number of degrees.
A wider f gives greater PSD confidence in terms of smoothing the spectral components.

The number of degrees N dof is defined as

N dof = 2BT (17.15)

The variable B is the bandwidth of an ideal rectangular filter. It is equal to the frequency step f for an
ideal rectangular filter. The B T product is unity and the number of degrees-of-freedom is always equal
to two for a given record assuming a rectangular filter. This number may be increased for a given record

185
by subdividing it into a set of small segments, where each segment has two degrees. The total number of
degrees is then equal to twice the number of segments. This subdivision approach tends to smooth the
PSD. The penalty, however, is that the frequency resolution widens as the record is subdivided which
could smear narrow peaks. These tradeoffs are shown by example in Table 17.2.

Table 17.2. Time History with 4096 Samples Taken over 16 Seconds, Rectangular Filter

Period of Frequency
Number of dof per
Number of Each Record Resolution
Records Record Total dof
Time Samples Ti Bi = 1 / Ti
= 2Bi TI
per Record (sec) (Hz)
NR

1 4096 16 0.0625 2 2
2 2048 8 0.125 2 4
4 1024 4 0.25 2 8
8 512 2 0.5 2 16

16 256 1 1 2 32
32 128 0.5 2 2 64

64 64 0.25 4 2 128

Additional rows could be added to the above table but are omitted for brevity. Selecting the optimum
processing strategy for a given time history is a matter of engineering judgment. Continuing with the
table example, a frequency step of 4 Hz with 128 dof might be a good choice for demonstrating that the
broadband base input for a shaker table test complied with the specified PSD’s tolerance bands. On the
other hand, a narrow frequency resolution from the first or second row may be needed if the purpose of
analyzing the PSD it to identify lightly-damped, low-frequency resonant peaks.

Recall that a rectangular filter leaves the time history data unmodified. Leakage error may be a concern,
however, as discussed in Section 16.3. A Hanning window can be applied to each data segment to reduce
leakage, but it is a non-ideal filter which reduces the number of degrees-of-freedom. The key to applying
the Hanning window is to subdivide the data into overlapping segments. Nearly 90% of the degrees-of-
freedom are recovered with a 50% overlap, as shown in Figure 17.4.

186
Figure 17.4. White Noise and Subdivision into Overlapped Segments with Hanning Windows

187
17.5 PSD CALCULATION EXAMPLES

17.5.1 White Noise PSD

Recall the white noise time history in Figure 16.25. A PSD is calculated for this signal with Hanning window
and 50% overlap for each of two frequency resolution cases. The PSDs are shown in Figure 17.5.

Figure 17.5. White Noise PSDs

The blue curve has a finer frequency resolution and more variation than the red curve. The red curve is
smoother due to its wider frequency resolution and higher dof number. The red curve effectively averages
the PSD points into wider bands. Note that pure white noise would have a flat PSD curve. The time history
in Figure 13.5 departs from this ideal because it was taken over a short, 10 second duration. It was also
band-limited via low-pass filtering such that it has a roll-off beginning at about 85 Hz. Filtering is covered
in Section 21.

188
17.5.2 Taurus Auto PSD

Figure 17.6. 2003 Ford Taurus Auto

The author mounted an accelerometer in his auto’s console as shown in Figure 17.6. The vehicle was
driven on a highway at 65 mph. An accelerometer time history from this road test is given in Figure 17.7.
The PSD is shown in Figure 17.8, with a narrow frequency resolution to aid in spectral peak identification.
The fundamental mode damping is calculated via the half-power bandwidth method in Figure 17.9.

189
Figure 17.7. Taurus Time History

Figure 17.8. Taurus PSD, f=0.3 Hz, 6 dof

190
Figure 17.9. Taurus PSD, Half-Power Bandwidth

f 1.5 Hz
Q= =  1.5 (33% damping) (17.16)
 f (1.9 − 0.88) Hz

The Taurus auto spring-mass Frequency is 1.5 Hz with 33% damping, which is typical damping for shock
absorbers. Common automobile natural frequencies are given in Table 17.3.

Table 17.3. Automobile Spring-Mass Frequencies

Fundamental
Vehicle
Frequency
Passenger Car 1 to 1.5 Hz
Sports Car 2 to 2.5 Hz
Hummer 4.5 Hz

191
The spectral peaks in Figure 17.8 at 14.6 and 29.1 Hz can be explained as tire imbalance frequencies per
the following calculations steps.

a. Assume 25 inch tire outer diameter at 65 mph.

b. Circumference =  (25 inch) = 78.5 inch
c. 65 mph = 1144 in/sec
d. Speed/Circumference = (1144 in/sec) / 78.5 in = 14.6 Hz
e. 2X harmonic = 29.1 Hz

17.6 PSD TIME HISTORY SYNTHESIS

Table 17.4. Time History Synthesis Steps to Meet PSD Specification

Step Description
1 Generate a white noise time history
2 Take the FFT
3 Scale the FFT amplitude per the PSD for each frequency
4 The time history is the inverse FFT
Use integration, polynomial trend removal, fade in and out, and differentiation so
5
that corresponding mean velocity and mean displacement are both zero
Scale the time history so that its GRMS value matches the specification’s overall
6
GRMS value
Take a PSD of the synthesized time history to verify that it matches the PSD
7
specification

The method for synthesizing a time history to satisfy a PSD are shown in Table 17.4 . These steps are used
to synthesize an acceleration time history to satisfy the Navmat P-9492 specification. The resulting
broadband random time history is shown in Figure 17.10, with its normal distribution histogram. It began
as white noise but was modified such that the final time history is no longer white noise. Its final,
corresponding PSD is shaped and defined over a finite frequency domain, as shown in Figure 17.11. White
noise would have a flat PSD in contrast.

192
Figure 17.10. Acceleration Time History and Histogram

Table 17.5. Acceleration Time History Statistics

Parameter Value Parameter Value

Mean 0G Kurtosis 3.01
Std Dev 6.06 G Crest Factor 5.375
RMS 6.06 G Maximum 32.56 G
Skewness 0.011 Minimum -30.93 G

193
Figure 17.11. PSD Verification

The synthesized time history satisfies the specification well within the tolerance bands. The frequency
step is 2.84 Hz with 1024 statistical degrees-of-freedom.

194
Figure 17.12. Synthesized Time History, Acceleration, Velocity and Displacement

Each of the three response time histories has a stable oscillation about its zero baseline. Each also has a
brief fade in and out, which could be seen in a close-up view of the start and finish.

195
Now assume an SDOF system with (fn=200 Hz, Q=10) subjected to the Navmat P-9492 PSD base input,
using the model in Figure 13.12. Solve for the acceleration response in the time domain using the
synthesized base input. The numerical engine is the Smallwood ramp invariant digital recursive
relationship. The response time history is shown in Figure 17.13.

Figure 17.13. SDOF Response to Synthesized Base Input

Table 17.6. Acceleration Response Time History Statistics

Parameter Value Parameter Value

Mean 0G Kurtosis 3.02
Std Dev 11.3 G Crest Factor 4.548
RMS 11.3 G Maximum 51.36 G
Skewness 0.003 Minimum -50.95 G

196
The overall response of 11.3 GRMS is very close to the theoretical values from the frequency domain
methods which will be shown in Sections 17.7 and 17.8. The theoretical crest factor for this case per
equation (14.18) is 4.71, slightly higher than the 4.548 value in the table.

The time domain synthesis method could whimsically be referred to as a Rube Goldberg approach, after
the famous inventor and cartoonist.

Figure 17.14. Input and Response Acceleration including Response Velocity, Close-up View

The base input is broadband random vibration. The response acceleration and velocity time histories are
narrowband random. The SDOF system prefers to oscillate at its natural frequency. The positive slope
zero cross rate is 199.6 Hz.

197
Figure 17.15. Input and Response PSDs for Navmat P-9492 Synthesis

The response PSD tracks the input at the low frequency end with nearly unity gain. The resonant response
occurs at and near the 200 Hz natural frequency. The energy above √2 times the natural frequency is
attenuated. Compare Figure 17.15 with the same set of PSD curves derived from the frequency domain
method which will be shown in Figure 17.18.

198
Figure 17.16. Power Transmissibility for Input and Response PSDs

The peak power transmissibility reaches nearly 100 G2/ G2 at 200 Hz which is equivalent to Q2 where
Q=10. But this is a special case of an SDOF system subjected to base excitation. A more robust method
for estimating the Q value, as if this were experimental data, is to use the half-power bandwidth method
shown in Figure 17.17.

199
Figure 17.17. Power Transmissibility with Half-Power Bandwidth Points

The half-power bandwidth method yields the following amplification factor.

f 199.4 Hz
Q= = = 9.63  10 (17.17)
f ( 209.6 − 188.9 ) Hz

200
17.7 MILES EQUATION FOR BASE EXCITATION

Consider an SDOF system subjected to base excitation where the input PSD is white noise over the
frequency domain from 0 to infinity Hz, using the model in Figure 13.12. The resulting acceleration
x GRMS can be determined by Miles equation from Reference [11]. The constant power spectral density
amplitude is represented by A with unit of (G2/ Hz).


x GRMS ( f n , Q ) =   f n Q A (17.18)
 
2

An equivalent form is

  1 
x GRMS ( f n ,  ) =     f n A (17.19)
 2   2 

The Miles equation is widely used due to its simplicity, but its assumption of white noise over an infinite
domain does not exist in physical reality. A rule-of-thumb states that it may be used as long as the input
PSD is flat within one-octave on either side of the natural frequency. In practice, this rule is used with
some compromise.

Assume an SDOF system with (fn=200 Hz, Q=10) subjected to the Navmat P-9492 PSD base input. The
overall response is


( )
x GRMS ( f n , Q ) =   ( 200Hz )(10 ) 0.04 G 2 / Hz = 11.2 GRMS
 
2
(17.20)

A more robust method for performing the response calculation is given in Section 17.8.

17.8 GENERAL METHOD FOR BASE EXCITATION

Real-world PSD specifications are shaped and have lower and upper frequency limits. Miles’ equation
cannot account for a PSD with ramps and plateaus. The problem is exacerbated if the PSD is from

201
narrowband measured data with peaks and dips. These practical characteristics can be readily accounted
for by applying the power transmissibility function to the base input PSD and then by doing a point-by-
point multiplication calculate the response PSD. The overall response is then the square root of the area
under the response PSD curve.

Recall the SDOF transmissibility function from equation (13.47).

1 + ( 2 )
2
H () = ,  = f / fn
( )
(17.21)
2 2
+ ( 2 )
2
1− 

The power transmissibility is equal the transmissibility squared. The response PSD is calculated from the
power transmissibility to the input PSD.

 
1 + ( 2 )
2
 
X PSD ( f ) =  YPSD ( f ) (17.22)
( )
2
 1−  + ( 2 ) 
2 2
 

 
1 + ( 2i )
2
N  
x GRMS ( f n , Q ) =   YPSD ( f i ) f i , i = f i f (17.23)
( )
2 n
i =1  1 −  + ( 2i ) 
2 2
 i 

Equation (17.23) is cumbersome but is readily implement via a software program. It does not appear to
have a name in the literature but is referred to as the “general method” in this document.

202
Figure 17.18. General Method SDOF Response to Navmat P-9492 Base Input PSD

Repeat the example of a SDOF system with (fn=200 Hz, Q=10) subjected to the Navmat P-9492 PSD base
input using the general method. The overall response in Figure 17.18 agrees with the Miles result. It is
also very close to the time domain overall response in Table 17.6. The general method would yield a more
accurate result if the natural frequency fell on either of the input PSD slopes.

Compare Figure 17.18 with the same set of PSD curves derived from the time domain synthesis in Figure
17.15.

203
18 VIBRATION RESPONSE SPECTRUM

18.1 VRS BASICS

The vibration response spectrum overall level is calculated using equation (17.23) by varying the natural
frequency for a given Q value. It can be plotted in terms of the overall 1-sigma, 3-sigma, and peak values
as a function of natural frequency. The GRMS and 1-sigma values are taken to be equal because the mean
is assumed to be zero. The expected peak value, or crest factor, is calculated from equation (14.15) and
is the only of three response metrics that depends on duration.

The vibration response spectrum is demonstrated in the following example. First calculate the SDOF
response to the Navmat P-9492 PSD base input for natural frequencies at 100, 200, and 300 Hz, with Q=10,
using the model in Figure 13.12.

Figure 18.1. SDOF Response to Navmat P-9492 Base Input PSD for Three Natural Frequencies

204
Figure 18.2. Vibration Response Spectrum of Navmat P-9492 Base Input PSD

The calculation of the response PSDs in Figure 18.1 could be repeated for a larger set of natural
frequencies. This process can be automated via a software program. The resulting VRS is shown in Figure
18.2. The Figure 18.1 legend’s three coordinates are shown explicitly in the plot. A duration of 180
seconds is assumed, but this only affects the peak curve. The peak curve shows that peaks up to 4 and
5-sigma may occur over the duration.

205
18.2 VRS COMPARISON OF TWO PSDS

Figure 18.3. PSD Comparison

The following example is actual case history. A launch vehicle avionics component was tested to the red
curve in Figure 18.3. The blue curve represents the flight data from an accelerometer mounted adjacent
to the component. The flight data has three peaks that exceeded the test level. Was the component
under-tested? The natural frequency of the component is unknown but could have been at any of the
three exceedance peaks. Assume that each curve represents the same duration for simplicity.

206
Figure 18.4. Vibration Response Spectrum Comparison

The test level was 8.6 dB higher that the flight data in terms of overall level. The vibration response
spectrum comparison shows that the component would have a higher response to the test level than to
the flight level across the entire natural frequency spectrum from 20 to 2000 Hz. The component test
level enveloped the flight data with margin to spare.

207
18.3 PSD ENVELOPE DERIVATION VIA VRS

The vibration response spectrum can also be used an indirect method for deriving a PSD envelope for a
measured environment. The PSD envelope will have a few coordinates in log-log format as a
simplification. The goal is to derive the least PSD possible that will still envelop the measured PSD in terms
of the respective VRS curves. This allows for some peak clipping as demonstrated in the following example
for an actual case history. There are some industry rules-of-thumb for peak clipping, but the VRS is a more
thorough approach.

Figure 18.5. Suborbital Launch Vehicle Flight PSD and Envelope

The flight data from a launch vehicle internal bulkhead location is shown in Figure 18.5. The spectral peaks
were unexpected. A realization was later made that a small access cover had been prematurely jettisoned
at liftoff. The result was that cylindrical module housing the bulkhead became a large whistle. This was
an anomaly which could be corrected for the next flight in this vehicle series. Nevertheless, it presented
an opportunity to develop the enveloping method. The following example proceeds assuming the whistle
condition would be repeated. The proposed envelope is shown with the flight data in Figure 18.5. The
VRS justification for the envelope is given in Figure 18.6. Assume the flight and envelope duration are the

208
same for simplicity. The envelope derivation process uses trial-and-error with random number
generation, along with some convergence routines.

Figure 18.6. Vibration Response Spectrum Comparison

209
19 SHOCK TESTING & ANALYSIS

19.1 SHOCK ISOLATOR PHOTO GALLERY

Figure 19.1. Titan II Missile Silo, Launch Control Room

The control room is mounted underground via huge isolation springs. A typical spring is shown in the
background. The purpose is to isolate the control room from mechanical shock and vibration in the event
of a nuclear strike above the launch site. The springs allow 18 inches of relative displacement. The control
room could thus carry out a retaliatory strike, as ordered by the U.S. president. This site is located south
of Tucson, Arizona. It has been decommissioned and is now a museum.

210
19.2 SIMPLE DROP SHOCK

Figure 19.2. iPhone 6 Accidental Drop

The first iPhone 6 was purchased in Perth, Australia on September 19, 2014. The event was covered by a
live TV report. The buyer mishandled the phone as he unboxed it. The phone survived the drop onto the
ground but may have had some fatigue or fracture damage.

Portable electronic devices (PEDs) are expected to survive multiple drops. Most original equipment
suppliers specifying between 30 and 50 drops. Recommended test methods are given in Reference [29].

x
m
x
m
k
k

Figure 19.3. Drop Shock Analytical Model

The accidental drop shock of a component onto a hard surface is difficult to model accurately. The item
may undergo rigid-body translation and rotation in each of its six degrees-of-freedom during freefall. The
item may have a nonlinear response with plastic deformation during impact. It may or may not bounce.
Furthermore, a box-shaped object may strike the ground on any of its corners, edges or faces.

211
A very simple method, as a first approximation, is to assume that the object is a linear, undamped, single-
degree-of-freedom subjected to initial velocity excitation as it strikes the ground and remains attached to
it via its spring. The object then undergoes free vibration in his configuration. The initial velocity is
calculated using a familiar physics formula where the change in kinetic energy is equal to the change in
potential energy due to gravity.

Assume that the object is dropped from rest. The initial velocity as it strikes the ground is

x ( 0 ) = 2g h (19.1)

The equation of motion is

m x + kx = 0 (19.2)

Or equivalently

x + n2 x = 0 (19.3)

The peak displacement is

2g h (19.4)
x (t) =
n

The peak velocity is equal to the initial velocity in equation (19.1).

The peak acceleration is

x ( t ) = n 2g h (19.5)

An example is shown in Table 11 for three natural frequency cases.

Table 19.1. Peak Response Values for 36 inch Drop Height

Natural Freq (Hz) Displacement (in) Velocity (in/sec) Acceleration (G)

200 0.133 167 543
600 0.044 167 1630

1000 0.027 167 2710

212
19.3 CLASSICAL SHOCK

Figure 19.4. Idealized Classical Pulse Examples

Classical pulses are the simplest base excitation pulses. They are deterministic and can be represented by
simple mathematical functions. They are typically one-sided. An SDOF system’s response to a classical
pulse can be solved for exactly using Laplace transforms.

Four classical pulse types are shown in Figure 19.4. Other types include rectangular and trapezoidal
pulses. These pulses do not necessarily represent real field environments, but they are still used
throughout industry to test equipment ruggedness for convenience.

Shock tests are performed on military equipment [27] to:

a. Provide a degree of confidence that materiel can physically and functionally withstand the
relatively infrequent, non-repetitive shocks encountered in handling, transportation, and service
environments. This may include an assessment of the overall materiel system integrity for safety
purposes in any one or of the handling, transportation, and service environments.
b. Determine the materiel's fragility level, in order that packaging may be designed to protect the
materiel's physical and functional integrity.
c. Test the strength of devices that attach materiel to platforms that can crash.

213
Potential equipment failure modes due to shock excitation include:

a. Materiel failure resulting from increased or decreased friction between parts, or general
interference between parts.
b. Changes in materiel dielectric strength, loss of insulation resistance, variations in magnetic and
electrostatic field strength.
c. Materiel electronic circuit card malfunction, electronic circuit card damage, and electronic
connector failure. (On occasion, circuit card contaminants having the potential to cause short
circuit may be dislodged under materiel response to shock.)
d. Permanent mechanical deformation of the materiel resulting from overstress of materiel
structural and nonstructural members.
e. Collapse of mechanical elements of the materiel resulting from the ultimate strength of the
component being exceeded.
f. Accelerated fatiguing of materials (low cycle fatigue).
g. Potential piezoelectric activity of materials, and materiel failure resulting from cracks in fracturing
crystals, ceramics, epoxies, or glass envelopes.

214
Figure 19.5. Drop Shock Test Machine, Initial Velocity Excitation

Classical pulse shock testing has traditionally been performed on a drop tower. The component is
mounted on a platform which is raised to a certain height. The platform is then released and travels
downward to the base, which has pneumatic pistons to control the impact of the platform against the
base. In addition, the platform and base both have cushions for the model shown. The pulse type,
amplitude, and duration are determined by the initial height, cushions, and the pressure in the pistons.
This is a textbook example of case where the initial potential energy of the raised platform and test item
are converted to kinetic energy. The final velocity of the freefall becomes the initial velocity of the shock
excitation.

215
Figure 19.6. A 50 G, 11 msec, Terminal Sawtooth Pulse for Shaker Shock Test

Classical pulse shock testing can sometimes be performed on shaker tables, but there some constraints.
The net velocity and net displacement must each be zero. Also, the acceleration, velocity and
displacement peaks must each be within the shaker table stroke limits. Pre and post pulses are added to
classical pulses to meet these requirements. A hypothetical terminal sawtooth pulse suitable for shaker
shock testing is shown in Figure 19.6.

216
19.4 HALF-SINE SHOCK EXAMPLE

Consider a single-of-freedom is subjected to a 50 G, 11 msec half-sine pulse applied as base excitation per
Figure 13.12. Set the amplification factor to Q=10. Allow the natural frequency to be an independent
variable. Solve for the absolute response acceleration. The equation of motion is for the relative
displacement z is

Asin (  t / T ) , for t  T
z + 2 n z + n 2 z =  (19.6)
 0 , for t > T

The primary response occurs during the half-sine pulse input. The residual response occurs during the
quiet period thereafter. The total response is the combination of primary and residual. The exact
response for a given time can be calculated via a Laplace transform. Note that the quiet period solution
is free vibration with initial velocity and displacement excitation.

Figure 19.7. Fourier Transform Half-Sine Shock

A common misunderstanding is to regard the half-sine shock pulse as having a discrete frequency which
would be the case if it were extended to a full-sine pulse. The Fourier transform in Figure 19.7 shows that
the half-sine pulse has a continuum of spectral content beginning at zero and then rolling-off as the
frequency increases. There are also certain frequencies where the magnitude drops to zero. The

217
magnitude represents the acceleration, but the absolute magnitude depends on the total duration include
the quiet period after the pulse is finished. The post-pulse duration was 10 seconds in the above example,
for a total of 10.011 seconds.

Figure 19.8. SDOF System Response to Half-Sine Shock, 10 Hz

The response in Figure 19.8 has an absolute peak value less than the peak input. This is an isolation case.
The response positive and negative peaks occur after the base input pulse is over.

218
Figure 19.9. SDOF System Response to Half-Sine Shock, 75 Hz

The response in Figure 19.9 has an absolute peak value that is 1.65 times the peak input. This is resonant
amplification case. The absolute peak response occurs during the base input pulse.

219
Figure 19.10. SDOF System Response to Half-Sine Shock, 400 Hz

The response converges to the base input as the natural frequency becomes increasingly high. This
becomes a unity gain case. The system is considered as hard-mounted.

The peak results from the three cases are shown in Table 19.2. The calculation can be repeated for a
family of natural frequencies. The peak acceleration results can then be plotted as a shock response
spectrum (SRS) as shown in Figure 19.11. The peak relative displacement values can likewise be plotted
as an SRS as shown in Figure 19.12.

Figure 19.11. Acceleration SRS, 50 G, 11 msec Half-Sine Pulse

The two curves in Figure 19.11 contain the coordinates in Table 19.2. The initial slope of each SRS curve
is 6 dB/octave indicating a constant velocity line. The curves also indicate that the peak response can be
lowered by decreasing the natural frequency. A low natural frequency could be achieved for a piece of
equipment by mounting it via soft, elastomeric isolator bushings or grommets. But a lower natural
frequency leads to a higher relative displacement as shown in Figure 19.12.

220
Table 19.2. Summary of Peak Response, 50 G, 11 msec, Half-Sine Base Input

Natural Absolute Value of

Peak Positive (G)
Frequency (Hz) Peak Negative (G)
10 20.3 17.3

75 82.5 65

400 51.8 4.38

Figure 19.12. Relative Displacement SRS, 50 G, 11 msec Half-Sine Pulse

The curves in Figure 19.12 could be used for designing isolator mounts for a component. The isolators
must be able to take up the relative displacement without bottoming or topping out. There must also be
enough clearance and sway space around the component.

221
Figure 19.13. Acceleration SRS, 50 G, 11 msec Terminal Sawtooth Pulse

The positive and negative SRS curves are reasonably close for the terminal sawtooth pulse shown in Figure
19.13. In contrast, the positive and negative SRS curves for the half-sine pulse in Figure 19.11 diverge as
the natural frequency increases above 80 Hz. The terminal sawtooth pulse is thus usually preferred over
the half-sine pulse for classical shock testing.

Another means of visualizing the SRS concept is given in Figure 19.14.

222
Figure 19.14. Half-Sine Shock Applied as Base Input to Independent SDOF Systems

spring_mass_srs.bmp

The systems are arranged in order of ascending natural frequency from left to right and subjected to a
common half-sine base input. The Soft-mounted system on the left has high spring relative deflection,
but its mass remains nearly stationary. The Hard-mounted system on the right has low spring relative

223
deflection, and its mass tracks the input with near unity gain. The Middle system ultimately has high
deflection for both its mass and spring. The peak positive and negative responses of each system are
plotted as a function of natural frequency in the shock response spectrum.

19.5 RESPONSE TO ARBITRARY EXCITATION

Recall from Section 13.3 that the response of a single-degree-of-freedom system to base excitation can
be expressed in terms of a second order, ordinary differential equation for the relative displacement z
for a base acceleration y .

z + 2 n z + n 2 z = − y (19.7)

Equation (19.7) can be solved via Laplace transforms if the base acceleration is deterministic such as a half-
sine pulse. A convolution integral is needed if the excitation varies arbitrary with time. The convolution
integral is computationally inefficient, however. An alternative is to use the Smallwood ramp invariant
digital recursive filtering relationship [18], [20].

Again, the recursive filtering algorithm is fast and is the numerical engine used in almost all shock response
spectrum software. It is also accurate assuming that the data has a sufficiently high sample rate and is
free from aliasing. One limitation is that it requires a constant time step.

The equation for the absolute acceleration is

xi= 2exp  −n t  cos  d t  x i −1

− exp  −2n t  x i −2

  1 
 

+ 1−   exp ( −n t ) sin ( d t )  y i
  d t 
 



   1  
+  2exp ( −n t )  − cos ( d t ) +   sin ( d t )   y i −1
 

   d t  


  1  

+  exp ( −2n t ) −   exp ( −n t ) sin ( d t )  y i −2

  d t  

(19.8)

224
The damped natural frequency is

 d =  n 1- 2 (19.9)

The digital recursive filtering relationship for relative displacement is omitted for brevity but is available
in Reference [20].

The relationship in equation (19.8) is recursive because the response at the current time depends on the
two previous responses, which are the first two terms on the righthand side of the equation. This is a
feedback loop in terms of control theory. The relationship is filtering because the energy at and near the
natural frequency is amplified whereas higher frequency energy above √2 times the natural frequency is
attenuated. See Figure 13.14.

19.5.1 Shock Response Spectrum Test Specification Objective

Figure 19.15. Mid-Field Pyrotechnic Shock Time History

225
Consider the measured mid-field shock time history in Figure 19.15 as taken from Reference [30]. The
time history would be essentially impossible to reproduce in a test lab given that it is a high-frequency,
high-amplitude complex, oscillation pulse. The aerospace practice instead is to derive an SRS to represent
the damage potential of the shock event. The test conductor may then use an alternate pulse to satisfy
the SRS specification within reasonable tolerance bands. This is an indirect method of achieving the test
goal. There are some limitations to this approach. One is that the test item is assumed to be linear.
Another is that it behaves as a single-degree-of-freedom system. Nevertheless, this method is used in
aerospace, military and earthquake engineering fields, for both analysis and testing.

Figure 19.16. Mid-Field Shock Response Spectrum & Envelope

The customary approach is to draw a conservative envelope over the measured SRS. The ramp-plateau
format is the most common, although there are variations. The enveloping process shown in the above
vibration is very conservative in the mid frequency domain. An additional dB uncertainty factor may be
needed to develop the envelope into a test specification, given that the SRS envelope is derived from a
single time history. Industry standards, such as Reference [31], give guidelines for the dB factor.

226
19.6 SEISMIC SHOCK

19.6.1 Seismic Waveforms

Figure 19.17. P-Wave

The primary wave, or P-wave, is a body wave that can propagate through the Earth’s core. This wave can
also travel through water. The P-wave is also a sound wave. It thus has longitudinal motion. Note that
the P-wave is the fastest of the four waveforms.

Figure 19.18. S-Wave

227
The secondary wave, or S-wave, is a shear wave. It is a type of body wave. The S-wave produces an
amplitude disturbance that is at right angles to the direction of propagation. Note that water cannot
withstand a shear force. S-waves thus do not propagate in water.

Figure 19.19. Love Wave

Love waves are shearing horizontal waves. The motion of a Love wave is similar to the motion of a
secondary wave except that Love wave only travel along the surface of the Earth. Love waves do not
propagate in water.

Figure 19.20. Rayleigh Wave, Retrograde

Rayleigh waves produce retrograde elliptical motion. The ground motion is thus both horizontal and
vertical. The motion of Rayleigh waves is similar to the motion of ocean waves in

Figure 19.21 except that ocean waves are prograde. Rayleigh waves resulting from airborne acoustical
sources may either be prograde or retrograde per Reference [32]. In some cases, the motion may begin

228
as prograde and then switch to retrograde. Airborne acoustical sources include above ground explosions
and rocket liftoff events.

Wave Propagation

Figure 19.21. Ocean Surface Wave Particle Motion, Prograde

The Love and Rayleigh waves are both surface waves. These are the two seismic waveforms which can
cause the most damage to building, bridges and other structures.

As an aside, seismic and volcanic activity at the ocean floor generates a water-borne longitudinal wave
called a T-wave, or tertiary wave. These waves propagate in the ocean’s SOFAR channel, which is
centered on the depth where the cumulative effect of temperature and water pressure combine to create
the region of minimum sound speed in the water column. SOFAR is short for “Sound Fixing and Ranging
channel.” These T-waves may be converted to ground-borne P or S-waves when they reach the shore.
Acoustic waves travel at 1500 m/s in the ocean whereas seismic P and S-waves travel at velocities from
2000 to 7000 m/s in the crust.

19.6.2 Seismic Response Spectrum Method

Professors Theodore von Kármán and Maurice Biot were very active in the early 1930s in the theoretical
dynamics aspects of what would later become known as the response spectrum method in earthquake
engineering. Biot proposed that rather than being concerned with the shape of the input time history,
engineers should instead use a method describing the response of systems to those shock pulses. The
emphasis should instead be on the effect, as represented by the response of a series of single-degree-of-
freedom oscillators, similar to that previously shown for the case of a half-sine input in Figure 19.14.

Practical use of the response spectrum method had to wait until the 1970s due to the intricacy of the
response calculation for complex, oscillating pulses which required digital computers. Time was also
needed to establish and publicize databases of strong motion acceleration time histories.

The response spectrum method was adopted for pyrotechnic shock in the aerospace industry and
renamed as shock response spectrum.

229
19.6.3 El Centro Earthquake

Figure 19.22. El Centro, Imperial Valley Earthquake Damage

Nine people were killed by the May 1940 Imperial Valley earthquake. At Imperial, 80 percent of the
buildings were damaged to some degree. In the business district of Brawley, all structures were damaged,
and about 50 percent had to be condemned. The shock caused 40 miles of surface faulting on the Imperial
Fault, part of the San Andreas system in southern California. Total damage has been estimated at about
$6 million. The magnitude was 7.1. The was the first major earthquake for which strong motion
acceleration data was obtained that could be used for engineering purposes.

230
Figure 19.23. El Centro Earthquake, Triaxial Time History

The highest excitation is in the North-South axis, parallel to the ground.

231
Figure 19.24. El Centro Earthquake North-South SRS, Three Damping Cases

The acceleration levels reached 1.5 G for the 1% damping curve. Recall that large civil engineering
structures can have nonlinear damping. The damping values tend to increase as the base excitation levels
increase as shown for the Transamerica Title Building in Section 11.6.

232
Figure 19.25. El Centro Earthquake Tripartite SRS

Seismic SRS curves are often plotted in tripartite format which displays relative displacement, pseudo
velocity and acceleration responses all on the same graph. The pseudo velocity PV is calculated from the
relative displacement Z as

PV =  n Z (19.10)

Stress can be calculated from pseudo velocity using the methods in Section 20. The acceleration curve
might be the most important design metric for equipment mounted inside a building. The relative
displacement might be the most important concern for analyzing the foundational strength. The curves
also show design tradeoffs. Lowering the building’s natural frequency below 1 Hz reduces the acceleration
response but increases the relative displacement. Note that Q=10 is the same as 5% damping.

233
19.6.4 Golden Gate Bridge

Figure 19.26. Golden Gate Suspension Bridge, San Francisco, California

In addition to traffic loading, the Golden Gate Bridge must withstand the following environments:

1. Earthquakes, primarily originating on the San Andreas and Hayward faults

2. Winds of up to 70 miles per hour
3. Strong ocean currents

The Golden Gate Bridge has performed well in all earthquakes to date, including the 1989 Loma Prieta
Earthquake. Several phases of seismic retrofitting have been performed since the initial construction. The
bridge’s fundamental mode is a transverse mode with a natural frequency of 0.055 Hz, with a period of
18.2 seconds

Note that California Department of Transportation (CALTRANS) standards require bridges to withstand an
equivalent static earthquake force (EQ) of 2.0 G. This level is plausibly derived as a conservative envelope
of the El Centro SRS curves in Figure 19.24.

234
19.6.5 Vandenberg AFB

Figure 19.27. Rocket Launch from Vandenberg AFB, California

The Vandenberg launch site is near the San Andreas fault system. The vehicle is mounted on the pad as
a tall cantilever beam. The vehicle must be analyzed to verify that it can withstand a major seismic event.
The vehicle may be mounted to the pad for only two weeks prior to launch. The odds of an earthquake
occurring to that time window are miniscule. But the launch vehicle and payload together may cost well
over $1 billion. The risk thus necessitates the analysis. Areas of concern are the loads imparted at the
launch vehicle’s joints and to the payload.

235
Figure 19.28. NASA SRS Curves for Launch Vehicles at Vandenberg

SRS curves are given for three damping cases. The curves are taken from Reference [30]. The vehicle
would typically be analyzed as a multi-degree-of-freedom system via a finite element model. Each SRS
curve could be applied to the model using a modal combination method. An alternative is to synthesize
a time history to satisfy a selected SRS curve. The time history could then be applied to the model via a
modal transient analysis.

236
19.6.6 Seismic Testing

Figure 19.29. Electrical Power Generator Seismic Shock Test

The diesel generator is mounted onto a platform at the top of a shaker table which is located below the
ground floor. This could be an emergency power generator for a hospital in an active seismic zone. A
video clip of the test is available on YouTube at: https://youtu.be/5uSqI7kSYE4

237
19.7 PYROTECHNIC SHOCK

19.7.1 Flight Events

Figure 19.30. Stage Separation Ground Test, Linear Shaped Charge

The plasma jet cuts the metal inducing severe mechanical shock energy, but the smoke and fire would not
occur in the near-vacuum of space.

Launch vehicle avionics components must be designed and tested to withstand pyrotechnic shock from
various stage, fairing and payload separation events that are initiated by explosive devices. Solid rocket
motor ignition is another source of pyrotechnic shock. The source shock energy can reach thousands of
acceleration Gs at frequencies up to and beyond 100 kHz. The corresponding velocities can reach a few
hundred in/sec, well above the severity thresholds in the empirical rules-of-thumb in Section 20.4.

Empirical source shock levels for a variety of devices are given in References [30], [33], [34]. These levels
are intended only as preliminary estimates for new launch vehicle designs. The estimates should be
replaced by ground test data once the launch vehicle hardware is built and tested.

238
19.7.2 Pyrotechnic Device Photo Gallery

Figure 19.31. Metal Clad Linear Shaped Charge

The chevron focuses a pyrotechnic plasma jet at the launch vehicle’s separation plane material. Severe
shock levels are generated as a byproduct.

Figure 19.32. Frangible Joint

A frangible joint may be used for stage or fairing separation. The key components of a frangible joint:

• Mild Detonating Fuse (MDF)

• Explosive confinement tube
• Separable structural element
• Initiation manifolds
• Attachment hardware

The hot gas pressure from the MDF detonation cause the internal tube to expand and fracture the joint.

239
Figure 19.33. Frangible Nuts, Hold Down Posts

The purpose of the nuts was to hold the SRBs in place against wind and ground-borne excitation. The nuts
were separated just before liftoff.

240
(image courtesy of European space agency)

Figure 19.34. Clamp Band with Pyrotechnic Bolt-Cutters

Clamp bands are often used for payload separation from launch vehicle adapters. They are also commonly
used for stage separation in suborbital launch vehicles, similar to the one in Figure 5.2. A pyrotechnic
bolt-cutter uses an explosive charge to drive a chisel blade to cut the band segments’ connecting bolt.
The cutters produce some shock energy, but much of the shock is due to the sudden release of strain
energy in the preloaded clamp band. This action can excite a ring mode in the radial axis. Recall Section
8.3.

The total clamp band release shock is significantly less than linear shaped charge and frangible joint. Note
that an analysis must be performed to verify that no gapping will occur in the between the band and the
joint as the vehicle undergoes bending mode vibration during powered flight.

241
19.8 PYROTECHNIC SHOCK DATA

19.8.1 Initial SRS Slopes

Figure 19.35. Expected Pyrotechnic SRS Initial Slope Limits

Near-field pyrotechnic shock can be difficult to measure accurately. The accelerometer data may have a
baseline shift or spurious low frequency transient. This error could be a result of the accelerometer’s own
natural frequency being excited or to some other problem. Aerospace pyrotechnic shock SRS
specifications usually begin at 100 Hz due to the difficult in accurately measuring low-amplitude, low-
frequency shock, while simultaneously measuring high-amplitude, high-frequency shock.

There are several methods for checking whether the data is acceptable. One is to the check the initial
slopes of both the positive and negative SRS curves. Each should have an overall trend of 6 to 12
dB/octave. The actual slopes may have local peaks and dips due to low frequency resonances. A 12
dB/octave slope represents constant displacement and zero net velocity change. A 6 dB/octave slope
indicates constant velocity. Recall the slope formulas in Section 17.1.2.

A second method for checking data accuracy is to verify that the positive and negative SRS curves are
within about 3 dB of each other across the entire natural frequency domain. A third method is to integrate

242
the acceleration time history to velocity. The velocity time history should oscillate in a stable manner
about the zero baselines.

These verification goals are challenging to meet with near-field shock measurements of high-energy
source shock, such as that from linear shaped charge and frangible joints. In practice, some high-pass
filtering or spurious trend removal may be necessary. There is no one right way to perform this “data
surgery.” It is a matter of engineering judgment.

19.8.2 Re-entry Vehicle Separation, Flight Data

Figure 19.36. Re-Entry Vehicle, Separation Shock, Near-Field Measurement, Time History

The source device was linear shaped charge.

243
Figure 19.37. Re-Entry Vehicle, Separation Shock, Near-Field Measurement, SRS

The SRS reached 20,040 G at 2400 Hz, which is an extreme level per the rule-of-thumb Section 20.4.

244
19.9 WATER IMPACT SHOCK

Figure 19.38. Solid Rocket Booster Water Impact Shock

Each of the two Space Shuttle Boosters was recovered and refurbished after every flight. Each booster
contained sensitive avionics components which underwent shock testing for the water impact event.

245
Figure 19.39. Solid Rocket Booster Water Impact Shock, Time History

The data is from the STS-6 Mission. The accelerometer was mounted at the forward end of the booster
adjacent to a large IEA avionics box. This was the worst-case shock event for this component.

246
Figure 19.40. Solid Rocket Booster Water Impact Shock, Tripartite SRS

The maximum acceleration response is 257 G at 85 Hz. The maximum pseudo velocity response is 201
in/sec at 76 Hz, which is severe per the rule-of-thumb Section 20.4.

247
19.10 SHOCK RESPONSE SPECTRUM SYNTHESIS

19.10.1 Synthesis Objectives

Consider an SRS specification for a component, subsystem or a large structure. The article should be
tested, if possible, to verify that it can withstand the shock environment. Certain shock tests can be
performed on a shaker table, like the generator test in Figure 19.29. This requires synthesizing an
acceleration time history to satisfy the SRS. The net velocity and net displacement must each be zero for
this test. These requirements can be met with a certain type of wavelet series, where the individual
wavelets may be nonorthogonal. The resulting wavelet time history should meet the SRS within
reasonable tolerance bands, but it may not “resemble” the expected time history which is a limitation of
this method. Another requirement is that the peak acceleration, velocity and displacement values must
be within the shaker table’s capabilities.

An innovative method for meeting the SRS with a wavelet series that resembles one or more measured
shock time histories is given in Reference [35].

A synthesized time history can also be used for modal transient analysis. This could be done for articles
which are too large or heavy for shaker tables. This analysis can also be done on small components prior
to shock testing to determine whether they will pass the test. Or the analysis could be done in support of
isolator mounting design. There is again a need for the synthesized time history to have net velocity and
net displacement which are each zero, to maintain numerical stability in stress calculations from relative
displacement values.

Damped-sines can be used for modal transient analysis where the goal is to meet the SRS with a time
history that plausibly resembles the expected field shock event. But damped-sines do not meet the
desired zero net velocity and displacement goals. The workaround is to first synthesize a damped-sine
series to meet the SRS and then reconstruct it via a wavelet series, in Rube Goldberg fashion.

Wavelet and damped-sine synthesis are shown in the following examples.

248
19.10.2 Wavelet Synthesis

19.10.2.1 Wavelet Equation

A wavelet is a sine function modulated by another sine function. The equation for an individual wavelet
is

 0 , for t  t dm


  2 f m   
A m sin  (t − t dm ) sin 2 f m (t − t dm ) , for t dm  t   t dm + N m 
Wm ( t ) =   Nm   2f m 


  N 
 0 , for t   t dm + m 
  2f m 

(19.11)

where

Wm ( t ) = acceleration of wavelet m at time t

Am = wavelet acceleration amplitude

fm = wavelet frequency
Nm = number of half-sines
t dm = wavelet time delay

Note that N m must be an odd integer greater than or equal to 3. This is required so that the net velocity
and net displacement will each be zero. The total acceleration at time t for a set of n wavelets is
n
x(t) =  Wm (t) (19.12)
m=1

249
Figure 19.41. Sample Wavelet

A sample, individual wavelet is shown in Figure 19.41. This wavelet was a component of a previous
analysis for an aerospace project.

19.10.2.2 Wavelet Synthesis Example

Consider the specification: MIL-STD-810E, Method 516.4, Crash Hazard for Ground Equipment.

Table 19.3. SRS Q=10, Crash Hazard Specification

Natural Peak Accel

Frequency (Hz) (G)
10 9.4
80 75
2000 75

Synthesize a series of wavelets as a base input time history for a shaker shock test to meet the Crash
Hazard SRS. The goals are:

1. Satisfy the SRS specification

2. Minimize the displacement, velocity and acceleration of the base input

250
The synthesis steps are shown in the following table.

Table 19.4. Wavelet Synthesis Steps

Step Description
Generate a random amplitude, delay, and half-sine number for each wavelet.
1
Constrain the half-sine number to be odd. These parameters form a wavelet table.
2 Synthesize an acceleration time history from the wavelet table.
3 Calculate the shock response spectrum of the synthesis.
Compare the shock response spectrum of the synthesis to the specification. Form a
4
scale factor for each frequency.
5 Scale the wavelet amplitudes.
6 Generate a revised acceleration time history.
Repeat steps 3 through 6 until the SRS error is minimized or an iteration limit is
7
reached.
Calculate the final shock response spectrum error.
Also calculate the peak acceleration values.
8
Integrate the signal to obtain velocity, and then again to obtain displacement.
Calculate the peak velocity and displacement values.
9 Repeat steps 1 through 8 many times.
Choose the waveform which gives the lowest combination of SRS error, acceleration,
10
velocity and displacement.

The resulting time history and SRS are shown in Figure 19.42 and Figure 19.43, respectively.

251
Figure 19.42. Crash Hazard Time History Synthesis

The acceleration time history has a reverse sine sweep character. It is an efficient and optimized
waveform for a shaker shock test, and it satisfies the SRS as shown in the next figure. A drawback is that
it does not resemble an actual crash shock time history.

252
Figure 19.43. Crash Hazard SRS

The positive and negative curves are from the synthesized waveform. The tolerance bands are set at +3
dB.

253
19.10.3 Damped-Sine Synthesis

19.10.3.1 Damped-Sine Equation

The equation for an individual damped sinusoid is


 0, for t  t dn
Wn (t) = 

A exp  −n  n ( t − t dn ) sin  n ( t − t dn ) , for t  t dn
 n
(19.13)

where

Wn (t) = acceleration of damped sinusoid n at time t

An = acceleration amplitude
n = angular frequency
n = damping value
t dn = time delay

The total acceleration at time t for a set of damped sinusoids is

N
x(t) =  Wn (t) (19.14)
n =1

19.10.3.2 Damped-Sine Example

Consider the following specification which could represent a stage separation shock level as some location
in a launch vehicle. A modal transient finite element analysis is to be performed on a component to verify
that the component will pass its eventual shock test. The immediate task is to synthesize a time history
to satisfy the SRS. The time history should “resemble” an actual pyrotechnic shock pulse.

Note that pyrotechnic SRS specifications typically begin at 100 Hz. The author’s rule-of-thumb is to
extrapolate the specification down to 10 Hz in case there are any component modes between 10 and 100

254
Hz. This guideline is also to approximate the actual shock event which should have an initial ramp
somewhere between 6 and 12 dB/octave.

Table 19.5. SRS Q=10, Stage Separation Shock

Natural Peak Accel

Frequency (Hz) (G)
10 10
2000 2000
10,000 2000

The specification has an initial slope of 6 dB/octave. The synthesis steps are shown in the following
table.

255
Table 19.6. Damped-Sine Synthesis Steps

Step Description
Generate random values for the following for each damped sinusoid: amplitude,
1 damping ratio and delay. The natural frequencies are taken in one-twelfth octave
steps.
2 Synthesize an acceleration time history from the randomly generated parameters.
3 Calculate the shock response spectrum of the synthesis.
Compare the shock response spectrum of the synthesis to the specification. Form a
4
scale factor for each frequency.
5 Scale the amplitudes of the damped sine components.

6 Generate a revised acceleration time history.

7 Repeat steps 3 through 6 as the inner loop until the SRS error diverges.
8 Repeat steps 1 through 7 as the outer loop until an iteration limit is reached.

9 Choose the waveform which meets the specified SRS with the least error.
Perform wavelet reconstruction of the acceleration time history so that velocity and
10
displacement will each have net values of zero.

The resulting time history and SRS are shown in Figure 19.44 and Figure 19.45, respectively.

256
Figure 19.44. Stage Separation Time History Synthesis

The acceleration time history somewhat resembles a mid or far-field pyrotechnic shock. The velocity and
displacement time histories each have a stable oscillation about their respective baselines.

257
Figure 19.45. Stage Separation SRS

The positive and negative curves are from the damped-sine synthesis.

258
19.10.4 Modal Transient Finite Element Analysis for Uniform Base Excitation

Consider a rectangular plate mounted to a base at each of its four corners. The plate is to be subject to
uniform seismic excitation. There are two methods to apply the base excitation in a finite element
analysis, as shown in Figure 19.46 and Figure 19.47.

Figure 19.46. Direct Enforce Acceleration Method

The direct enforcement method is computationally intensive, requiring matrix transformations and a
matrix inversion as shown in Reference [36].

Rigid links
F

Figure 19.47. Seismic Mass Method

The seismic mass is chosen to be several orders of magnitude higher in mass than the plate. An equivalent
force is calculated and apply to the seismic mass to excite the desired acceleration at each of the plate’s

259
corners. This method adds a degree-of-freedom to the plate system resulting in a rigid-body mode at zero
frequency. But the remaining natural frequencies and mode shapes should be the same as if the plate
were mounted normally to its joining structure. The author’s experience is that the seismic mass method
is faster and more accurate and reliable than the direct enforced method. This method was introduced
in Section 13.3.3.

19.10.5 Shock Fields

19.10.5.1 Near-Field

The near-field environment is dominated by direct stress wave propagation from the source. Peak
accelerations in excess of 5000 G occur in the time domain with a frequency content extending beyond
100 kHz. The near-field usually includes structural locations within approximately 15 cm of the source for
severe devices such as linear shaped charge and frangible joint. No shock-sensitive hardware should be
mounted where it would be exposed to a near-field environment.

19.10.5.2 Mid-Field

The mid-field environment is characterized by a combination of wave propagation and structural

resonances. Peak accelerations may occur between 1000 and 5000 G, with substantial spectral content
above 10 kHz. The mid-field for intense sources usually includes structural locations between
approximately 15 and 60 cm of the source, unless there are intervening structural discontinuities.

19.10.5.3 Far-Field

The far-field environment is dominated by structural resonances. The peak accelerations tend to fall
below 2000 G, and most of the spectral content below 10 kHz. The far-field distances occur outside the
mid-field. The typical far-field SRS has a knee frequency corresponding to the dominant frequency
response. The knee frequency is the frequency at which the initial ramp reaches the plateau in the log-
log SRS plot.

19.10.5.4 Joint & Distance Attenuation

The source shock energy is attenuated by intervening material and joints as it propagates from the near-
field to the far-field. Empirical distance and joint attenuation factors for the SRS reduction are given in
References [30], [33], [34] . A typical attenuation curve from [33] is given in Figure 19.48, which assumes
an input source shock SRS consisting of a ramp and plateau in log-log format. Such curves should be used
with caution given that the attenuation is highly dependent on damping and structural details. The curves

260
can be used for preliminary estimates for new launch vehicle designs, but ground separation tests are still
needed for the actual launch vehicle hardware. These tests are needed to measure the source shock as
well as the levels at key component mounting locations.

Figure 19.48. Shock Response Spectrum Versus Distance from Pyrotechnic Shock Source

261
19.10.6 Shock Mitigation

Figure 19.49. Sensor Electronics, Wire Rope Isolators (NASA/JPL photo)

The source shock energy is attenuated as it propagates to avionics mounting locations through the launch
vehicle’s material and joints. The input shock to a component can be mitigated by mounting the
component as far away from the source device as possible. Another effective attenuation technique is to
mount the component via elastomeric bushings or wire rope isolators. The NASA Mars Science
Laboratory Sensor Support Electronics unit is mounted on vibration isolators in Figure 19.49.

262
Figure 19.50. SCUD-B Avionics Component Isolation

The bushings are made from some type of rubber or elastomeric compound. The bushings provide
damping, but their main benefit is to lower the natural frequency of the system. The isolators thus
attenuate the shock and vibration energy which flows from the instrument shelf into the avionics
component.

263
19.10.7 Pyrotechnic Shock Testing Methods

Figure 19.51. Near-Field Shock Simulation using a Plate

The test component is mounted on other side of plate. The source device is a textile explosive cord with
a core load of 50 grains/ft (PETN explosive). Up to 50 ft of Detonating Cord has been used for some high
G tests. The maximum frequency of shock energy is unknown so analog anti-aliasing filters are needed
for the accelerometer measurements per the guidelines in Section 15.

Note that a component may be mounted in the mid-field shock region of a launch vehicle. The SRS test
level derivation for this zone may include a significant dB factor for uncertainty or as a qualification
margin. This conservatism may require a near-field-type shock test for a component that is actually in a
mid-field zone. This situation can also occur for components mounted in far-fields.

264
Figure 19.52. JPL Tunable Shock Beam (photos courtesy of NASA JPL)

The NASA/JPL Environmental Test Laboratory developed and built a tunable beam shock test bench based
on a design from Sandia National Laboratory many years ago. The excitation is provided by a projectile
driven by gas pressure. The beam is used to achieve SRS specifications, typically consisting of a ramp and
a plateau in log-log format. The intersection between these two lines is referred to as the “knee
frequency.” The beam span can be varied to meet a given knee frequency. The high frequency shock
response is controlled by damping material.

265
19.10.8 Shock Failure Modes

Figure 19.53. Sensitive Electronic Parts

Pyrotechnic shock can cause crystal oscillators to shatter. Large components such as DC-DC converters
can detached from circuit boards. In addition, mechanical relays can experience chatter or transfer.

Figure 19.54. Shock Test Case History, Adhesive & Solder Joint Failure

The image shows adhesive failure and rupture of solder joint after a stringent shock test. A large
deflection of the PCB resulting from an insufficient support/reinforcement of the PCB combined with high
shock loads can lead to these failures. Staking is needed for parts weighing more than 5 grams.

266
Figure 19.55. Shock Test Case History, Lead Wire Failure

The image shows a sheared lead between solder joint and winding of coil. The lacing cord was insufficient
by itself. The lacing should be augmented by staking with an adhesive.

20 PSEUDO VELOCITY & THE STRESS-VELOCITY RELATIONSHIP

20.1 PSEUDO VELOCITY

The Shock Response Spectrum (SRS) models the peak response of a single-degree-of-freedom (SDOF)
system to a base acceleration, where the system's natural frequency is an independent variable. The
absolute acceleration and the relative displacement of the SDOF system can be readily calculated. The
relative velocity can also be calculated, but the common practice is to approximate the relative velocity
in terms of the "pseudo velocity."

The preferred method for calculating the peak pseudo velocity is to multiply the peak relative
displacement by the angular natural frequency as shown in equation (19.10). An alternate method is to
divide the peak acceleration by this frequency.

267
20.2 STRESS-VELOCITY RELATIONSHIP

Howard A. Gaberson (1931-2013) championed the idea that dynamic stress correlated better with pseudo
velocity than with either acceleration or relative displacement. Gaberson was a shock and vibration
specialist who was with the U.S. Navy Civil Engineering Laboratory and later the Facilities Engineering
Service Center from 1968 to 2000, mostly conducting dynamics research. He specialized in shock and
vibration signal analysis and published more than 100 papers and articles. Gaberson pointed out that the
kinetic energy in a mechanical system is proportional to velocity squared. The pseudo velocity is thus a
measure of the stored peak energy in the system at a particular frequency and has a direct relationship
to the survival or failure of the system. The following equations are taken from Reference [37]:

Consider the stress-velocity relationship for a traveling wave in an infinite rod.

Figure 20.1. Propagating Wave in an Infinite Rod

268
The stress  is proportional to the velocity as

(x, t) = − c v ( x, t ) (20.1)

The velocity depends on natural frequency, but the stress-velocity relationship does not. The lack of a
frequency term in equation (20.1) was another reason that Gaberson promoted the stress-velocity
relationship.

Figure 20.2. Propagating Wave in a Finite Rod

The principle can be extended to a finite rod with common boundary conditions where the peak modal
stress   n  is
max

n  ma x = c  Vn  ma x (20.2)

The index n represents the mode number. The peak stress is mainly due to the fundamental mode in
most cases. But higher modes can contribute significant stress in certain cases, such as pyrotechnic shock
response.

269
20.3 BEAM STRESS

Consider the bending vibration of a simply-supported beam. The modal stress due to bending is
proportional to the modal velocity as follows.

n  max = ĉ E A I  Vn  max (20.3)

The distance from the beam’s neutral axis to the outer fiber is ĉ .

Note this equation applies to other boundary condition cases per Reference [38].

Equation (20.3) can be simplified as follows:

n  ma x = k̂ c  Vn  ma x (20.4)

Values for the k̂ constant for typical cross-sections are:

Cross-section k̂
Solid Circular 2
Rectangular 3

The method can also be used for complex structures, but with some uncertainty. The stress-velocity
formula is

n  max = K c  Vn  max (20.5)

The constant of proportionality K depends on geometry, etc.

Bateman [39] claimed the following range for complex equipment 4 < K < 8. Gaberson [40] gave a larger
range with 1 < K < 10.

Equation (20.5) could thus be used to “bound the problem” by assuming a high and low K value.

Note that the stress-velocity equations do not require that peak stress and peak velocity occur at the same
location. Any stress concentration factor must be applied separately. The equations are useful as quick

270
estimates prior to any advanced modeling such as finite element analysis. The equations can also be used
to “sanity check” any finite element results.

20.4 SHOCK VELOCITY THRESHOLD RULE-OF-THUMB

An empirical rule-of-thumb in MIL-STD-810E [27] and SMC-TR-06-11 [41] states that a shock response
spectrum is considered severe only if one of its components exceeds the level

Threshold = [ 0.8 (G/Hz) * Natural Frequency (Hz) ] (20.6)

For example, the severity threshold at 100 Hz would be 80 G.

This rule is effectively a velocity criterion. MIL-STD-810E states that it is based on unpublished
observations that military-quality equipment does not tend to exhibit shock failures below a shock
response spectrum velocity of 100 inches/sec (254 cm/sec).

Equation (20.6) actually corresponds to 50 inches/sec. It thus has a built-in 6 dB margin of conservatism.

The 100 inches/sec threshold is defined in part by the observation that the severe velocities which cause
yield point stresses in mild steel beams turn out to be about 130 inches/sec.

20.5 VELOCITY THRESHOLD FROM YIELD STRESS

The velocity threshold Vth can be calculated from the yields stress as follow for a beam with rectangular
cross section.

 yield
Vth = (20.7)
K c

The parameters for mild steel are

 yield 33 ksi
ρ 0.00075 lbf sec^2/in^4

c 199,000 inches/sec

271
 33, 000 lbf / in 2
Vth = = 128 in / sec
( )
3 0.00075 lbf sec 2 / in 4 (199, 000 in / sec )
(20.8)

Gaberson [40] gave the limits in Table 20.1. Severe Velocities, Fundamental Limits to Modal Velocities in
Structures.

Table 20.1. Severe Velocities, Fundamental Limits to Modal Velocities in Structures

Rod Beam Plate

yield ρ
Material E (psi) Vmax Vmax Vmax
(psi) (lbm/in3)
(ips) (ips) (ips)
Douglas Fir 1.92e+06 6450 0.021 633 366 316

Aluminum
10.0e+06 35,000 0.098 695 402 347
6061-T6
Magnesium
6.5e+06 38,000 0.065 1015 586 507
AZ80A-T5
33,000 226 130 113
Structural Steel,
29e+06 0.283
High Strength
100,000 685 394 342

20.6 MORSE CHART

272
Figure 20.3. Morse Chart, Shock Response Spectra, Damage Thresholds

The Morse Chart in Figure 20.3 was derived for aerospace components subjected to shock environments
[42]. The curves are defined by the following formulas.

Table 20.2. Morse Chart Severity Thresholds

Threshold Formula
300 ips [ 4.8 (G/Hz) * Natural Frequency (Hz) ]

100 ips [ 1.6 (G/Hz) * Natural Frequency (Hz) ]

50 ips [ 0.8 (G/Hz) * Natural Frequency (Hz) ]

273
21 DIGITAL FILTERING

21.1 FILTERING INTRODUCTION

Filtering is a tool for resolving signals. It can be performed on either analog or digital signals. Analog anti-
aliasing filtering was covered in Section 15.4. The present section is limited to digital filtering with an
emphasis on the Butterworth filter, which is implemented in the time domain as a digital recursive filtering
relationship.

21.2 HIGH-PASS & LOW-PASS FILTERS

A high-pass filter is a filter which allows the high-frequency energy to pass through. It is thus used to
remove low-frequency energy from a signal.

A low-pass filter is a filter which allows the low-frequency energy to pass through. It is thus used to
remove high-frequency energy from a signal.

A band-pass filter may be constructed by using a high-pass filter and low-pass filter in series

21.3 BUTTERWORTH FILTER CHARACTERISTICS

274
 BUTTERWORTH LOWPASS FILTER L=ORDER c = 1 rad/sec

1.1
L=6
1.0 L=4
L=2
0.9

0.8

0.7
MAGNITUDE

0.6

0.5

0.4

0.3

0.2

0.1

0
0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0

FREQUENCY (rad/sec)

Figure 21.1. Butterworth Filter, Three Orders, Transfer Function Magnitudes

Butterworth_three_orders.emf

A Butterworth filter is one of several common infinite impulse response (IIR) filters. Other filters in this
group include Bessel and Chebyshev filters. These filters are classified as feedback filters. The
Butterworth filter can be used either for high-pass, low-pass, or band-pass filtering. It is characterized by
its cut-off frequency. The cut-off frequency is the frequency at which the corresponding transfer function
magnitude is –3 dB, equivalent to 0.707. The transfer function curves in Figure 21.1 each pass through
the same cut-off frequency point.

A Butterworth filter is also characterized by its order, with the sixth-order preferred in this document. A
property of Butterworth filters is that the transfer magnitude is –3 dB at the cut-off frequency regardless
of the order. Other filter types, such as Bessel, do not share this characteristic.

Consider a lowpass, sixth-order Butterworth filter with a cut-off frequency of 100 Hz. The corresponding
transfer function magnitude is given in Figure 21.2.

275
Figure 21.2. Butterworth Low-Pass Filter, 100 Hz, Transfer Function

Note that the curve in Figure 21.2 has a gradual roll-off beginning at about 70 Hz. Ideally, the transfer
function would have a rectangular shape, with a corner at (100 Hz, 1). This ideal is never realized in
practice due to stability concerns in the time domain. Thus, a compromise is usually required to select the
cut-off frequency. The transfer function could also be represented in terms of a complex function, with
real and imaginary components. A transfer function magnitude plot for a sixth-order Butterworth filter
with a cut-off frequency of 100 Hz as shown in Figure 21.3.

276
Figure 21.3. Butterworth High-Pass Filter, 100 Hz, Transfer Function

277
The curves in the previous figures suggests that filtering could be achieved as follows:

1. Take the Fourier transform of the input time history

2. Multiply the Fourier transform by the filter transfer function, in complex form
3. Take the inverse Fourier transform of the product

The above frequency domain method is valid. Nevertheless, the filtering algorithm is usually implemented
in the time domain for computational efficiency, maintaining stability, avoid leakage error, etc.

Time domain
Input x Response y
equivalent of
H()

Figure 21.4. Filter Block Diagram

The filter transfer function can be represented by H ( ) as shown in Figure 21.4.

The filtering transfer function is implemented as a digital recursive filtering relationship in the time
domain. The response y k at time index k is

 L   L 
y k =   b n x k − n  −   a n y k −n  (21.1)
n =0  n =0 

The filter order is L . The filter coefficients a n and b n can be derived using the method is Reference [43].
The equation is recursive because the output at any time depends on the output at previous times.

278
21.4 PHASE CORRECTION

x Time Time domain

Time Time domain y
Reversal equivalent of equivalent of
Reversal
H() H()

Figure 21.5. Refiltering for Phase Correction

An ideal filter should provide linear phase response. This is particularly desirable if shock response spectra
calculations are required. Butterworth filters, however, do not have a linear phase response. Other IIR
filters share this problem. A number of methods are available, however, to correct the phase response.
One is the refiltering method in Reference [43]. An important note about refiltering is that it reduces the
transfer function magnitude at the cut-off frequency to –6 dB, as shown in example in Figure 21.6.

Figure 21.6. Butterworth Filter 6th Order Low-Pass Filter with Refiltering, Transfer Function

279
21.5 PSD CALCULATION VIA DIGITAL FILTER
Recall the synthesized time history from Figure 17.10 which satisfied the Navmat P-9492 PSD specification.
Calculate the power spectral density by applying band-pass filtering over successive bands. This method
is useful as a means of better understanding the PSD function. The filtering is performed using a
Butterworth sixth-order filter without refiltering for phase correction. The bands are full octave. The
results are shown in Table 21.1.

Table 21.1. Band-pass Filter Results, Octave Bands

fl (Hz) fc (Hz) fu (Hz) Δf (Hz) GRMS GRMS^2 GRMS^2/Hz

14 20 28 14 0.36 0.13 0.0091
28 40 57 29 0.80 0.64 0.0226
57 80 113 56 1.47 2.17 0.0384
113 160 226 113 2.16 4.67 0.0413
226 320 453 227 2.92 8.54 0.0377
453 640 905 452 3.11 9.70 0.0214
905 1280 1810 905 3.02 9.15 0.0101
1810 2560 3620 1810 1.31 1.71 0.0009

The band center frequency is fc. The lower and upper band limits are fl and fc, respectively. The
bandwidth is Δf which is the difference between the band limits. A GRMS value is calculated for each
band. Then the GRMS value is squared and divided by the bandwidth. The last column is the PSD in
GRMS^2/Hz. The time histories for the two rows in the table are shown in Figure 21.7. The GRMS^2/Hz
points are plotted along with Navmat PSD in Figure 21.8. Recall that GRMS^2/Hz is abbreviated as G^2/Hz.

280
Figure 21.7. Synthesized Time History with Two Filtering Cases, Close-up View

The bands are in full octave format.

281
Figure 21.8. PSD Comparison

The PSD points calculated via band-pass filtering track the specification fairly well. The dropout for the
last point is not a concern because the bandwidth extended from 1810 to 3620 Hz, but the specification
stopped at 2000 Hz.

282
21.6 SEISMIC DATA FILTERING

The following example shows how filtering can aid in analyzing seismic data from a seismometer. Filtering
is used to find onset of P-wave in seismic time history from Solomon Island earthquake, Magnitude 6.8,
October 8, 2004. The measured data is from a homemade seismometer in Mesa, Arizona.

Figure 21.9. Author’s Homemade Lehman Seismometer

The boom is a horizontal pendulum. It has a period of 14.2 seconds, equivalent to a natural frequency of
0.071 Hz. A sensor at the free end measures the displacement. The boom length is 64 inch. The total frame
height is 35 inch. The boom has a knife edge that pivots against a bolt head in the lower cross-beam of
the frame.

The boom is suspended from the frame by a wire cable. The cable is attached to the top cross-beam of
the frame. The other end of the cable is attached to the boom, about two-thirds of the distance from the
pivot to the free end of the boom. The pivot point is offset from the top cable attachment point. Thus,
the boom oscillates as if it were a "swinging gate."

The plate supporting the frame has three adjustable mounting feet. The feet can be adjusted to tune the
pendulum to the desired natural frequency. Furthermore, the wire cable has a turnbuckle which is used

283
to adjust height of the free end of the boom. The detached frame in the center of the figure is used for
assembly and to limit the displacement during tuning.

Figure 21.10. Non-contact Relative Displacement Sensor

The classic sensor for Lehman seismometers is a magnetic coil attached to the boom. The coil moves
between the poles of a magnet, thus inducing a voltage proportional to the velocity. The signal requires
tremendous amplification. The author has better success with the inductive displacement sensor.

284
Figure 21.11. Ballast Mass Partially Submerged in Oil

The damping method is shown in this figure.

Figure 21.12. Pivot at the End of the Boom

The chisel blade butted up against a chrome-plated bolt head.

285
Figure 21.13. Seismometer Free Vibration Response

The seismometer was given an initial displacement and then allowed to vibrate freely. The period was
14.2 seconds, with 9.8% damping. The corresponding frequency is 0.0705 Hz.

286
Figure 21.14. Measured Solomon Earthquake Time History, Lehman Seismometer Solomon_1.emf

The trace in Figure 21.14 is the displacement time history of the Solomon Islands earthquake on October
8, 2004 as measured by the Lehman seismometer via the inductive sensor. The seismometer was located
at home in Mesa, Arizona. The data was acquired by a Nicolet Vison system. The Nicolet sample rate was
set to 50 samples per second with its lowpass filter set to 5 Hz.

A Krohn-Hite filter, model 3343, was used to high-pass filter the analog displacement signal at 0.03 Hz
prior to its input to the Nicolet system. It also provided a 20 dB gain.

The time is referenced to the earthquake occurrence using the USGS data. The plot’s Y-axis is labeled as
relative displacement because it is the response of the boom relative to the ground. Further calculation
would be required to estimate the true ground motion. The time history shows that the Earth is
remarkably reverberant. The oscillations last well over one hour. The phase components are P primary
wave, S secondary or shear wave, and LQ Love wave. Recall the seismic wave diagrams in Section 19.6.1.

The P-wave is indiscernible against the background microseismic noise. Nevertheless, it can be extracted
by additional filtering, as shown in the next figure.

287
Figure 21.15. High-pass Filtered Seismic Waveform Showing Onset of P-wave

The P-wave arrives at in Mesa, Arizona about 760 seconds after it was generated near the Solomon
Islands. The S-wave arrives around 1400 seconds. There are also some intermediate waveforms which
begin as P-waves and transform into other types through reflection and refraction. Both the P and S-
waves are body waves which travel faster than surface waves. The frequency content of body and surface
waves is shown in Table 21.2, as taken from Reference [1].

Table 21.2. Seismic Waveform Periods and Frequencies

Wave Type Period Natural Frequency

(sec) (Hz)
Body 0.01 to 50 0.02 to 100
Surface 10 to 350 0.003 to 0.1

288
21.7 INTEGRATION AND DIFFERENTIATION OF ACCELEROMETER SIGNALS

There are occasional needs to integrate accelerometer time histories to velocity. This may be done to
evaluate the accuracy of pyrotechnic shock data where the velocity is expected to have a stable oscillation
about its zero baseline. Accelerometer data may also be double-integrated to displacement, but the
accuracy depends on the frequency response of the accelerometers. Some accelerometer types such as
variable capacitance and servo motor designs can accurately measure acceleration down to zero
frequency. The common piezoelectric accelerometer may have a practical lower limit of a few Hertz,
however. Some type of high-pass filter is needed in this case. The high-pass filtering may be performed
using an analog filter in the accelerometer’s signal conditioner. Or the filtering may be performed digitally
as a post-processing step.

Another integration scenario is the case where an acceleration time history is synthesized for use in a
modal transient analysis. The corresponding velocity and displacements should each be stable, especially
if dynamic stresses are to be calculated. Consider the simple case of the synthesized white noise
acceleration time history in Figure 21.16.

The integration steps are performed using the trapezoidal method which has sufficient accuracy. The
resulting velocity displacement each has a drift. These effects arise in part because the integration does
not account for initial velocity and initial displacement, each of which are either unknown or unspecified.
The next task is to correct the acceleration time history so that its velocity and displacement will each be
stable. There is no one right way to do this, but an example is shown in Figure 21.17. The correction steps
were:

1. Remove mean from acceleration

2. High-pass filter acceleration at 6 Hz, Butterworth 6th order with refiltering for phase
correction
3. Apply fade in and out to acceleration using 2% of the total duration at each end
4. Integrate to velocity
5. Remove mean from velocity
6. Apply fade in and out to velocity using 2% of the total duration at each end
7. Integrate to displacement
8. Perform first-order trend removal on displacement
9. Apply fade in and out to displacement using 2% of the total duration at each end
10. Differentiate displacement to velocity
11. Differentiate velocity to acceleration Integrate corrected acceleration to velocity
12. Integrate velocity to displacement
13. Verify that the resulting velocity and displacement each have stable oscillations about their
respective zero baselines and that each begins and ends at zero

289
Figure 21.16. White Noise Acceleration Time History with Corresponding Velocity and Displacement

290
Figure 21.17. Acceleration Time History Correction and Verification

Nine subplots appear as follows:

Row 1, Col 1 - processed acceleration prior to integration

Row 2, Col 1 - integrated velocity
Row 3, Col 1 - double integrated displacement

Row 3, Col 2 - double integrated displacement (repeated)

Row 2, Col 2 - differentiated velocity from double integrated displacement
Row 1, Col 2 - double differentiated acceleration from double integrated displacement

291
 Row 1, Col 3 - double differentiated acceleration (repeated)
Row 2, Col 3 - integrated velocity from double differentiated acceleration
Row 3, Col 3 - double integrated displacement from double differentiated acceleration

The corrected acceleration is the time history shown in Row 1, Col 2 and repeated in Row 1, Col 3.

Figure 21.18. Fade In and Out Windows

Note that the fade in and out windows are similar to those shown in Figure 21.19.

292
21.8 BAND LIMITED WHITE NOISE

Figure 21.20. White Noise Segment Before and After Low-Pass Filtering

White noise is an idealized concept. True white noise would have a flat PSD level from zero to infinity
frequency units. Digital white noise is inherently frequency limited by its sample rate. Furthermore,
white noise may be subject to low-pass filtering to smooth it prior to its use as an excitation function, as
shown in

Figure 21.20. The top plot is a white noise sampled at 10,000 samples/sec and has a choppy appearance.
Aliasing would be suspected if the top plot were measured data.

The bottom plot is the same time history low-pass filtered at 2000 Hz which smooths the data. Either
time history should be suitable for either applied force or base excitation of an SDOF system with a natural
frequency of 200 Hz, for example. Recall that the SDOF system itself behaves a mechanical filter. But the
low-pass filtered version could have some potential numerical accuracy benefits, although an
investigation of this hypothesis is left as future work.

293
21.9 LAUNCH VEHICLE COUPLED-LOADS ANALYSIS FILTERING

Figure 21.21. Sources of Aero-Mechanical Loads

Coupled Loads Analysis (CLA) predicts payload & launch vehicle responses due to major dynamic and
quasi-static loading events. This analysis is performed prior to launch. It can also be performed as post-
flight data reconstruction using flight accelerometer data. Flight accelerometer data is low-pass filtered
for coupled-loads analyses. The cut-off frequency varies by launch vehicle, payload, key events, etc. The
primary sources of these low frequency loads are:

1. Pre-launch events: ground winds, seismic loads

2. Liftoff: engine/motor thrust build-up, ignition overpressure, pad release
3. Air-loads: buffet, gust, static-elastic
4. Liquid engine ignitions and shutdowns

Reference [44] gives the following frequency guidelines for CLA. The low-frequency dynamic response of
the launch vehicle/payload system to transient flight events is typically analyzed from 0 Hz to 100 Hz. The
response frequency domain can be up to 150 Hz for some small launch vehicles.

294
22 THE GREAT SPACECRAFT BASE INPUT VIBRATION TEST DEBATE

22.1 DEBATE INTRODUCTION

Figure 22.1. Shaker Table Vibration Test

Engineers prepare the MESSENGER spacecraft for a vibration test on a shaker table at The Johns Hopkins
University Applied Physics Laboratory, Laurel, Maryland. The purpose of this test is to verify that the
spacecraft can withstand the launch vehicle environment during liftoff and ascent.

295
Flight configured spacecraft are subjected to base input vibration tests for certain programs. The
spacecraft are often one-of-a-kind, so the vibration test is effectively a proto-qualification test covering
both design and workmanship verification. The tests may be sinusoidal or random. The sine vibration is
typically at low frequencies, below 100 Hz. The author’s colleagues are divided on whether these
spacecraft system-level tests are prudent and effective.

22.2 ARGUMENTS AGAINST SPACECRAFT VIBRATION TESTING

The following assertions are made by A.M. Kabe and E. Perl from The Aerospace Corporation in informal
working groups meetings.

1. Vibration tables cannot replicate the impedance of the launch vehicle interface, nor the
interaction that occurs between the launch vehicle and spacecraft when they are a coupled
system; hence, the modes of vibration will not be the same as in flight.
2. Only translational motions are applied at the base, one axis at a time, whereas during flight,
the launch vehicle/spacecraft system will vibrate simultaneously in all six degrees of freedom
at each mass point and at each interface point between the launch vehicle and spacecraft.
3. The total acceleration load during powered flight also depends on the spacecraft rigid-body
acceleration which a shaker cannot replicate.
4. Derivation of a “base input” environment from a few accelerometer locations at the
launch vehicle/spacecraft interface will generally lead to an over prediction of the motions
at the interface, since local deformations are mapped on the assumption that the interface
acts as a rigid plane.
5. The use of (response)/(base motion) ratios to extract damping, a common practice, is not a
valid approach for multi-degree-of-freedom systems it fails to account for the mode
participation factor.
6. The test article may not include the actual spacecraft launch vehicle adapter or the
propellant mass in the tanks because of safety and contamination concerns.
7. The test requirement forces the spacecraft organization to design its system to not only
survive the launch environment, but also to survive an artificial test that more often than
not produces overly conservative loads in many parts of the structure while not adequately
testing others.
8. The test can pose unnecessary risk of damaging flight hardware late in the program.

Note that A.M. Kabe advocates acoustic reverberant chamber testing of spacecraft as a workmanship
screen, as an alternative to base shake testing. He also favors shaker table vibration testing on a
component or subsystem level.

296
22.3 ARGUMENTS FOR SPACECRAFT VIBRATION TESTING

The following justification points are made by NASA engineers Daniel Kaufman, Scott Gordon, Steve
Hendricks and Dennis Kern.

1. Essentially all current launch vehicle organizations (Delta, Atlas, Taurus, Pegasus, Ariane, HII,
Proton, Long March, Falcon, etc.) specify and require or strongly recommend a spacecraft
sine or random vibration test.
2. Testing is also required by NASA documents, such as References [45] and [46].
3. Insurance companies require vibration tests on all commercial communications satellites.
4. Some test facilities have the capability to perform simultaneous multi-axis vibration testing,
as needed.
5. The vibration test provides qualification for tertiary/ancillary hardware that would not
otherwise be tested. This includes: Cable harnesses, bellows, connectors, actuators,
plumbing lines, wave guides, brackets, dampers, shades and shields, articulation &
deployment mechanisms, shunt heaters, louvers, purge equipment, hinges and restraints,
blankets & supports.
6. The test provides an opportunity to determine the structural linearity in the operational
vibration range of response. Note that linearity is typically assumed in coupled loads analysis
(CLA).
7. Force limiting reasonably accounts for the interaction with the base motion, and has been
effectively employed in spacecraft vibration testing, thus reducing the potential for an over-
test at the spacecraft’s natural frequencies in the test configuration. The force limiting takes
into account the CLA response levels.
8. Force gauges under the spacecraft provide a very accurate method of measuring and limiting
to the CLA loads during the vibration test for mid to high apparent mass modes.
9. Numerous case histories have shown that vibration testing is effective for uncovering design
or workmanship flaws which would have otherwise caused mission degradation or failure.

22.4 SPACECRAFT TESTING EXPERIENCES

Note that NASA/GSFC typically uses sine vibration testing, whereas JPL tends to use random.
A few examples from sine testing at GSFC are:
▪ TRMM: During Observatory sine testing, found that the NASDA supplied PAF clamp band had
insufficient tension and gapped during the test. As a result, the clamp band tension was increased
for flight.
▪ GOES had a workmanship problem involving a missing or loose bolt which caused structural failure
of a mission-critical antenna. It was detected during the lateral sine test.

297
 ▪ NOAA-K experienced IMU saturation during sine sweep testing. Because the spacecraft IMU
provides guidance information for the Titan II launch vehicle during ascent, IMU saturation during
launch would have resulted in a mission failure. Changes were made and launch vehicle restraints
were implemented to resolve the problem, including wind restrictions at launch and a commanded
first stage shutdown vs. fuel depletion.
▪ TDRS-H: During the sine vibration test, the first two modes for the Space Ground Link antenna (SGL)
were lower than predicted by the model. The first mode dropped from 15 Hz to 11 Hz and the
second mode dropped from 33 Hz to 25 Hz. It turned out that the mathematical model of this
“simple” antenna was wrong and therefore the Verification Loads Cycle had to be rerun.

A few examples from random vibration testing at JPL are:

▪ Cassini: Experienced an RTG electrical short to its spacecraft mount in system random vibration
test. Significant degradation in spacecraft electrical power could have resulted. Spacecraft mount
was redesigned.
▪ CloudSat: Cloud Profiling Radar waveguide failure in spacecraft random vibration test due to
apparent poor workmanship of adhesive bonding. Possible loss of science data averted.
▪ MER 1: Fundamental modes of the Rover in spacecraft random vibration test were 20% greater
than predicted in all three axes. (Fixed base modal test had been performed on Rover, Lander, and
Cruise Stage separately; FE models were then combined. Estimated stiffness of Lander attachment
to Rover was too low.) FE model was updated just in time for the verification CLA cycle. Vibration
test also revealed improper torque of bolts on some tanks in low level runs. Bolts were properly
torqued and test completed successfully.
▪ MSL Rover: experienced several motor encoder screws backed out of at least one of the Rover
actuators during Rover random vibration test. The actuators are used throughout Rover and the
issue was unlikely to have otherwise been found before launch, which could have been a serious
threat to the mission.
JPL prefers random vibration because it easier to control, particularly with respect to force limiting.
Note that JPL tested the SMAP spacecraft to the following workmanship PSD:
20 Hz to 250 Hz, 0.01 G^2/Hz
The overall level was 1.5 GRMS. The duration was one minute. The vibration was applied in the vertical
axis and in one “45 degree” lateral axis.

22.5 SINE SWEEP CONTROL

Sine sweep vibration is more difficult to control than random especially for the case of lightly-damped
modes. The sweep rate, compression factor, and tracking filter must be selected with great care.
INVAP experienced control issues during sine vibration testing of the ARSAT-1 structural test model. But
note that dynamic simulator were used for the test, which can have much less damping than the actual
flight hardware.

298
23 TACOMA NARROWS BRIDGE FAILURE

23.1 BRIDGE INTRODUCTION

Figure 23.1. Torsional Mode of the Tacoma Narrows Bridge

The Tacoma Narrows Bridge failure is an example of self-excited vibration of such importance that it
warrants its own section in this document. The original Tacoma Narrows Bridge was opened to traffic on
July 1, 1940. It was located in Washington State, near Puget Sound.

The Tacoma Narrows Bridge was the third-longest suspension bridge in the United States at the time, with
a length of 5939 feet including approaches. Its two supporting towers were 425 feet high. The towers
were 2800 feet apart.

Prior to this time, most bridge designs were based on trusses, arches, and cantilevers to support heavy
freight trains. Automobiles were obviously much lighter. Suspension bridges were both more elegant
and economical than railway bridges. The suspension design thus became favored for automobile traffic.
Unfortunately, engineers did not fully understand the forces acting upon bridges. Neither did they
understand the response of the suspension bridge design to these poorly understood forces.

Furthermore, the Tacoma Narrows Bridge was built with shallow plate girders instead of the deep
stiffening trusses of railway bridges. Note that the wind can pass through trusses. Plate girders, on the
other hand, present an obstacle to the wind. As a result of its design, the Tacoma Narrows Bridge
experienced rolling undulations which were driven by the wind. It thus acquired the nickname "Galloping
Gertie."

299
23.2 BRIDGE COLLAPSE

Figure 23.2. Tacoma Narrows Bridge Collapse

Strong winds caused the bridge to collapse on November 7, 1940. Initially, 35 mile per hour winds excited
the bridge's transverse vibration mode, with an amplitude of 1.5 feet. This motion lasted 3 hours. The
wind then increased to 42 miles per hour. In addition, a support cable at mid-span snapped, resulting in
an unbalanced loading condition. The bridge response thus changed to a 0.2 Hz torsional vibration mode,
with an amplitude up to 28 feet. The torsional mode is shown in Figure 23.1.

The torsional mode shape was such that the bridge was effectively divided into two halves. The two halves
vibrated out-of-phase with one another. In other words, one half rotated clockwise, while the other
rotated counter-clockwise. The two half spans then alternate polarities.

One explanation of this is the "law of minimum energy." A suspension bridge may either twist as a whole,
or divide into half spans with opposite rotations. Nature prefers the two half-span option since this
requires less wind energy. The dividing line between the two half-spans is called the "nodal line." Ideally,
no rotation occurs along this line. The bridge collapsed during the excitation of this torsional mode.
Specifically, a 600 foot length of the center span broke loose from the suspenders and fell a distance of
190 feet into the cold waters below. The failure is shown in Figure 23.2.

300
23.3 FAILURE THEORIES

23.3.1 Candidates

The fundamental weakness of the Tacoma Narrows Bridge was its extreme flexibility, both vertically and
in torsion. This weakness was due to the shallowness of the stiffening girders and the narrowness of the
roadway, relative to its span length. Engineers still debate the exact cause of its collapse, however. Three
theories are:

1. Random turbulence
2. Periodic vortex shedding
3. Aerodynamic instability (negative damping)

These theories are taken from Reference [47]. Aerodynamic instability is the leading candidate.

23.3.2 Random Turbulence

An early theory was that the wind pressure simply excited the natural frequencies of the bridge. This
condition is called "resonance." The problem with this theory is that resonance is a very precise
phenomenon, requiring the driving force frequency to be at, or near, one of the system's natural
frequencies in order to produce large oscillations. The turbulent wind pressure, however, would have
varied randomly with time. Thus, turbulence would seem unlikely to have driven the observed steady
oscillation of the bridge.

23.3.3 Vortex Shedding

Theodore von Karman, a famous aeronautical engineer, was convinced that vortex shedding drove the
bridge oscillations. A diagram of vortex shedding around a spherical body is shown in Figure 23.3. Von
Karman showed that blunt bodies such as bridge decks could also shed periodic vortices in their wakes.

A problem with this theory is that the natural vortex shedding frequency was calculated to be 1 Hz, as
shown in Section 23.4. This frequency is also called the "Strouhal frequency." The torsional mode
frequency, however, was 0.2 Hz. This frequency was observed by Professor F. B. Farquharson, who
witnessed the collapse of the bridge. The calculated vortex shedding frequency was five times higher than
the torsional frequency. It was thus too high to have excited the torsional mode frequency.

In addition to "von Karman" vortex shedding, a flutter-like pattern of vortices may have formed at a
frequency coincident with the torsional oscillation mode. Whether these flutter vortices were a cause or
an effect of the twisting motion is unclear.

301
Figure 23.3. Vortex Shedding around a Spherical Body

23.3.4 Aerodynamic Instability

Aerodynamic instability is a self-excited vibration. In this case, the alternating force that sustains the
motion is created or controlled by the motion itself. The alternating force disappears when the motion
disappears. This phenomenon is also modeled as free vibration with negative damping.

Airfoil flutter and transmission line galloping are related examples of this instability. Further explanations
of instability are given in References [48], [49] and [50].

The following scenario shows how aerodynamic instability may have caused the Tacoma Narrows Bridge
to fail. For simplicity, consider the motion of only one span half. Assume that the wind direction was not
perfectly horizontal, perhaps striking the bridge span from below, as shown in Figure 23.4.

302
 Wind

Figure 23.4. Span at Initial Rest Position

Thus, the bridge is initially at an angle-of-attack with respect to the wind. Aerodynamic lift is generated
because the pressure below the span is greater than the pressure above. This lift force effectively places
a torque, or moment, on the bridge. The span then begins to twist clockwise as show in Figure 23.5.
Specifically, the windward edge rotates upward while the leeward edge rotates downward.

Wind

Figure 23.5. Span Rotates Clockwise

303
 Wind

Figure 23.6. Span Rotates Counter-Clockwise

The span has rotational stiffness, however. Thus, elastic strain energy builds up as the span rotates.
Eventually, the stiffness moment overcomes the moment from the lift force. The span then reverses its
course, now rotating counter-clockwise

The span's angular momentum will not allow it to simply return to its initial rest position, however. The
reason is that there is little or no energy dissipation mechanism. Thus, the span overshoots its initial rest
position. In fact, it overshoots to the extent that the wind now strikes the span from above as shown in
Figure 23.6. The wind's lift force now effectively places a counter-clockwise moment on the span.

Once again, strain energy builds up in the span material. Eventually, the stiffness moment exceeds the
moment from the wind's lift force. The span thus reverse course, now rotating clockwise. Again, it
overshoots its rest position. The cycle of oscillation begins anew from the position shown in Figure 23.4,
except that the span now has rotational velocity as it passes through the original rest position.

The cycles of oscillation continue in a repetitive manner.

Note that the wind force varies as a function of the span angle during the cycle. The wind force may also
vary with the angular velocity. The wind force is not a function of time, however.

Eventually, one or two failure modes occurs. One possibility is that the span experiences fatigue failure
due to an excessive number of stress reversals. The other is that the angular displacement increased in
an unstable manner until the material is stressed beyond its yield point, and then beyond its ultimate
stress limit.

These two failure modes can be interrelated. For example, accumulated fatigue effectively lowers the
yield and ultimate stress limits. Regardless, the bridge collapses.

As a final note, the aerodynamic instability oscillation is not a resonant oscillation since the wind does not
have a forcing frequency at, or near, the bridge's torsional mode frequency. Some physics and
engineering textbooks mistakenly cite the Tacoma Narrows Bridge as an example of resonance. This

304
problem is discussed in Reference [51]. Nevertheless, the bridge's collapse remains the most famous
structural failure due to vibration.

23.4 BRIDGE’S STROUHAL FREQUENCY

The original Tacoma Narrows Bridge collapsed in 1940. It experienced severe torsional oscillations driven
by a 42 mile per hour wind. The fundamental weakness of the Tacoma Narrows Bridge was its extreme
flexibility, both vertically and in torsion. This weakness was due to the shallowness of the stiffening
girders and the narrowness of the roadway, relative to its span length. A cross-section of the bridge is
shown in Figure 23.7.

39 ft

8 ft

Figure 23.7. Cross-section of Tacoma Narrows Bridge Span

The Strouhal number for the bridge cross-section is S = 0.11 according to Reference [51]. Furthermore,
the characteristic dimension is the girder height 8 ft, which may be substituted for the diameter. The
Strouhal frequency for a 42 mph (61.6 ft/sec) wind is thus

fs = S U / D (23.1)

fs = ( 0.11) ( 61.1 ft / sec ) / ( 8 ft ) = 0.84 Hz (23.2)

The Strouhal frequency is rounded to 1 Hz in Reference [51]. On the other hand, the bridge's torsional
oscillation was observed to be 0.2 Hz. The Strouhal frequency was thus at least two octaves greater than
the torsional mode frequency. Billah and Scanlan thus argue in Reference [51] that the Karman vortex
shedding could not have driven the torsional mode oscillation. They argue instead that an aerodynamic
instability resulted in a self-excited oscillation, which caused the failure. Several decades after the failure,
the matter is not completely settled.

23.5 REPLACEMENT BRIDGE

305
Figure 23.8. New Tacoma Narrows Bridge

A new Tacoma Narrows Bridge was built in 1950, as shown in

Figure 23.8. The new bridge had truss-girders which allowed the winds to pass through. It also had
increased torsional stiffness because it was thicker and wider. Furthermore, wind tunnel testing was
performed to verify the design of the new bridge prior to its construction. A similar bridge was later built
alongside the new bridge in the figure to ease traffic flow, with one-way traffic on each bridge.

306
24 BIBLIOGRAPHY

[1] T. Lay and T. Wallace, Modern Global Seismology, New York: Academic Press, 1995.

[2] ISO 2631-1, Evaluation of human exposure to whole-body vibration, Part 1: General requirements,
International Organization for Standardization, 1997.

[3] T. Irvine, "Free Vibration of a Single-Degree-of-Freedom System, Revision B," Vibrationdata, 2005.

[4] T. Irvine, "Free Vibration of a Two-Degree-of-Freedom System Subjected to Initial Velocity and
Displacement, Revision C," Vibrationdata, 2011.

[5] T. Irvine, "Longitudinal Natural Frequencies of Rods and Response to Initial Conditions, Revision
D," Vibrationdata, 2010.

[6] M. Bai and W. Thatcher, High G Pyrotechnic Shock Simulation using Metal-to-Metal Impact,
Washington DC: The Shock and Vibration Bulletin, Naval Research Laboratory, September 1979.

[7] R. Blevins, Formulas for Natural Frequency and Mode Shapes, Malabar, Florida: R. Krieger, 1979.

[8] T. Irvine, "Bending Frequencies of Beams, Rods and Pipes, Revision T," Vibrationdata, 2018.

[9] A. Leissa, Vibration of Plates, NASA SP-160, Washington D.C.: National Aeronautics and Space
Administration, 1969.

[10] E. Baker, Structural Analysis of Shells, Malabar, Florida: Krieger, 1972.

[11] D. S. Steinberg, Vibration Analysis for Electronic Equipment, Third Edition, New York: Wiley, 2000.

[12] R. Marshall, L. Phan and M. Celebi, "Full-Scale Measurement of Building Response to Ambient
Vibration and the Loma Prieta Earthquake," Fifth U.S. National Conference on Earthquake
Engineering, Earthquake Awareness and Mitigation Across the Nation, vol. 11, 1994.

[13] M. Collins, Liftoff: The Story of America's Adventure in Space, New York City: Grove Press, 1989.

[14] M. Griffin, Handbook of Human Vibration, London: Academic Press, 1990.

[15] T. Irvine, "Modal Transient Analysis of a System Subjected to an Applied Force, Revision L,"
Vibrationdata, 2013.

[16] J. Emmerson, "“On Shaky Ground”," Letters, Smithsonian Air & Space, May 2001.

[17] T. Irvine, "The Generalized Coordinate Method for Discrete Systems, Revision F," Vibrationdata,
2012.

307
[18] D. O. Smallwood, "An Improved Recursive Formula for Calculating Shock Response Spectra," Shock
and Vibration Bulletin, No. 51, May 1981.

[19] T. Irvine, " Two-Stage Isolation for Harmonic Base Excitation, Revision C," Vibrationdata, 2015.

[20] T. Irvine, "Derivation of the Filter Coefficients for the Ramp Invariant Method as Applied to Base
Excitation of a Single-degree-of-Freedom System, Revision D," Vibrationdata, 2018.

[21] T. Irvine, "Modal Transient Vibration Response of a Cantilever Beam Subjected to Base Excitation,"
Vibrationdata, 2016.

[22] T. Scharton and D. Pankow, "Extreme Peaks in Random Vibration Testing," in Aerospace/JPL,
Spacecraft and Launch Vehicle Dynamic Environments Workshop, Hawthorne, California, 2006.

[23] K. Ahlin, "Comparison of Test Specifications and Measured Field Data," Sound & Vibration,
September 2006.

[24] H. Himelblau, A. Piersol and e. al, IES Handbook for Dynamic Data Acquisition and Analysis, Mount
Prospect, Illinois, 1994.

[25] R. Randall, Frequency Analysis, 3rd edition, Bruel & Kjaer, 1987.

[26] F. Harris, Trigonometric Transforms, Technical Publication DSP-005, San Diego, CA: Scientific-
Atlanta.

[27] MIL-STD-810E, Environmental Engineering Considerations and Laboratory Tests, United States
Department of Defense, 1989.

[28] NAVMAT P-9492, Navy Manufacturing Screening Program, Washington D.C.: Department of the
Navy, Naval Material Command, 1979.

[29] JEDEC JESD22-B111A, Board Level Drop Test Method of Components for Handheld Electronic
Products Standard, JEDEC Solid State Technology Association, 2016.

[30] NASA-HDBK-7005, Dynamic Environmental Criteria, Washington DC: National Aeronautics and
Space Administration, 2001.

[31] "Mechanical Shock Design and Verification Handbook, ECSS-E-HB-32-25A," European Cooperation
for Space Standardization, Noordwijk, The Netherlands, July 2015.

[32] D. Albert and J. Orcutt, "Observations of Low-frequency Acoustic-to-seismic Coupling in the

Summer and Winter," Journal of the Acoustical Society of America , vol. (86), July 1989.

[33] W. Kacena, M. McGrath and A. Rader, "Aerospace Systems Pyrotechnic Shock Data, Vol. VI., NASA
CR 116406," Goddard Space Flight Center, Maryland, 1970.

[34] European Cooperation for Space Standardization, "Mechanical Shock Design and Verification
Handbook," Noordwijk, The Netherlands, July 2015.

308
[35] R. Ferebee and T. Irvine, NASA/TM-2008-215253, An Alternative Method of Specifying Shock Test
Criteria, 2008.

[36] T. Irvine, "Modal Transient Analysis of a Multi-degree-of-freedom System with Enforced Motion,
Revision C," Vibrationdata, 2011.

[37] R. Chalmers and H. Gaberson, "Modal Velocity as a Criterion of Shock Severity," Shock and
Vibration Bulletin, Naval Research Lab, December 1969.

[38] T. Irvine, "Shock and Vibration Stress as a Function of Velocity, Revision G," Vibrationdata, 2013.

[39] A. Piersol and T. Paez, Harris’ Shock and Vibration Handbook, 6th Edition; V. Bateman, Shock
Testing Machines, Chapter 27, New York City: McGraw-Hill, 2010.

[40] H. Gaberson, "Shock Severity Estimation," Sound & Vibration, January 2012.

[41] SMC-TR-06-11 AEROSPACE REPORT NO. TR-2004(8583)-1 Rev. A, Test Requirements for Launch,
Upper-Stage, and Space Vehicles, El Segundo, CA: Aerospace Corporation.

[42] R. Morse, Presentation at Spacecraft & Launch Vehicle Dynamics Environments Workshop
Program, Aerospace Corp., El Segundo, CA, June 2000.

[43] S. Sterns and R. David, Signal Processing Algorithms in Fortran and C, New Jersey: , Prentice Hall,
1993.

[44] ECSS-E-HB-32-26A, Spacecraft Mechanical Loads Analysis Handbook, Noordwijk, The Netherlands:
European Cooperation for Spacecraft Standardization (ECSS), 2013.

[45] NASA-STD-7002B, Payload Test Requirements, Washington D.C.: National Aeronautics and Space
Administration, 2018.

[46] GSFC-STD-7000 REV A, General Environmental Verification Statndard (GEVS), Greenbelt, MD: NASA
Goddard Space Flight Center, 2013.

[47] J. Koughan, The Collapse of the Tacoma Narrows Bridge, Evaluation of Competing Theories of its
Demise, and the Effects of the Disaster of Succeeding Bridge Designs, The University of Texas at
Austin, 1996.

[48] D. Hartog, Mechanical Vibrations, New York: McGraw-Hill, 1940.

[49] M. Levy and M. Salvadori, Why Buildings Fall Down, New York: Norton, 1992.

[50] H. Bachmann, Vibration Problems in Structures, Berlin: Birkhauser Verlag, 1995.

[51] K. Billah and R. Scanlan, "Resonance, Tacoma Narrows Bridge Failure, and Undergraduate Physics,
Textbooks," American Journal of Physics, 1991.

309
310