Numerical analysis of cracked beams using nonlinear output frequency response functions

Available online at www.sciencedirect.com

Computers and Structures 86 (2008) 1809–1818 www.elsevier.com/locate/compstruc

Numerical analysis of cracked beams using nonlinear output frequency response functions Z.K. Peng a,*, Z.Q. Lang a, F.L. Chu b a

Department of Automatic Control and Systems Engineering, University of Sheffield, Mappin Street, Sheffield, S1 3JD, UK b Department of Precision Instruments, Tsinghua University, Beijing 100084, PR China Received 2 April 2007; accepted 21 January 2008 Available online 17 March 2008

Abstract Nonlinear output frequency response functions (NOFRFs) are a new concept recently developed for the analysis of nonlinear systems in the frequency domain. Based on this concept, the output frequency behaviors of nonlinear systems can be expressed using a number of one-dimensional functions similar to the approach used in the traditional frequency response analysis of linear systems. In this paper, the NOFRFs are first employed to analyze a typical representation for cracked structures, a single degree of freedom system (SDOF) bilinear model, to explain the occurrence of the nonlinear phenomena when a cracked structure is subjected to sinusoidal excitations, including the generation of super-harmonic components and sub-resonances. Then the NOFRF concept is used to analyze the crack induced nonlinear response of a beam represented by a finite model. The results show that higher order NOFRFs are extremely sensitive to the appearance of cracks in the beam, and can therefore be used as crack damage indicators to indicate the existences and the sizes of cracks. This research study establishes an important basis for the application of the NOFRF concept in fault diagnosis of mechanical structures. Ó 2008 Elsevier Ltd. All rights reserved. Keywords: Nonlinear vibration; Crack detection; Volterra series; Fault detection; Bilinear oscillator

1. Introduction Fatigue cracks are a potential source of catastrophic failure in civil structures or mechanical machines. To avoid the failure caused by cracks, many researchers have performed extensive investigations to develop structural integrity monitoring techniques. Most of the techniques are based on vibration measurement and analysis because, in most cases, vibration based methods can offer an effective and convenient way to detect fatigue cracks in structures. Generally, vibration based methods can be classified into two categories: the linear approach and the nonlinear approach. By linear approaches the presences of cracks in a target object are detected through the monitoring of changes in resonant frequencies [1] or in mode shapes [2] *

Corresponding author. E-mail addresses: [email protected] (Z.K. Peng), z.lang@ sheffield.ac.uk (Z.Q. Lang), chufl@mail.tsinghua.edu.cn (F.L. Chu). 0045-7949/$ - see front matter Ó 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.compstruc.2008.01.011

or in damping factors [3]. However, some researches have shown the linear detection procedures do not always come up to practical requirement because of low sensitivity to defects. For instance, in Ref. [4], the numerical results show that the presence of a crack, which makes up about 10–20% of the undamaged cross-section area, reduces the natural frequencies of a beam only 0.6–1.9%. Increasingly over recent years, attentions have been focused on the investigations of the nonlinear effects caused by the presence of cracks and associated applications to the problem of crack detections. When a cracked object is subjected to a single harmonic input, main distinctive features of such a vibration system are the appearance of super-harmonic components and subharmonic resonances [4–8]. In the literatures [4–8], these phenomena are termed as the nonlinear effects. In Ref. [4], Bovsunovsky and Surace have claimed that the nonlinear effects are more sensitive to the presence of a crack than the changes of natural frequencies, mode shapes, and they studied the effect of the damping on the

1810

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

nonlinear effects. Also based on the nonlinear effect of subharmonic resonances, Tsyfansky and Beresnevich [5] have developed a new approach for the detection of fatigue cracks in flexible, geometrically nonlinear beam-type structural elements. Sundermeyer and Weaver [6] have studied the forced responses of a bilinear model subjected to a excitation with two frequencies, and based which they have further exploited the weakly nonlinear character of a cracked beam for the purpose of determining crack location. In Ref. [7], Saavedra and Cuitino have studied dynamic behaviors of different multi-beams systems containing a transverse crack both theoretically and experimentally. Some results about the nonlinear effects that would be useful for crack detection have been given. Darpe et al. [8] studied the transient response and breathing behavior of a cracked Jeffcott rotor passing through its critical speed and subharmonic resonances. Their work has further confirmed the existence of subharmonic resonances for a cracked mechanical structure. In summary, the nonlinear effect analysis based methods would be much more sensitive to presences of cracks than the linear vibration based methods. The nonlinear effects of super-harmonic components and subharmonic resonances caused by cracks have been observed in numerical studies as well as experimental studies and, during the last decades, great attention has been paid by research scientists to the diagnosis of cracks in vibrating structures and in rotating machinery by using the crack–associated nonlinear effects. Excellent reviews about the crack subject study can be found in [9–12]. In [9], Imam et al. have developed a technology for on-line rotor crack detection and monitoring. This technique is based on the analytical modeling of the dynamics of the system where a 3-D finite element crack model and a nonlinear rotor dynamic code have also been developed to accurately model a cracked rotor system. The early studies about the dynamics of cracked rotors have been reviewed in Ref. [10] by Wauer, where various modeling of the cracked part of structures and different detection procedures to diagnose fracture damage have been summarized. Later, Dimarogonas [11] has contributed a more comprehensive review about the cracked object studies, in which over 300 papers about the various cracked objects including beam, rotor, plates, shell and blades, et al. have been presented, different methods and theories related to the crack have been briefly introduced and, moreover, a list of topics of interest for further research have been suggested. This really provided an excellent overall view to the scope of the crack study. In 2004, Sabnavis and his colleagues [12] have made a review of the literature on cracked shaft detection and diagnostics published after 1990. However few efforts have been dedicated to theoretically explain why and how the nonlinear effects occur in a cracked object. It is well known that the nonlinearity caused by the opening and closing of a crack can be well modeled as a piecewise-linear model, known as a bilinear oscillator model, and therefore many studies about cracks [13–17] have adopted the bilinear oscillator model. Based

