7
Response of SDOF Systems
to Periodic Excitation:
Frequency-Domain Analysis
In Chapter 4 you studied the response of SDOF systems to harmonic excitation and
became familiar with important concepts such as resonance. You also learned how to
simplify the analysis of viscous-damped systems by using complex frequency response.
These concepts from Chapter 4 are now extended to determine the response of SDOF
systems to periodic excitation.
Upon completion of this chapter you should be able to:
•
•
•
•
•
•
Determine the Fourier series representation of a periodic function using the real
form of the Fourier series.
Determine the Fourier series representation of a periodic function using the complex form of the Fourier series.
Determine the steady-state response of an SDOF system to periodic excitation
using either the real form or the complex form.
Use the basic definition to determine the Fourier transform of a transient function
p(t).
Show that the impulse-response function for a linear system and the frequencyresponse function for the system form a Fourier transform pair.
Apply the FFT algorithm to compute the Fourier transform of a periodic function
p(t), and plot the magnitude and phase of the Fourier transform computed.
7.1 RESPONSE TO PERIODIC EXCITATION: REAL
FOURIER SERIES
Forces acting on structures are frequently periodic, or can be approximated closely by
periodic forces. For example, the forces exerted on an automobile traveling at constant
speed over certain roadway surfaces can be considered to be periodic. Figure 7.1 shows
a periodic function with period T1 , that is,
p(t + T1 ) = p(t)
(7.1)
A periodic function can be separated into its harmonic components by means of a
Fourier series expansion. In this section we consider real Fourier series. In Section 7.2
184
7.1 Response to Periodic Excitation: Real Fourier Series
185
p(t )
t
T1
Figure 7.1
Periodic function with period T1 .
complex Fourier series are introduced. The complex form is very useful when combined
with the complex frequency-response function of Chapter 4 to study the steady-state
response of damped systems.
7.1.1
Real Fourier Series
The periodic function p(t) may be separated into its harmonic components by means of
a Fourier series expansion. The real Fourier series expansion of p(t) may be defined as
p(t) = a0 +
∞
!
n=1
an cos n.1 t +
∞
!
bn sin n.1 t
(7.2)
n=1
where
.1 =
2π
T1
(7.3)
is the fundamental frequency (in rad/s), and an and bn are the coefficients of the nth
harmonic. The coefficients a0 , an , and bn are related to p(t) by the equations
a0 =
1
T1
2
an =
T1
bn =
2
T1
%
τ +T1
τ
%
τ +T1
p(t) dt = average value of p(t)
p(t) cos n.1 t dt,
n = 1, 2, . . .
p(t) sin n.1 t dt,
n = 1, 2, . . .
τ
%
τ
τ +T1
(7.4)
where τ is an arbitrary time.
Although theoretically, a Fourier series representation of p(t) may require an infinite number of terms, in actual practice p(t) can generally be approximated with
186
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
sufficient accuracy by a relatively small number of terms. Example 7.1 illustrates the
Fourier series representation of a square wave.
Example 7.1 (a) Determine expressions for the coefficients of a real Fourier series
representation of the square wave shown in Fig. 1. Write the Fourier series representation of p(t). (b) Plot truncated series employing, respectively one, two, and three terms
of the Fourier series.
p(t )
p0
t
p0
T1
Figure 1 Square wave.
SOLUTION (a) The integrals in Eqs. 7.4 can be evaluated over the period −T1 /2 <
t < T1 /2 and thus written in the form
%
1 T1 /2
p(t) dt
(1)
a0 =
T1 −T1 /2
%
2 T1 /2
an =
p(t) cos n.1 t dt
(2)
T1 −T1 /2
%
2 T1 /2
bn =
p(t) sin n.1 t dt
(3)
T1 −T1 /2
where
p(t) =
<
−p0
p0
− T1 /2 ≤ t < 0
0 ≤ t < T1 /2
(4)
Substituting Eqs. 4 into Eqs. 1 and 2, we get
a0 = an = 0
Ans. (a) (5)
This results from the fact that p(t) is an odd function of t [i.e., p(t) = −p(−t)],
whereas a0 and an are coefficients of even terms in the Fourier series. The coefficient
for the odd terms is
%
4 p0 T1 /2
sin n.1 t dt
(6)
bn =
T1 0
so
9T1 /2
9
4p0 −1
bn =
cos n.1 t 99
(7)
T1 n .1
0
7.1 Response to Periodic Excitation: Real Fourier Series
But .1 T1 = 2π, so
bn = −
2p0
(cos nπ − 1)
nπ
187
(8)
or
4p0
,
n = 1, 3, 5, . . .
nπ
The Fourier series representation of the square wave is thus
4p0 ! 1
sin n.1 t
p(t) =
π
n
bn =
Ans. (a) (9)
Ans. (a) (10)
n=1,3,...
(b) The plots in Fig. 2 show the contributions over one period of the first three
nonzero terms of the Fourier series representation of the square wave in Fig. 1.
p (t )
p0
t
(a)
T1
p0
(b)
(c)
(d)
t
t
t
Figure 2 Fourier-series representations of a square wave: (a) square wave; (b) sin(.1 t) term;
(c) two-term representation; (d ) three-term representation.
7.1.2
Steady-State Response to Periodic Excitation
Having determined the response of SDOF systems to harmonic excitation in Chapter 4,
and having determined how to represent a periodic function in terms of its harmonic
components, we can now determine the response of SDOF systems to periodic excitation. In Example 7.2, an undamped SDOF system is subjected to the square-wave
excitation of Example 7.1.
188
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
Example 7.2 The undamped SDOF system in Fig. 1 is subjected to a square-wave
excitation p(t) like the one in Example 7.1. Determine a Fourier series expression for
the steady-state response of the system if ωn = 6.1 .
u (t )
k
m
p (t )
Figure 1
Spring–mass oscillator.
SOLUTION From Eq. 4.9, the steady-state response of an undamped SDOF system
to cosine excitation p0 cos .t is
p0 /k
cos .t
1 − r2
u=
(1)
√
where r = ./ωn = . m/k.
From Example 7.1 we can write the square-wave excitation p(t) in the form
p(t) =
where
∞
!
Pn sin n.1 t
(2)
n=1

 4p0
