Comparison of vibration and wave propagation approaches applied to assess damage influence on the behavior of Euler–Bernoulli beams

Computers and Structures 89 (2011) 1820–1828

Contents lists available at ScienceDirect

Computers and Structures journal homepage: www.elsevier.com/locate/compstruc

Comparison of vibration and wave propagation approaches applied to assess damage influence on the behavior of Euler–Bernoulli beams R.A. Tenenbaum ⇑, L.T. Stutz, K.M. Fernandes Dynamics, Acoustics and Vibration Lab., Graduate Program of Computer Modeling, IPRJ, State University of Rio de Janeiro, Brazil

a r t i c l e

i n f o

Article history: Received 10 November 2008 Accepted 4 June 2010 Available online 5 November 2010 Keywords: Structural damage Structural vibration Wave propagation Flexibility matrix Euler–Bernoulli beams

a b s t r a c t In the present work, the damage influence on the vibrational behavior and on the wave propagation issues of a slender Euler–Bernoulli beam is investigated. In the vibration framework, the damage is assessed in the modal space by considering changes in the reduced flexibility matrix of the structure, which is only related to the measured degrees of freedom and may be accurately estimated from a few of the lower frequency modes in a modal test. In the wave propagation framework, on the other hand, the presence of damage is perceived in the time domain by an early echo output generated by the inhomogeneity within the beam in an impact test. A slender aluminium beam with distinct imposed damage scenarios is considered for numerical analysis. The modal properties required to estimate the reduced flexibility matrix of the beam and its impulse response are assumed to be available from modal and impact tests, respectively. In the numerical analysis performed, both approaches showed sensitivities to damage that enable them to be applied for damage identification purposes. Ó 2010 Civil Comp Ltd and Elsevier Ltd. All rights reserved.

1. Introduction Damage identification is an essential issue for determining safety reliability and remaining lifetime of aerospace, civil and mechanical structures. The technological and scientific challenges posed by damage identification problems yielded a great research activity, within the engineering community, on this subject [1,2]. Different nondestructive damage identification approaches are proposed in the literature [3,1]. These ones, encompassing deterministic or statistical perspectives, consider different types of data (modal parameters, time series, frequency responses), several forms of excitations and experimental setups, distinct mathematical formulations and numerical algorithms for solving the corresponding inverse problem. Most of the proposed approaches are built on the vibrational behavior of the structure, more specifically on the traditional modal analysis. Many vibration based damage identification approaches presented in the literature are built on the general framework of Finite Element Model (FEM) updating methods [4,2]. These methods are intended for identifying structural damage through determination of changes in some parameters of a FEM of the structure. Hence, the damage identification problem may be cast as a minimization one and a set of parameters, supposed to describe the damage scenario, is sought in order to minimize an arbitrarily defined error function. This error may be defined, for instance, as the difference ⇑ Corresponding author. E-mail address: [email protected] (R.A. Tenenbaum).

between some outputs or matrices of the FEM of the undamaged structure and the corresponding ones obtained from a modal testing on the damaged structure [5,6]. The basic idea of the vibration based damage identification approach is that the modal properties of the structure (frequencies, mode-shapes and modal damping) are functions of the physical ones (mass, stiffness and damping) and changes in these last due to damage will be reflected in the modal characteristics, which can be measured and used to infer the damage. Although these methods were proven to succeed in many damage identification problems, the high frequency effects of small defects, such as cracks, may be slightly or even not reflected in the required modal properties of the structure, making the damage identification a very difficult task. This great difficulty is not present in damage identification approaches built on wave propagation, since these approaches are highly sensitive to changes in local dynamic impedance [7–9], such as that caused by small defects. Although these approaches are much more uncommon in the literature than the vibration ones, they have the advantage of being fast techniques with good accuracy [10]. Wave propagation has also been considered as an useful model for damage assessment of structural elements, using Lamb waves [11,12]. Some recent works report its application in hysteretic materials [13] and in nonlinear plates [14]. Succeeded applications are reported in the fields of geophysics [15,16], medical ultrasonics [17], fatigue testing [18,19], non-destructive testing [20], and in the evaluation of the integrity of drilled piles [21–23].

0045-7949/$ - see front matter Ó 2010 Civil Comp Ltd and Elsevier Ltd. All rights reserved. doi:10.1016/j.compstruc.2010.10.006

1821

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828

In the problem under concern, in the wave propagation approach, the damage is considered as an inhomogeneity of the beam, and the basic idea is as follows. In a wave propagation test, an incident longitudinal stress is applied at the boundary of the structure and, consequently, a progressive wave propagates along it. This progressive wave is reflected whenever it encounters a local change of impedance (for instance, a boundary condition or a damage), generating a regressive wave. Finally, the generated echo is measured at the sensor locations and used to infer about the position and shape of the damage. The main goal of this research is to develop damage identification techniques based on vibration and wave propagation approaches and compare them in terms of practice, accuracy and robustness. Although having the inverse problem in mind, the concern here is to comparatively solve the analysis problem, i.e. to assess the influence of the damage in the vibration and wave propagation behavior, and study the sensitivity of each method to structural damage. 2. Problem modeling

In the vibration framework, in order to describe the damage influence on the vibrational behavior of a system, its structural

