10 FATIGUE
ANALYSIS
FREQU ENCY
IN T H E
DOMA! N
Yu NG-L! LEE DAIMLERCHRYSLER
10. 1
INTRODUCTION
Fatigue analysis is often performed in the time domain, in which all input loading and output stress or strain response are time-based signals. In this case, the response time history can be calculated either by the static stress analysis (or inertia relief method) (Lee et al., 1995) or by the modal transient analysis (Crescimanno and Cavallo, 1999; Vellaichamy, 2002). In the static stress analysis, stress can be obtained by superimposing all stress influences from the applied loads at every time step. In the modal transient analysis, stress can be calculated by using the modal stress and the modal coordinates. The modal transient analysis is recommended if resonant fatigue is of primary concern; that is, the loading frequencies appear closer to the system natural frequencies. In some situations, however, response stress and input loading are preferably expressed as frequency-based signals, usually in the form of a power spectral density (PSD) plot. In this case, a system function (a characteristic of the structural system) is required to relate an input PSD of loading to the output PSD of response. The PSD techniques for dynamic analysis have been widely used in the offshore engineering field. PSD represents the energy of the time signal at different frequencies and is another way of denoting the loading signal in time domain. The Fast Fourier Transform (FFT) of a time signal can be used to obtain the PSD of the loading, whereas the Inverse Fourier Transform (IFT) can be used to transform the frequency-based signal to the time-based loading. The transform of loading history between 369
370
F A T I G U E A N A L Y S I S IN T H E F R E Q U E N C Y D O M A I N
the time domain and frequency domain is subject to certain requirements, as per which the signal must be stationary, random, and Gaussian (normal). This chapter has two objectives. First, the fundamental theories of random vibration, which can be found elsewhere (Newland, 1984; Bendat and Piersol, 1986; Wirsching et al., 1995), are summarized in Sections 10.1-10.5 for the benefit of readers. Second, in Sections 10.6 and 10.7, the methods of predicting fatigue damage from an output PSD of stress response are introduced.
10.2
A RANDOM SAMPLE TIME HISTORY
A system produces a certain response under excitation. If the excitation or the response motion, X(t), is unpredictable, the system is in random vibration because the exact value of X(t) cannot be precisely predicted in advance. It can only be described probabilistically. The probability density function (PDF) of a time history X(t) can be obtained by calculating its statistical properties. Figure 10.1 shows a sample time history for a random process, x(t), during the time interval 7", where X(t) exits between the values of x and x + d x for a total time of (dtl +dt2 +dt3 +dt4). The probability that x<-X(t)<-x + d x is therefore given by
P[x<_X(t)<_x + dx] -
dtl + dt2 + dt3 + dt4
(lO.2.1)
If the duration T is long enough, the PDFJx(x) is given by
x(o
x+dx
.
.
.
.
.
.
.
.
.
.
.
.
_
_
X
A
A
//y , /VD H
FIGURE 10. 1 Probabilitydensity function (PDF) for a random process X(t).
37 1
A RANDOM S A M P L E T I M E HISTORY
k
E clt~ fx(x) = P[x <- x(t) -< x + dx] - i=l T
(10.2.2)
Equation 10.2.2 is correct if the time duration T goes to infinity, and implies that the sample time history continues forever. Measurement of the time segments, ~ dti, for the P D F fx(x) in Equation 10.2.2 is very cumbersome.
Alternatively, the P D F can be determined by the fraction of the total number of samples in the band between x and x + dx. This can be done by digitizing the time history at a certain sampling rate in the time interval T. As shown in Figure 10.2, fx(x) is then given by
fx(x) - P[x <- X(t) <- x + dx] =
# sampleband 6 = # sample v 69
(10.2.3)
With a given P D F fx(x), some statistical properties of the random process X(t) can be obtained. The mean,/~x, and the variance, azx, of the process can be calculated as +~
~x =
l
T
1J
xfx(x)dx ~- -~ X(t)dt
- ')c +~c
(10.2.4)
0 T
a 2 = J [ x - l~xl2fx (x)dx ~ -~1 l [ X ( t ) - laxl2dt -~
(10.2.5)
0
x+dx X
t
T FIGURE
1 0.2
P D F for a digitized r a n d o m process X(t).
I
FATIGUE ANALYSISIN THE FREQUENCYDOMAIN
372
tx(X) 0.05 [ Px
O.O4
0.03
1
0.02
0.01
I -20
~~0 0 FIGURE
1 0.3
i
1
20
40
~
IX 60
P D F for a Gaussian process.
where T is the record length. When #x = O, ax is the root mean square (RMS) of the random process X(t). The RMS is a measure of the amplitude of the process. The process X ( t ) is called Gaussian if its PDF f x ( x ) follows the bell-shape distribution of Figure 10.3. The PDF is given by
f x ( x ) - ~-~--------~exp - ~
-~x #x
- e c < x < +oc
(10.2.6)
where /~x and a x are the mean and standard deviation of the process, respectively.
1 0.3
A RANDOM
PROCESS
A collection of an infinite number of sample time histories, such as Xl(t), X2(t), and X3(t) makes up the random process X(t) as shown in Figure 10.4. Each time history can be generated from a separate experiment. For example, the front axle shaft torque of a vehicle is of interest. It is necessary to learn how the torque varies during long-distance driving. A large number of time histories for various road and weather conditions must be recorded. In engineering, the ensemble of a sufficiently large number of sample time histories approximates the infinite collection representing a random process. Instead of being measured along a single sample, the ensemble statistical properties are determined across the ensemble, as shown in Figure 10.5. For a Gaussian random process, the ensemble probability density at each time instant and at any two time units must be Gaussian. A random process is said to be stationary if the probability distributions for the ensemble remain
A
RANDOM
~ _;, _u _~ L .t
--
1
373
PROCESS
LJFT Fz(Ibs) .
Front
5
1
5|
0
=..~LF:o=t~rj I[:'.., /
~
_,
r,,il~i
I
I
I
I
I
I
I
,, ..
I
-
I J
~" .
- - ji.~l, S 9 0
-m,-
~.,Lt
Front giFT Fz(Ibs)
.
,--m
0oo
.
.
-
"/
~ooo
I
I
WFT F z ( I b = )
9 s,e,
...do.t
/
~o~o I
Lt F r o n t
""
/
~(x,(,);
__i_.. 1 AA/-,i
p I
rand^m
~ o-
lJi,J
.
s~v.~
~ ./t U ) I
3ooo I I G o n d I
"(t)) ,,
.
-r
.
.
.
'
~-~m-m--l-r
i
.
~opo
.
.
.
" ~ ~ -
. . . . . .
-- .
/-~
rIndom~
.
.
.
.
.
5ooo
.
. *random3
,x,(o,' i i
"
""
I
_
./ _
I -
~ obo
~ ~ ~ ao
1 0.4
FIGURE
=~cl@
~o
Random process---ensemble of random sample time histories.
h
t2
,
I
,,(1(0 xdt)" k x3(t)
x,(,)/
LvW
L/
V J V "~
i
v
vtk/
i
v~v ~
I-W W FIGURE 1 0 . 5
~1
v~
~'~'
v ,,wvy-
lllustration of ensemble statistical data.
the same (stationary) for each time instant. This implies that the ensemble mean, standard deviation, variance, and mean square are time invariant. A stationary process is called ergodic if the statistical properties along any single sample are the same as the properties taken across the ensemble. That means each sample time history completely represents the ensemble. Note that if a random process is ergodic, it must be stationary. However, the converse is not true: a stationary process is not necessary ergodic.
374
FATIGUE
ANALYSIS
IN T H E F R E Q U E N C Y
DOMAIN
The autocorrelation function of a random process is the mean value of the product X(tl)X(t2), and is denoted by E[X(tl)X(t2)]. It can be approximated by the average value of the product X(tl)X(t2), which can be obtained by sampling the random variable x at times tl and t2. For a stationary random process, the value of E[X(tl)X(t2)] is time invariant. However, it depends only on the time difference r = ]t2 - tl 1. Therefore, the autocorrelation function of X ( t ) , denoted by R(r), is actually a function of r. It is expressed as (10.3.1)
R(~) = E[X(tl)X(t2)]
Also, because X ( t ) is stationary, its mean and standard deviation are independent of t. Thus, E[X(tl)] = E[X(t2)] = #x
(10.3.2)
ax(tl) = ox(t2) = ax
(10.3.3)
The correlation coefficient, p, for X ( t l ) and X(t2) is given by P-
R x ( r ) - #2x a2x
(10.3.4)
If p = +1, there is perfect correlation between X ( t l ) and X(t2), and if p = 0, there is no correlation. Because the value of p lies between - 1 and 1, it has /~r - a2 -< Rx(r)-< #2 + a~,
(10.3.5)
Figure 10.6 shows the properties of the autocorrelation function R(r) of a stationary random process X ( t ) . When the time interval approaches infinity, the random variables at times t~ and t2 are not correlated. That means R x ( r --* oc) ~- #x
and
p --, 0
(10.3.6)
Rx(0
/dx + O'x
f
r 2
2
/4x - O'x
FIGU RE 10. 6
Autocorrelation function Rx(r) of a stationary random process.
375
FOURIER ANALYSIS
When the time interval is zero,
Rx(v - O) - p2 + cr2 _ E[X21
(10.3.7)
and R(0) becomes the mean square value of the process. Also, for a stationary random process, R(T) is an even function, i.e.,
Rx( - r)
Rx(r) = 1 0.4
(10.3.8)
FOURIER ANALYSIS
Any periodic time history can be represented by the summation of a series of sinusoidal waves of various amplitude, frequency, and phase. If X(t) is a periodic function of time with a period T, X(t) can be expressed by an infinite trigonometric series of the following form:
X(t) - Ao + Z
Ak c o s - - f - + Bk sin
(10.4.1)
k=l
where
r/2 A o - ~ I X(t)dt -T/2 1
T/2
_2 Ak -- -~
J X(t) cos ~2r&t ~ dt -T/2 T/2
Bk -- ~2
I X(t) sin T2r&t dt -T/2
The Fourier series can be also expressed by using complex coefficients as O(3
X(t)- Z
Ckei2nkt/T
(10.4.2)
k=l
where the complex coefficients Ck are given by T/2
Ck --
X(t)e-i2rtkt/rdt
(10.4.3)
-T/2
The Fourier transform can be considered as the limit of the Fourier series of X(t) as T approaches infinity. This can be illustrated as follows by rewriting Equation 10.4.2 with infinite T:
376
FATIGUE ANALYSIS
X(t)-
1
lim ~ T---,oc k = - ~
IN T H E F R E Q U E N C Y
X(t)e-iZ'~t/rdt ei2nkt/T
DOMAIN
(10.4.4)
- T/2
If the frequency of the kth harmonic, wk, in radians per second, is 2~zk T and the spacing between adjacent periodic functions, Aw, is wk -
2rt
A ~ -- --T
(10.4.5)
(10.4.6)
E q u a t i o n 10.4.4 becomes
~A ~
X(t) - lim
X( t)e- ika~t d t
T--,cx~ k = -
e ika~
(10.4.7)
- T/2
As T goes to infinity, the frequency spacing, A~, becomes infinitesimally small, denoted by dw, and the sum becomes an integral. As a result, Equation 10.4.7 can be expressed by the well-known Fourier transform pair X(t) and X(~o): OO
l I
X(co) - ~
X(t)e-i~'dt
(10.4.8)
--OG OO 1
X(t) -
I X(c~
(10.4.9)
-- (X)
The function X(a~) is the forward Fourier transform of X(t), and X(t) is the inverse Fourier transform of X(aJ). The Fourier transform exists if the following conditions are met: 1. The integral of the absolute function exists, i.e., J'_~ IX(t)ldt < ~ . 2. Any discontinuities are finite.
1 O.5
SPECTRAL
DENSITY
The Fourier transform of a stationary random process X(t) usually does not exist because the condition O(3
I IX(t)ldt < c~ --00
(10.5.1)
377
SPECTRAL DENSITY
is not met. However, the Fourier transform of the autocorrelation function Rx(r) always exists. If the stationary random process X(t) is adjusted (or normalized) to a zero-mean value, i.e.,
Rx(r --, oc) - 0
(10.5.2)
i
(10.5.3)
the condition
lRx(r)ldt < ~
--OC
is met. In this case, the forward and inverse Fourier transforms of Rx(r) are given by
'J
7Y,=
Sx(~) - ~
Rx(r)-
i
(~0.5.4)
Rx(r)e-i~*dr
Sx(a2)ei";rda;
(lo.5.5)
where Sx(~) is the spectral density of the normalized stationary random process X(t). If r - 0, Equation 10.5.5 reduces to E [ X 2] -
Rx(O)- i Sx(a;)da;- a2
(10.5.6)
This means that the square root of the area under a spectral density plot Sx(~) is the RMS of a normalized stationary random process. Sx(~) is also called mean square spectral density and is illustrated in Figure 10.7.
sx(c0) / 9
v
0
m
FIGURE 1 0.7 Relationship between the spectral density and root mean square of a normalized stationary r a n d o m process.
378
FATIGUE ANALYSIS
IN T H E F R E Q U E N C Y
DOMAIN
The idea of a negative frequency has been introduced for mathematical completeness. However, it has no physical meaning. It is common practice to consider the frequency from zero to infinity and to have the frequency expressed in hertz (cycle/second) rather than radians/second. Therefore, the two-sided spectral density, Sx(co), can be transformed into an equivalent one-sided spectral density W x ( f ) as follows: E[xZ]- a2-
I Wx(f)df
(10.5.7)
0
where
W x ( f ) - 4rcSx(cO)
(10.5.8)
is the PSD and u)
f = 2n
(10.5.9)
The following spectral density relationships exist for first and second derivatives of a stationary random process X(t):
a22 =
S f ( ( w ) = 032 S x ( o 3 )
(1 0.5.10)
Wit(f) = (2n)2f 2 W x ( f )
(10.5.11)
S2(co)dco-oc
~2Sx(w)dco = (2r02 f2 W x ( f ) d f -oo
0
Sjd(w) = ~o4Sx(co)
(10.5.13)
Wjd(f) = (2r04f 4 W x ( f )
(10.5.14)
~} = [ Sji(co)dw = -c~
(10.5.12)
co4Sx(w)dw _ (2r04 f4 W x ( f ) d f -oo
(10.5.15)
0
A random process is called a narrow-band process if its spectral density has only a narrow band of frequencies. In contrast, a broad-band process is a process whose spectral density covers a broad band of frequencies. White noise is a broad-band process. Figure 10.8 shows examples of a narrow-band process, a wide-band process, and white noise. The PSD function is usually presented on a log-log scale. An octave (oct) is a doubling of frequency. The increase in octaves from fl to f2 is
1 oct = ~ l o g 7
IOgl0 2
f2 J!
(10.5.16)
37 9
SPECTRAL DENSITY
S t r a in(uE)
nBrrow
- I. 804
1
~
3
Strm in(uE)
~'~" ~ ~ ~ -2.53
u ide
~,, ,, ,, ,, ~ ~ ~
M,, ,, ,,, ~ ~ ~ . ~ ~ ,, ,, ,,~
,!"" "" ~ w'w~w""" "~ wT~ w""" " , w ~ ~ ~ l
~
s***,,,~
uh i teno i se
S t r a tn ( u E )
_,.~
~ 'Vtr'~'~llr~~~r,'~~~l~"p~ ~ , . ~ ~l,r~( ~,r,~,vv,
(a) I~M5 Power(uE^2.
Hz ^-I)
nerrow
,b I~I"15 P o u e r ( u E ^ 2 .
8b
6b
wide
H~^-I )
0'''II I O0
~ ~ ~ @MS
(b)
o
.,b
20 Pouer(uE^2.
20
~6o
m,-~:m,
60
"
sb
l-lz^-I )
~6o
HII.
uhj L e n o i ~ e
-o
~o
8o
.j?o
FIGURE 10.8 (a) Narrow-band, wide-band, and white noise random processes, (b) power spectral densities (PSDs) of narrow-band, wide-band, and white noise random processes.
A bel is the common logarithm of the ratio of two measurements of power. A decibel (db) is one-tenth of a bel and is defined by
wx(/~)
db - 10 lOgl0 Wx([i)
(10.5.17)
A doubling of PSD corresponds to an increase of approximately 3 db. If the PSD doubles for each doubling of frequency, the spectrum increases
380
FATIGUE
ANALYSIS
IN THE
FREQUENCY
DOMAIN
at 3 db/octave or the spectrum has a positive roll-off rate of 3 db/octave. Figure 10.9 shows a random vibration test profile with a constant positive and negative roll-off rate of 6 db/octave. It is noted that on a log-log scale the constant decibel per octave lines are straight. The negative roll-off rate of 6 db/octave indicates a power decrease of 4:1 in one octave from 1000 to 2000 Hz.
Example 10.1 Calculate the RMS acceleration for the Wx(f) profile shown in Figure 10.9.
Solution. The calculation of the RMS acceleration may appear relatively simple, at first, as being the RMS sum of the three areas under the PSD curve in Figure 10.9. But because the profile is drawn on a log scale, the slopes that appear to be straight lines are actually exponential curves. Thus, we need to replot using linear graph paper as in Figure 10.10 and subsequently calculate the area under the exponential curves. Mathematically, the exponential curve can be expressed in the following form:
(]0.5.18)
W x ( f ) - A • (f)S
where A and B are the exponential parameters. B is the slope of the roll-off line on log-log scales and is defined as
B
logWx(f2) wxqi )
m
logs ~ jl
1 x ( lOlog Wx(f2))) lo wxoc~ log 2 x
('
log2
log
(10.5.19) j~/
1
x(the roll-off rate in db/octave)
=
3-
2
3
A N
6 db/octave
+ 6 db/octave C O
"o
E O O L_
3: o a.
0.1
i
1000
100
10
10000
Frequency (Hz)
FIGURE
1 0.9
Example of a random vibration test profile on log log scales.
38 1
S P E CT RA L DENSITY
. . . . . . . . . . .
1
i
. . . . . . . . .
~:
-
~" 0.8
.
0.6
i
.
.
.
.
. . . . . . .
:o,i i i
.........
"
~ "0 0
500
1000
1500
2000
Frequency (Hz) FIGURE 10. 1 0
Example of a random vibration test profile on linear scales.
A is the exponential coefficient and is obtained as
A-
Wx~l) = (fl) 8
Wxq2~
(10.5.20)
(f2) 8
The area under the roll-off curve is therefore J~ P
Area-
J~ .]9 W x ( f ) d f - j A x ( J ) ' B d 9J - ~A B+I Yi .li
[(f2)8+ i
- ( A ) e+']
For the first segment 1-2, that the following is known: Ji - 50 Hz,
W x ( J i ) - 0.25 g2/Hz
and J~-
100Hz,
Wx(f2)-
1.00gZ/Hz.
B and A are therefore equal to log B
~
1.0
0.25 = 2 100 log 50
and m
m
wx(A) (A) 8
0.25 _ 10-4 502
(10.5.21)
FATIGUE ANALYSIS IN THE FREQUENCY DOMAIN
382
The area under the segment 1-2 is then equal to 10-4 Areal_2 - B A+____I[(fz)B+l --(fl)B+l] -- 2 + 1 [1002+1 -- 502+1]- 29"282 For the second segment 2-3, the calculation of the rectangular area is Area2_3 = 1.0 • (1000 - 100) = 900 g2
For the third segment 3-4, the following is known: Ji = 1000 Hz, Wx(f~) = 1.00 g2/Hz and f2 = 2000 Hz, Wx(f2) = 0.25 g2/Hz. By using Equations 10.5.19 and 10.5.20, B and A are equal to - 2 and 106, respectively. The area under the segment 3-4 is A [ B+I )B+l] _ 106 [2000-2+1 _ 1000-2+l] Area3_4 - B +----1 (f2) -(fl -2 + 1
_ 500 g2
Thus, the total area is equal to Areal_2 + Area2_3 + Area3_4 = 29.2 + 900 + 500 = 1429.2g 2 and the R M S acceleration is equal to v/1429.2 = 37.8 g.
10.6
LEVEL
CROSSING RATE OF NARROW-BAND RANDOM PROCESSES
For a continuous and differentiable stationary process, X(t), the expected number of positively sloped crossing (upcrossing) in an infinitesimal interval is only dependent on dt. We have
E[N.+ (dt)] = oo+ dt
(l 0.6.1)
where o,,+ is the expected rate of upcrossing per time unit. If A denotes the event that any random sample from X(t) has an upcrossing x - a in an infinitesimal time interval dt, the probability of Event A, P(A), is given
P(A) = oo.dt Equation 10.6.2 allows us to express exist, we must have
a - X(t) < X(t) < a
Oa +
(10.6.2)
in terms of P(A). For Event A to
and
~'(t) > 0
L E V E L C R O S S I N G R A T E OF N A R R O W - B A N D
383
RANDOM PROCESSES
Combining these two conditions, P(A) can be written as P(A) = e ( a - X(t) < X(t) < a f3 X(t) > 0)
(10.6.3)
These conditons define a triangle area in the X(t) - ~'(t) plane, as shown in Figure 10.11. The probability of Event A is calculated by integrating the joint PDF of X(t) and ~'(t) over this region, i.e., P(A)
i fxx(U, v)dudv 0 a-vdt
(10.6.4)
Substituting Equation 10.6.4 into Equation 10.6.2 leads to the following expression of the level upcrossing rate for a stationary random process:
(10.6.5)
l)a+ -- i I)fX'~ (a, o)do o If X(t) is Gaussian, the expected upcrossing rate of x - a is 1 e2 exp
(10.6.6)
The expected rate of zero upcrossings E[0 +] is found by letting a Equation 10.6.6: 1 tr2 E[0+] = 2r~ trx
0 in
(10.6.7)
By using Equations 10.5.7 and 10.5.12, the expected rate of zero upcrossing is
h.._
v
X
a
FIGURE
10. 1 1
The region where Event A occurs.
384
F A T I G U E A N A L Y S I S IN T H E F R E Q U E N C Y D O M A I N
O0
f f2 W x ( f ) d f 00(3
El0 +] _
(10.6.8)
,[ Wx(f)df
0
The expected rate of peak crossing, E[P], is found from a similar analysis of the velocity process X(t). The rate of zero downcrossing of the velocity process corresponds to the occurrence of a peak in ~'(t). The result for a Gaussian process is 1 a~,"
E[P] - 2rt ox
(10.6.9)
By using Equations 10.5.12 and 10.5.15, we have OG
f f4 Wx(f)df o ~[f2 Wx(f)df
E[P]-
OC
(10.6 10)
0
A narrow-band process is smooth and harmonic. For every peak there is a corresponding zero upcrossing, meaning El0 +] is equal to E[P]. However, the wide-band process is more irregular. A measure of this irregularity is the ratio of the zero upcrossing rate to the peak crossing rate. The ratio is known as the irregularity factor, ~, i.e., E[O+] -
E[p]
(10.6.11)
Figure 10.12 shows a simple example of how to calculate the irregularity factor. Note that when ~, ~ 0, there is an infinite number of peaks for every zero upcrossing. This is considered a wide band random process. The vaue of ~, = 1 corresponds to one peak per one zero upcrossing and it represents a narrow-band random process. Alternatively, a narrow- or wide-band process can be judged by the width of its spectrum. For this reason, the spectral width parameter, 2, is introduced as 2 = V/1 - y2
(10.6.12)
Note that 2 ~ 0 represents a narrow-band random process. If Mj is the jth moment of a one-sided PSD function (Figure 10.13) defined as OG
Mj - If i Ws~(f )df o
(10.6.13)
L E V E L C R O S S I N G R A T E OF N A R R O W - B A N D
385
RANDOM PROCESSES
E [0 §
4
E[P] 7
/ Ix
d
v
t
T FIGURE 1 0 . 1 2
v
Calculation of the irregularity factor 7.
%
..-t, 9"-(Hz)
-.I' ~ d, FIGURE 1 0 . 1 3
Moments from a one-side PSD.
the rate of zero crossings E[0 +] and rate o f peaks
E[P] are
given by (10.6.14)
E[P]- V/--~
(10.6.15)
MoM4
(10.6.16)
7
-
386
F A T I G U E A N A L Y S I S IN T H E F R E Q U E N C Y D O M A I N
/ 2-
M2
~/]
10.7
(10.6.17)
MoM4
V
MODELS OF FATIGUE DAMAGE UNDER NARROWBAND RANDOM PROCESSES
In this section, a linear S - N curve in log-log coordinates and the linear damage accumulation rule are employed. It is assumed that all stress cycles have the same mean value. If a mean stress exists, it can be accounted for by various mean stress correction models (see Chapter 4). Variable amplitude loading is simulated by a sequence of blocks of constant amplitudes. Damage D is defined as k
D - Z N fi=l ,
ni
(10.7.1)
i
where n; is the total number of cycles in the ith block of constant-stress amplitude S,~,i, Nf, i is the number of cycles to failure under S~,i, and k is the total number of blocks. Failure occurs when D > 1. Recall that the relationship between constant stress amplitude S~,~ (or constant stress range S~,~) and the fatigue life Nf, i (see Chapter 4) has the following expression: m
Sa, i : S;(2NT, i) b
(10.7.2)
Nf , , = -~
(10.7.3)
or
where S) is the fatigue strength coefficient and b is the fatigue strength exponent. If m = - 1 / b S - N curve is given by
and A = 0.5 x Nf, i
-
-
A x Sa,m
, an alternative form of the (10.7.4)
This is a very convenient formula to predict Nf, i for a given S,,~. The cycle counting histogram (also called discrete stress spectrum) for a narrow-band stress process S(t) can be established by either performing the rainflow cycle counting technique (see Chapter 3) or by the counting the number of peaks ni in the window AS around a stress level. Suppose k that the total number of peaks counted in the S(t) process is denoted by ~ n;. The probability fi that a stress amplitude Sa = Sa, i will occur is i=1
M O D E L S OF F A T I G U E D A M A G E U N D E R N A R R O W - B A N D
387
RANDOM PROCESSES
ni
(10.7.5)
J~= k
~-~ ni
i=1
Thus,f. is the PDF of the random variable S~. In this case, the total fatigue damage can be written as k
k
k fi ~"~ ?/i D - / ~ 1 nN~, - ~ i=fl, 9= i .= i
(10.7.6)
By using the linear S-N model, the expression for fatigue damage is
k ~-~ ni k
D
-i=lA ZfiS a',i
(10.7.7)
i=l
Also, the expected value of S m is
k
s7
(10.7.8)
i=1 Therefore, k ~ni
D = ~=!. E(S~) A
(10.7.9)
k where the total count of cycles ~ n i is equal to the rate of upcrossing i=1 multiplying the total time period T ( = E[0 +] x T). Assume that the PDF of stress amplitude S~ can be treated as a continuous random variable as shown in Figure 10.14. The expected value of Sam is O(3
E(S~) = J ~fs,(s~)ds,
(10.7.10)
0
Even though any statistical model of S~ can be employed, it is common to use the Weibull distribution with the following CDF"
Fs,(Sa) - 1 - exp -
(10.7.11)
where ~ and /3 are the scale parameter (characteristic life) and shape parameter (Weibull slope), respectively. For the Weibull distribution,
388
F A T I G U E A N A L Y S I S IN T H E F R E Q U E N C Y D O M A I N
Sa oo c 0 "0 .m .Q 0 Q.
fsa(Sa)
~
a
~
~
stress amplitude FIGURE 1 0 . 1 4
Continuous PDF of stress amplitude.
where F(.) is the gamma function. In the special case in which fl = 2, the Weibull distribution reduces to the Rayleigh distribution. This is an important case, because Rayleigh is the distribution of peaks or ranges or amplitude in a stationary narrow-band Gaussian process that has an RMS value of as. Also, it can be shown that = 2x/~s
(10.7.13)
Therefore, if S(t) is stationary narrow-band Gaussian and the stress amplitudes follow the Rayleigh distribution,
+ l)
(10
where as - x/--M-0. In such a case, therefore, the expected total fatigue damage DNB of a zeromean stationary narrow-band Gaussian stress process over a time interval can be written as k
~-~ ni D U B - /:1 E ( s m ) _ E[0 +] • T A A
10.8
MODELS
1-' m
+1
OF FATIGUE DAMAGE UNDER RANDOM PROCESSES
(10.7.15)
WIDE-BAND
Based on the rainflow counting method, a model for predicting fatigue damage under stationary wide band Gaussian stress process has been
M O D E L S OF F A T I G U E
DAMAGE
UNDER WIDE-BAND
RANDOM
389
PROCESSES
proposed by Wirsching and Light (1980). By using the narrow-band approach as a starting point, the general expression for the damage DWB,Wirschingover a time interval ~ is (10.8.1 )
DWB, Wirsching = ~ W DNB
where DNB is the fatigue damage under a narrow-band random process and ~w is the rainflow correction factor. ~w is an empirical factor derived from extensive Monte Carlo simulations that include a variety of spectral density functions. It is expressed as follows:
~s = as + [1 - as](1 -
2) bw
(10.8.2)
where
as =
0.926 - 0.033m
bs =
1.587m - 2.323
Note that m is the slope of the S-N curve defined in Equation 10.7.4, and 2 is the spectral width parameter defined in Equations 10.6.12 or 10.6.17. Oritz and Chen (1987) also derived another similar expression for fatigue damage, Dwa, Ortiz, under wide-band stresses as Dwu, Oritz =
{oDNB
(10.8.3)
where
~o
1, / M2Mk
--'TVa/~M-k+2and
k-
2.0 m
(10.8.4)
The irregularity factor ), is defined in Equation 10.6.11 or 10.6.16. Instead of using the damage correction factor from the narrow-band stresses to the wide-band stresses, Dirlik (1985) developed an empirical closed-form expression for the P D F of rainflow amplitude, Jso(sa), based on extensive Monte Carlo simulations of the stress amplitudes. Dirlik's solutions were sucessfully verified by Bishop (1988) in theory. Dirlik's damage model for a time period of ~: is as follows: O(3
DWB, Dirlik --
E[Plz l .W~fs~ A
(10.8.5)
0 9
DI
fs..(s")-2v/-M~ Q
e~ •
D2 x Z
-z__E2xs2
+ ~ e2R: 2x/-MTR 2
D3 x Z
z2
+ ~ x e : • 2v/Mo
where Z
z
1 2x/-Mo
M2 7=v/MoM4
M1
Xm--~o
M~M ~
(10.8.6)
390
FATIGUE
D1 -
2(Xm - 72) 1 +
7 - Xm -- DI2
R =
~,?+
D3=
ANALYSIS
1 -
7-
1-DI-D2
D I +
Q=
IN T H E F R E Q U E N C Y
DOMAIN
1 - 7 - DI -F DI2
D2
-
D~
1 -
R
1 . 2 5 ( 7 - D3 - D2 x R) D1
Bishop a n d S h e r r a t t (1989) a n d Bishop (1994) c o n c l u d e d t h a t the Dirlik s o l u t i o n gives the better results as c o m p a r e d with the c o r r e s p o n d i n g time d o m a i n results.
Example 10.7.1. A h o t - r o l l e d c o m p o n e n t m a d e o f S A E
1008 steel is subjected to r a n d o m l o a d i n g process. T h e fatigue strength coefficient (Sf) a n d the fatigue s t r e n g t h e x p o n e n t , (b) are 1297 M P a a n d - 0 . 1 8 , respectively T h e stress response at a critical location is calculated in terms o f the P S D in F i g u r e I 0.15. T h e P S D has two frequencies o f 1 a n d 10 Hz, c o r r e s p o n d i n g to 10,000 a n d 2500 M P a 2 / H z , respectively. D e t e r m i n e the fatigue d a m a g e o f this c o m p o n e n t , using the e q u a t i o n s for w i d e - b a n d stresses. N o t e t h a t a crest factor is defined as the ratio o f p e a k ( m a x i m u m ) to R M S values. A sine wave has a crest factor o f v/2 = 1.414.
Solution. W e will calculate the fatigue d a m a g e o f the c o m p o n e n t in time d o m a i n first a n d use it as baseline i n f o r m a t i o n to c o m p a r e with the p r e d i c t e d d a m a g e based on the P S D . T h e original stress time history can be o b t a i n e d by a d d i n g t w o sine waves, o n e for each block in the P S D . F o r a sine wave, the a m p l i t u d e o f each is calculated f r o m 1.414 times the R M S value, i.e., the area o f each P S D block. T h e stress a m p l i t u d e So, l o f the sine wave at 1 Hz is
Sa, 1
x/10000 x 1 x 1.414 = 1 4 1 . 4 M P a
=
12000---------!
:
I
I
!
:
!
I
!
f
i
~
;
~
10000-
"~" -1- 8000 t'l;I
a..
.....~ I
+
I
i
'_. ~ t ] i ' 1 : i ~ ! i ....+++.~_[---+I:a~........+_.... i....4. . . . . . . +-.-_-+~.+-~-L----J,........ ~-.++......... ;...........++..... + i ; I + ; i ', I + + + + + + + + +I ' + ++ ~ I +, .1 ...... ;+ --'..+ --:- ' ,-[ i'+ '9 . . . .+ . . . 9. + " +' ' + ' I
l i l + J
6000-
............
a o} IX. 4 0 0 0 -
i ..... I ....
i
+~+....
++................... +
..... .... .....
I
+--I
+
+
+-+--i
I
i
I
I
',
I
;,
+ 'I ........I
.I++I++,
:
...... T--+
.......... f
i ....
+,......... ,+
....... i +............
~t ++ ....
t,
+,,+
+
",...... r+................ r +
+.......... .
+ + + , +
++
........
+.......... + ...... ~ ....
+
t
i ....
+.....
ll-+ -r +..... +
'
'
:
+~
:++.....
1~,....
+
+ ............
+
?;
..... t....
+-+ .... :+ ....... ~ ..... T .... + .... i .....
:
'
+
.... I --~+- -I.... "-- ,+........... ,~...................+.... ~.... ,++-t ....... <.- -+........ +-.... .-.... - ,---4 .............-+..... 2OOO
Iii+ . . . . . .
'.
,
t
FIGURE
1
i ..... + ........
10. 1 5
2
+
i
+
2 ++
0 0
:
3
"
,--t-+-i; 4
:
I ++
i
+
'
!
+
+ ++
.... ++-+-+--++ .... +-4
5 6 7 Frequency (Hz)
+
"
'
+
~
'
'
+,,:+.,+.z I..-L..I
,
I
' +
.... T ..... +.....
8
9
10
11
12
PSD of the stress response of a component made of SAE 1008 steel.
M O D E L S OF F A T I G U E D A M A G E U N D E R W I D E - B A N D
39 1
RANDOM PROCESSES
The stress amplitude Sa, 2 of the second sine wave at 10 Hz is S,,2 = v/2500 x 1 • 1.414 = 7 0 . 7 M P a The material S - N curve for vibration fatigue usually follows the following expression:
Nf, i = A x S a ,
m
where
m = -1/b = -1/(A - 0.5 x
0.18) = 5.56
- 0.5 x 12975.56 - 1.02 x l017 M P a
With this S - N curve, we determine the fatigue life for each sine wave as follows:
Nf, l - 1.02
x
1017 • 141.4 . 5 . 5 6 - 1.13 X 105 cycles
N7,2 - 1.02 x l017 x 70.7 -5s6 - 5.32 x l06 cycles In a 1-sec time interval, the sine waves at 1 and l0 Hz represent 1 cycle (i.e., nl = 1) and 10 cycles (i.e., n2 = 10) of reversed loading, respectively. The linear d a m a g e calculation for this time interval gives Dtime
nl n• -~-~= Nf, l Nf, 2
=
1 10 -1.08• 1.13 • l0 s + 5.32 x 106 -
O_ 5
This corresponds to a fatigue life of 93,000 seconds (= 1/1.08 • 10-5). We then calculate the fatigue d a m a g e to the c o m p o n e n t based on the given PSD. The aforementioned frequency d o m a i n methods will be used for the d a m a g e estimation. The j t h m o m e n t of the PSD can be easily calculated as O(3
Mj
-
| f Ws~
lJ •
10000 • l +
10J • 2500 • 1
,3
0
M0=
1~
10000z 1+10 ~215215
1 = 12500
M2 = 12 • 10000 z 1 + 102 • 2500 • 1 = 260000 M4 = 14 • 10000 • 1 + 104 • 2500 • 1 = 25010000 F r o m this we can c o m p u t e El0+] - ~ ~ 0 2 -
E[P] -
M~_
V/262~176 - 4.56 ,/25010000 V 260000 - 9.81
392
FATIGUE ANALYSIS
.
E[0 +] . . .
4.56 . 9.81
E[P]
IN T H E F R E Q U E N C Y
DOMAIN
0.465
2 -- V/1 - 72 -- V/I - 0.4652 -- 0.885
1. W i r s c h i n g a n d L i g h t M e t h o d F a t i g u e d a m a g e DNB o f a z e r o - m e a n s t a t i o n a r y n a r r o w - b a n d stress p r o c e s s o v e r a t i m e i n t e r v a l r - 1 sec c a n be w r i t t e n as DNB-
Gaussian
E[0+]x r ( 2V/~o0)m(2) 6 + 1) 4.56 x 1 ( x/'2 x 12500 )5561-'(~5"5 A
F
DNB
--
1.02 x 1017
DNB
--
0.000345
+
1
T h e r a i n f l o w c o r r e c t i o n f a c t o r is c a l c u l a t e d as aw - 0.926 - 0 . 0 3 3 m - 0.926 - 0.033 x 5.56 - 0.743 bw-
1.587m-
2.323 -
1.587 x 5 . 5 6 -
2.323 - 6.501
~w - a w + [1 - aw](1 - 2) bw - 0.743 + (1 - 0.743) x (1 - 0.885) 65~ - 0.743 F i n a l l y , the f a t i g u e d a m a g e DwB, Wirsching is c o m p u t e d as DWB, Wirsching = ~wDNB -- 0.743 X 0.000345 -- 0 . 0 0 0 2 5 6 / s e c o n d
T h i s c o r r e s p o n d s to a f a t i g u e life o f 3900 s e c o n d s ( = 1 / 0 . 0 0 0 2 5 6 ) , w h i c h is v e r y c o n s e r v a t i v e as c o m p a r e d w i t h the b a s e l i n e fatigue life (93,000 seconds). 2. O r i t z a n d C h e n M e t h o d DNB - - 0 . 0 0 0 3 4 5 is the s a m e as t h a t c a l c u l a t e d p r e v i o u s l y . C a l c u l a t i o n o f the r a i n f l o w c o r r e c t i o n f a c t o r ~'o is r e q u i r e d . G i v e n the slope o f the S - N curve, m - 5.56, 2.0 2.0 k . . . . . m 5.56 Mk -- 1~
Mk+2-
0.3597
x 10000 x 1 + 10 0.3597 • 2500 • 1 -- 15723
10.3597+2 X 10000 X 1 + 10 0.3597+2 • 2500 • 1 -- 582338
1 /M2Mk
1
,/ 60000
~o - ~ V M o M k + 2 - 0.465 V i 2-if0--0 x ~582-338 T h e f a t i g u e d a m a g e DwB, Lutes is c o m p u t e d as
=0.432
M O D E L S OF F A T I G U E D A M A G E U N D E R W I D E - B A N D
DwB, Oritz -- U,o D N B
393
RANDOM PROCESSES
- - 0.432 x 0.000345 -
0.00015/second
This c o r r e s p o n d s to a fatigue life of 6 7 0 0 s e c o n d s ( = 1/0.00015), which is very conservative as c o m p a r e d with the baseline fatigue life (93,000 seconds). 3. Dirlik M e t h o d It is required to d e t e r m i n e the following p a r a m e t e r s for the P D F o f stress a m p l i t u d e s that have been rainflow cycle c o u n t e d . Z =
1 2 x/--M~
Xm-
1
=
= 0.004472
2x/12500
MI ~ / ~ 2 _ 35000 ~ / 260000 M00 V M4 i-~0-0 g 2~i0-0-00 - 0.2855 - ;2) 2(0.2855 - 0.4652) 1 + ,,2 = 1 + 0.4652 = 0.1074
2(Xm
D1 = R-
7 -
Xm -
D~
_-
1-7-D~+D~
0.465-0.2855-0.10742 1-0.465-0.1074+0.10742
--0.3825
1 - 7 - D~ + D~ 1 - 0.465 - 0.1074 + 0.10742 _ 0.7111 1- R l - 0.3825 D3 - 1 - DI - D2 - 1 - 0.1074 - 0.7111 - 0.1815 D
-)
--
z
1.25(7 - D3 - D2 x R)
Q-
Dl
--
1.25(0.465 - 0.1815 - 0.7111 x 0.3825)
=
0.1074
=0.1339
Substituting these values into E q u a t i o n 10.8.6 gives the Dirlik's P S D of stress a m p l i t u d e s as follows: fs',,(sa)
- - 0 . 0 0 3 5 8 7 e -~176
+ 9.72 x 10-Se -684x10 5,.2, + 3.63 x 10-6e - l ~ 1 7 6 ss~
T h e fatigue d a m a g e in a given time period of 1 second is
DwB, Dirlik - -
E[P]r J ,-,,,,r A
_
~
0
9.82 x 1 ,5.56,~. (s,)d&, 1.02 x 1017 ' aa .IS,, 0
oc ill
-- 9.627 x 10 --17 / aa'5"56"./s,,(sa)dsa /
o
The previous e q u a t i o n can be a p p r o x i m a t e d by the discrete i n t e r g r a t i o n scheme listed next:
394
F A T I G U E A N A L Y S I S IN T H E F R E Q U E N C Y D O M A I N OQ
,',", 5 56 -0 0334s,~,i -(0.0035~51Sa~ i e " -+-
DwB, Dirlik ~ 9.627 X 10 -17 i=1
9.72 x 1"-SUSa, i5.56e-6.84xlO5~,i)msa..~ -
lO-6~5"56e-l'OxlO-S~")Asa'a,i
O(3
9.627
x
10-17(~-~ 3.63
x
i=1
where As a
sa, i = - - ~
+ (j -
1)Asa
The numerical integration technique leads to DwB, Dirlik 1.268 x 10-5/second. This corresponds to a fatigue life of 79,000 seconds, which correlates better to the baseline fatigue life (93,000 seconds). Therefore, the Dirlik method is the preferred method for fatigue damage calculation based on the PSD and has been widely adopted by many commercial fatigue software packages. =
REFERENCES Bendat, J. S., and Piersol, A. G., Random Data, Analysis and Measurement Procedures, 2nd ed., Wiley Interscience, New York, 1986. Bishop, N. W. M., Spectral methods for estimating the integrity of structural components subjected to random loading. In Handbook of Fatigue Crack Propagation in Metallic Structures, Vol. 2, Carpinteri, A. (Eds.), Elsevier, Dordrecht, 1994, pp. 1685-1720. Bishop, N. W. M., The Use of Frequency Domain Parameters to Predict Structural Fatigue, Ph.D. Thesis, Warwick University, 1988. Bishop, N. W. M., and Sherratt, F., Fatigue life prediction from power spectral density data. Part 2: Recent Development, Environmental Engineering, Vol. 2, Nos. 1 and 2, 1989, pp. 5-10. Crescimanno, M. R., and Cavallo, P., On duty simulation of a trimmed body under dynamics loads: modal superposition approach to evaluate fatigue life, SAE Paper 1999-01-3150, 1999. Dirlik, T., Application of Computers in Fatigue Analysis, Ph.D. Thesis, Warwick University, 1985. Lee, Y., Raymond, M. N., and Villaire, M. A., Durability design process of a vehicle suspension component, Journal of Testing and Evaluation, Vol. 23, No. 5, 1995, pp. 354-363. Newland, D. E., An Introduction to Random Vibrations and Spectral Analysis, 2nd ed., Longman, New York, 1984. Oritz, K., and Chen, N. K., Fatigue damage prediction for stationary wide-band stresses. Presented at the 5th International Conference on the Applications of Statistics and Probability in Civil Engineering, Vancouver, Canada, 1987. Vellaichamy, S., Transient dynamic fatigue analysis using inertia relief approach with modal resonance augmentation, SAE Paper 2002-01-3119, 2002. Wirsching, P. H., and Light, M. C., Fatigue under wide band random stresses, ASCE Journal of Structural Division, Vol. 106, 1980, pp. 1593--1607. Wirsching, P. H., Paez, T. L., and Oritiz, K., Random Vibrations." Theory and Practice, Wiley, New York, 1995.