nπ
Pn =

0
n = odd
(3)
n = even
Hence, the steady-state response can be written in the form
u(t) =
∞
!
Un sin n.1 t
(4)
n=1
From Eq. 1, the coefficients Un of the steady-state-response series have the form
Un =
Pn
k (1 − rn2 )
(5)
n.1
ωn
(6)
where
rn =
Finally, combining Eqs. 3, 4, and 5, we get the following Fourier series expression for
the steady-state response:
u(t) =
4p0
kπ
!
n=1,3,...
1
sin n.1 t
n [1 − (n/6)2 ]
since, for the present problem, ωn = 6.1 .
Ans. (7)
7.2 Response to Periodic Excitation: Complex Fourier Series
189
pPn
4p0
1
wn = 6Ω1
Ωn
Ω1
3Ω1
5Ω1
7Ω1
9Ω1
7Ω1
9Ω1
pkUn
4p0
1
Ω1
3Ω1
Ωn
5Ω1
Figure 7.2 Excitation and response spectra based on Example 7.2.
It is very convenient to visualize periodic functions in terms of their spectra, that
is, plots of the amplitude of each harmonic component versus frequency. The spectra of
the periodic excitation p(t) and resulting steady-state response u(t) from Examples 7.1
and 7.2 are plotted in Fig. 7.2. From Eqs. 3 and 5 of Example 7.2, the following
nondimensional forms are obtained:

 4
n = odd
Pn
= nπ
(7.5a)

p0
n
=
even
0

4

n = odd
Un
nπ[1
−
(n/6)2 ]
=
(7.5b)
p0 /k 
n = even
0
and
.n = n.1
(7.5c)
Note that for the particular frequency ratio ωn = 6.1 , some of the Fourier components of the excitation are at frequencies below resonance, whereas others are above
resonance. The response would have a very large Fourier component if n.1 were to
fall close to ωn for some value of n.
7.2 RESPONSE TO PERIODIC EXCITATION: COMPLEX
FOURIER SERIES
In Section 4.3, the complex frequency-response function H (.) was introduced as a
convenient means of representing the response of a viscous-damped system to harmonic
190
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
excitation. A complex Fourier series representation of periodic excitation and response
functions will also prove to be very useful.
7.2.1
Complex Fourier Series
Let the periodic function p(t) be separated into its harmonic components by means of
a complex Fourier-series expansion, given by1
p(t) =
∞
!
P n (.)ei(n.1 t)
(7.6)
n=−∞
where the fundamental frequency .1 (rad/s) is related to the period of the function by
.1 T1 = 2π. The bar over Pn symbolizes the fact that the coefficients of the series may
be complex, even though the series as a whole represents a real function of time. To
evaluate these components, note that
:
% τ +T1
0
n &= m
i(n .1 t) −i(m .1 t)
e
e
dt =
(7.7)
T1
n=m
τ
So, multiplying Eq. 7.6 by e−i(m .1 t) and integrating over one period, we get
Pn =
1
T1
%
τ +T1
p(t)e−i(n.1 t) dt,
τ
n = 0, ±1, . . .
(7.8)
Note that
∗
P −n = P n = complex conjugate of P n
(7.9)
and that
1
P0 =
T1
%
τ +T1
τ
p(t) dt = average value of p(t)
(7.10)
[Actually, P 0 is real-valued since p(t) is real-valued.]
Example 7.3 (a) Show that if p(t) is a real-valued function, the right-hand side of
Eq. 7.6 will turn out to be real-valued, as it should be. (b) Show that if p(t) is an odd
function, P n (.) is purely imaginary and the coefficients P n and P −n are related by the
equation
P −n = −P n
SOLUTION (a) From Euler’s formula (Eq. 3.14),
e±iθ = cos θ ± i sin θ
1 By
and e0 = 1
(1)
including negative n terms in this series, the Fourier series can represent a real function of time, p(t).
7.2 Response to Periodic Excitation: Complex Fourier Series
191
Therefore, Eq. 7.6 can be expanded into the following form:
∞
!
p(t) = P 0 +
+
n=1
∞
!
P n (cos n.1 t + i sin n.1 t)
(2)
P −n (cos n.1 t − i sin n.1 t)
n=1
Since P n is assumed to be complex, it can be expressed in terms of its real and imaginary
components as
P n = -(P n ) + i .(P n )
(3)
Then, from Eq. 7.9,
∗
P −n = P n = -(P n ) − i .(P n )
(4)
Combining Eqs. 2 through 4, we get
p(t) = P 0 + 2
∞
!
n=1
[-(P n ) cos n.1 t − .(P n ) sin n.1 t]
which is real.
Q.E.D.
(b) Making use of Euler’s formula, Eq. 1, we can write Eq. 7.8 as
% τ +T1
1
Pn =
p(t) (cos n.1 t − i sin n.1 t) dt
T1 τ
Since p(t) is said to be an odd function of t, the cosine term drops out and
%
−i τ +T1
P +n =
p(t) sin n.1 t dt
T1 τ
(5)
(6)
(7)
Then
P −n
−i
=
T1
i
=
T1
%
%
τ +T1
p(t) sin(−n.1 t) dt
τ
τ +T1
(8)
p(t) sin n.1 t dt
τ
Hence, from Eq. 7, P n is purely imaginary, and from Eqs. 7 and 8, P −n = −P n .
Q.E.D.
By comparing Eq. 5 of Example 7.3 with Eq. 7.2, you will notice that the coefficients of the real Fourier series and the coefficients of the complex Fourier series are
related by the following equations:
ao = P 0 ,
an = 2-(P n )
bn = −2.(P n )
(7.11)
192
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
Example 7.4 (a) Determine an expression for the Fourier coefficients P n of the complex Fourier series representation for the square wave of Example 7.1. (b) Sketch
spectra of -(P n ), .(P n ), and |P n |.
SOLUTION (a) Equation 7.8 can be evaluated over the period 0 < t < T1 , giving
%
%
1 T1 /2
1 T1
Pn =
(p0 )e−i(n.1 t) dt +
(−p0 )e−i(n.1 t) dt
(1)
T1 0
T1 T1 /2
9 '
−p0 & −i(n.1 t) 99T1 /2
−i(n.1 t) 9T1
e
(2)
Pn =
−
e
0
T1 /2
in.1 T1
But .1 T1 = 2π, so
<
+1
n = even
e−i(n.1 T1 /2) = e−i(nπ) =
−1
n = odd
(3)
e−i(n.1 T1 ) = e−i(2nπ) = 1
Therefore, the complex Fourier coefficients are given by
ip0
[2e−i(nπ) − 1 − e−i(2nπ) ]
Pn =
2πn
or