36bi þ 24bq þ 36bj

IðxÞ ¼ bðxÞI0

ð3Þ

and, hence, the damage field reads as

bðxÞ ¼

IðxÞ : I0

ð4Þ

Considering the finite bidimensional Euler–Bernoulli beam element depicted in Fig. 1 and assuming that the damage field b(x), within this element, is interpolated by three nodal values (bi, bq, bj) and classical Lagrangian piecewise linear shape functions, according to Eq. (2), the element stiffness matrix reads as

  3 Ke ¼ E0 I0 =8le A;

2.1. Vibration approach

2

Element Model order, which depends only on the discretization of the displacement field. Besides, from Eq. (2), one may conclude that the damage may be interpreted as a change in the Young modulus E(x) or in the area moment of inertia I(x) or in both, through the bend stiffness E(x)I(x). For example, if one considers that the Young modulus E0 is uniform along the beam and the damage is only due to changes in its geometric properties, the area moment of inertia of the damaged structure, in the spatial domain, is given by

2le ð13bi þ 6bq þ 5bj Þ

ð5Þ

with the matrix A given by

36bi  24bq  36bj

2le ð5bi þ 6bq þ 13bj Þ

3

6 7 2 6 2le ð13bi þ 6bq þ 5bj Þ l2e ð19bi þ 10bq þ 3bj Þ 2le ð13bi þ 6bq þ 5bj Þ le ð7bi þ 2bq þ 7bj Þ 7 6 7; 6 36b  24b  36b 2le ð13b þ 6b þ 5b Þ 36bi þ 24bq þ 36bj 2le ð5bi þ 6bq þ 13bj Þ 7 i q j i q j 4 5 2le ð5bi þ 6bq þ 13bj Þ

2

le ð7bi þ 2bq þ 7bj Þ

2le ð5bi þ 6bq þ 13bj Þ

integrity is considered to be continuously described by a scalar attribute b(x), named cohesion parameter, defined over the entire spatial domain x 2 (0, l). This parameter is related with the links among material points and can be interpreted as a measure of the local cohesion state of the material and, therefore, is supposed to describe the state of damage within a mechanical system. The cohesion parameter possesses as extreme values:

b ¼ 1 ðmaterial properties preservedÞ; b ¼ 0 ðlocal ruptureÞ:

ð1Þ

The modeling considers that only the elastic terms contained in its basic equations are affected by the damage field [5,6]. Therefore, in the context of an Euler–Bernoulli beam, the stiffness matrix obtained by a spatial discretization using the FEM reads as

Kðbh Þ ¼

Z

l

bh ðxÞE0 I0 HT ðxÞHðxÞdx;

ð2Þ

0

where x is the structural elastic domain, l is the beam length, E0 and I0 are, respectively, the nominal Young modulus and area moment of inertia, H is the standard discretized differential operator, bh(x) is the spatial discretization of the damage field b(x), and T means transpose. It is worth noting that the capability of the present method to represent the damage is greatly supported by the discretization bh of the damage field, which is not necessarily coincident with the one used for the displacement field. Hence, different spatial discretizations for the damage and displacement fields may be defined. This favors a more detailed description of the damage within an element without requiring a higher Finite

ð6Þ

2

le ð3bi þ 10bq þ 19bj Þ

where le is the length of the element, bi and bj are, respectively, the nodal cohesion parameters associated with nodes i and j at the ends of the finite element, and bq is the one associated with the middle node q. In Fig. 1, u and h are, respectively, the translational and rotational degrees of freedom of the finite beam element. In the present work, the influence of the damage on the dynamic behavior of a system will be assessed by considering changes in its flexibility matrix [2,24,25]. The flexibility matrix of a system is defined as the inverse of the stiffness one and, therefore, it is a function of the cohesion parameters describing the damage scenario. The flexibility matrix G can be directly expressed as a function of the system modal parameters as [24]

G¼

n X 1 i¼1

x2i

/ðiÞ  /ðiÞ ;

ð7Þ

where n is the number of degrees of freedom (DOF) of the FEM, xi and /(i) are the ith system natural frequency and mass-normalized mode-shape, respectively, and  denotes tensor product. It is worth mentioning that the flexibility matrix does not depend on the damping characteristics of the structure, which is generally difficult

ui

uj j

i q i

j

le Fig. 1. Bidimensional Euler–Bernoulli beam element.

1822

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828

to model and, therefore, represents an important source of uncertainty in the results. As one can see from Eq. (7), the higher the frequency, the smaller the modal contribution to the flexibility matrix. Hence, a good estimate of the flexibility matrix can be obtained from a few of the lower frequency modes. However, in a modal testing, the number of DOF at which the mode shapes are sampled is typically much smaller than that of the FEM. Therefore, if the mode-shapes are sampled at only m DOF, the reduced flexibility matrix is obtained as

Gmm ¼

n X 1 i¼1

x2i

ðiÞ ðiÞ /m  /m ;

ð8Þ

1

cðxÞ ¼ ðEðxÞ=qðxÞÞ2

ðiÞ /m

where is the ith mass-normalized mode shape sampled at a subset of m DOF. Considering Eq. (8), an estimate GE for the flexibility matrix Gmm may be obtained from a modal testing as