on a simple bilinear oscillator, in this paper the nonlinear effects caused are interpreted by cracks using a new concept known as nonlinear output frequency response functions (NOFRFs) [18,19]. The concept of NOFRFs in general allows the analysis of nonlinear systems to be implemented in a manner, which is similar to the analysis of the linear system frequency response. This provides a great insight into how nonlinear phenomena including the generation of new frequency components and the occurrence of the sub-resonances in practices. Moreover, in this study the NOFRFs are used to analyze a cracked cantilever beam described by a finite element model. The results show that the NOFRFs are a very sensitive indicator of the presence and size of cracks. This research study establishes an important basis for the application of the NOFRF concept in fault diagnosis of mechanical structures. 2. Nonlinear output frequency response functions The definition of NOFRFs is based on the Volterra series theory of nonlinear systems. The Volterra series extends the well-known convolution integral description for linear systems to a series of multi-dimensional convolution integrals, which can be used to represent a wide class of nonlinear systems [20]. Consider the class of nonlinear systems which are stable at zero equilibrium and which can be described in the neighbourhood of the equilibrium by the Volterra series Z 1 N Z 1 n X Y yðtÞ ¼  hn ðs1 ; . . . ; sn Þ uðt  si Þdsi ð1Þ n¼1

1

1

i¼1

where yðtÞ and uðtÞ are the output and input of the system, hðs1 ; . . . ; sn Þ is the nth order Volterra kernel, and N denotes the maximum order of the system nonlinearity. Lang and Billings [20] derived an expression for the output frequency response of this class of nonlinear systems to a general input. The result is 8 N P > > Y ðjxÞ ¼ Y n ðjxÞ for 8x > > > n¼1 > < pffiffi R 1= n ð2Þ H n ðjx1 ; . . . ; jxn Þ Y n ðjxÞ ¼ ð2pÞ n1 x1 þ...þxn ¼x > > > n > Q > >  U ðjxi Þ drnx : i¼1

This expression reveals how nonlinear mechanisms operate on the input spectra to produce the system output frequency response. In Eq. (2), Y ðjxÞ is the spectrum of the system output, Y n ðjxÞ represents the nth order output frequency response of the system, Z 1 Z 1 ... hn ðs1 ; . . . ; sn Þeðx1 s1 þþxn sn Þj H n ðjx1 ; . . . ; jxn Þ ¼ 1

1

 ds1 . . . dsn

ð3Þ

is the nth order generalised frequency response function (GFRF) [20], and

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

Z

H n ðjx1 ; . . . ; jxn Þ

x1 þ...þxn ¼x

n Y

un(t)

U ðjxi Þ drnx

i¼1

Qn denotes the integration of H n ðjx1 ; . . . ; jxn Þ i¼1 U ðjxi Þ over the n-dimensional hyper-plane x1 þ    þ xn ¼ x. Eq. (2) is a natural extension of the well-known linear relationship Y ðjxÞ ¼ H ðjxÞU ðjxÞ, where H ðjxÞ is the frequency response function, to the nonlinear case. For linear systems, the possible output frequencies are the same as the frequencies in the input. For nonlinear systems described by Eq. (1), however, the relationship between the input and output frequencies is more complicated. Given the frequency range of an input, the output frequencies of system (1) can be determined using the explicit expression derived by Lang and Billings in [20]. Based on the above results for the output frequency response of nonlinear systems, a new concept known as the nonlinear output frequency response function (NOFRF) was recently introduced by Lang and Billings [18]. The NOFRF is defined as R Qn H n ðjx1 ; . . . ; jxn Þ i¼1 U ðjxi Þ drnx x1 þþxn ¼x R Qn Gn ðjxÞ ¼ ð4Þ i¼1 U ðjxi Þ drnx x1 þþxn ¼x under the condition U n ðjxÞ 6¼ 0. Actually, U n ðjxÞ is the Fourier Transform of un ðtÞ, i.e., Z n Y U ðjxi Þ drnx U n ðjxÞ ¼ x1 þþxn ¼x i¼1 Z 1 1 n jxt

¼ pffiffiffiffiffiffi 2p

u ðtÞe

ð5Þ

dt

1

Notice that Gn ðjxÞ is valid over the frequency range of U n ðjxÞ, which can be determined using the algorithm in [20]. By introducing the NOFRFs Gn ðjxÞ, n ¼ 1; . . . ; N , Eq. (2) can be written as Y ðjxÞ ¼

N X

Y n ðjxÞ ¼

n¼1

N X

Gn ðjxÞU n ðjxÞ

ð6Þ

n¼1

which is similar to the description of the output frequency response for linear systems. For a linear system, the relationship between Y ðjxÞ and U ðjxÞ can be illustrated as in Fig. 1. Similarly, the nonlinear system input and output relationship of Eq. (6) can be illustrated as in Fig. 2. The NOFRFs reflect a combined contribution of the system and the input to the system output frequency response behaviour. It can be seen from Eq. (4) that Gn ðjxÞ depends not only on H n ðn ¼ 1; . . . ; N Þ but also on the input U ðjxÞ. According to Eq. (4), the NOFRF Gn ðjxÞ is a weighted sum of H n ðjx1 ; . . . ; jxn Þ over x1 þ    þ xn ¼ x with the weights depending on the test