n = even
0
ip0
−i(nπ)
(2e
Pn =
− 2) = −2ip0

2nπ
n = odd
nπ
(b) See Fig. 1.
(4)
Ans. (1) (5)
(Pn)
n
−5
2p0
5p
−5
−3
−1
1
2p0
p
(Pn)
3
5
3
5
2p0
3p
1
−3
−1
−2p0
3p
n
−2p0
5p
−2p0
p
|Pn|
2p0
p
−5
Figure 1
−3
−1
1
2p0
3p
2p0
5p
3
5
n
Various spectra for the square wave of Example 7.1.
7.2 Response to Periodic Excitation: Complex Fourier Series
7.2.2
193
Complex Frequency Response
The topic of complex frequency response was introduced in Section 4.3. Here we extend
that discussion to cover the complex representation for steady-state response of an SDOF
system to periodic excitation. From Eqs. 4.30 and 4.33, the steady-state response can
be written in complex form as
u(t) = U (.)ei .t = H (.)p0 ei.t
(7.12)
where the frequency-response function H is given by2
H (.) ≡ H u/p (.) =
1/k
[1 − (./ωn ) ] + i [2ζ(./ωn )]
2
(7.13)
where . is the forcing frequency and ωn is the undamped natural frequency of the
SDOF system.
When the excitation is periodic, we can use the complex Fourier series representation of Eq. 7.6, repeated here:
p(t) =
∞
!
P n ei(n.1 t)
(7.6)
n=−∞
Then the steady-state response can be written as
u(t) =
∞
!
U n ei(n.1 t)
(7.14)
n=−∞
Noting from Eq. 7.12 that for harmonic excitation U = H p0 , we see that the corresponding expression for periodic response is
U n = H n P n = |H n ||P n | ei(αHn +αPn )
(7.15)
where for a viscous-damped SDOF system,3
H n ≡ H (n.1 ) =
1/k
[1 − (n .1 /ωn ) ] + i [2ζ(n .1 /ωn )]
2
(7.16)
Example 7.5 illustrates the use of complex Fourier series in determining the steadystate response of a SDOF system subjected to periodic excitation. The method would
be even more beneficial if the system were a damped system.
Example 7.5 (a) Repeat Example 7.2 by determining an expression for the Fourier
coefficients U n for the undamped SDOF system subjected to square-wave excitation
with ωn = 6.1 . (b) Sketch magnitude and phase spectra, |U n | and αUn .
2
Here we use the dimensional form of H rather than the nondimensional form given in Eq. 4.33.
The subscript n of the undamped natural frequency ωn should not be confused with the Fourier series
index n.
3
194
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
SOLUTION (a) From Eq. 7.15,
U n = H n P n = |H n ||P n | ei(αHn +αPn )
(1)
and from Eq. 7.16, with ζ = 0,
Hn =
1/k
1/k
=
1 − (n.1 /ωn )2
1 − (n/6)2
(2)
From Eq. 5 of Example 7.4, the Fourier coefficients of the square wave are

0
n = even
P n = −2ip0

n = odd
nπ
Hence, the nth term of the frequency response is


0
U n = H nP n =
−i(2p0 )


