- Email: [email protected]
A numerical method for nonlinear complex modes with application to active–passive damped sandwich structures
Engineering Structures 31 (2009) 284–291
Contents lists available at ScienceDirect
Engineering Structures journal homepage: www.elsevier.com/locate/engstruct
A numerical method for nonlinear complex modes with application to active–passive damped sandwich structures Hakim Boudaoud a,b , Salim Belouettar a,∗ , El Mostafa Daya b , Michel Potier-Ferry b a
Centre de Recherche Public Henri Tudor, 29, Avenue John F. Kennedy, L-1855 Luxembourg, Luxembourg
b
LPMM, UMR CNRS 7554, I.S.G.M.P., Université Paul Verlaine-Metz, Ile du Saulcy, F-57045 Metz Cedex 01, France
article
info
Article history: Received 24 January 2008 Received in revised form 16 June 2008 Accepted 5 August 2008 Available online 23 September 2008 Keywords: Sandwich structures Viscoelasticity Piezoelectricity Active Passive Hybrid Finite element method Vibrations Nonlinear eigenvalues Perturbation technique Homotopy Asymptotic numerical method
a b s t r a c t In this paper, a numerical method is proposed for determining complex vibrations modes of sandwich structures with piezoelectric and viscoelastic layers. Based on homotopy and asymptotic numerical techniques, this method leads to the damping properties calculation (loss factor and natural frequency per mode) of the hybrid sandwich structure. The numerical results of the loss factor and natural frequency are compared to those obtained from analytical beam model and from numerical studies the modal strain energy method. © 2008 Elsevier Ltd. All rights reserved.
1. Introduction In many industrial and defence applications, noise and vibration are important problems. Control of sound and vibration has been and is always the subject of a lot of research in recent years, and examples of applications are now numerous [1,2]. The three most common classifications of vibration control are: (1) passive; (2) active and (3) hybrid passive–active control. Passive control involves the use of reactive or resistive devices that either load the transmission path of the disturbing vibration or absorb vibratory energy. Active control also loads the transmission path but achieves this loading through the use of actuators that generally require external energy. In recent years, researchers have been directed to simultaneous use of piezoelectric and viscoelastic materials, in the so called ‘‘hybrid vibration control’’. A large part of studies highlighted structure with an hybrid arrangement called ACLD (SCLD) for active (smart) constrained layer damping. In many theoretical and experimental investigations the active
∗
Corresponding author. Tel.: +352 54 55 80 530; fax: +352 42 59 91 333. E-mail address: [email protected] (S. Belouettar).
0141-0296/$ – see front matter © 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.engstruct.2008.08.008
system increases the passive transverse shear deformation of the viscoelastic constrained layer and appears to be an effective means to control vibration [3,4]. To provide reliable and effective damping treatments, other interesting configurations are studied [5,6]. Typically, these structures are five layers sandwich ones (Fig. 1) in which a thin and soft viscoelastic layer is confined between identical elastic and stiff layers and two thin bonded piezoelectric layers. Much investigation has been devoted to the linear dynamic analysis of these structures. Their stiffness matrix is a complex one and depends nonlinearly on the vibration frequency. So, the vibration analysis leads to a non-linear complex eigenvalue problem. Solving this problem yields, first complex modes that can be slightly different from undamped modes, second complex eigenvalues whose real and imaginary parts are related to the loss factors and to the viscous eigen frequencies. From an engineering point of view, the most relevant quantity is the loss factor that is associated with each mode. A lot of methods have been presented to predict these structural damping properties and a detailed review of these works can be found in Reference [7]. Previously, there was no direct method available for solving such a problem efficiently in the case of large-scale structures. Often the stiffness matrix is assumed to be constant in a frequency
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
285
Fig. 1. Five layered sandwich structure with viscoelastic middle layer and two piezoelectric layers.
the asymptotic method of reference by [11], but more general since the introduction of the continuation procedure permits us to apply it to any case, and not only for small damping. The proposed investigation is an extension of a previous work of [6] for viscoelastically damped sandwich structures and the work of [24] for the nonlinear vibrations analysis of viscoelastically damped sandwich shells. 2. Material modelling 2.1. Piezoelectric constitutive equations
Fig. 2. Detailed Sketch of the considered five layered sandwich structure.
The electroelastic response of a piezoelectric body of volume Ω and regular boundary surface S, is governed by the mechanical and electrostatic equilibrium equations,
band, leading to eigenvalue problems with constant complex matrices. Sometimes, this complex eigenvalue problem is solved by the QR method [8], but this does not work for structures with a large number of degrees of freedom and can only used when the [K (ω)] is a constant complex matrix. Another way consists of computing the response to harmonic forcing for many frequencies, to yield resonance peaks from which loss factors and eigenfrequencies can be determined [9]: this can be applied to any structural problem, but the computational cost is very high. Approximate techniques such as the modal strain energy method [9] can also be used, but only to estimate the loss factor. An asymptotic technique was proposed in [10], but the procedure was limited to a small order of truncation and did not use a continuation algorithm, so this method cannot be used for large damping or for a strongly varying stiffness modulus. Recently, the algorithm proposed by Chen et al. [11] has proved much more satisfactory, because the nonlinear complex eigenvalue problem is solved by coupling a first order perturbation, an iterative algorithm and a reduced basis technique. In this way, the problem can be solved accurately, even for large-scale structures. Initially, analytical methods were developed to yield approximate loss factors and natural frequencies of sandwich beams or plates with simple boundary conditions [12–18]. In practice, it is necessary to design sandwich structures with a complex geometry and generic boundary conditions by using finite element simulations. Finite elements are still being used for analysis of damped sandwich structures [9–11,19–23] with the help of the aforementioned computational methods. This paper aims to address this deficiency and propose an efficient numerical method able to solve accurately the nonlinear complex eigenvalue problem presented above with a moderate computational cost. The proposed method is convenient for large-scale problems and applicable with any constitutive relation. The presented numerical methodology is similar to
σij = Cijkl εkl − ekij Ek Dk = ekij εij + kl El
(1)
where σ and are the second rank stress and strain tensors, E and D are, respectively, the electric field and the electric displacement vectors, ε is the second rank dielectric tensor, e is the third rank coupling tensor and C is the fourth rank elasticity tensor. The superscripts E and e indicate that the elasticity and permittivity constants are determined under conditions of zero or constant electric field and strain, respectively. We assume that no body electric charges are applied to the piezoelectric layers. The electric potential layer-wise assumption implies that in each of the layers the electric potential scalar field is represented by an independent function φ . A linear type variation is chosen to represent the electric potential in each piezoelectric layer, which implies that two electric potential nodes are needed for each lamina. The electric field-electric potential relation is the following: Ei = −(∇ φ)i .
(2)
2.2. Viscoelastic material constitutive modelling A general linear viscoelastic constitutive law is assumed for the core. Classically, such a law involves convolution product, [25]. For instance the relation between the axial stress and strain can be written in the form of a convolution product:
Z
t
Y (t − τ ) −∞
∂ ∂ (x, τ )dτ = Y ∗ ∂τ ∂t
(3)
where Y is the relaxation function of the viscoelastic material. Considering that the structure is subjected to harmonic vibrations, the viscoelastic behaviour can be described by a complex Young modulus which is depending nonlinearly on the frequency and temperature [25].
286
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
3. Finite element formulation From the application of the principle of virtual work, the following equilibrium equation is obtained
δ Wext + δ Ec − δ Eh = 0
(4)
where δ Wext denotes the virtual work of the external loading, δ Ec represents the kinematical energy and δ Eh is the internal electromechanical energy. Ec and Eh are given as:
Z 1X 2 2 2 ˙ ˙ ˙ d Ωk + W + V U ρ E = k k k k c 2 k Ωk k = s, e1 , Zc , e2 , a Z X 1X E = 1 k k σij ij dΩk − Eil Dli dΩl h 2
2
Ωk
k
Ωl
l
(5) l = s , a.
We presume that the piezoelectric layers are thin enough to consider the hypothesis of constant electric field through the thickness of each piezoelectric lamina and nile in the other directions as valid. We also consider as in [26] the central to layer to be conductive with a uniform potential fixed to zero. Thus, the Eq. (2) takes the following form: E3 = −
(6)
hp
e i
= Bu u (7) i where Bu is the product of the differential operating matrix i i S
n o {Fu (x, y, z , t )} = F˜u (x, y, z ) eiωt {u (x, y, z , t )} = {U (x, y, z )} eiωt {φ (x, y, z , t )} = {Φ (x, y, z )} eiωt
(11)
√
where ω is the frequency and i = −1. Substituting Eq. (11) in Eq. (9) and considering separately the sensor and actuator potentials, the global dynamic finite element equations can be rewritten as:
"K (ω) uu Kφ u s Kφ u a ( ) Fu = 0
Ku φ a 0 Kφφ a
Kuφ s Kφφ s 0
#
Muu 0 0
" −ω
2
0 0 0
0 0 0
#! (
U
)
Φs Φa
(12)
0
φ
where φ denotes the superficial electrical potential and hp is the thickness of the piezoelectric layer. Introducing the strain-displacement relation in terms of the nodal displacement lead, in each i layer, to the discretized expression of the mechanical strain:
i
The matrices Bi depend only on the shape functions of the considered finite element and Ci are the matrices obtained from the strain stress laws. [Muu ] is assumed to be constant in this analysis. In order to estimate the damping properties (frequency and loss factor per vibration mode) of the investigated adaptive structure when subjected to harmonic loads, a frequency response analysis is performed. The applied force, the displacement field and the electrical potential are expressed as:
relating S to the shape function matrix Nu . Similarly we get the same relations for the piezoelectric layer:
= −∇ Φ k = −∇ Nφ k φ e = − Bφ k φ e (8) k where Bφ is the product of the differential operating matrix relating the electrical field E k and the shape function matrices k
T
where Kuφ = Kφ u . From the previous equation, the generated potential on the sensor is given by:
−1 {Φs } = − Kφφ s Kφ u s {U } .
(13)
In the proposed analysis, a velocity and a proportional feedback control laws are considered. According to these laws, the controlled voltage for each layer is expressed as:
φ a = gd φ s + gv φ˙ s .
(14)
gd and gv are the feedback control gain and the velocity feedback gain. The Eq. (14) permits to write:
k E
Nφ . Substituting Eqs. (7) and (8) into Eq. (4) one obtains the equilibrium equation for an element, which after taking into account the contribution of all elements in the domain, yields the equilibrium equations of the system obtained after performing the usual finite element assembling techniques.
[K ] {u} + Kuφ {φ} + [Muu ] {¨u} = {Fu } uu T Kuφ {u} + Kφφ {φ} = {0} .
(9)
The quantities [Muu ], [Kuu ] are the mass and the stiffness matrices used in the structural analysis, Kuφ is electrical
mechanical coupled stiffness matrix, Kφφ is the electrical stiffness matrix, {Fu } is the force vector of the mechanical loads applied to the structure.
XZ t [Muu ] = ρ Nu i Nu i dΩ vi i X Z t Bu i C i Bu i dΩ [Kuu ] = v i i Z k X k t t k K = Bu [e] Bφ dΩ uφ k Zvk k t k X Kφφ = [ ] Bφ k dΩ Bφ n
vk
−1 {Φa } = − (gd + iωgv ) Kφφ s Kφ u s {U } .
(15)
Substituting the Eqs. (13) and (15) into Eq. (12), we get the dynamic finite element equations with condensed electrical degrees of freedom. [Kreal (ω)] + i Kcomp (ω) − ω2 [Muu ] {U } = {Fu }
(16)
where the complex nodal vibration eigenmode is denoted by U and with the following definition for the real and imaginary parts:
] = Re {[Kuu (ω)]} [Kreal (ω) −1 − Kuφ s + gd Kuφ a Kφφ s Kφ u s −1 Kcomp (ω) = Im {[Kuu (ω)]} − ωgv Kuφ a Kφφ s Kφ u s for simplicity the real and imaginary part of the complex stiffness matrix Re {R[K uu (ω) ]}I and Im {[Kuu (ω)]} will be respectively denoted Kuu and Kuu . 4. Application of the asymptotic numerical method 4.1. Homotopy technique
i = e 1 , e 2 , s , a, c k = s , a.
(10)
The object of this first stage is to transform the nonlinear eigenvalue problem using the homotopy technique [27]. This allows us to introduce artificially a parameter ε in Eq. (16) as:
[Kreal (ε)] + iε Kcomp − ω2 [Muu ] {U } = 0 0≤ε≤1
(17)
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
287
I {U0 }t −gd [Kactuator ]t + i Kuu − p0 gv [Kactuator ]t {U0 } C1 = {U0 }t [Muu ] {U0 } k=j−1 j −1 t k=P t I t P pk [Kactuator ]t {U0 } Kuu − gd [Kactuator ]t − igv (Uj−k−1 ) − Ck Uj−k [Muu ] + i Uj−1 k=0 k=1 Cj = t { } { } ] [M U U 0 uu 0 j = 2, 3, . . . Box I.
with
(
−1 R [Kreal (ε)] = Kuu − Kuφ s + ε gd Kuφ a Kφφ s Kφ u s I −1 Kcomp = Kuu − ωgv Kuφ a Kφφ s Kφ u s .
One can easily note that the initial problem (Eq. (16)) corresponds to ε = 1. The system corresponding to ε = 0 is the vibration problem of an undamped structure which can be written in the following form
R Kuu − Kuφ
s
Kφφ
−1 s
Kφ u
s
− ω2 [Muu ] {U0 } = 0.
(18)
This equation represents a real eigenvalue problem that can be solved by classical methods such as subspace iteration or the Lanczos method, and we suppose that we have computed the required eigenpairs of (18). To get the solutions of the nonlinear eigenvalue problem (16), we try to solve the Eq. (17) deduced by the homotopy transformation. In fact, the Eq. (17) does not define a well-posed problem, because the unknown can be multiplied by an arbitrary constant. To obtain a well-posed problem, it is sufficient to require an additional condition, for instance
{U0 }t {U − U0 } = {0} .
(19)
The Eqs. (16) and (18) define a well-posed path following problem ε → {{U (ε)} , ω(ε)}, having to be solved by an appropriate method. We could use an incremental-iterative method like Newton–Raphson or the modified Newton method, but we preferred an asymptotic-numerical method, permitting us to restrict the number of matrices to be inverted. 4.2. Perturbation technique At this stage, we describe the asymptotic-numerical method being applied to solve the perturbed eigenvalue problem (16) and (18). This technique has been extensively studied for other nonlinear problems, such as nonlinear elasticity, Navier–Stokes equations and structural problems involving unilateral contact, plasticity or viscoplasticity [28]. The asymptotic numerical method associates a perturbation technique with a classical discretization procedure, which is generally the finite element method. It was first proposed by [29] and we refer to [30] for a complete reference list. Within asymptotic-numerical techniques, one searches for the unknowns {U } and ω in form of truncated integro-power series with respect to the parameter ε :
N X { } U = ε j Uj j =0 N X ω= ε j pj j =0 N X 2 ω = ε j Cj j=0 0≤ε≤1
(20)
where N is the truncation order of the series and ω the complex eigenfrequency. In this analysis, the modulus E (w) is assumed to be dependent on frequency and after insertion of the series (20) in the Eq. (17), one gets a set of recurrent linear problems that coming from the identification of the like powers of ε : The first order problem:
[A] {U0} = 0 R [A] = Kuu − [Ksensor ] − C0 [Muu ] (21) 2 C0 = ω0 −1 where [Ksensor ] = Kuφ s Kφφ s Kφ u s and [Kactuator ] = Kuφ a −1
Kφφ s Kφ u s are respectively the added static stiffness due to sensor and actuator. The second order problem:
[A] {U1 } = {F1 } I {F } = g [K ] − i Kuu − p0 gv [Kactuator ] 1+ C [Md ]actuator 1 uu ) {U0 }.
(22)
The (j + 1)th order problem:
[A] Uj = Fj k=j X I Fj = Ck [Muu ] Uj−k − i Kuu − gd [Kactuator ] Uj−1 k=1 kX =j −1 pk [Kactuator ] Uj−k−1 . +igv
(23)
k=0
Notice that for the first order, we recover the undamped eigenvalue problem and we propose to compute the required eigenpairs by any classical method for solving a symmetric eigenvalue problem. The problems for the higher orders are nonhomogeneous linear ones having the same matrix. Therefore, since the method combines the classical and asymptotic-numerical approaches, we can compute many terms of the series by inverting only one matrix per frequency. In the present case, the method has an additional advantage because the vectors in Box I are complex, while the matrix to be inverted is real. Note also that the matrix is singular and its kernel is generated by the undamped mode. Then the problems (Box I) have a solution if and only their right-hand sides j = 1 to N satisfy the following solvability condition:
t Fj
{U0 } = 0.
(24)
This condition permits us to compute the term Cj in the expansion (Eq. (20)) of the eigenvalue as the equation in Box I. We note that this number depends only on the previously computed vectors {Un } and Cn obtained in the previous order. So, we can compute C1 by Box I, next C2 by solving Box I, next C2 again by Box I and so on. If the solvability condition is satisfied, the linear equations (22) and (23) are well-posed and they have a unique solution if one takes into account the additional condition (25), that leading to:
{U0 }t Uj = 0 j = 1, . . . , N .
(25)
288
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
−1
Re {[Kuu (ω)]} + Kuφ s + ε Gd Kuφ a Kφφ s Kφ u s + iε [C (ω)] − ω2 [Muu ] {U (ε)}
R (ε) = k(Re {[Kuu (ω)]} {U (ε)})k N X {U } = Uj j =0 v u N uX Cj ω = it j =0
Box II.
One can rearrange the equations Eqs. (23) and (25) by introducing a Lagrange multiplier k, and it will be easy to check that the solution of (23) and (25) satisfies the following matrix equation:
[A]
{U0 }t
{U0 } 0
Uj k
=
Fj 0
.
(26)
The matrix in [A] has an inverse and we can factorize it by the Crout method. This permits us to solve many linear problems (22) and (23) with a single matrix factorization. In this way, we get solutions of the perturbed problem (17) in the form:
N X {U (ε)} = ε j Uj v j =0 u N uX ω = i t ε j Cj .
followed in classical matrix form: σ C C12 C13 0 0 11 1 C C C 0 0 σ 12 22 23 2 σ 3 C13 C23 C33 0 0 0 0 0 C 0 σ 44 4 0 0 0 0 C σ5 = 55 0 0 0 0 0 σ6 0 0 0 e15 0 D1 D 0 0 0 e 0 2 24 D3
e31
e32
e33
0
0
0
0
0
0
0
0
0
0
0
0
0
−e24
0
−e15
0
C66
0
0
0
11
0
0
0
22
0
0
0
1 2 3 γ23 × γ13 γ12 E 1 E2
(27)
j=0
−e31 −e32 −e33 0 0 0 0 0
33
(28)
E3
This is only an estimate of the solution of (17) because the order of truncation N is necessarily finite (generally, N = 20 is chosen) and also because can be outside the radius of convergence of the series. The latter has to be greater than or equal to 1, because we have to solve the problem that corresponds to the case ε = 1. If the estimate is good enough, i.e. the relative residual Box II is smaller than a tolerance (such as 10−8 . . . 10−6 ), the solution of problem is obtained from (27) by putting ε = 1. Otherwise we have to improve the method by replacing the polynomial approximation by a rational one. This technique, named Padé approximants, allow us to increase the radius of convergence of the series. More details of this technique may be found in Cochelin et al. [31] and Elhage-Hussein et al. [32]. If the the radius of convergence of the series still fewer than one we have to use a continuation method [6] or iterative algorithms [24].
5. Numerical applications The present numerical procedure has been applied to the free vibrations analyses of a 5-layered simply supported sandwich beam, Figs. 1 and 2. The materials composing the elastic layers and the viscoelastic layers are characterised by the values of the Young modulus and the Poisson ratio: Ee = 45.54 GPa, ν = 0.33 for the elastic material and Ev = 7.25 MPa, ν = 0.45 for the viscoelastic core. The elastic, viscoelastic and piezoelectric materials have respectively the following densities: ρe = 2040 kg/m3 , ρv = 1200 kg/m3 and ρp = 7500 kg/m3 . The constitutive equations of the PZT4 orthotropic piezoelectric layers may be summarised as
where:
C11 = C22 = 139.0 GPa; C33 = 115.3 GPa C = C = 74 . 3 GPa ; C = 77.8 GPa 13 23 12 C = C = 25.6 GPa; C66 = 30.6 GPa 44 55 2 e33 = 15.1 C/m2 ; e31 = e32 = −5.21 C/m ; 2 e15 = e24 = 12.7 C−/8m 11 = 22 = 1.31e F/m; 33 = 1.15e−8 F/m.
(29)
In the application of the ANM, the truncation is taken as N = 20 and the relative residual R () is smaller than 10−4 . 5.1. Case of five layered sandwich beam We propose to analyse the modal properties of a five-layered simply supported beam of length L = 500 mm and high h = 10 mm. The thicknesses of the various layers composing the beam are the same as those used in [33] i.e. he = 69h/150, hv = h/15 and hp = h/150 respectively for the elastic, viscoelastic and piezoelectric layers. The viscoelastic properties of the central layer are introduced by a complex modulus that is assumed to be constant (see Table 1). Two tests are presented that mainly differ from one another by the control law i.e. a direct proportional feedback one used is the first test and velocity feedback control law is used in the second one. 5.1.1. Proportional direct control law The free vibration of the beam is investigated in the case of direct proportional feedback (gv = 0). Several values of the proportional feedback control gain (gd ) and different viscoelastic loss factors ηvisco values were considered (see Tables 1 and 2). This test has been investigated quite extensively in [33] from
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
289
Table 1 Proportional feedback control gd = −1, gv = 0: Natural frequencies and loss factors of simply supported sandwich beam for for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
95.2 96.5 96.5 95.5 96.8 96.8 96.2 97.5 97.5 97.4 98.7 98.7 99.0 100.2 100.3 100.8 102.0 102.1
0 0 0 0.0740 0.0741 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
298.9 303.1 303.1 299.1 303.5 303.3 299.7 304.6 304.0 300.7 306.6 305.5 302.2 309.5 306.5 304.0 313.4 308.4
0 0 0 0.057 0.057 0.057 0.113 0.114 0.114 0.168 0.169 0.168 0.221 0.223 0.221 0.272 0.275 0.275
612.5 622.2 622.2 612.6 622.6 622.3 613. 623.9 622.3 613.6 626.1 623.3 614.5 629.2 624.2 615.6 633.4 625.4
0 0 0 0.035 0.035 0.035 0.07 0.07 0.07 0.105 0.105 0.105 0.14 0.139 0.14 0.174 0.173 0.174
Table 2 Proportional feedback control gd = −20, gv = 0: Natural frequencies and loss factors of simply supported sandwich beam for for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
89.8 91.0 91.0 90.1 91.3 91.3 90.8 92.0 92.0 91.9 93.1 93.1 93.4 94.6 94.7 95.1 96.3 96.4
0 0 0 0.0740 0.0741 0.074 0.143 0.144 0.143 0.204 0.205 0.204 0.254 0.256 0.254 0.293 0.296 0.293
280.8 284.8 284.8 281.0 285.2 285.1 281.6 286.35 285.7 282.7 288.3 286.8 284.1 291.2 288.2 286.0 295.2 290.1
0 0 0 0.059 0.059 0.059 0.117 0.117 0.117 0.173 0.174 0.174 0.227 0.230 0.227 0.279 0.282 0.279
573.4 582.6 582.6 573.5 583.0 582.8 573.9 584.4 583.2 574.6 586.6 583.8 575.5 589.8 584.7 576.6 594.1 585.9
0 0 0 0.037 0.037 0.037 0.074 0.074 0.074 0.11 0.11 0.11 0.146 0.145 0.146 0.181 0.18 0.181
Table 3 Proportional feedback control: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco and control gain gd (gv = 0)
ηvisco
0
0.2
0.4
0.6
0.8
1
gd = −1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
gd = −20
gd = −150
f1 in (Hz)
η1
f1 in (Hz)
η1
f1 in (Hz)
η1
95.2 96.5 96.5 95.5 96.8 96.8 96.2 97.5 97.5 97.4 98.7 98.7 99.0 100.2 100.3 100.8 102.0 102.1
0 0 0 0.0740 0.0741 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
89.8 91.0 91.0 90.1 91.3 91.3 90.8 92.0 92.0 91.9 93.1 93.2 93.4 94.6 94.7 95.1 96.3 96.4
0 0 0 0.074 0.0741 0.074 0.143 0.144 0.143 0.204 0.205 0.204 0.254 0.256 0.254 0.279 0.282 0.279
44.8 45.5 45.6 45.1 45.8 45.8 45.8 46.5 46.5 46.8 47.6 47.6 48.1 49.0 48.9 49.6 51.3 50.4
0 0 0 0.12 0.12 0.119 0.23 0.23 0.23 0.32 0.31 0.31 0.38 0.39 0.376 0.42 0.41 0.42
theoretical points of view. In [33], the analytical results are validated by using a 2D finite element model on Abaqus. The FE mesh of the sandwich beam is composed of 240 CPE8R: 8-node plane strain elements with a reduced integration. A steady-state linear dynamic analysis procedure [34] is used to monitor the
linear response of the beam structure subjected to a continuous harmonic excitation. In Table 3 are presented the first modal property (frequency and damping) obtained using the CSME, the analytical model [33] and those using the ANM method.In all analysed situations, we observe a very good agreement between the results of the ANM
290
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
Table 4 Calculation time needed with CSME method and new ANM method considering different number of d.o.f. d.o.f.
Time using CSME
Time using ANM
600 1000 1500
2 min 3 s 6 min 33 s 21 min 48 s
46 s 3 min 1 s 7 min 45 s
using three methods: the analytical approach by [33], the CSME and the asymptotic numerical method. Beyond the agreement between the results of the proposed method and those obtained analytically or by the CSME, the results of the ANM are valid for higher values of damping in the core or higher values of the control gains and for the higher modes (Tables 5–7). 6. Conclusion
method and those obtained analytically [33], those obtained by using the CSME (see tables). However, when the damping in the core becomes more significant or when the control increases, the values obtained with the ANM are in agreement with the analytical values whereas those obtained with the CSME tend not to be in agreement. In a second part we consider the necessary time to get the four first modal properties using the CSME and the ANM method. In Table 4 we observe, for 3 different meshes, that in each considered case the ANM method is less time-consuming than the CSME method. 5.1.2. Velocity feedback control law The second example has been designed to establish the ability of the method to analyse adaptive structures with a velocity feedback control. As for the previous case, different values of feedback control gain, gv , and damping were considered. In Table 5 are presented the results of the frequency and damping obtained by
In this paper, an asymptotic numerical method has been developed for solving an eigenvalue problem for a five layers hybrid damped sandwich structure. The method consists of introducing an artificially modified problem that depends on a parameter , and we seek eigenvalues and eigenmodes in the form of integro-power series expansions with respect to . The asymptotic numerical analysis permits us to decompose the perturbed problem into a sequence of linear problems. The first one is a real eigenvalue problem that can be solved by classical procedures. The other problems are non-homogeneous with the same linear operator. In this way, the method permits us to considerably simplify the analysis of the nonlinear eigenvalue problem. The method presented has three advantages, as compared with the classical techniques that are commonly used to estimate the loss factor of a viscoelastic structure. First, it solves the complex nonlinear eigenvalue problem accurately, so that the
Table 5 Velocity feedback control gd = 0, gv = 0: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
95.5 96.8 96.8 95.8 97.0 97.1 96.5 97.8 97.8 97.8 99.0 99.0 99.3 100.5 100.6 101.1 102.3 102.4
0 0 0 0.074 0.074 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
299.8 304.1 304.1 300.0 304.6 304.3 300.6 306.4 304.9 301.7 309.4 306.0 303.1 313.7 307.4 304.9 319.6 309.3
0 0 0 0.057 0.057 0.057 0.113 0.113 0.113 0.168 0.167 0.168 0.221 0.217 0.221 0.271 0.264 0.271
614.5 624.2 624.2 614.6 624.6 624.3 615.0 625.9 624.7 615.6 628.1 625.4 616.5 631.2 626.2 617.6 635.4 627.4
0 0 0 0.035 0.035 0.035 0.07 0.07 0.07 0.105 0.105 0.105 0.140 0.139 0.139 0.174 0.173 0.173
Table 6 Velocity feedback control gd = 0, gv = −0.01: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
95.5 96.8 96.8 95.6 96.9 96.8 96.1 97.5 97.4 97.1 98.5 98.4 98.5 100. 99.8 100.1 101.7 101.5
0.036 0.036 0.036 0.11 0.11 0.11 0.18 0.181 0.18 0.243 0.243 0.243 0.295 0.296 0.296 0.337 0.338 0.337
299.5 305.5 302.9 298.4 306.7 301.8 297.8 309.1 301.1 297.5 312.8 300.9 297.7 318.0 301.1 298.0 325.0 301.7
0.117 0.119 0.119 0.175 0.179 0.176 0.233 0.241 0.234 0.291 0.302 0.292 0.347 0.36 0.348 0.402 0.422 0.402
610.2 642.0 612.1 608.7 643.5 610.5 605.7 647.4 607.6 603.0 652.3 604.9 600.5 658.3 602.5 598.0 665.7 600.4
0.28 0.286 0.278 0.30 0.31 0.30 0.33 0.355 0.33 0.37 0.404 0.373 0.412 0.455 0.412 0.45 0.51 0.45
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
291
Table 7 Velocity feedback control: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco and control gain gv (gd = 0)
ηvisco
0
0.2
0.4
0.6
0.8
1
gv = −0.01
gv = 0 Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
gv = −0.05
f1 in (Hz)
η1
f1 in (Hz)
η1
f1 in (Hz)
η1
95.5 96.8 96.8 95.8 97.0 97.0 96.5 97.8 97.8 97.8 99.0 99.0 99.3 100.5 100.6 101.1 102.3 102.4
0 0 0 0.074 0.074 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
95.5 96.8 96.8 95.6 96.9 96.8 96.1 97.5 97.4 97.1 98.5 98.4 98.5 100 99.8 100.1 101.7 101.5
0.036 0.036 0.036 0.11 0.11 0.11 0.18 0.185 0.181 0.243 0.243 0.243 0.295 0.296 0.296 0.337 0.338 0.337
95.2 96.5 95.8 94.6 96.1 95.1 94.2 96.1 94.9 94.5 96.6 95.2 95.2 97.6 95.8 96.2 98.8 96.9
0.180 0.181 0.18 0.258 0.258 0.258 0.334 0.332 0.333 0.404 0.4 0.401 0.464 0.458 0.461 0.514 0.505 0.51
computational technique is still valid for large damping and higher values of control gain. Secondly it is able to account for a strong variation of the modulus with frequency, as is often the case for polymers. Thirdly, the method remains efficient for problems involving many degrees of freedom, because it is only necessary to solve a classical eigenvalue problem and to decompose a few matrices. Nevertheless, the present method can probably be significantly improved. Indeed, there are many ways to set up a homotopy transformation and, using the one presented in [35], one could avoid decomposing complex matrices. Furthermore, the power series do not yield the best possible representation of a solution path and the algorithm becomes much more efficient if one introduces Padé approximants [35,30]. Acknowledgements This work has been done in the framework of the FNR-MAT08-01 ADYMA project of the the European ERA-NET Materials (MATERA) and CASSEM (Contact N: FP6-NMP3-CT-2005-013517). The authors acknowledge the financial support of the Luxembourgish National Research Fund (www.fnr.lu) and the European Commission. References [1] Nakra BC. Vibration control with viscoelastic materials. Shock Vib Dig 1976;8: 3–12. [2] Hu H, Belouettar S, Daya EM, Potier-Ferry M. Review and assessment of various theories for modeling sandwich composites. Compos Struct 2008; 84(3):282–92. [3] Vasques CMA, Mace BR, Gardonio P, Dias Rodrigues J. Arbitrary active constrained layer damping treatments on beams: Finite element modelling and experimental validation. Comput Struct 2006;84:1384–401. [4] Vasques CMA, Dias Rodrigues J. Combined feedback/feedforward active control of vibration of beams with acld treatments: Numerical simulation. Comput Struct 2008;86:292–306. [5] Agnes GS, Napolitano K. Active constrained layer viscoelastic damping. In: 34th AIAA/ASME/ASCE/AHS/ASC Struct., structural dyn., and mater. conf., AIAA. 1993. p. 3499–506. [6] Daya EM, Potier-Ferry M. A numerical method for nonlinear eigenvalue problems: Application to vibrations of viscoelastic structures. Comput Struct 2001;79:533–41. [7] Ganapathi M, Patel BP, Boisse P, Polit O. Flexural loss factors of sandwich and laminated composite beams using linear and non linear dynamic analysis. Compos Part B: Eng 1999;30:245–56. [8] Bathe KJ. Finite element procedures in engineering analysis. NJ: Prentice-Hall; 1982. [9] Soni ML. Finite element analysis of viscoelastically damped sandwich structures. Shock Vib Bull 1981;55(1):97–109. [10] Ma BA, He JF. Finite element analysis of viscoelastically damped sandwich plates. J Sound Vibration 1992;52:107–23.
[11] Chen WF. Constitutive equations for engineering materials, vol. 2: Plasticity and modelling. Elsevier; 1994. [12] Kerwin EM. Damping of flexural waves by a constrained viscoelastic layer. J Acoust Soc Am 1959;31(7):952–62. [13] DiTaranto RA, Blasingame W. Composite damping of vibrating sandwich beams. J Eng Ind 1967;89(B):633–8. [14] Mead DJ, Markus S. The forced vibration of three-layer damped sandwich beam with arbitrary boundary conditions. J Sound Vibration 1969;10: 163–175. [15] Yan MJ, Dowell EH. Governing equations for vibrating constrained-layer damping sandwich plates and beams. J Appl Mech 1972;94:1041–7. [16] Oravsky V, Markus S, Simkova O. New approximate method of finding the loss factors of sandwich cantilever. J Sound Vibration 1974;33:335–52. [17] Rao DK. Frequency and loss factor of sandwich beams under various boundary conditions. J Mech Eng Sci 1978;20(5):271–82. [18] Sadek EA. Dynamic optimisation of a sandwich beam. Comput Struct 1984; 1984(4):605–15. [19] Rikards R, Chate A, Barkanov E. Finite element analysis of damping the vibrations of laminated composites. Comput Struct 1993;47(6):1005–15. [20] Johnson CD, Kienholz DA, Rogers LC. Finite element prediction of damping in beams with constrained viscoelastic layer. Shock Vib Bull 1981;51(1):71–81. [21] LU YP, Killian JW, Everstine GC. Vibrations of three layered damped sandwich plate composites. J Sound Vibration 1979;64(1):63–71. [22] Sainsbury MG, Zhang QJ. The galerkin element method applied to the vibration of damped sandwich beams. Comput Struct 1999;71:239–56. [23] Ramesh TC, Ganesan N. Finite element analysis of conical shells with a constrained viscoelastic layer. J Sound Vibration 1994;171(5):577–601. [24] Duigou L, Daya EM, Potier-Ferry M. Iterative algorithms for non-linear eigenvalue problems. Application to vibrations of viscoelastic shells. Comput Methods Appl Mech Eng 2003;192:1323–35. [25] Ferry JD. Viscoelastic properties of polymers. New-York: Wiley; 1980. [26] Gopinathan SV, Varandan VV, Varandan VK. A review and critique of theories for piezoelectric laminates. Smart Mater Struct 2000;9:24–48. [27] Cochelin B, Damil N, Potier-Ferry M. The Asymptotic-Numerical Method: An efficient pertubation technique for non-linear structural mechanics. Revue Européenne des éléments finis 1994;3(2):281–97. [28] Cochelin B. Mthode asymptotique-numrique pour le calcul non-linaire gomtrique des structures lastiques. 1994. [29] Thompson JMT, Walker AC. The non-linear perturbation analysis of discrete structure systems. Internat J Solids Struct 1968;4:219–59. [30] Najah A, Cochelin B, Damil N, Potier-Ferry M. A critical review of Asymptotic Numerical Methods. Arch Comput Methods Eng 1998;5:31–50. [31] Cochelin B, Damil N, Potier-Ferry M. Asymptotic-Numerical Methods and Padé approximants for non-linear elastic structures. Internat J Numer Methods Eng 1994;37:1187–213. [32] Elhage-Hussein A, Potier-Ferry M, Damil N. A numerical continuation method based on padé approximants. Internat J Solids Struct 2000;37:6981–7001. [33] Boudaoud H, Belouettar S, Daya EM, Potier-Ferry M. A shell element for active–passive vibration control of composite structures with piezoelectric and viscoelastic layers. Mech Adv Mater Struct 2008;15(3):208–19. [34] Abaqus Standard. Version 5.8. Pawtucket (RI): Hibbitt, Karlsson & Sorensen, Inc.; 1998. [35] Damil N, Potier-Ferry M, Najah A, Chari R, Lahmam H. An iterative method based upon padé approximants. Commun Numer Methods Eng 1999;15: 701–8.
Contents lists available at ScienceDirect
Engineering Structures journal homepage: www.elsevier.com/locate/engstruct
A numerical method for nonlinear complex modes with application to active–passive damped sandwich structures Hakim Boudaoud a,b , Salim Belouettar a,∗ , El Mostafa Daya b , Michel Potier-Ferry b a
Centre de Recherche Public Henri Tudor, 29, Avenue John F. Kennedy, L-1855 Luxembourg, Luxembourg
b
LPMM, UMR CNRS 7554, I.S.G.M.P., Université Paul Verlaine-Metz, Ile du Saulcy, F-57045 Metz Cedex 01, France
article
info
Article history: Received 24 January 2008 Received in revised form 16 June 2008 Accepted 5 August 2008 Available online 23 September 2008 Keywords: Sandwich structures Viscoelasticity Piezoelectricity Active Passive Hybrid Finite element method Vibrations Nonlinear eigenvalues Perturbation technique Homotopy Asymptotic numerical method
a b s t r a c t In this paper, a numerical method is proposed for determining complex vibrations modes of sandwich structures with piezoelectric and viscoelastic layers. Based on homotopy and asymptotic numerical techniques, this method leads to the damping properties calculation (loss factor and natural frequency per mode) of the hybrid sandwich structure. The numerical results of the loss factor and natural frequency are compared to those obtained from analytical beam model and from numerical studies the modal strain energy method. © 2008 Elsevier Ltd. All rights reserved.
1. Introduction In many industrial and defence applications, noise and vibration are important problems. Control of sound and vibration has been and is always the subject of a lot of research in recent years, and examples of applications are now numerous [1,2]. The three most common classifications of vibration control are: (1) passive; (2) active and (3) hybrid passive–active control. Passive control involves the use of reactive or resistive devices that either load the transmission path of the disturbing vibration or absorb vibratory energy. Active control also loads the transmission path but achieves this loading through the use of actuators that generally require external energy. In recent years, researchers have been directed to simultaneous use of piezoelectric and viscoelastic materials, in the so called ‘‘hybrid vibration control’’. A large part of studies highlighted structure with an hybrid arrangement called ACLD (SCLD) for active (smart) constrained layer damping. In many theoretical and experimental investigations the active
∗
Corresponding author. Tel.: +352 54 55 80 530; fax: +352 42 59 91 333. E-mail address: [email protected] (S. Belouettar).
0141-0296/$ – see front matter © 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.engstruct.2008.08.008
system increases the passive transverse shear deformation of the viscoelastic constrained layer and appears to be an effective means to control vibration [3,4]. To provide reliable and effective damping treatments, other interesting configurations are studied [5,6]. Typically, these structures are five layers sandwich ones (Fig. 1) in which a thin and soft viscoelastic layer is confined between identical elastic and stiff layers and two thin bonded piezoelectric layers. Much investigation has been devoted to the linear dynamic analysis of these structures. Their stiffness matrix is a complex one and depends nonlinearly on the vibration frequency. So, the vibration analysis leads to a non-linear complex eigenvalue problem. Solving this problem yields, first complex modes that can be slightly different from undamped modes, second complex eigenvalues whose real and imaginary parts are related to the loss factors and to the viscous eigen frequencies. From an engineering point of view, the most relevant quantity is the loss factor that is associated with each mode. A lot of methods have been presented to predict these structural damping properties and a detailed review of these works can be found in Reference [7]. Previously, there was no direct method available for solving such a problem efficiently in the case of large-scale structures. Often the stiffness matrix is assumed to be constant in a frequency
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
285
Fig. 1. Five layered sandwich structure with viscoelastic middle layer and two piezoelectric layers.
the asymptotic method of reference by [11], but more general since the introduction of the continuation procedure permits us to apply it to any case, and not only for small damping. The proposed investigation is an extension of a previous work of [6] for viscoelastically damped sandwich structures and the work of [24] for the nonlinear vibrations analysis of viscoelastically damped sandwich shells. 2. Material modelling 2.1. Piezoelectric constitutive equations
Fig. 2. Detailed Sketch of the considered five layered sandwich structure.
The electroelastic response of a piezoelectric body of volume Ω and regular boundary surface S, is governed by the mechanical and electrostatic equilibrium equations,
band, leading to eigenvalue problems with constant complex matrices. Sometimes, this complex eigenvalue problem is solved by the QR method [8], but this does not work for structures with a large number of degrees of freedom and can only used when the [K (ω)] is a constant complex matrix. Another way consists of computing the response to harmonic forcing for many frequencies, to yield resonance peaks from which loss factors and eigenfrequencies can be determined [9]: this can be applied to any structural problem, but the computational cost is very high. Approximate techniques such as the modal strain energy method [9] can also be used, but only to estimate the loss factor. An asymptotic technique was proposed in [10], but the procedure was limited to a small order of truncation and did not use a continuation algorithm, so this method cannot be used for large damping or for a strongly varying stiffness modulus. Recently, the algorithm proposed by Chen et al. [11] has proved much more satisfactory, because the nonlinear complex eigenvalue problem is solved by coupling a first order perturbation, an iterative algorithm and a reduced basis technique. In this way, the problem can be solved accurately, even for large-scale structures. Initially, analytical methods were developed to yield approximate loss factors and natural frequencies of sandwich beams or plates with simple boundary conditions [12–18]. In practice, it is necessary to design sandwich structures with a complex geometry and generic boundary conditions by using finite element simulations. Finite elements are still being used for analysis of damped sandwich structures [9–11,19–23] with the help of the aforementioned computational methods. This paper aims to address this deficiency and propose an efficient numerical method able to solve accurately the nonlinear complex eigenvalue problem presented above with a moderate computational cost. The proposed method is convenient for large-scale problems and applicable with any constitutive relation. The presented numerical methodology is similar to
σij = Cijkl εkl − ekij Ek Dk = ekij εij + kl El
(1)
where σ and are the second rank stress and strain tensors, E and D are, respectively, the electric field and the electric displacement vectors, ε is the second rank dielectric tensor, e is the third rank coupling tensor and C is the fourth rank elasticity tensor. The superscripts E and e indicate that the elasticity and permittivity constants are determined under conditions of zero or constant electric field and strain, respectively. We assume that no body electric charges are applied to the piezoelectric layers. The electric potential layer-wise assumption implies that in each of the layers the electric potential scalar field is represented by an independent function φ . A linear type variation is chosen to represent the electric potential in each piezoelectric layer, which implies that two electric potential nodes are needed for each lamina. The electric field-electric potential relation is the following: Ei = −(∇ φ)i .
(2)
2.2. Viscoelastic material constitutive modelling A general linear viscoelastic constitutive law is assumed for the core. Classically, such a law involves convolution product, [25]. For instance the relation between the axial stress and strain can be written in the form of a convolution product:
Z
t
Y (t − τ ) −∞
∂ ∂ (x, τ )dτ = Y ∗ ∂τ ∂t
(3)
where Y is the relaxation function of the viscoelastic material. Considering that the structure is subjected to harmonic vibrations, the viscoelastic behaviour can be described by a complex Young modulus which is depending nonlinearly on the frequency and temperature [25].
286
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
3. Finite element formulation From the application of the principle of virtual work, the following equilibrium equation is obtained
δ Wext + δ Ec − δ Eh = 0
(4)
where δ Wext denotes the virtual work of the external loading, δ Ec represents the kinematical energy and δ Eh is the internal electromechanical energy. Ec and Eh are given as:
Z 1X 2 2 2 ˙ ˙ ˙ d Ωk + W + V U ρ E = k k k k c 2 k Ωk k = s, e1 , Zc , e2 , a Z X 1X E = 1 k k σij ij dΩk − Eil Dli dΩl h 2
2
Ωk
k
Ωl
l
(5) l = s , a.
We presume that the piezoelectric layers are thin enough to consider the hypothesis of constant electric field through the thickness of each piezoelectric lamina and nile in the other directions as valid. We also consider as in [26] the central to layer to be conductive with a uniform potential fixed to zero. Thus, the Eq. (2) takes the following form: E3 = −
(6)
hp
e i
= Bu u (7) i where Bu is the product of the differential operating matrix i i S
n o {Fu (x, y, z , t )} = F˜u (x, y, z ) eiωt {u (x, y, z , t )} = {U (x, y, z )} eiωt {φ (x, y, z , t )} = {Φ (x, y, z )} eiωt
(11)
√
where ω is the frequency and i = −1. Substituting Eq. (11) in Eq. (9) and considering separately the sensor and actuator potentials, the global dynamic finite element equations can be rewritten as:
"K (ω) uu Kφ u s Kφ u a ( ) Fu = 0
Ku φ a 0 Kφφ a
Kuφ s Kφφ s 0
#
Muu 0 0
" −ω
2
0 0 0
0 0 0
#! (
U
)
Φs Φa
(12)
0
φ
where φ denotes the superficial electrical potential and hp is the thickness of the piezoelectric layer. Introducing the strain-displacement relation in terms of the nodal displacement lead, in each i layer, to the discretized expression of the mechanical strain:
i
The matrices Bi depend only on the shape functions of the considered finite element and Ci are the matrices obtained from the strain stress laws. [Muu ] is assumed to be constant in this analysis. In order to estimate the damping properties (frequency and loss factor per vibration mode) of the investigated adaptive structure when subjected to harmonic loads, a frequency response analysis is performed. The applied force, the displacement field and the electrical potential are expressed as:
relating S to the shape function matrix Nu . Similarly we get the same relations for the piezoelectric layer:
= −∇ Φ k = −∇ Nφ k φ e = − Bφ k φ e (8) k where Bφ is the product of the differential operating matrix relating the electrical field E k and the shape function matrices k
T
where Kuφ = Kφ u . From the previous equation, the generated potential on the sensor is given by:
−1 {Φs } = − Kφφ s Kφ u s {U } .
(13)
In the proposed analysis, a velocity and a proportional feedback control laws are considered. According to these laws, the controlled voltage for each layer is expressed as:
φ a = gd φ s + gv φ˙ s .
(14)
gd and gv are the feedback control gain and the velocity feedback gain. The Eq. (14) permits to write:
k E
Nφ . Substituting Eqs. (7) and (8) into Eq. (4) one obtains the equilibrium equation for an element, which after taking into account the contribution of all elements in the domain, yields the equilibrium equations of the system obtained after performing the usual finite element assembling techniques.
[K ] {u} + Kuφ {φ} + [Muu ] {¨u} = {Fu } uu T Kuφ {u} + Kφφ {φ} = {0} .
(9)
The quantities [Muu ], [Kuu ] are the mass and the stiffness matrices used in the structural analysis, Kuφ is electrical
mechanical coupled stiffness matrix, Kφφ is the electrical stiffness matrix, {Fu } is the force vector of the mechanical loads applied to the structure.
XZ t [Muu ] = ρ Nu i Nu i dΩ vi i X Z t Bu i C i Bu i dΩ [Kuu ] = v i i Z k X k t t k K = Bu [e] Bφ dΩ uφ k Zvk k t k X Kφφ = [ ] Bφ k dΩ Bφ n
vk
−1 {Φa } = − (gd + iωgv ) Kφφ s Kφ u s {U } .
(15)
Substituting the Eqs. (13) and (15) into Eq. (12), we get the dynamic finite element equations with condensed electrical degrees of freedom. [Kreal (ω)] + i Kcomp (ω) − ω2 [Muu ] {U } = {Fu }
(16)
where the complex nodal vibration eigenmode is denoted by U and with the following definition for the real and imaginary parts:
] = Re {[Kuu (ω)]} [Kreal (ω) −1 − Kuφ s + gd Kuφ a Kφφ s Kφ u s −1 Kcomp (ω) = Im {[Kuu (ω)]} − ωgv Kuφ a Kφφ s Kφ u s for simplicity the real and imaginary part of the complex stiffness matrix Re {R[K uu (ω) ]}I and Im {[Kuu (ω)]} will be respectively denoted Kuu and Kuu . 4. Application of the asymptotic numerical method 4.1. Homotopy technique
i = e 1 , e 2 , s , a, c k = s , a.
(10)
The object of this first stage is to transform the nonlinear eigenvalue problem using the homotopy technique [27]. This allows us to introduce artificially a parameter ε in Eq. (16) as:
[Kreal (ε)] + iε Kcomp − ω2 [Muu ] {U } = 0 0≤ε≤1
(17)
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
287
I {U0 }t −gd [Kactuator ]t + i Kuu − p0 gv [Kactuator ]t {U0 } C1 = {U0 }t [Muu ] {U0 } k=j−1 j −1 t k=P t I t P pk [Kactuator ]t {U0 } Kuu − gd [Kactuator ]t − igv (Uj−k−1 ) − Ck Uj−k [Muu ] + i Uj−1 k=0 k=1 Cj = t { } { } ] [M U U 0 uu 0 j = 2, 3, . . . Box I.
with
(
−1 R [Kreal (ε)] = Kuu − Kuφ s + ε gd Kuφ a Kφφ s Kφ u s I −1 Kcomp = Kuu − ωgv Kuφ a Kφφ s Kφ u s .
One can easily note that the initial problem (Eq. (16)) corresponds to ε = 1. The system corresponding to ε = 0 is the vibration problem of an undamped structure which can be written in the following form
R Kuu − Kuφ
s
Kφφ
−1 s
Kφ u
s
− ω2 [Muu ] {U0 } = 0.
(18)
This equation represents a real eigenvalue problem that can be solved by classical methods such as subspace iteration or the Lanczos method, and we suppose that we have computed the required eigenpairs of (18). To get the solutions of the nonlinear eigenvalue problem (16), we try to solve the Eq. (17) deduced by the homotopy transformation. In fact, the Eq. (17) does not define a well-posed problem, because the unknown can be multiplied by an arbitrary constant. To obtain a well-posed problem, it is sufficient to require an additional condition, for instance
{U0 }t {U − U0 } = {0} .
(19)
The Eqs. (16) and (18) define a well-posed path following problem ε → {{U (ε)} , ω(ε)}, having to be solved by an appropriate method. We could use an incremental-iterative method like Newton–Raphson or the modified Newton method, but we preferred an asymptotic-numerical method, permitting us to restrict the number of matrices to be inverted. 4.2. Perturbation technique At this stage, we describe the asymptotic-numerical method being applied to solve the perturbed eigenvalue problem (16) and (18). This technique has been extensively studied for other nonlinear problems, such as nonlinear elasticity, Navier–Stokes equations and structural problems involving unilateral contact, plasticity or viscoplasticity [28]. The asymptotic numerical method associates a perturbation technique with a classical discretization procedure, which is generally the finite element method. It was first proposed by [29] and we refer to [30] for a complete reference list. Within asymptotic-numerical techniques, one searches for the unknowns {U } and ω in form of truncated integro-power series with respect to the parameter ε :
N X { } U = ε j Uj j =0 N X ω= ε j pj j =0 N X 2 ω = ε j Cj j=0 0≤ε≤1
(20)
where N is the truncation order of the series and ω the complex eigenfrequency. In this analysis, the modulus E (w) is assumed to be dependent on frequency and after insertion of the series (20) in the Eq. (17), one gets a set of recurrent linear problems that coming from the identification of the like powers of ε : The first order problem:
[A] {U0} = 0 R [A] = Kuu − [Ksensor ] − C0 [Muu ] (21) 2 C0 = ω0 −1 where [Ksensor ] = Kuφ s Kφφ s Kφ u s and [Kactuator ] = Kuφ a −1
Kφφ s Kφ u s are respectively the added static stiffness due to sensor and actuator. The second order problem:
[A] {U1 } = {F1 } I {F } = g [K ] − i Kuu − p0 gv [Kactuator ] 1+ C [Md ]actuator 1 uu ) {U0 }.
(22)
The (j + 1)th order problem:
[A] Uj = Fj k=j X I Fj = Ck [Muu ] Uj−k − i Kuu − gd [Kactuator ] Uj−1 k=1 kX =j −1 pk [Kactuator ] Uj−k−1 . +igv
(23)
k=0
Notice that for the first order, we recover the undamped eigenvalue problem and we propose to compute the required eigenpairs by any classical method for solving a symmetric eigenvalue problem. The problems for the higher orders are nonhomogeneous linear ones having the same matrix. Therefore, since the method combines the classical and asymptotic-numerical approaches, we can compute many terms of the series by inverting only one matrix per frequency. In the present case, the method has an additional advantage because the vectors in Box I are complex, while the matrix to be inverted is real. Note also that the matrix is singular and its kernel is generated by the undamped mode. Then the problems (Box I) have a solution if and only their right-hand sides j = 1 to N satisfy the following solvability condition:
t Fj
{U0 } = 0.
(24)
This condition permits us to compute the term Cj in the expansion (Eq. (20)) of the eigenvalue as the equation in Box I. We note that this number depends only on the previously computed vectors {Un } and Cn obtained in the previous order. So, we can compute C1 by Box I, next C2 by solving Box I, next C2 again by Box I and so on. If the solvability condition is satisfied, the linear equations (22) and (23) are well-posed and they have a unique solution if one takes into account the additional condition (25), that leading to:
{U0 }t Uj = 0 j = 1, . . . , N .
(25)
288
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
−1
Re {[Kuu (ω)]} + Kuφ s + ε Gd Kuφ a Kφφ s Kφ u s + iε [C (ω)] − ω2 [Muu ] {U (ε)}
R (ε) = k(Re {[Kuu (ω)]} {U (ε)})k N X {U } = Uj j =0 v u N uX Cj ω = it j =0
Box II.
One can rearrange the equations Eqs. (23) and (25) by introducing a Lagrange multiplier k, and it will be easy to check that the solution of (23) and (25) satisfies the following matrix equation:
[A]
{U0 }t
{U0 } 0
Uj k
=
Fj 0
.
(26)
The matrix in [A] has an inverse and we can factorize it by the Crout method. This permits us to solve many linear problems (22) and (23) with a single matrix factorization. In this way, we get solutions of the perturbed problem (17) in the form:
N X {U (ε)} = ε j Uj v j =0 u N uX ω = i t ε j Cj .
followed in classical matrix form: σ C C12 C13 0 0 11 1 C C C 0 0 σ 12 22 23 2 σ 3 C13 C23 C33 0 0 0 0 0 C 0 σ 44 4 0 0 0 0 C σ5 = 55 0 0 0 0 0 σ6 0 0 0 e15 0 D1 D 0 0 0 e 0 2 24 D3
e31
e32
e33
0
0
0
0
0
0
0
0
0
0
0
0
0
−e24
0
−e15
0
C66
0
0
0
11
0
0
0
22
0
0
0
1 2 3 γ23 × γ13 γ12 E 1 E2
(27)
j=0
−e31 −e32 −e33 0 0 0 0 0
33
(28)
E3
This is only an estimate of the solution of (17) because the order of truncation N is necessarily finite (generally, N = 20 is chosen) and also because can be outside the radius of convergence of the series. The latter has to be greater than or equal to 1, because we have to solve the problem that corresponds to the case ε = 1. If the estimate is good enough, i.e. the relative residual Box II is smaller than a tolerance (such as 10−8 . . . 10−6 ), the solution of problem is obtained from (27) by putting ε = 1. Otherwise we have to improve the method by replacing the polynomial approximation by a rational one. This technique, named Padé approximants, allow us to increase the radius of convergence of the series. More details of this technique may be found in Cochelin et al. [31] and Elhage-Hussein et al. [32]. If the the radius of convergence of the series still fewer than one we have to use a continuation method [6] or iterative algorithms [24].
5. Numerical applications The present numerical procedure has been applied to the free vibrations analyses of a 5-layered simply supported sandwich beam, Figs. 1 and 2. The materials composing the elastic layers and the viscoelastic layers are characterised by the values of the Young modulus and the Poisson ratio: Ee = 45.54 GPa, ν = 0.33 for the elastic material and Ev = 7.25 MPa, ν = 0.45 for the viscoelastic core. The elastic, viscoelastic and piezoelectric materials have respectively the following densities: ρe = 2040 kg/m3 , ρv = 1200 kg/m3 and ρp = 7500 kg/m3 . The constitutive equations of the PZT4 orthotropic piezoelectric layers may be summarised as
where:
C11 = C22 = 139.0 GPa; C33 = 115.3 GPa C = C = 74 . 3 GPa ; C = 77.8 GPa 13 23 12 C = C = 25.6 GPa; C66 = 30.6 GPa 44 55 2 e33 = 15.1 C/m2 ; e31 = e32 = −5.21 C/m ; 2 e15 = e24 = 12.7 C−/8m 11 = 22 = 1.31e F/m; 33 = 1.15e−8 F/m.
(29)
In the application of the ANM, the truncation is taken as N = 20 and the relative residual R () is smaller than 10−4 . 5.1. Case of five layered sandwich beam We propose to analyse the modal properties of a five-layered simply supported beam of length L = 500 mm and high h = 10 mm. The thicknesses of the various layers composing the beam are the same as those used in [33] i.e. he = 69h/150, hv = h/15 and hp = h/150 respectively for the elastic, viscoelastic and piezoelectric layers. The viscoelastic properties of the central layer are introduced by a complex modulus that is assumed to be constant (see Table 1). Two tests are presented that mainly differ from one another by the control law i.e. a direct proportional feedback one used is the first test and velocity feedback control law is used in the second one. 5.1.1. Proportional direct control law The free vibration of the beam is investigated in the case of direct proportional feedback (gv = 0). Several values of the proportional feedback control gain (gd ) and different viscoelastic loss factors ηvisco values were considered (see Tables 1 and 2). This test has been investigated quite extensively in [33] from
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
289
Table 1 Proportional feedback control gd = −1, gv = 0: Natural frequencies and loss factors of simply supported sandwich beam for for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
95.2 96.5 96.5 95.5 96.8 96.8 96.2 97.5 97.5 97.4 98.7 98.7 99.0 100.2 100.3 100.8 102.0 102.1
0 0 0 0.0740 0.0741 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
298.9 303.1 303.1 299.1 303.5 303.3 299.7 304.6 304.0 300.7 306.6 305.5 302.2 309.5 306.5 304.0 313.4 308.4
0 0 0 0.057 0.057 0.057 0.113 0.114 0.114 0.168 0.169 0.168 0.221 0.223 0.221 0.272 0.275 0.275
612.5 622.2 622.2 612.6 622.6 622.3 613. 623.9 622.3 613.6 626.1 623.3 614.5 629.2 624.2 615.6 633.4 625.4
0 0 0 0.035 0.035 0.035 0.07 0.07 0.07 0.105 0.105 0.105 0.14 0.139 0.14 0.174 0.173 0.174
Table 2 Proportional feedback control gd = −20, gv = 0: Natural frequencies and loss factors of simply supported sandwich beam for for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
89.8 91.0 91.0 90.1 91.3 91.3 90.8 92.0 92.0 91.9 93.1 93.1 93.4 94.6 94.7 95.1 96.3 96.4
0 0 0 0.0740 0.0741 0.074 0.143 0.144 0.143 0.204 0.205 0.204 0.254 0.256 0.254 0.293 0.296 0.293
280.8 284.8 284.8 281.0 285.2 285.1 281.6 286.35 285.7 282.7 288.3 286.8 284.1 291.2 288.2 286.0 295.2 290.1
0 0 0 0.059 0.059 0.059 0.117 0.117 0.117 0.173 0.174 0.174 0.227 0.230 0.227 0.279 0.282 0.279
573.4 582.6 582.6 573.5 583.0 582.8 573.9 584.4 583.2 574.6 586.6 583.8 575.5 589.8 584.7 576.6 594.1 585.9
0 0 0 0.037 0.037 0.037 0.074 0.074 0.074 0.11 0.11 0.11 0.146 0.145 0.146 0.181 0.18 0.181
Table 3 Proportional feedback control: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco and control gain gd (gv = 0)
ηvisco
0
0.2
0.4
0.6
0.8
1
gd = −1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
gd = −20
gd = −150
f1 in (Hz)
η1
f1 in (Hz)
η1
f1 in (Hz)
η1
95.2 96.5 96.5 95.5 96.8 96.8 96.2 97.5 97.5 97.4 98.7 98.7 99.0 100.2 100.3 100.8 102.0 102.1
0 0 0 0.0740 0.0741 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
89.8 91.0 91.0 90.1 91.3 91.3 90.8 92.0 92.0 91.9 93.1 93.2 93.4 94.6 94.7 95.1 96.3 96.4
0 0 0 0.074 0.0741 0.074 0.143 0.144 0.143 0.204 0.205 0.204 0.254 0.256 0.254 0.279 0.282 0.279
44.8 45.5 45.6 45.1 45.8 45.8 45.8 46.5 46.5 46.8 47.6 47.6 48.1 49.0 48.9 49.6 51.3 50.4
0 0 0 0.12 0.12 0.119 0.23 0.23 0.23 0.32 0.31 0.31 0.38 0.39 0.376 0.42 0.41 0.42
theoretical points of view. In [33], the analytical results are validated by using a 2D finite element model on Abaqus. The FE mesh of the sandwich beam is composed of 240 CPE8R: 8-node plane strain elements with a reduced integration. A steady-state linear dynamic analysis procedure [34] is used to monitor the
linear response of the beam structure subjected to a continuous harmonic excitation. In Table 3 are presented the first modal property (frequency and damping) obtained using the CSME, the analytical model [33] and those using the ANM method.In all analysed situations, we observe a very good agreement between the results of the ANM
290
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
Table 4 Calculation time needed with CSME method and new ANM method considering different number of d.o.f. d.o.f.
Time using CSME
Time using ANM
600 1000 1500
2 min 3 s 6 min 33 s 21 min 48 s
46 s 3 min 1 s 7 min 45 s
using three methods: the analytical approach by [33], the CSME and the asymptotic numerical method. Beyond the agreement between the results of the proposed method and those obtained analytically or by the CSME, the results of the ANM are valid for higher values of damping in the core or higher values of the control gains and for the higher modes (Tables 5–7). 6. Conclusion
method and those obtained analytically [33], those obtained by using the CSME (see tables). However, when the damping in the core becomes more significant or when the control increases, the values obtained with the ANM are in agreement with the analytical values whereas those obtained with the CSME tend not to be in agreement. In a second part we consider the necessary time to get the four first modal properties using the CSME and the ANM method. In Table 4 we observe, for 3 different meshes, that in each considered case the ANM method is less time-consuming than the CSME method. 5.1.2. Velocity feedback control law The second example has been designed to establish the ability of the method to analyse adaptive structures with a velocity feedback control. As for the previous case, different values of feedback control gain, gv , and damping were considered. In Table 5 are presented the results of the frequency and damping obtained by
In this paper, an asymptotic numerical method has been developed for solving an eigenvalue problem for a five layers hybrid damped sandwich structure. The method consists of introducing an artificially modified problem that depends on a parameter , and we seek eigenvalues and eigenmodes in the form of integro-power series expansions with respect to . The asymptotic numerical analysis permits us to decompose the perturbed problem into a sequence of linear problems. The first one is a real eigenvalue problem that can be solved by classical procedures. The other problems are non-homogeneous with the same linear operator. In this way, the method permits us to considerably simplify the analysis of the nonlinear eigenvalue problem. The method presented has three advantages, as compared with the classical techniques that are commonly used to estimate the loss factor of a viscoelastic structure. First, it solves the complex nonlinear eigenvalue problem accurately, so that the
Table 5 Velocity feedback control gd = 0, gv = 0: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
95.5 96.8 96.8 95.8 97.0 97.1 96.5 97.8 97.8 97.8 99.0 99.0 99.3 100.5 100.6 101.1 102.3 102.4
0 0 0 0.074 0.074 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
299.8 304.1 304.1 300.0 304.6 304.3 300.6 306.4 304.9 301.7 309.4 306.0 303.1 313.7 307.4 304.9 319.6 309.3
0 0 0 0.057 0.057 0.057 0.113 0.113 0.113 0.168 0.167 0.168 0.221 0.217 0.221 0.271 0.264 0.271
614.5 624.2 624.2 614.6 624.6 624.3 615.0 625.9 624.7 615.6 628.1 625.4 616.5 631.2 626.2 617.6 635.4 627.4
0 0 0 0.035 0.035 0.035 0.07 0.07 0.07 0.105 0.105 0.105 0.140 0.139 0.139 0.174 0.173 0.173
Table 6 Velocity feedback control gd = 0, gv = −0.01: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco
ηvisco 0
0.2
0.4
0.6
0.8
1
Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
f1 in (Hz)
η1
f2 in (Hz)
η2
f3 in (Hz)
η3
95.5 96.8 96.8 95.6 96.9 96.8 96.1 97.5 97.4 97.1 98.5 98.4 98.5 100. 99.8 100.1 101.7 101.5
0.036 0.036 0.036 0.11 0.11 0.11 0.18 0.181 0.18 0.243 0.243 0.243 0.295 0.296 0.296 0.337 0.338 0.337
299.5 305.5 302.9 298.4 306.7 301.8 297.8 309.1 301.1 297.5 312.8 300.9 297.7 318.0 301.1 298.0 325.0 301.7
0.117 0.119 0.119 0.175 0.179 0.176 0.233 0.241 0.234 0.291 0.302 0.292 0.347 0.36 0.348 0.402 0.422 0.402
610.2 642.0 612.1 608.7 643.5 610.5 605.7 647.4 607.6 603.0 652.3 604.9 600.5 658.3 602.5 598.0 665.7 600.4
0.28 0.286 0.278 0.30 0.31 0.30 0.33 0.355 0.33 0.37 0.404 0.373 0.412 0.455 0.412 0.45 0.51 0.45
H. Boudaoud et al. / Engineering Structures 31 (2009) 284–291
291
Table 7 Velocity feedback control: Natural frequencies and loss factors of simply supported sandwich beam for various values of ηvisco and control gain gv (gd = 0)
ηvisco
0
0.2
0.4
0.6
0.8
1
gv = −0.01
gv = 0 Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM Anal CSME ANM
gv = −0.05
f1 in (Hz)
η1
f1 in (Hz)
η1
f1 in (Hz)
η1
95.5 96.8 96.8 95.8 97.0 97.0 96.5 97.8 97.8 97.8 99.0 99.0 99.3 100.5 100.6 101.1 102.3 102.4
0 0 0 0.074 0.074 0.074 0.144 0.144 0.144 0.205 0.206 0.205 0.256 0.258 0.256 0.296 0.299 0.296
95.5 96.8 96.8 95.6 96.9 96.8 96.1 97.5 97.4 97.1 98.5 98.4 98.5 100 99.8 100.1 101.7 101.5
0.036 0.036 0.036 0.11 0.11 0.11 0.18 0.185 0.181 0.243 0.243 0.243 0.295 0.296 0.296 0.337 0.338 0.337
95.2 96.5 95.8 94.6 96.1 95.1 94.2 96.1 94.9 94.5 96.6 95.2 95.2 97.6 95.8 96.2 98.8 96.9
0.180 0.181 0.18 0.258 0.258 0.258 0.334 0.332 0.333 0.404 0.4 0.401 0.464 0.458 0.461 0.514 0.505 0.51
computational technique is still valid for large damping and higher values of control gain. Secondly it is able to account for a strong variation of the modulus with frequency, as is often the case for polymers. Thirdly, the method remains efficient for problems involving many degrees of freedom, because it is only necessary to solve a classical eigenvalue problem and to decompose a few matrices. Nevertheless, the present method can probably be significantly improved. Indeed, there are many ways to set up a homotopy transformation and, using the one presented in [35], one could avoid decomposing complex matrices. Furthermore, the power series do not yield the best possible representation of a solution path and the algorithm becomes much more efficient if one introduces Padé approximants [35,30]. Acknowledgements This work has been done in the framework of the FNR-MAT08-01 ADYMA project of the the European ERA-NET Materials (MATERA) and CASSEM (Contact N: FP6-NMP3-CT-2005-013517). The authors acknowledge the financial support of the Luxembourgish National Research Fund (www.fnr.lu) and the European Commission. References [1] Nakra BC. Vibration control with viscoelastic materials. Shock Vib Dig 1976;8: 3–12. [2] Hu H, Belouettar S, Daya EM, Potier-Ferry M. Review and assessment of various theories for modeling sandwich composites. Compos Struct 2008; 84(3):282–92. [3] Vasques CMA, Mace BR, Gardonio P, Dias Rodrigues J. Arbitrary active constrained layer damping treatments on beams: Finite element modelling and experimental validation. Comput Struct 2006;84:1384–401. [4] Vasques CMA, Dias Rodrigues J. Combined feedback/feedforward active control of vibration of beams with acld treatments: Numerical simulation. Comput Struct 2008;86:292–306. [5] Agnes GS, Napolitano K. Active constrained layer viscoelastic damping. In: 34th AIAA/ASME/ASCE/AHS/ASC Struct., structural dyn., and mater. conf., AIAA. 1993. p. 3499–506. [6] Daya EM, Potier-Ferry M. A numerical method for nonlinear eigenvalue problems: Application to vibrations of viscoelastic structures. Comput Struct 2001;79:533–41. [7] Ganapathi M, Patel BP, Boisse P, Polit O. Flexural loss factors of sandwich and laminated composite beams using linear and non linear dynamic analysis. Compos Part B: Eng 1999;30:245–56. [8] Bathe KJ. Finite element procedures in engineering analysis. NJ: Prentice-Hall; 1982. [9] Soni ML. Finite element analysis of viscoelastically damped sandwich structures. Shock Vib Bull 1981;55(1):97–109. [10] Ma BA, He JF. Finite element analysis of viscoelastically damped sandwich plates. J Sound Vibration 1992;52:107–23.
[11] Chen WF. Constitutive equations for engineering materials, vol. 2: Plasticity and modelling. Elsevier; 1994. [12] Kerwin EM. Damping of flexural waves by a constrained viscoelastic layer. J Acoust Soc Am 1959;31(7):952–62. [13] DiTaranto RA, Blasingame W. Composite damping of vibrating sandwich beams. J Eng Ind 1967;89(B):633–8. [14] Mead DJ, Markus S. The forced vibration of three-layer damped sandwich beam with arbitrary boundary conditions. J Sound Vibration 1969;10: 163–175. [15] Yan MJ, Dowell EH. Governing equations for vibrating constrained-layer damping sandwich plates and beams. J Appl Mech 1972;94:1041–7. [16] Oravsky V, Markus S, Simkova O. New approximate method of finding the loss factors of sandwich cantilever. J Sound Vibration 1974;33:335–52. [17] Rao DK. Frequency and loss factor of sandwich beams under various boundary conditions. J Mech Eng Sci 1978;20(5):271–82. [18] Sadek EA. Dynamic optimisation of a sandwich beam. Comput Struct 1984; 1984(4):605–15. [19] Rikards R, Chate A, Barkanov E. Finite element analysis of damping the vibrations of laminated composites. Comput Struct 1993;47(6):1005–15. [20] Johnson CD, Kienholz DA, Rogers LC. Finite element prediction of damping in beams with constrained viscoelastic layer. Shock Vib Bull 1981;51(1):71–81. [21] LU YP, Killian JW, Everstine GC. Vibrations of three layered damped sandwich plate composites. J Sound Vibration 1979;64(1):63–71. [22] Sainsbury MG, Zhang QJ. The galerkin element method applied to the vibration of damped sandwich beams. Comput Struct 1999;71:239–56. [23] Ramesh TC, Ganesan N. Finite element analysis of conical shells with a constrained viscoelastic layer. J Sound Vibration 1994;171(5):577–601. [24] Duigou L, Daya EM, Potier-Ferry M. Iterative algorithms for non-linear eigenvalue problems. Application to vibrations of viscoelastic shells. Comput Methods Appl Mech Eng 2003;192:1323–35. [25] Ferry JD. Viscoelastic properties of polymers. New-York: Wiley; 1980. [26] Gopinathan SV, Varandan VV, Varandan VK. A review and critique of theories for piezoelectric laminates. Smart Mater Struct 2000;9:24–48. [27] Cochelin B, Damil N, Potier-Ferry M. The Asymptotic-Numerical Method: An efficient pertubation technique for non-linear structural mechanics. Revue Européenne des éléments finis 1994;3(2):281–97. [28] Cochelin B. Mthode asymptotique-numrique pour le calcul non-linaire gomtrique des structures lastiques. 1994. [29] Thompson JMT, Walker AC. The non-linear perturbation analysis of discrete structure systems. Internat J Solids Struct 1968;4:219–59. [30] Najah A, Cochelin B, Damil N, Potier-Ferry M. A critical review of Asymptotic Numerical Methods. Arch Comput Methods Eng 1998;5:31–50. [31] Cochelin B, Damil N, Potier-Ferry M. Asymptotic-Numerical Methods and Padé approximants for non-linear elastic structures. Internat J Numer Methods Eng 1994;37:1187–213. [32] Elhage-Hussein A, Potier-Ferry M, Damil N. A numerical continuation method based on padé approximants. Internat J Solids Struct 2000;37:6981–7001. [33] Boudaoud H, Belouettar S, Daya EM, Potier-Ferry M. A shell element for active–passive vibration control of composite structures with piezoelectric and viscoelastic layers. Mech Adv Mater Struct 2008;15(3):208–19. [34] Abaqus Standard. Version 5.8. Pawtucket (RI): Hibbitt, Karlsson & Sorensen, Inc.; 1998. [35] Damil N, Potier-Ferry M, Najah A, Chari R, Lahmam H. An iterative method based upon padé approximants. Commun Numer Methods Eng 1999;15: 701–8.