GE ¼

nE X 1

ðiÞ /E 2 i;E

x

i¼1



ðiÞ /E ;

ð9Þ

where nE is the number of modes obtained from a modal test and xi,E and /EðiÞ are, respectively, the ith experimental natural frequency and the mass-normalized mode-shape, sampled at the same subset of DOF, as its analytical counterpart /ðiÞ m. Aiming at obtaining a relation between the analytical flexibility matrix Gmm and a m  m matrix which contains information about all the stiffness properties of the structure, the original stiffness matrix K should be reorganized according to the Guyan’s reduction method as follows [26]. Partitioning the n DOF of the system into m measured DOF and o omitted ones, the ith mode shape casts as

" /ðiÞ ¼

ðiÞ /m /oðiÞ



Kmo

KTmo

Koo

 :

ð11Þ

It can be shown [27], that the analytical flexibility matrix Gmm is equal to the inverse of the statically condensed system stiffness matrix K with respect to the same set of DOF, viz.

Gmm ¼ K

h

¼ Kmm 

T Kmo K1 oo Kmo

i1

:

ð12Þ

ð14Þ

where r(x, t) is the longitudinal stress field, u(x, t) is the longitudinal displacement field, the independent variable t stands for the time, and the subscript denotes partial derivative. The linear momentum balance, in the other hand, is given as

ðrðx; tÞAðxÞÞx ¼ qðxÞAðxÞv t ðx; tÞ;

ð15Þ

where vt(x, t) is the time partial derivative of the particle velocity field in the x-direction. Eqs. (14) and (15), when combined, provide, after some algebraic manipulation, a second-order differential equation for the longitudinal stress, which takes the form [30]



rtt ðx; tÞ ¼ c2 ðxÞ rxx ðx; tÞ þ 

ð10Þ

Kmm

1

rðx; tÞ ¼ EðxÞux ðx; tÞ;

þ

:

ð13Þ

is also x-dependent. The linear constitutive equation is given as

#

According to this partition, the stiffness matrix can be partitioned as

K¼

will be sought. The starting point will be the classical characteristic variables formulation [29], leading to a compact form for the wave propagation equation. This equation and its algebraic solution for noiseless data can be applied as a first approach to model the Euler–Bernoulli beam with an inclusion or obstacle in it, as a crosssection variation, a crack, a bubble, or whatever else characterizing a damage scenario. Let us consider an essentially one-dimensional solid linear elastic beam with its mechanical properties varying along the space coordinate x: the cross-section area A(x); the density q(x); and the elastic modulus E(x). As a consequence, the wave propagation speed,

 0  A ðxÞ q0 ðxÞ  rx ðx; tÞ AðxÞ qðxÞ !

0 A0 ðxÞ qðxÞrðx; tÞ qðxÞAðxÞ

¼ 0;

ð16Þ

where the prime denotes total derivative with respect to the function argument (x, in this case). In the following, the arguments of the functions and parameters will be dropped out. Since the hyperbolic differential Eq. (16) has no analytical closed-form solution, another approach, based on the characteristic variables will be adopted. The idea is to substitute the second-order equation by a system of two first-order equations. Eqs. (14) and (15) can be expressed as:

rt ¼ E v x ; rx þ ðA0 =AÞr ¼ qv t :

ð17Þ

Now, a new independent variable, s, called the travel time coordinate, is defined as

Hence, as it has been shown in Eq. (12), the analytical flexibility matrix Gmm, only related to the measured DOF, contains information about all the stiffness properties of the system and, hence, about all the structural damage scenario. Therefore, the damage state of a structure may be inferred from the difference between the analytical flexibility matrix Gmm of the nominal structure, without damage, and the experimental one, GE, obtained from a modal testing on the supposed damaged structure.

where g is a dummy integration variable and the wave propagation speed is defined as in Eq. (13). Introducing the medium acoustic (mechanical) impedance,

2.2. Wave propagation approach

zðxÞ ¼ qðxÞcðxÞ;

The scattering of acoustic waves by one-dimensional inhomogeneous media can be applied as an alternative approach to the damage identification problem. Since the damage may be considered as an inhomogeneity, the complete nonuniform wave equation can be taken into account [28]. The time-domain solution for the acoustic scattering problem in one-dimensional media, probed with an arbitrarily shaped pulse

sðxÞ ¼

Z

x

0

dg ; cðgÞ

ð18Þ

ð19Þ

one can define the characteristic variables as:

r ¼ t þ s; s ¼ t  s:

ð20Þ

Note that in the plane (s, t) the characteristic curves become straight lines r = constant and s = constant, making angles of ± p/ 4 rad with the coordinate axes, as seen in Fig. 2.

1823

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828

t s

progressive wave component. Analogously, the echo observed at x = 0, due to the damage, will be the output signal g(t) = V(s, s), a regressive wave component. Assuming now the beam under study as a Goupillaud layered medium [32], i.e., a sectionally homogeneous one with equal travel time layers of Ds = Dt, the discrete incoming pulse is

r

F j ¼ f ð2ðj  1ÞDtÞ and the discrete outgoing echo will be

~

Gj ¼ gð2jDtÞ;

Fig. 2. Characteristic curves in the plane s, t.