u(t)

FFT

U( jω)

H( jω)=G1( jω)

Y( jω)

Fig. 1. The output frequency response of a linear system.

u2(t) u(t)

FFT

FFT FFT

UN ( jω) U2 ( jω) U1 ( jω)

1811

GN ( jω) G2( jω) G1( jω)

YN ( jω) Y2 ( jω) Y1 ( jω)

Y( jω)

Fig. 2. The output frequency response of a nonlinear system.

input. Therefore Gn ðjxÞ can be used as an alternative representation of the dynamical properties described by H n . The most important property of the NOFRF Gn ðjxÞ is that it is one-dimensional, and thus allows the analysis of nonlinear systems to be implemented in a convenient manner similar to the analysis of linear systems. Moreover, there is an effective algorithm [18] available which allows the estimation of the NOFRFs to be implemented directly using system input output data. 3. Finite element modelling of a cracked beam The beam under consideration is a solid-spar homogeneous viscoelastic beam shown in Fig. 3. The equation for the beam bending vibrations excited by the force uðtÞ is   o2 yðx; tÞ oyðx; tÞ o2 o2 yðx; tÞ qA þ ðEI þ c  EIðx; x ; d ÞÞ 0 c c ot2 ot ox2 ox2 ¼ uðtÞdðx  LÞ

ð7Þ

where dðxÞ is the Dirac delta function, and q is the mass density of the beam, A ¼ w  d is the cross-sectional area, c is the damping of the beam, E is the Young’s modulus of material, I 0 is the moment of inertia of the beam cross-section, Iðx; xc ; d c Þ is the reduction of the moment of inertia caused by the open crack of depth d c located at xc , and, in this study, an exponential function suggested by Christides and Barr [21] is used, as follows: Iðx; xc ; d c Þ ¼

I0 1 þ C expð2ajx  xc j=dÞ

ð8Þ

where C ¼ ðI 0  I c Þ=I c , and I 0 ¼ wd 3 =12 and 3 I c ¼ wðd  d c Þ =12 are the second moment of the areas of the undamaged section and damaged section respectively. According to the exponentially decay distribution (8), the stiffness reduction distribution will reach maximum at the crack position and decrease to zero outside the crack influence zone. This distribution is quite similar to the additional strain energy distribution caused by the presence of crack according to the fracture theory. Moreover, the stress decay rate used in exponentially decay distribution can be determined by either experimental means or numerical approaches. Shen and Pierre [22] employed the finite element method to determine the decay rate where the singularity crack elements were used, and the obtained rate 0.667 is

1812

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

y

xc

u(t)

w

dc

d

x L Fig. 3. The beam with a crack.

very close to the result experimentally tested by Christides and Barr [21]. To obtain the finite element model of the cracked beam, the two-noded beam element shown in Fig. 4 is employed. Adopting the Hermitian shape functions, the stiffness matrix of the two-noded beam element without a crack is obtained using the standard integration based on the variation in flexural rigidity as Z l T ½K e  ¼ ½BðxÞ EI½BðxÞdx ð9Þ 0



 where ½BðxÞ ¼ ðxÞ , and H 1 ðxÞ, H 2 ðxÞ, H 3 ðxÞ and H 4 ðxÞ are the Hermitian shape functions defined as H 001 ðxÞ

2

H 002 ðxÞ

H 003 ðxÞ

H 004

3

3x 2x þ 3 2 l l 2x2 x3 þ 2 H 2 ðxÞ ¼ x  l l 3x2 2x3 H 3 ðxÞ ¼ 2  3 l l 2x2 x3 þ 2 H 4 ðxÞ ¼  l l H 1 ðxÞ ¼ 1 

2l2

6l

k 24

k 14

k 44

where k 11

"  2 # 12EðI 0  I c Þ 2l3c 2n ¼ þ 3lc 2  1 l2 l4 l

k 12

   12EðI 0  I c Þ l3c 7n 6n2 ¼ þ lc 2  þ 2 l l2 l l3

k 14

   12EðI 0  I c Þ l3c 5n 6n2 ¼ þ lc 1  þ 2 l l2 l l3 "  2 # 12EðI 0  I c Þ 3l3c 3n 2 ¼ þ 2lc l l2 l2

k 24 ¼

k 44

   12EðI 0  I c Þ 3l3c 9n 9n2 þ þ 2l 2  c l l2 l2 l2

"  2 # 12EðI 0  I c Þ 3l3c 3n 1 ¼ þ 2lc l l2 l2

4l2

Assume the stiffness reduction caused by an open crack falls within a single element, and then the stiffness matrix ½K ec  of the cracked element can be written as ½K ec  ¼ ½K e   ½K c 

y υi

k 14

k 22

Assume the beam rigidity EI is constant given by EI 0 within the element, and then the element stiffness matrix is 2 3 12 6l 12 6l 4l2 6l 2l2 7 EI 0 6 6 6l 7 ½K e  ¼ 3 6 ð10Þ 7 l 4 12 6l 12 6l 5 6l

where ½K c  is the reduction in the stiffness matrix due the crack. According to Sinha et al. [14], the matrix ½K c  is 2 3 k 11 k 12 k 11 k 14 6 k k 22 k 12 k 24 7 6 12 7 ð12Þ ½K c  ¼ 6 7 4 k 11 k 12 k 11 k 14 5

ð11Þ