nkπ[1 − (n/6)2 ]
n = even
(3)
Ans. (a) (4)
n = odd
(b) From Eq. 4 we can evaluate expressions for the magnitude and phase angle.
|U n | =
2p0 /kπ
|n [1 − (n/6)2 ]|
(5)
Because the nonzero terms in Eq. 4 are pure imaginary, the phase angles will all be
±π/2. Hence,
:
−π/2
n = +1, +3, +5, −7, −9, . . .
αUn =
(6)
π/2
n = −1, −3, −5, +7, +9, . . .
Note that sign changes occur as the resonance frequency is passed at ωn = 6.1 , that
is, at n = 6.
mn =
pk |Un|
zp0
1
n
−9
−7
−5
−3
−1
1
3
5
7
9
an
p/2
n
−wn = −6Ω1
−p/2
wn = 6Ω1
Figure 1 Response spectrum (magnitude and phase).
7.3
Response to Nonperiodic Excitation: Fourier Integral
195
For sketching purposes, evaluate the nondimensionalized magnitude
µn ≡
(π/2) |U n |
|1/n |
=
,
p0 /k
|[1 − (n/6)2 ]|
n = ±1, ±3, . . .
(7)
The values to be plotted are
36
= 1.029,
35
36
=
= 0.396,
91
4
= 0.444,
9
4
=
= 0.089
45
µ1 = µ−1 =
µ3 = µ−3 =
µ7 = µ−7
µ9 = µ−9
µ5 = µ−5 =
36
= 0.655
55
The results are shown in Fig. 1.
Compare the sketch of µn in Example 7.5 with Fig. 7.2b and note that the complex
coefficients µn in Example 7.5 have an amplitude that is half that of the corresponding
real coefficients, which are plotted in Fig. 7.2b. The contribution of the −n terms in
the complex Fourier series accounts for this difference.
It is convenient to think of problems such as Example 7.5 in terms of transformations from the time domain to the frequency domain (spectrum), and from the frequency
domain to the time domain. Figure 7.3 illustrates this.
7.3 RESPONSE TO NONPERIODIC EXCITATION:
FOURIER INTEGRAL
In previous sections you have seen that a periodic function can be represented by a
Fourier series, as in Eq. 7.2 or 7.6. When the function to be represented is not periodic,
it can be represented by a Fourier integral. In developing expressions for the Fourier
integral transform pair, it will be convenient to employ the complex Fourier series,
letting the period T1 approach infinity.
7.3.1
Fourier Integral Transforms
Equation 7.6, which defines the complex Fourier series, is repeated here:
p(t) =
∞
!
P n ei(n.1 t)
Time
domain
p (t )
(7.6)
n=−∞
Frequency
domain
Eq. 7.8
Pn (Ω)
Eq. 7.15
u (t )
Figure 7.3
Eq. 7.14
Un (Ω)
Solution of a periodic response problem by transformation to the frequency domain.
196
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
where, as given by Eq. 7.8, the Fourier coefficients P n are related to p(t) by
% τ +T1
1
Pn =
p(t)e−i(n .1 t) dt,
n = 0, ±1, . . .
(7.8)
T1 τ
provided that the integral exists.
In letting T1 → ∞, it will be convenient to introduce the following notation4 :
.1 = /.,
.n = n .1
P (.n ) = T1 P n =
2π
Pn
/.
(7.17)
(7.18)
Then Eq. 7.6 can be written in the form
p(t) =
∞
1 !
P (.n )ei(.n t) /.
2π n=−∞
(7.19)
where, from Eqs. 7.8 and Eq. 7.18,
P (.n ) =
%
T1 /2
p(t)e−i(.n t) dt
(7.20)
−T1 /2
The limits of integration on Eq. 7.20 have been taken as shown so that when T1 → ∞,
the entire time history of p(t) will be included regardless of the specific form of p(t).
As T1 → ∞, .n becomes the continuous variable ., and /. becomes the differential d.. Then Eqs. 7.20 and 7.19, respectively, can be written as
% ∞
P (.) =
p(t)e−i.t dt
(7.21)
−∞
1
p(t) =
2π
%
∞
P (.)ei.t d.
(7.22)
−∞
Equations 7.21 and 7.22 are called a Fourier transform pair. P (.) is known as the
Fourier transform of p(t); and p(t) is called the inverse Fourier transform of P (.).
The representation of p(t) by its Fourier transform requires the existence of the integral
in Eq. 7.21. Conditions that must be satisfied for the Fourier transform integral to exist
are discussed in texts on integral transforms (e.g., Ref. [7.1]). These conditions are met
by most physically realizable functions representing forces, displacements, and so on.
Finally, the Fourier transform pair can be written in a more symmetric form if
written in terms of the frequency f = ./2π. Then
P (f ) ≡ F [p(t)] =
−1
%
p(t) ≡ F [P (f )] =
4 Note
that Eq. 7.18 introduces a scaling factor, T1 .
∞
p(t)e−i(2πf t) dt
(7.23)
−∞
%
∞
−∞
P (f )ei(2πf t) df
(7.24)
7.3
Response to Nonperiodic Excitation: Fourier Integral
197
Example 7.6 illustrates the use of straightforward time-domain integration to determine the Fourier transform of a rectangular pulse that is symmetric about t = 0. Note
that the resulting Fourier transform is a real function of frequency.
Example 7.6 Let p(t) be the rectangular pulse defined by

t < −T
0


−T ≤ t ≤ T
p(t) = p0