In the new plane (r, s), Eq. (17) become then, after some algebraic manipulation:

 0 cA rr  zv r þ r ¼ 0; 2A  0 cA rs þ zv s  r ¼ 0: 2A

Gj ¼

j X

Rn þ

n2 X

! Q pn

F jnþ1 ;

ð29Þ

p¼1

n¼1

" Q pn

¼ Rnp

X l Q pn1  Rnp1 Rn1 þ Q n1 Rnp1 l¼1 p1

n ¼ 1; 2; . . . ; N; p ¼ 1; 2; . . . ; n  2:

!# ; ð30Þ

In Eq. (30), Ri stands for the reflection coefficient at the ith layer of the medium, defined as

ð22Þ Ri ¼

Z i  Z i1 ; Z i þ Z i1

i ¼ 1; 2; . . . ; N:

ð31Þ

The mathematical procedure, in the wave propagation approach, consists then in the following steps. Given the medium, with A0 and Z0, it is discretized into N elements. Then, the reflection coefficients are computed by Eq. (31). In the sequence, the polynomials Q pn are calculated from Eq. (30). Finally, the echo output is computed from Eq. (29). 3. Damage influence

ð23Þ

_ will Then, the derivative of Z(s) with respect to s, denoted as Z, be

_ sÞ ¼ dZ ¼ cðzðxÞAðxÞÞ0 : Zð ds

ð24Þ

The system of differential equations describing the wave propagation, Eq. (21), combined with Eqs. (22)–(24), yields:

Z_ U ¼ 0; 4Z Z_ V ¼ 0: Vs  4Z Ur þ

ð25Þ

Eqs. (25) are a compact and uncoupled pair of differential equations that describes the wave propagation phenomenon in a more convenient way. To integrate it, boundary conditions in the (r, s) plane must be furnished, corresponding to the physical situation under consideration. Let us, for instance, consider the probing of a semi-infinite medium, x P 0, by a pulse excitation at x = 0. Assuming also the Sommerfeld radiation hypothesis [31], the boundary conditions can be stated as:

Uðs; sÞ ¼ FðsÞ ¼ f ðtÞ;

ð28Þ

where the polynomials Q pn have the general recursive formula

where A0 and z0 are, respectively, the nominal cross-section area and the acoustic impedance of the beam. The physical meaning of this variable transform is that U(r, s) and V(r, s) represent, respectively, the progressive and regressive waves travelling along the r and s directions. With the aim at obtaining a closed-form solution for Eq. (22), let us introduce the generalized acoustic impedance, given by

ZðsÞ ¼ zðxðsÞÞAðxðsÞÞ:

j ¼ 1; 2; . . . ; N;

where NDt is the total time interval under consideration. It can be shown that Eq. (25), with the boundary conditions given in Eq. (26), has, after the discretization given in Eqs. (27) and (28), the following algebraic solution for the echo [28]

ð21Þ

One of the advantages of the system of Eq. (21), in plane (r, s), is that it is uncoupled. Now, a change of dependent variables will be performed:

sffiffiffiffiffiffiffiffi 1 z0 A Uðr; sÞ ¼ ðrðr; sÞ  zv ðr; sÞÞ; 2 zA0 sffiffiffiffiffiffiffiffi 1 z0 A Vðr; sÞ ¼ ðrðr; sÞ þ zv ðr; sÞÞ; 2 zA0

ð27Þ

Vðr; 0Þ ¼ 0;

ð26Þ

where f(t) is the incident longitudinal stress being applied at the boundary r = s (x = s = 0) and the second equation ensures that there is no disturbance in s 6 0 (t 6 s). Note that f(t), being the longitudinal stress at the physical boundary x = 0, corresponds to U(s, s), a

In this paper, the damage influence on the vibrational behavior and on the wave propagation issues of a slender aluminium Euler– Bernoulli beam is considered. In the numerical analysis that follows, the nominal geometric properties of the beam are: length l0 = 1 m; width w0 = 30 mm; and height h0 = 10 mm (see Fig. 3), yielding the area moment of inertia I0 = 2.5  109 m4. Besides, the assumed nominal material properties are: Young modulus E0 = 7.1  1010 N/m2; and density q0 = 2.7  103 kg/m3, yielding the sound speed c0  5128 m/s. With the aim at simulating a damage state, it is firstly assumed that a triangular prismatic cut, with base length a = 25 mm and depth d = 5 mm, crossing its whole width, is taken away from the middle part of the beam. Fig. 3 depicts the homogeneous and the inhomogeneous (damaged) beam, along with a detail of the imposed damage scenario. Different damage states are further considered in a sensitivity analysis of the approaches. 3.1. Vibration approach In order to assess the influence of the damage on the vibrational behavior of the structure, a FEM comprising 40 bidimensional Euler–Bernoulli beam elements, with length le = 25 mm, is considered. Three nodal values (bi, bq, bj) for the cohesion parameter are used to linearly interpolate the damage field b(x) within a finite element (see Fig. 1), so that the element stiffness matrix Ke is given by Eq. (6). Besides, according to Eq. (4), the damage field b(x) along the beam is given by

1824

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828

1m a

30 mm

0.5m

d

10 mm

500 mm

(a) homogeneous and inhomogeneous beam