lc ¼ 1:5d, and n is the distance between the left node and the crack, as illustrated in Fig. 5. It is supposed that the crack does not affect the mass distribution of the beam. Therefore, the consistent mass matrix of the beam element can be formulated directly as

y

θi

υ i+1

υc

θ i+1

θc ζ

υc+1

θ c+1

x

x xi=0

xi+1=l

Fig. 4. Two-noded beam element.

xc=0

xc+1=l

Fig. 5. The cracked beam element.

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

e

½M  ¼

Z

l

qA½H ðxÞ 0 2 156 6 qAl 6 22l ¼ 6 420 4 54

T

13l

½H ðxÞdx 22l 4l2

54 13l

13l 3l2

156 22l

3 13l 3l2 7 7 7 22l 5

ð13Þ

4l2

where ½H ðxÞ ¼ f H 1 ðxÞ H 2 ðxÞ H 3 ðxÞ H 4 ðxÞ g. In the dynamic analysis, the system matrix is usually required to be inverted. From this aspect, a diagonalized mass matrix has a computational advantage. In this study, a diagonalized mass matrix is adopted, which is developed from the consistent mass matrix using the approach introduced in [23] and is given below 2 3 39 0 0 0 2 0 07 qAl 6 6 0 l 7 ½M e  ¼ ð14Þ 6 7 78 4 0 0 39 0 5 0

0

0

l2

Moreover, although the damping is definitely affected by the presence of crack, as widely used in the cracked beam analysis [13,17], the damping is still assumed to be proportional to the stiffness and not affected by the crack. Then, after assembling the element matrices, the equations of the motion of a cracked beam can be written as the following bilinear form M€ q þ Dq_ þ Kq ¼ F ðtÞ if crack is close ð15Þ M€ q þ Dq_ þ ðK  K c Þq ¼ F ðtÞ if crack is open where M, K, D are the mass, stiffness and damping matrices of the beam, and K c is the stiffness reduction caused by the crack, and q is the displacement vector of the beam. The status of the crack is determined by the curvature c at the crack location, which is given by [13] as c ¼ ½BðnÞf tc

hc

tcþ1

hcþ1 g

T

ð16Þ

Normally, the crack is regarded to be close when c > c0 , vice versa, and c0 is the threshold value of the bend curvature. In this study, the crack is assumed to open and close when the curvature changes sign and c0 ¼ 0 as the weight of the beam has not been taken into consideration.

1813

Method [17], a popular tool widely used in the nonlinear system analysis, can predict the amplitudes of the superharmonic components and the occurrences of the sub-resonances, this method is essentially an empirical and posterior method which is based on the assumption that the responses of the cracked beam consist of only harmonic components, and cannot provide theoretical explanation to the generation of super-harmonics and sub-resonances. In this section, the two well-known nonlinear phenomena will first be explained by using the NOFRFs on a SDOF bilinear oscillator, then the sensitivity of the NOFRFs to the crack ratio defined as a ¼ d c =d will be numerically studied. 4.1. Super-harmonic components and sub-resonances When a beam vibrates mainly in its first mode and a crack does not affect the mass of the beam, then its response may be approximated using a single degree of freedom bilinear oscillator. Although this situation rarely occurs in practices, as the SDOF bilinear oscillator and Eq. (15) has the same characteristics [13], it is reasonable to use the SDOF bilinear oscillator as a demonstration. The motion equation for a cracked object described by a bilinear oscillator is given by [13] m€y þ c_y þ aky ¼ uðtÞ y P 0 ð17Þ m€y þ c_y þ ky ¼ uðtÞ y < 0 where m and c are the object mass and damping coefficient respectively; yðtÞ is the displacement; k is the stiffness; a is known as the stiffness ratio ð0 6 a 6 1Þ. uðtÞ is the external force exciting the model. Define SðyÞ as the restoring force of a bilinear oscillator as follows: aky y P 0 SðyÞ ¼ ð18Þ ky y<0 Obviously SðyÞ is a piecewise-linear continuous function of displacement y as illustrated in Fig. 6. In mathematics, the Weierstrass Approximation Theorem [24] guarantees that any continuous function on a closed and bounded interval can be uniformly approximated on that interval by a polynomial to any degree of accuracy. Since the restoring force SðyÞ is a continuous function of displacement y, it can be

4. NOFRF analysis of nonlinear effects induced by cracks The bilinear form motion Eq. (15) of a cracked beam implies that a cracked beam is a nonlinear system. There are two well-known nonlinear phenomena observed on cracked beams: one is when the cracked beam is subjected to a harmonic input, super-harmonic components can appear in the output response, and the other is the sub-resonance where if the excitation frequency happens to be half of any eigenfrequency of a cracked object, vibrations often become significant. Unfortunately, the conventional Frequency Response Function cannot explain the two phenomena. In addition, although the Harmonic Balance

a

b

S(y) Order = 4

S(y)

y

y

ky Fig. 6. The restoring force of a bilinear oscillator and its polynomial approximation ða ¼ 0:8Þ.

1814

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

well approximated by a polynomial. Fig. 6 gives the result of using a fourth order polynomial to approximate SðyÞ where the stiffness ratio a is taken as 0.8. It can be seen that a fourth order polynomial can fit SðyÞ very well. If using a polynomial P ðyÞ to represent for SðyÞ and ignoring the tiny approximation error, the SDOF model Eq. (17) can be rewritten as m€y þ c_y þ P ðyÞ ¼ uðtÞ

P ðyÞ ¼

ki ky i

ð20Þ

i¼1