t >T
0
(a) Determine the Fourier transform of this rectangular pulse. Express the transform as
a function of the frequency variable .. (b) Plot the Fourier transform.
SOLUTION (a) From Eq. 7.21,
% ∞
%
−i.t
P (.) =
p(t)e
dt =
−∞
Therefore,
P (.) =
T
p0 e−i.t dt
(1)
−T
p0 −i.T
− ei.T )
(e
−i .
(2)
which can also be written in terms of the sinc function (sin θ )/θ .
P (.) = 2 p0 T
sin .T
.T
Ans. (a) (3)
(b) The Fourier transform P (.) is therefore a real function. It can be plotted
versus the frequency variable . (Fig. 1) and compared with the discrete Fourier series
of Example 7.4.
P (Ω)
2p0T
Ω
0
−4p/T
−2p/T
2p/T
4p/T
Figure 1 Fourier transform of a symmetric rectangular pulse.
Reference [7.1] contains a table of a number of Fourier transform pairs.
198
7.3.2
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
Frequency-Response Functions
In Eq. 7.14 we found that the response of an SDOF system to periodic excitation can
be expressed in the form
u(t) =
∞
!
U n ei(n.1 t)
(7.14)
n=−∞
where, from Eq. 7.15,
U n = H nP n
(7.15)
Following the procedure of Eqs. 7.17 through 7.24, we obtain the following Fourier
transform pair for the response of an SDOF system:
U (f ) =
u(t) =
where
%
∞
u(t)e−i(2πf t) dt
(7.25)
U (f )ei(2πf t) df
(7.26)
−∞
%
∞
−∞
U (f ) = H (f )P (f )
(7.27)
which is the product of the system frequency-response function, H (f ), and the Fourier
transform of the excitation, P (f ). Therefore, the response can be expressed by the
following inverse Fourier transform:
% ∞
H (f )P (f )ei(2πf t) df
(7.28)
u(t) =
−∞
In some cases a table of Fourier transform pairs can be used to evaluate this inverse
transform.[7.1,7.2] However, evaluation of this definite integral generally involves contour
integration in the complex plane, which is beyond the scope of this book.
7.3.3
Parameter Identification
Of great importance (as demonstrated in Chapter 18) is the fact that Eq. 7.27 can be
written symbolically in the form
H (f ) =
U (f )
P (f )
(7.29)
The system frequency-response function, H (f ), can thus be obtained from Fourier
transforms of the measured time histories of excitation p(t) and response u(t). Important
system parameters (e.g., the undamped natural frequency ωn and damping factor ζ of
the system) can then be extracted from this system frequency-response function. For
example, for a viscous-damped SDOF system, H (f ) is given by (Eq. 7.13)
1/k
H (f ) =
(7.30)
[1 − (f/fn )2 ] + i(2ζf/fn )
Such parameter identification procedures are discussed in Chapter 18 and in modal
analysis references such as Refs. [7.4] and [7.5].
7.4 Relationship Between Complex Frequency Response and Unit Impulse Response
199
In Section 7.5 we describe how the Fourier integral can be approximated by the
discrete Fourier transform (DFT), which, in turn, can be evaluated numerically by use
of the fast Fourier transform (FFT) algorithm.
7.4 RELATIONSHIP BETWEEN COMPLEX FREQUENCY
RESPONSE AND UNIT IMPULSE RESPONSE
The complex frequency response, H (.) or H (f ), of a linear system describes its
response characteristics in the frequency domain, and the unit impulse response h(t)
describes the system’s response in the time domain. We now show that the unit impulse
function and the system frequency-response function form a Fourier transform pair.
From Eq. 7.23, the Fourier transform of the unit impulse excitation is
% ∞
P (f ) =
p(t)e−i(2πf t) dt = 1
(7.31)
−∞
Therefore, the response to this unit impulse is, from Eq. 7.28,
h(t) =
%
∞
H (f )ei(2πf t) df
(7.32)
−∞
But from Eq. 7.24, the equation that defines the inverse Fourier transform, this unit
impulse response function is just the expression for the inverse Fourier transform of
the system frequency-response function H (f ). Conversely, it follows that the system
frequency-response function H (f ) is the Fourier transform of the unit impulse response
function h(t), that is,
H (f ) ≡ F [h(t)] =
%
∞
h(t)e−i(2πf t) dt
(7.33)
−∞
This Fourier transform pair relationship between a system’s impulse-response function
in the time domain and its frequency-response function in the frequency domain is
illustrated in Fig. 7.4. For example, for a viscous-damped SDOF system, H (f ) is
given by Eq. 7.30 and h(t) by Eq. 5.30, both repeated here.
H (f ) =
h(t) =
Time Domain
h(t)
Impulse Response
<
1/k
[1 − (f/fn )2 ] + i(2ζf/fn )
(7.30)
1 −ζ ωn t
e
sin ωd t,
m ωd
(5.30)
F −→
←− F −1
0
t >0
Frequency Domain
H (f )
Frequency-Response Function
Figure 7.4 Relationship between impulse response and frequency response.
200
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
It is left as an exercise for the reader to show that these SDOF system functions satisfy
Eq. 7.33.
In Chapter 18 it will be noted that some parameter-estimation algorithms work
directly with the frequency-response function H (f ); these are called frequency-domain
algorithms. Other parameter-estimation algorithms make use of the impulse-response
function h(t); these are called time-domain algorithms.
7.5 DISCRETE FOURIER TRANSFORM AND FAST
FOURIER TRANSFORM
Although the Fourier integral techniques discussed in Section 7.3 provide a means for
determining the transient response of a system, numerical implementation of the Fourier
integral became a practical reality only with the publication of the Cooley–Tukey
algorithm for the fast Fourier transform in 1965.[7.3] Since that date, the FFT has virtually
led to a revolution in many areas of technology, including the area of vibration testing
(e.g., see Chapter 18 and Refs. [7.4] and [7.5]).
Two steps are involved in the numerical evaluation of Fourier transforms.[7.1] First,
discrete Fourier transforms (DFTs), which correspond to Eqs. 7.23 and 7.24, are derived.
Then an efficient numerical algorithm, the fast Fourier transform (FFT), is used to the
compute the DFTs.
7.5.1
Discrete Fourier Transform Pair
For numerical treatment of the Fourier transform, it is necessary first to define a discrete Fourier transform pair corresponding to the continuous Fourier transform pair of
Eqs. 7.23 and 7.24. To be transformed, a continuous function must first be sampled at
discrete time intervals /t. Second, due to computer memory and execution-time limitations, only a finite number N of these sampled values can be utilized. Figure 7.5a
illustrates a sampled waveform with N = 16 samples taken at /t = 0.25-sec intervals.
1.5
5
Mag(P )
p(t )
1
0.5
0
0
−0.5
−1
0
1
2
t (sec)
(a)
3
4
5
−5
0
4
8
n
12
16
(b)
Figure 7.5 Sampled waveform (e−t with /t = 0.25 sec, N = 16): (a) discrete-time representation; (b) discretefrequency representation.
7.5
Discrete Fourier Transform and Fast Fourier Transform
201
The effects of this sampling and truncation are to approximate the continuous signal by a periodic signal of period T1 = N /t, sampled at times tm = m /t, m =
0, 1, . . . , (N − 1).
The total sample time is T1 , so the fundamental-frequency sinusoid that fits within
this sample time has a period T1 . Therefore, the frequency interval of the discrete Fourier
transform is
/f =
1
1
=
T1
N /t
(7.34)
The DFT consists of N samples extending up to a maximum frequency of N /f .
Specific features of sampling and truncation that are illustrated in Fig. 7.5 are discussed following Example 7.7. Sampling and truncation are discussed in greater detail
in Section 18.4.
Since a period T1 consists of N samples, the integral of Eq. 7.23 is replaced by the
following finite sum:
P (fn ) =
N−1
!
p(tm )e−i2π(m/t)(n/f ) /t
(7.35)
m=0
Finally, the discrete Fourier transform (DFT) can be written
P (fn ) = /t
N−1
!
p(tm )e−i(2πm(n/N)) ,
m=0
n = 0, 1, . . . , N − 1
(7.36)
The inverse DFT can be obtained from Eq. 7.24 in a similar manner. Thus,
p(tm ) =
N−1
!
P (fn )ei2π(m/t)(n/f ) /f
(7.37)
n=0
This expression for the inverse Fourier transform (IFT) can be written as
p(tm ) =
N −1
1 !
P (fn )ei(2πm(n/N)) ,
N /t
n=0
m = 0, 1, . . . , N − 1
(7.38)
Equations 7.36 and 7.38 define a discrete Fourier transform pair that is consistent with
the continuous Fourier transform.
The discrete Fourier transform approximates the continuous Fourier transform at
discrete frequencies fn . The accuracy of a DFT representation depends on the sampling
interval /t and the number of samples, or block size, N . In Section 18.4 we discuss
these effects in much greater detail.
Reference [7.1] presents a graphical derivation and a theoretical derivation of the
DFT pair. The resulting discrete Fourier transform (DFT), written in present notation, is
P (fn ) =
N−1
!
m=0
p(tm )e−i(2πm(n/N)) ,
n = 0, 1, . . . , N − 1
(7.39)
202
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
where fn = n /f . The corresponding inverse DFT is
p(tm ) =
N−1
1 !
P (fn )ei(2πm(n/N)) ,
N
n=0
m = 0, 1, . . . , N − 1
(7.40)
This form is the one that is utilized in computing FFTs.[7.6] However, a scale factor
/t is required to produce an equivalence between this DFT form and the continuous
Fourier transform.
7.5.2
Fast Fourier Transform Algorithm
The fast Fourier transform (FFT) is not a new type of transform but rather, an efficient
numerical algorithm for evaluating the DFT. Its importance lies in the fact that by
eliminating most of the repetition in the calculation of a DFT, it permits much more
rapid computation of the DFT.
Either Eq. 7.39 or 7.40 can be cast in the form
Am =
where
N−1
!
Bn WNmn ,
n=0
m = 0, 1, . . . , N − 1
WN = e−i(2π/N )
(7.41)
(7.42)
A measure of the amount of computation involved in Eq. 7.41 is the number of complex
products implied by the form of the equation and the range of m. It is clear that there are
N sums, each of which requires N complex products, or there are N 2 products required
for computing all of the Am ’s. By taking advantage of the cyclical nature of powers
of WN , the total computational effort can be drastically reduced. Figure 7.6 shows the
repetition cycle for W8mn . The number of complex products for the FFT algorithm is
given by (N/2) log2 N . For example, if N = 512, the number of FFT operations is less
than 1% of the corresponding number of DFT operations. Signal-processing software
products invariably provide for computation of the FFT.[7.6]
W 68 = W 814 = …
W 78 = W 815 = …
W 58 = W 813 = …
W 08 = W 88 = …
W 48 = W 812 = …
W 18 = W 89 = …
W 38 = W 811 = …
W 28 = W 810 = …
Figure 7.6 Cyclical nature of WNmn for N = 8.
7.5
203
FFT Computations
In Example 7.7 the FFT command in the Matlab computer program is used to compute the Fourier transform of a square wave, similar to the one treated analytically in
Example 7.4 in Section 7.2. Computation of FFTs of nonperiodic signals is discussed
in Section 18.4.
Example 7.7 Let p(t) be the square wave defined over one period by the function
:
2
0.0 ≤ t < 0.5 sec
p(t) =
0
0.5 sec ≤ t < 1.0 sec
represented by 16 samples over the finite interval 0 sec ≤ t < 1 sec. Use the FFT
command in Matlab to determine the Fourier transform of this square wave, and plot
the real and imaginary parts of the resulting Fourier transform.
SOLUTION Relative to an offset average value of 1.0, the square wave is antisymmetric in the time window 0.5 sec ≤ t ≤ 0.5 sec. The 16 discrete-time sample points
are shown in Fig. 1a. These 16 data points constitute the input to the FFT algorithm.
(See the Comments following this example.)
1.5
(P )
1
0.5
0
2
−0.5
0
2
4
6
8
6
8
f (Hz)
1
0
1
0
0.5
t (sec)
(a )
0.5
1
(P )
p (t )
7.5.3
Discrete Fourier Transform and Fast Fourier Transform
0
−0.5
−1
0
2
4
f (Hz)
(b)
Figure 1 (a) Time- and (b) frequency-domain representations of a square wave.
204
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
In Fig. 1b the results of a 16-point FFT are plotted versus the discrete frequency
fn for 0 ≤ fn ≤ fN/2 . The real part (top right) and the imaginary part(lower right) can
be compared with the plots in Fig. 1 of Example 7.4, noting that the square wave in
Example 7.4 is an antisymmetric function with an average value of zero, whereas the
square wave in this example is an antisymmetric function relative to an average value
of 1.0. Other differences are explained in the Comments that follow.
Comments Regarding Matlab FFT Input/Output
1. The FFT algorithm is executed with the command
y = fft (x, N )
where N is the block size and x is the input sequence, that is, the vector of
discrete-time samples of the input function. For structural dynamics applications,
x will be a vector of real numbers. It is desirable to let N be a power of 2.
2. The FFT treats the first N numbers in the input sequence as one period of a
periodic function. For example, the rectangular pulse in Fig. 1a of Example 7.7
will correspond to a square wave that is represented by a Fourier series whose
fundamental frequency is 1.0 sec, with the sampled value at t = 1 sec repeating the sampled value at t = 0 sec. This topic is treated in greater detail in
Section 18.4.
3. Where there is a step discontinuity in the input function, the sampled value is
taken as the average of the value prior to the jump and the value after the jump,
as shown in Fig. 1a of Example 7.7 and in Fig. 7.5a.
4. Because P (fn ) and P (f−n ) are complex conjugates of each other, the second
half of the Fourier coefficients are usually not plotted (e.g., as in Fig. 1b of
Example 7.7). Thus, the complex values in the output vector are sequenced by
Matlab as follows:
y(1) = P (f0 ) = average value of input (real)
y(2) = P (f+1 )
..
.
y(N/2 + 1) = P (f+N/2 )
y(N/2 + 2) = P (f−N/2+1 )
..
.
y(N ) = P (f−1 )
5. To treat the FFT output values as approximations to the complex Fourier coefficients given by Eq. 7.8, it is necessary to divide the output values y(n) by the
block size N .5
5
This scaling is necessitated by the scaling factor T1 introduced in Eq. 7.18, together with the additional
scaling factor of /t between Eqs. 7.36 and 7.39.
Problems
205
6. To convert the input sequence number m to input sample time tm in seconds,
it is necessary to use tm = m/t. Similarly, to convert the frequency scale of
the output from coefficient number n to frequency fn in hertz, one must use
fn = n/f = n/N /t.
REFERENCES
[7.1] E. O. Brigham, The Fast Fourier Transform, Prentice-Hall, Englewood Cliffs, NJ, 1974.
[7.2] R. A. Gabel and R. A. Roberts, Signals and Systems, 3rd ed., Wiley, New York, 1987.
[7.3] J. W. Cooley and J. W. Tukey, “An Algorithm for Machine Calculation of Complex
Fourier Series,” Math Computation, Vol. 19, 1965, pp. 297–301.
[7.4] D. J. Ewins, Modal Testing: Theory, Practice and Application, 2nd ed., Research Studies
Press, Baldock, Hertfordshire, England, 2000.
[7.5] N. M. M. Maia and J. M. M. Silva, ed., Theoretical and Experimental Modal Analysis,
Wiley, New York, 1997.
[7.6] Using MATLAB, Version 6, The MathWorks, Natick, MA, 2002, pp. 12-41 to 12-48.
PROBLEMS
p (t )
For problems whose number is preceded by a C,
you are to write a computer program and use it to
produce the plot(s) requested. Note: Matlab .m-files
for many of the plots in this book may be found on
the book’s website.
Half-sine waves
p0
t
T1/2
T1
Figure P7.2
Problem Set 7.1
7.1 Determine the real Fourier series for the periodic
excitation function in Fig. P7.1.
p (t )
p0
t
−T1/4
3T1/4
−p0
C 7.3 (a) Determine the real Fourier series for the
square wave shown in Fig. P7.3. (b) By comparing
your result from part (a) with the answer in Eq. 10 of
Example 7.1, what do you observe to be the effect(s) of
the phase shift of the square wave? (c) Run the Matlab .m-file sd2hw7 3.m, and print out the plots produced. Modify this .m-file to produce the input signal of
Example 7.1, and print out the plots produced.
p (t )
p0
T1/4
5T1/4
3T1/4
Figure P7.1
−p0
7.2 Determine the real Fourier series for the periodic
excitation function in Fig. P7.2.
Figure P7.3
t
206
Response of SDOF Systems to Periodic Excitation: Frequency-Domain Analysis
7.4 (a) Determine the real Fourier series for the periodic
excitation function in Fig. P7.4. (b) By comparing your
result from part (a) with the answer in Eq. 10 of Example
7.1, what do you observe to be the effect(s) of the upward
shift of the square wave?
7.8 (a) Determine the coefficients of the complex
Fourier series that corresponds to the excitation p(t)
in Fig. P7.3. (b) Sketch the real and imaginary parts
and the magnitude of the complex coefficients as in
Example 7.4.
7.9 (1) Determine the coefficients of the complex
Fourier series that corresponds to the excitation p(t)
in Fig. P7.4. (2) Sketch the real and imaginary parts
and the magnitude of the complex coefficients as in
Example 7.4.
p(t)
2p0
t
T0/2
T0/2
7.10 The undamped SDOF system of Fig. P7.10a is
subjected to the sawtooth base motion z(t) illustrated
in Fig. 7.10b. Let
z(t) =
Figure P7.4
7.5 The undamped SDOF system shown in Fig. P7.5
is subjected to the excitation p(t) given in Problem 7.1.
(a) Determine a Fourier series expression for the steadystate response of the system if ωn = 4.1 . (b) Sketch the
spectra for Pn and for Un , similar to those shown in
Fig. 7.2.
∞
!
Z n ei(n.1 t) ,
n=−∞
m
∞
!
W n ei(n.1 t)
n=−∞
where W n = H n Z n . (a) Determine H n . (b) Determine
Z n . (These will have the same form as P n of Problem
7.1.) (c) Determine W n , and sketch |W n |/Z versus
frequency order, as in Example 7.4.
z (t )
u (t )
k
k
w(t) =
w=u–z
m
p (t )
( a)
Figure P7.5
z (t )
7.6 The undamped SDOF system shown in Fig. P7.5
is subjected to the excitation p(t) given in Fig. P7.4.
(a) Determine a Fourier series expression for the steadystate response of the system if ωn = 6.1 . (b) Sketch
the spectra for Pn and for Un similar to those shown
in Fig. 7.2.
Z
t
−Z
T1
(b)
Figure P7.10
Problem Set 7.2
7.7 (a) Determine the coefficients of the complex
Fourier series that corresponds to the excitation p(t)
given in Fig. P7.1. (b) Sketch the real and imaginary
parts and the magnitude of the complex coefficients as
in Example 7.4.
Problem Set 7.3
7.11 Determine the continuous Fourier transform for the
rectangular pulse in Fig. P7.11. (Note that this involves
a time shifting of the pulse in Example 7.6.)
Problems
p (t )
p0
t
2T
Figure P7.11
7.12 Determine the real and imaginary parts of the
continuous Fourier transform of each of the following
functions:
Note: Absolute value of time t.
(a) p(t) = p0 e−α|t|

0,
t <0





1
0≤t ≤
(b) p(t) = p0 sin(2πf0 t),
2f
0



1

 0,
t>
2f0
207
compute a 16-point Fourier transform to verify the (real
and imaginary) results shown in Figs. 1b. (Note: You can
find the Matlab .m-file for Example 7.7 as sd2ex7 7.m
on the book’s website.) (b) What is the resulting value
of the coefficient (P 1 )16pt ? Compare this answer with
the value given in Eq. 5 of Example 7.4. (c) Repeat
part (a) to calculate a 32-point transform and a 64point transform. Are the values of (P 1 )16pt , (P 1 )32pt ,
and (P 1 )64pt converging to the value given in Eq. 5 of
Example 7.4?
C 7.15 Let p(t) be the periodic function defined by
the 16-point discrete-time sawtooth wave shown in
Fig. P7.15. (a) Use the FFT command in Matlab to
determine the 16-point discrete Fourier transform (DFT)
of this sawtooth wave. (b) Plot the real and imaginary
parts of this DFT transform.
1
0.5
7.13 Using the unit impulse response function h(t) of
Eq. 5.30 and the definition in Eq. 7.33 of the continuous
Fourier transform H (f ), show that the Fourier transform
of the unit impulse response function is the frequency
response function H (f ) given by Eq. 7.30. That is,
using Eq. 7.33, show that h(t) and H (f ) form a Fourier
transform pair, as illustrated in Fig. 7.4.
p (t )
Problem Set 7.4
0
−0.5
−1
0
0.25
0.5
Time
Problem Set 7.5
C 7.14 (a) Using the input sequence given in Fig. 1a
of Example 7.7, use the Matlab FFT command to
Figure P7.15
0.75
1