(b) Detail of the imposed damage

Fig. 3. Euler–Bernoulli beam and the imposed damage scenario.

bðxÞ ¼ ðhðxÞ=h0 Þ3 ;

ð32Þ

where h(x) and h0 are, respectively, the height of the beam in the damaged and undamaged cases. Fig. 4 depicts the nominal flexibility matrix of the undamaged beam accounting for all modal content of the FEM of structure, i.e., the eighty natural frequencies and mass-normalized mode shapes sampled at the eighty DOF. In order to illustrate the damage influence on the modal properties of the beam, Fig. 5 depicts the error between the flexibility matrix of the damaged structure, considering the imposed damaged scenario, showed in Fig. 3, and the nominal flexibility matrix depicted in Fig. 4. From Fig. 5, one can clearly note the damage influence on the vibrational behavior of the beam reflected on

the flexibility matrix. Besides, changes in the magnitude of the elements of the error matrix may indicate the damaged region within the beam. In practice, however, the error matrix depicted in Fig. 5 is not available for analysis, since only a few of the lower frequency modes may be obtained from a modal test and, besides, the mode-shapes are only sampled at a reduced subset of DOF. Therefore, only an estimate for the reduced flexibility matrix, only related to the measured DOF, may be computed from Eq. (9). In this work, it is assumed that the modal properties required for computing an estimate of the flexibility matrix are available from a modal test. The mode-shapes are assumed to be sampled at only 4 vertical DOF by sensors located at positions: 0.2, 0.4, 0.6 and 0.8 m. Fig. 6 depicts the reduced flexibility Gmm of the damaged beam, i.e., the flexibility matrix only related to the measured DOF and that accounts for all frequency range of the model (the eighty modes, in this case). The relative error in the estimation of the reduced flexibility matrix from a few lower frequency modes may be computed as

Emm ¼

Fig. 4. Flexibility matrix of the undamaged beam.

Gmm  GE ; Gmm

ð33Þ

where Gmm and GE are defined, respectively, in Eqs. (8) and (9). Fig. 7 depicts the relative estimation error given in Eq. (33) when only the first four natural frequencies and mode-shapes, instead of the eighty-ones, are taken into account in the computation of the reduced flexibility matrix of the damaged beam. As one can see from Fig. 7, a good estimate GE of the flexibility matrix may be obtained from a few of the lower frequency modes. Fig. 8 depicts the influence of the damage on the reduced flexibility matrix through the difference between the flexibilities GE of the damaged beam, obtained from a modal test, and the analytical

−4

x 10

5

G

mm

( m/N )

10

0

1 2

4

3 Row

3 2

4 1

Fig. 5. Error between the flexibility matrices of the damaged and undamaged beam.

Column

Fig. 6. Reduced flexibility matrix of the damaged beam.

1825

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828 0.04

Relative Error in G

E

0.03 0.02

0.02

0.01

0.01

0

0

−0.01

−0.01 1

−0.02

2

4 3

3 Row

4

2 1

−0.03

Column

Fig. 7. Relative error in the estimation of the reduced flexibility matrix of the damaged beam.

−0.04 0

50

100

150

200 t(µs)

250

300

350

400

Fig. 9. Echo produced by the damage (impulse response).

As can be seen from Figs. 9 and 10, there is a noticeable effect on the output wave signal generated by the damage.

−6

x 10

10

4. Sensitivity analysis

5

In this section, the sensitivity of the vibration and wave propagation approaches to damage is assessed by considering different damage scenarios. Four cases are considered, as depicted in Fig. 11.

0

G

E

− G

mm

( m/N )

15

−5 1 2 Row

4 3

3 4

2 1

Column

Fig. 8. Error between the reduced flexibility matrices of the damaged and undamaged beam.