N is the order of the approximating polynomial, and ki k, i ¼ 1; . . . ; N are the polynomial coefficients. The model described by Eq. (19) is the polynomial-type nonlinear system where the term k1 ky represents the linear part and the other high order terms represent the nonlinear part. For the bilinear oscillator model, the polynomial coefficients are determined by the stiffness ratio a. Table 1 shows the results of using a fourth order polynomial to approximate the bilinear oscillator with different stiffness ratios. It is known from Table 1 that all coefficients, apart from k1 , will increase with a decrease of a. This means that the nonlinear strength of the bilinear oscillator will increase with the decrease of a. It is worth to noting that except for k1 , the values of k2 . . . and kN also depend on the range of x which the polynomial approximation is defined. In the case shown in Table 1, this range of x is [1, 1]. For the free undamped vibration of bilinear oscillators, its effective natural frequency xB is given by [14] xB ¼ 2x0 x1 =ðx0 þ x1 Þ where pffiffiffiffiffiffiffiffiffi x0 ¼ k=m and

a

xL =x0

xB =x0

jxB  xL j=xB (%)

1.00 0.95 0.90 0.85 0.80

1.0000 0.9874 0.9747 0.9618 0.9487

1.0000 0.9872 0.9737 0.9594 0.9443

0.0000 0.0203 0.1027 0.2501 0.4660

ð19Þ

where N X

Table 2 Comparison between xL and xB

fies using a polynomial-type nonlinear model to describe a bilinear oscillator. For the polynomial-type nonlinear system (19), the relationship between the displacement yðtÞ and the external force uðtÞ can be described by the Volterra series Eq. (1) in the time domain and by Eq. (6) using NOFRFs in the frequency domain. When the input uðtÞ is a sinusoidal force uðtÞ ¼ A cosðxF t þ bÞ

Peng et al. [19] have shown that the frequency components of the nth order output Y n ðjxÞ of the nonlinear system can be determined as Xn ¼ fðn þ 2kÞxF ; k ¼ 0; 1; . . . ; ng

ð26Þ

and the frequency components of the system output Y ðjxÞ can be determined as X¼