reduced flexibility Gmm of the undamaged beam, i.e., the one obtained from the FEM of the undamaged beam. Once again, one can clearly observe the damage influence on the structure’s vibrational behavior. 3.2. Wave propagation approach To test the beam with a wave pulse, an excitation f(t) = d(t), where d(t) stands for the Dirac delta, is longitudinally applied at its free end x = 0. For the homogeneous beam, of course, only the boundary at x = l0 generates an echo. Therefore, if one restrain the time interval to pick up the echo to NDt 2 [0, 2l0/c], no signal at all will be heard at the sensor position (x = 0). On the other hand, the inhomogeneous (damaged) beam generates an output signal inside the time interval NDt. The medium was discretized with 1000 points, corresponding to a time interval Dt = 1.95  107 s. Hence, the damage length was then discretized with 25 points. Fig. 9 shows the result of probing the considered damaged beam with an impulse f(t) = d(t), obtained with Eqs. (26)–(31). Fig. 10(a) shows a detail of the output signal, observed in the interval (185, 215) ls. It is worth noting that the first part of the echo is negative, corresponding to a decrease in the generalized acoustic impedance, while the second part presents the inverse behavior, due to an impedance increase. In Fig. 10(b) a zoom on the output signal is shown, in the interval (200, 215) ls. It can be observed the tail of the echo (after 200 ls), corresponding to the reverberation inside the obstacle. Theoretically, the tail length is infinite, but, of course, only the main part of it can be observed.

 Case 1: Different values of length a and depth d of a prismatic triangular cut located at the middle of the beam, as depicted in Fig. 11(a), are studied.  Case 2: Two damaged regions with the same geometric appearance (triangular cuts). The first is centered at 400 mm, with length a = 25 mm and depth d = 5 mm. The second is centered at 600 mm, with length a = 12.5 mm and depth d = 2.5 mm.  Case 3: Two triangular cuts close to each other so that they superpose somehow. The first is centered at 500 mm, with length a = 20 mm and depth d = 5 mm. The second is centered at 513 mm, with length a = 10 mm and depth d = 2.5 mm.  Case 4: A circular hole with radius r = 2.5 mm located at the middle of the beam, crossing its entire width. 4.1. Vibration approach In order to quantify the sensitivity of the flexibility matrix to different damage scenarios, the relative Frobenious norm of the difference between the damaged and undamaged reduced flexibility matrices is considered. The plot of the error matrix for the different cases considered is omitted here because the similar appearance between them does not give any visual insight about the distinct situations. Firstly, variations in the damage length a and depth d of a prismatic cut are separately considered (Case 1). Table 1 shows the sensitivity of the flexibility matrix to changes in the depth d, with length a = 25 mm, centered at the middle of the beam. From Table 1 one can clearly notice an increase in the relative norm with the increase in depth d. According to Eq. (32), the damage field b is proportional to the cube of the beam height, with corroborates with the relatively greater norms for greater values of depth d. Table 2 presents the sensitivity of the flexibility matrix to changes in the length a of the damage with depth d = 50 mm, centered in the middle of the beam. It is worth noting that the second entrance in Tables 1 and 2 are the same, since the damage has the same length (a = 25 mm) and depth (d = 5 mm). Comparing Tables 1 and 2, one may clearly conclude that the vibration approach

1826

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828

0.04

15

0.02

10

0

5

−0.02

0

−0.04 185

190

195

200 t(µs)

205

210

x 10−4

−5 200

215

205

210

215

t(µs)

(a) Detail of the echo

(b) Reverberant part of the echo

Fig. 10. Two different zooms on the output signal.

a

30 mm

25 mm

d

10 mm

12.5 mm 2.5 mm

5 mm

500 mm

400 mm

(a) Case 1

600 mm

(b) Case 2

10 mm 5 mm 2.5 mm

5 mm

5 mm 500 mm

513 mm

500 mm

(c) Case 3

(d) Case 4

Fig. 11. Four different damages imposed to the slender beam: (a) Case 1; (b) Case 2; (c) Case 3; and (d) Case 4.

Table 1 Sensitivity of the flexibility matrix to damage depth in Case 1.

kGE Gmm kF kGmm kF

d = 2.5 mm

d = 5.0 mm

d = 7.5 mm

0.0057

0.0198

0.0838

It is worth emphasizing that the damage identification procedures should be based on the difference between the flexibility matrices of the damaged and undamaged structure, as that depicted in Fig. 8. The results presented in Tables 1–3 show that the flexibility matrix is sensitive to damage and, therefore, can be considered for damage identification purposes. 4.2. Wave propagation approach

Table 2 Sensitivity of the flexibility matrix to damage length in Case 1.

kGE Gmm kF kGmm kF

a = 15 mm

a = 25 mm

a = 35 mm

0.0121

0.0198

0.0275

Table 3 Sensitivity of the flexibility matrix to distinct damages.

kGE Gmm kF kGmm kF

Case 2

Case 3

Case 4

0.0398

0.0122

0.0068

adopted here presents a greater sensitivity to damage depth than to its length. The sensitivity of the flexibility matrix to distinct damage scenarios is presented in Table 3, which includes Cases 2, 3 and 4. Note that in Case 2, where there are two damaged regions apart one from another, the norm of the flexibility matrix is considerably greater than that of Case 3, where two damaged regions are close to each other. The obtained results are coherent with the damage geometries, showing to be greater for greater damage intensities.

To examine now the sensitivity of the output signal (echo) to damage in the wave propagation approach, the same four cases depicted in Fig. 11 are considered. As described before, four distinct combinations of parameters a and d are studied in Case 1. It should be noted that the vertical scales are different. The outgoing echoes g(t) generated from Case 1 for different values of depth d and length a are depicted in Fig. 12. In Fig. 12(a), the obstacle has a depth d = 2.5 mm. In Fig. 12(b), the depth was changed to d = 7.5 mm. In the situation shown in Fig. 12(c), the depth is the original one, d = 5 mm, but the length was changed to a = 15 mm. Finally, in Fig. 12(d) the length was modified to a = 35 mm. Comparing Fig. 12(a) and (b), it can be observed that the change in amplitude of the output signal is noticeably greater when the depth d of the damage increases. On the other hand, comparing Fig. 12(c) and (d), one may clearly observe the difference in the support of the echo: It is wider for the greater damage length. It is worth noting that the smaller the parameter a the bigger are the module of the reflection coefficients (see Eqs. (31) and (23)) and, as a consequence, the echo presents higher values. The results depicted in Fig. 12 indicate that the echo signal has enough sensitivity to changes in damage parameters.

1827

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828 0.015

0.15

0.01

0.1 0.05

0.005

0

0 −0.05

−0.005

−0.1

−0.01

−0.15

−0.015 150

160

170

180

190

200 t(µs)

210

220

230

240

250

−0.2 150

160

170

(a) depth 2.5 mm

180

190

200 t(µs)

210

220

230

240

250

230

240

250

(b) depth 7.5 mm

0.1

0.03 0.02

0.05 0.01

0

0 −0.01

−0.05 −0.02

−0.1 150

160

170

180

190

200 t(µs)

210

220

230

240

250

−0.03 150

160

170

(c) length 15 mm

180

190

200 t(µs)

210

220

(d) length 35 mm

Fig. 12. Echo generated by four different damages, Case 1.

0.04

0.15 0.1

0.02 0.05 0

0

−0.05 −0.02 −0.1 −0.04 100

−0.15 130

160

190

220 t(µs)

250

280

310

340 −0.2 150

160

Fig. 13. Echo generated from Case 2.

170

180

190

200 t(µs)

210

220

230

240

250

Fig. 15. Echo generated from Case 4. 0.04 0.03 0.02 0.01 0 −0.01 −0.02 −0.03 −0.04 150

160

170

180

190

200 t(µs)

210

220

230

240

250

Fig. 14. Echo generated from Case 3.

Fig. 13 shows the double echo obtained for Case 2. This result pinpoints the presence of damage in two distinct regions. It is worth noting that there is a main echo, resulting from the first damage, around 156 ls, a second one around 234 ls, and a third and minor component of the echo around 312 ls, this last one resulting from the reverberation between the two damages.

The echo generated from the damage scenario composed of two triangular cuts close to each other, Case 3, depicted in Fig. 11(c), is presented in Fig. 14. The presence of two events in the damaged region may be perceived in the echo shape. Since the damage assessment built on wave propagation approach is established in the time domain, the first part of the echo is not dependent on the whole damage shape, as can be recognized comparing, for instance, Figs. 14 and 9. The system response for Case 4 is presented in Fig. 15. One may clearly note the difference between the echo general shape for this case and that of the previous ones. Indeed, the initial part has a strong derivative, resulting from the high value of the initial reflection coefficients; after that, however, the behavior changes to a smooth curve, as expected. 5. Conclusions The damage influence on the vibrational behavior and on the wave propagation issues of a slender Euler–Bernoulli beam was

1828

R.A. Tenenbaum et al. / Computers and Structures 89 (2011) 1820–1828

studied in this work. For this, an aluminium beam with distinct imposed damage scenarios was considered. In the vibration framework, the damage influence was assessed by considering changes on the flexibility matrix of the structure. It was shown that a good estimate of the reduced flexibility matrix of the damaged structure may be obtained from a standard modal test. The results show that the reduced flexibility matrix presents sensitivity to structural damage and, hence, may be considered for damage identification purposes. In the wave propagation framework, the damage is clearly perceived by measuring the generated echo. Both the damage location and its intensity provide information to the output signal. The main conclusion for now is that both approaches present a noticeable sensitivity of the results to damage, showing that they are promising ones for the synthesis of damage identification techniques. Of course, the actual capability of each technique to identify both the position and the intensity of the real damage can only be assessed by using the available identification methods— such as Levenberg–Marquardt, Particle Swarm Optimization, Simulated Annealing, or Genetic Algorithms methods, among others—to solve the inverse problem, for each approach. It is worth mentioning that, even for the case of simple structures, such as the Euler–Bernoulli beam considered in this paper, the vibrational approach based on the flexibility matrix presents further difficulties, not considered in the present work, related to the computation of the required modal parameters. For instance, the algorithm for extracting the normalized mode shapes from the system response and the number and position of sensors are important issues that may significantly affect the results. However, this approach has the advantage of being directly applicable to assessment of damage influence in more complex structures. The wave propagation approach, on the other hand, yielded a straightforward assessment of the damage influence on the behavior of the beam. Since this is a direct approach based on the time domain, significant economy of time and computational effort were obtained. In the case of one-dimensional wave propagation, a unique sensor suitably positioned on the structure may provide sufficient information for the identification of the damage position and intensity. However, if more complex structures are taken into account, the number and position of sensors are important issues that must be considered. In the wave propagation approach, the computational effort in solving the direct problem involves the calculus of the algebraic Eqs. (26)–(31), which is straightforward and costless. Considering a Dirac delta input, all the values of Fjn+1 in Eq. (29) vanish, except for F1. If any other kind of input is considered, Eq. (29) involves, of course, more terms. However, this will not significantly increase the computational effort. In the vibration approach, a relatively greater computational effort is required. The modal parameters: natural frequencies and mass-normalized mode shapes, used to compute the flexibility matrix of the structure are not directly obtained from standard modal tests. Therefore, these parameters must be obtained by post-processing the modal test data: frequency response functions or impulse response data. Acknowledgements The authors acknowledge the Brazilian National Council for Scientific and Technological Development, CNPq, and Rio de Janeiro’s Foundation for Research Support, FAPERJ, for their financial support to this research.