N [

Xn ¼ fkxF ; k ¼ N ; . . . ; 1; 0; 1; . . . ; N g

ð27Þ

n¼1

and the kth super-harmonic component of the system output can be expressed as

ð21Þ Y ðj kxF Þ ¼

x1 ¼

ð25Þ

½ðN kþ1Þ=2 X

Gkþ2ðn1Þ ðj kxF ÞAkþ2ðn1Þ ðj kxF Þ

n¼1

pffiffiffiffiffiffiffiffiffiffiffi ak=m

ð22Þ

ðk ¼ 0; 1; . . . ; N Þ

ð28Þ

where [] stands for taking the integer part, and

Therefore pffiffiffi rffiffiffiffi pffiffiffi 2 a k 2 a pffiffiffi pffiffiffi x0 xB ¼ ¼ ð1 þ aÞ m ð1 þ aÞ

ð23Þ

An ðjðn þ 2kÞxF Þ ¼

1 n! n jAj ejðnþ2kÞb n 2 k!ðn  kÞ! k

For the polynomial-type nonlinear system (19), the natural frequency of the linear part can be defined as pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffi ð24Þ xL ¼ k1 k=m ¼ k1 x0 Table 2 shows a comparison between xL and xB under different stiffness ratios. It can be seen that the xL is a good approximation of xB . To certain extent, this further justiTable 1 The polynomial approximation result for a bilinear oscillator a

k1

k2

k3

k4

1.00 0.95 0.90 0.85 0.80

1.0000 0.9750 0.9500 0.9250 0.9000

0.0000 0.0409 0.0818 0.1228 0.1637

0.0000 0.0000 0.0000 0.0000 0.0000

0.0000 0.0204 0.0407 0.0611 0.0814

ð29Þ

nk

zfflfflfflfflfflfflfflfflffl}|fflfflfflfflfflfflfflfflffl{ zfflfflfflfflfflfflfflfflfflfflfflffl}|fflfflfflfflfflfflfflfflfflfflfflffl{ Gn ðjðn þ 2kÞxF Þ ¼ H n ðjxF ; . . . ; jxF ; jxF ; . . . ; jxF Þ ð30Þ Notice H n ðjx1 ; . . . ; jxn Þ is a symmetric function, thus, for the case of sinusoidal inputs, Gn ðjxÞ over the nth order output frequency range Xn ¼ fðn þ 2kÞxF ; k ¼ 0; 1; . . . ; ng is equal to the H n ðjx1 ; . . . ; jxn Þ evaluated at x1 ¼    ¼ xk ¼ xF , xkþ1 ¼    ¼ xn ¼ xF , (k ¼ 0; . . . ; nÞ. Clearly Eq. (28) shows how the super-harmonic components are generated. For example, this equation shows that it is the 2nd, 4th . . . NOFRFs that shift the input energy from the input frequency to the second harmonic component, and it is the 3rd, 5th . . . NOFRFs transform the input energy from the input frequency to the third harmonic component. Therefore, from Eq. (28), the NOFRFs produce a straightforward explanation to how the super-har-

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

monic components are generated for bilinear oscillators subjected to a harmonic excitation. In addition, from Eq. (19), the GFRF up to 3th order can be calculated recursively using the algorithm by Billings and Peyton Jones [25] to produce the results below H 1 ðjxÞ ¼

1 mx2 þ j cx þ k1 k

H 2 ðjx1 ; jx2 Þ ¼ k2 kH 1 ðjx1 ÞH 1 ðjx2 ÞH 1 ðjx1 þ jx2 Þ

ð31Þ

cracked beam. Therefore, the nonlinear phenomena with cracked beams including the generation of super-harmonic components and the sub-resonance can be explained from the above NOFRF based analysis. Fig. 7a and b shows the amplitudes of the first harmonics and the third harmonics obtained by using Runge–Kutta method on a specific bilinear as follows:

ð32Þ

m ¼ 1 kg;

2k2 k ½H 1 ðjx1 ÞH 2 ðjx2 ; jx3 Þ 3 þ H 1 ðjx2 ÞH 2 ðjx1 ; jx3 Þ

H 3 ðjx1 ; jx2 ; jx3 Þ ¼ 

þ H 1 ðjx3 ÞH 2 ðjx1 ; jx2 ÞH 1 ðjx1 þ jx2 þ jx3 Þ

ð33Þ

Substituting above three equations into Eq. (30) follows: G1 ðjxF Þ ¼

mx2F

1 þ j cxF þ k1 k

ð34Þ

G2 ðj2xF Þ ¼ k2 kG21 ðjxF ÞG1 ðj2xF Þ

ð35Þ

2

G3 ðj3xF Þ ¼ 2ðk2 kÞ G31 ðjxF ÞG1 ðj2xF ÞG1 ðj3xF Þ

ð36Þ

Because jG1 ðjxF Þj reaches a maximum when xF ¼ xL , Eq. (35) obviously indicates that jG2 ðj2xF Þj can reach a maximum at xF ¼ xL and xF ¼ xL =2. Similarly, Eq. (36) indicates that jG3 ðj3xF Þj can reach a maximum at xF ¼ xL , xF ¼ xL =2 and xF ¼ xL =3. Therefore, xL and xL =2 can be regarded the resonant frequencies of jG2 ðj2xF Þj, and xL , xL =2 and xL =3 are the resonant frequencies of jG3 ðj3xF Þj. At the resonance frequencies, the energy transmission is efficient, and then considerable input signal energy may be transferred by the NOFRFs G2 ðj2xF Þ and G3 ðj3xF Þ from the driving frequency to the second and third order harmonic components respectively, and consequently the second and third order harmonic components can reach a maximum at these resonant frequencies. This is just the sub-resonance phenomenon observed on the

k ¼ 3:55  104 N s=m; c ¼ 23:5619 N=m

and the stiffness ratio a is 0.8. The amplitude of the sinusoidal external force is 1 and frequency xF changes within the range 0 6 xF 6 1:2x0 . Clearly, the second harmonic component reaches resonances around the frequencies 0:93x0 and 0:47x0 respectively which, as shown in Table 2, are the approximate values of xL and 1=2xL under a ¼ 0:8. Similarly, it can be seen that the third harmonic component reaches resonances around the frequencies xL , 1=2xL and 1=3xL . Therefore, the numerical results shown in Fig. 7 confirm above NOFRF based analysis results that the second and third order harmonic components can reach maxima at the frequencies xF ¼ xL , xF ¼ xL =2 and xF ¼ xL =3. Also these numerical results validate the explanations presented by NOFRFs to the nonlinear phenomena of cracked beams. 4.2. NOFRFs analysis of a cracked beam The cracked beam considered in this section is the endfixed beam shown in Fig. 3 whose length is 1 m long. The width of the beam is 0.02 m, the cross-sectional area is 0.004 m2, and the mass density is 1000 kg/m3. The elastic modulus is 100 GPa. The position of the single crack is xc ¼ 0:5 m, and the crack ratio a ¼ d c =d is varying from 0.1 to 0.4 in step of 0.1. The FE model of this beam is constructed using five Euler–Bernoulli beam elements, as shown in Fig. 8. Ignoring the fixed node, the FE model of this beam has 5 nodes and 10 degrees of freedom. Using

-5

x 10

x 10

-7

(0.47, 1.98e-005)

1.8

1815

(0.31, 6.85e-007) 6

1.6

5

Amplitude

Amplitude

1.4 1.2 1 0.8 (0.93, 5.25e-006)

0.6 0.4

(0.47, 4.17e-007)

4 3 2

(0.93, 9.09e-008)

1

0.2 0.2

0.4

0.6

0.8

ω F /ω 0

1

1.2

(a) the second harmonic component

0.2

0.4

0.6

0.8

1

1.2

ω F /ω 0 (b) the third harmonic component

Fig. 7. The second and third harmonic components obtained using Runge–Kutta method.

1816

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

1

2 Node 1

3 Node2

Node 3

4

5

Node 4

Node 5

Fig. 8. The FE model for the cracked beam.

this FE model, the modal frequencies of the cracked beam are calculated, and the results are shown in Table 3. It can be seen from Table 3 that the crack does not cause significant change in the modal frequencies. It implies that the modal frequencies are not sensitive to the presence of crack. The insensitivity of the modal frequencies to the crack is mainly due to the fact that the nonlinearity introduced by the crack is often weak, and therefore many of the common testing techniques will tend to linearize the response [26]. Sinusoidal forcing will tend to emphasize the nonlinearity, and crack detection methods based on detecting harmonics of the forced response have been proposed [27]. In this study, the sensitivity of the NOFRFs to the depth of the crack is studied. Using the algorithm in [18], the NOFRFs are estimated from the force responses of the cracked beam subjected to harmonic excitations with two level amplitudes A ¼ 1:0 and A ¼ 1:2, respectively. The frequency xF of the sinusoidal force uðtÞ varies from 1 to

Table 3 The natural frequencies (Hz) of the beam with one crack

1 2 3 4 5 6 7 8 9 10

Natural frequency (103) (Hz) d c =d ¼ 0:0

d c =d ¼ 0:1

d c =d ¼ 0:2

d c =d ¼ 0:3

d c =d ¼ 0:4

0.0100 0.0594 0.1580 0.2924 0.4365 0.9385 1.1151 1.3198 1.4839 1.5789

0.0099 0.0592 0.1579 0.2918 0.4355 0.9367 1.1141 1.3189 1.4836 1.5788

0.0098 0.0590 0.1578 0.2913 0.4343 0.9353 1.1133 1.3183 1.4833 1.5787

0.0097 0.0588 0.1578 0.2909 0.4338 0.9342 1.1127 1.3178 1.4831 1.5786

0.0097 0.0587 0.1577 0.2907 0.4332 0.9333 1.1123 1.3174 1.4829 1.5786

-4

4

x 10

Transient state

Steady state

2

Amplitude

Mode

200 Hz in step of 5 Hz. The simulation studies were conducted using a fourth order Runge–Kutta method to obtain the forced response of the system. It is worth noting here that, as the Volterra series is only valid for the description of the steady state responses of the nonlinear systems, therefore only the responses at the steady state are used to estimate the NOFRFs. A response passing from the transient state to the steady status is shown in Fig. 9, which is the time history of the transverse displacement of the node 5 of the cracked beam with crack ratio d c =d ¼ 0:4 and excitation frequency xF ¼ 20 Hz. Fig. 10 shows the NOFRFs of the transverse displacements of node 5 of the cracked beam with different crack ratios. It can be seen that G1 ðjxF Þ reaches maxima at frequencies 10 Hz, 60 Hz and 160 Hz which are the first three resonant frequencies in Table 3, clearly there is only tiny and negligible changes appearing on the NOFRF G1 ðjxF Þ with the increasing of the crack ratios, therefore the NOFRF G1 ðjxF Þ is not sensitive to the presence of crack. Fig. 10b shows that the NOFRF G2 ðj2xF Þ has four peaks at the frequencies 5 Hz, 30 Hz, 60 Hz and 150 Hz, and the resonances at the frequencies 5 Hz and 30 Hz are just the so-called sub-resonances. Similarly, it can be seen that the NOFRF G3 ðj3xF Þ reaches maxima at the frequencies 5 Hz, 30 Hz, 60 Hz, 95 Hz, 130 Hz and 150 Hz, among which the frequency of 95 Hz is about the one-third of the fourth natural frequency in Table 3. These analysis results are consistent with the analysis results of the sub-resonance using the bilinear model. In addition, as Fig. 10b and c show, the NOFRFs G2 ðj2xF Þ and G3 ðj3xF Þ are very sensitive to the presence of crack, especially at the resonant frequencies such as at the frequencies 5 Hz and 60 Hz where the G2 ðj2xF Þ and G3 ðj3xF Þ increases sharply with the increasing of the crack ratio. The results imply that the NOFRFs can be a sensitive indicator of the appearances and size of cracks in beams. Consequently, compared to traditional output spectrum based crack detection methods, the NOFRFs based method proposed in the present study should be more effective in practical applications. Actually, a Volterra series based method has been used to analyze the vibration of a cracked beam in [17] where the higher order transform function (HOTF) of a cracked cantilevered beam was estimated. The HOTF was defined

0 -2 -4

0

5

10

15

20

25

30

Time / sec Fig. 9. A response from transient state to steady state.

35

40

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

a -2

1

Freq (Hz)

10

-3

Amplitude

1817

10

2

1

10

2

60

3

160

-4

10

3 -5

10

20

40

60

80

100 120 140 160 180

200

Frequency (Hz)

b

-3

10

1

-4

10

-5

1

5

2

30

3

60

4

150

-4

1

Freq (Hz)

10

-5

10

10

-6

Amplitude

Amplitude

c

Freq (Hz)

3

2 -6

10

3

10

-7

10

1

5

2

30

3

60

4

95

5

130

6

150

2

4 -8

4

10

-7

10

6 5

-9

10 -8

10

50

100

150

200

50

Frequency (Hz)

100

150

200

Frequency (Hz)

Fig. 10. NOFRFs at different crack ratios [(a) G1 ðjxF Þ, (b) G2 ðj2xF Þ, (c) G3 ðj3xF Þ].

as the ratio between the output spectrum Y ðjxÞ and U n ðjxÞ for a particular n of interest and is also based on the Volterra series of nonlinear systems. Compared to NOFRFs, HOTF is not a theoretically well established concept although, under certain conditions, NOFRFs and HOTF can be related to each other. According to Eq. (28), the harmonic components of the frequency response of the cracked beam can be rewritten as

In such cases, the NOFRFs Gn ðjxÞ degenerates to the HOTF proposed in [17] where it was validated that the HOTFs are an accurate and extremely sensitive indicator of nonlinear behavior of a system and the evolution of HOTFs for the crack size is very useful for structural damage assessment. This observation and comparison with the HOTFs confirms the significance of applying the well defined and more general concept of NOFRFs to the detection of cracks in engineering structures.

Y ðj kxF Þ ¼ Gk ðj kxF ÞAk ðj kxF Þ þ

½ðN kþ1Þ=2 X 

5. Conclusions

Gkþ2ðn1Þ ðj kxF ÞAkþ2ðn1Þ ðj kxF Þ

n¼2

ðk ¼ 1; . . . ; N Þ

ð37Þ

If the super-harmonic Y ðj kxF Þ is mainly determined by Gk ðj kxF ÞAk ðj kxF Þ, and the contribution of the second term in Eq. (37) to the super-harmonic Y ðj kxF Þ can be ignored, then Y ðj kxF Þ can be approximated by Y ðj kxF Þ  Gk ðj kxF ÞAk ðj kxF Þ

ðk ¼ 1; . . . ; N Þ

ð38Þ

This paper presents an analysis of the dynamic behavior of beams with a closing crack in the frequency domain using the NOFRF concept recently developed by the authors. A simple SDOF bilinear oscillator was first employed to explain the well-known nonlinear phenomena occurring in cracked beams subjected to sinusoidal excitations, including the generation of super-harmonic components and the sub-resonances. Then the new NOFRF concept is used to analyze the crack induced nonlinear

1818

Z.K. Peng et al. / Computers and Structures 86 (2008) 1809–1818

response of the beam. The NOFRF based analysis shows that the high order NOFRFs are extremely sensitive to the appearance of crack in the beam, and can therefore be used as crack damage indicators to indicate the existence and size of cracks. This research study establishes an important basis for the application of the NOFRF concept in fault diagnosis of mechanical structures. Acknowledgement The authors gratefully acknowledge the support of the Engineering and Physical Science Research Council, UK and the Natural Science Foundation of China (No. 50425516) for this work. References [1] Chinchalkar S. Determination of crack location in beams using natural frequencies. J Sound Vib 2001;247:417–29. [2] Kim JT, Ryu YS, Cho HM, Stubbs N. Damage identification in beam-type structures: frequency-based method vs mode-shape-based method. Eng Struct 2003;25:57–67. [3] Panteliou SD, Chondros TG, Argyrakis VC, Dimarogonas AD. Damping factor as an indicator of crack severity. J Sound Vib 2001;241:235–45. [4] Bovsunovsky AP, Surace C. Considerations regarding superharmonic vibrations of a cracked beam and the variation in damping caused by the presence of the crack. J Sound Vib 2005;288:865–86. [5] Tsyfansky SL, Beresnevich VI. Detection of fatigue cracks in flexible geometrically non-linear bars by vibration monitoring. J Sound Vib 1998;213:159–68. [6] Sundermeyer JN, Weaver RL. On crack identification and characterization in a beam by nonlinear vibration analysis. J Sound Vib 1995;183:857–71. [7] Saavedra PN, Cuitino LA. Crack detection and vibration behavior of cracked beams. Comput Struct 2001;79:1451–9. [8] Darpe AK, Gupta K, Chawla A. Transient response and breathing behaviour of a cracked Jeffcott rotor. J Sound Vib 2004;272:207–43.

[9] Iman I, Azzaro SH, Bankert RJ. Development of an on-line rotor crack detection and monitoring system. J Vib Acoust Stress Reliab Des 1989;111:241–50. [10] Wauer J. On the dynamics of cracked rotors: a literature survey. Appl Mech Rev 1990;43:13–7. [11] Dimarogonas AD. Vibration of cracked structures: a state of the art review. Eng Fract Mech 1996;55:831–57. [12] Sabnavis G, Kirk RG, Kasarda M. Cracked shaft detection and diagnostics: a literature review. Shock Vib Digest 2004;36:287–96. [13] Friswell MI, Penny JET. Crack modeling for structural health monitoring. Struct Health Monitor 2002;1:139–48. [14] Sinha JK, Friswell MI, Edwards S. Simplified models for the location of cracks in beam structures using measured vibration data. J Sound Vib 2002;251:13–38. [15] Chu YC, Shen MHH. Analysis of forced bilinear oscillators and the application to cracked beam dynamics. AIAA J 1992;30:2512–9. [16] Chati M, Rand R, Mukherjee S. Modal analysis of a cracked beam. J Sound Vib 1997;207:249–70. [17] Ruotolo R, Surace C, Crespo P, Storer D. Harmonic analysis of the vibration of a cantilevered beam with closing crack. Comput Struct 1996;61:1057–74. [18] Lang ZQ, Billings SA. Energy transfer properties of nonlinear systems in the frequency domain. Int J Control 2005;78:354–62. [19] Peng ZK, Lang ZQ, Billings SA. Resonances and resonant frequencies for a class of nonlinear systems. J Sound Vib 2007;300:993–1014. [20] Lang ZQ, Billings SA. Output frequency characteristics of nonlinear system. Int J Control 1996;64:1049–67. [21] Christides S, Barr ADS. One-dimensional theory of cracked Bernoulli–Euler Beams. Int J Mech Sci 1984;26:639–48. [22] Shen MHH, Pierre C. Natural modes of Bernoulli–Euler beams with symmetric cracks. J Sound Vib 1990;138:115–34. [23] Kwon Young W, Bang Hyochoong. The finite element method using Matlab. 2nd ed. Boca Raton, FL: CRC Press, Inc.; 2000. [24] Jeffreys H, Jeffreys BS. Methods of mathematical physics. 3rd ed. Cambridge, England: Cambridge University Press; 1988. [25] Billings SA, Peyton Jones JC. Mapping nonlinear integro-differential equations into the frequency domain. Int J Control 1990;52:863–79. [26] Worden K, Tomlinson GR. Nonlinearity in structural dynamics: detection, identification and modeling. Bristol, UK: IoP Publishing; 2001. [27] Tsyfanskii SL, Magone MA, Ozhiganov VM. Using nonlinear effects to detect cracks in the rod elements of structures. The Sov J Nondestr Test 1985;21:224–9.