References [1] Santos JVA, Maia NMM, Soares CMM. Structural damage identification: a survey. In: Topping BHV, Papadrakakis M, editors. Trends in computational structures technology. Stirlingshire, Scotland: Saxe-Coburg Publications; 2008. [2] Doebling SW, Farrar CR, Prime MB, Sheritz DW. Damage identification and health monitoring of structural and mechanical systems from changes in their vibration characteristics: a literature review. Los Alamos National Lab., Rept. LA-13070-MS, Los Alamos, NM: 1996. [3] Carden EP, Fanning P. Vibration based condition monitoring: a review. Struct Health Monit 2004;3:355–77. [4] Stutz LT, Castello DA, Rochinha FA. Structural damage identification based on a finite element updating. In: Proceedings of the sixth US national congress on computational mechanics, Dearborn, Michigan: 2001. [5] Stutz LT, Castello DA, Rochinha FA. A flexibility-based continuum damage identification approach. J Sound Vib 2005;279:641–67. [6] Castello DA, Stutz LT, Rochinha FA. A structural defect identification based on a continuum damage model. Comput Struct 2002;80:417–36. [7] Gangadharan R, Roy Mahapatra D, Gopalakrishnan S, Murthy CRL, Bhat MR. On the sensitivity of elastic waves due to structural damages: time-frequency based indexing method. J Sound Vib 2009;320:915–41. [8] Grabowska J, Palacz M, Krawczuka M. Damage identification by wavelet analysis. Mech Syst Signal Process 2008;22:1623–35. [9] Fernandes KM, Stutz LT, Tenenbaum RA, Silva Neto AJ. Vibration and wave propagation approaches applied to assess damage influence on the behavior of Euler-Bernoulli beams—Part II: Inverse problem. In: Proceedings of the ninth international conference on computational structures technology, Paper 39, Athens: 2008. [10] Tenenbaum RA, Zindeluk M. A fast algorithm to solve the inverse scattering problem in layered media with arbitrary input. J Acoust Soc Am 1992;92(6):3371–8. [11] Zhongqing S, Lin Y, Ye L. Guided. Lamb waves for identification of damage in composite structures: a review. J Sound Vib 2006;295:753–80. [12] Pruell C, Kim JY, Qu J, Jacobs LJ. Evaluation of plasticity driven material damage using Lamb waves. Appl Phys Lett 2007;91:231911. [13] Vanaverbeke S, Abeele KVD. Two-dimensional modeling of wave propagation in materials with hysteretic nonlinearity. J Acoust Soc Am 2007;122:58–72. [14] Küchler S, Meurer T, Jacobs LJ, Qu J. Two-dimensional wave propagation in an elastic half-space with quadratic nonlinearity: a numerical study. J Acoust Soc Am 2009;125:1293–301. [15] Mendell JM, Ashrafi FH. A survey of approaches to solving inverse problems for lossless layered media systems. IEEE Trans Geosci Remote Sens 1980;GE18(4):320–30. [16] Robinson EA, Treitel S. Geophysical signal analysis. Englewood Cliffs: Prentice Hall; 1980. [17] Lefebre JP. La tomographie d’impédance acoustique. Trait Signal 1985;2(2):103–10. [18] Nagy PB. Fatigue damage assessment by nonlinear ultrasonic material characterization. Ultrasonics 1998;36:375–81. [19] Cantrell JH, Yost WT. Nonlinear ultrasonic characterization of fatigue microstructures. Int J Fatigue 2001;23:480–7. [20] Tenenbaum RA, Zindeluk M. Classical signal processing techniques applied to variable impedance identification by impulse testing. In: Proceedings of the 15th Internoise, vol. 2, 1986. p. 1241–46. [21] Nia SH, Charngc JJ, Lod KF. Low-strain integrity testing of drilled piles with high slenderness ratio. Comput Geotech 2006;33:283–93. [22] Ni SH, Lo KF, Lehmann L, Huang YH. Time-frequency analyses of pile-integrity testing using wavelet transform. Comput Geotech 2008;35:600–7. [23] Niederleithinger E, Taffe A. Early stage elastic wave velocity of concrete piles. Cem Concr Compos 2006;28(4):317–20. [24] Alvandi A, Cremona C. Assessment of vibration-based damage identification techniques. J Sound Vib 2006;292:179–202. [25] Pandey AK, Biswas M. Damage detection in structures using changes in flexibility. J Sound Vib 1994;169:3–17. [26] Guyan RJ. Reduction of stiffness and mass matrices. AIAA J 1965;3(2):380. [27] Alvin K. Method for determining minimum-order mass and stiffness matrices from modal test data. AIAA J 1995;33:128–45. [28] Tenenbaum RA, Zindeluk M. An exact solution for the one-dimensional elastic wave equation in layered media. J Acoust Soc Am 1992;92(6):3364–70. [29] Santosa F, Schwetlick H. The inversion of acoustical impedance profile by methods of characteristics. Geophys JR Astron Soc 1982;70:229–43. [30] Gordon P, Sanday SC. Elastic waves in heterogeneous bars of varying crosssection. J Franklin Inst 1977;303(2):129–45. [31] Pierce AD. Acoustics: an introduction to its physical principles and applications. Woodbury: Acoustical Society of America; 1989. [32] Hooshyar MA. Goupillaud layers and construction of wave speed and density of a layered acoustic medium. J Acoust Soc Am 1990;87:2310–8.