# Vibration of axisymmetric composite piezoelectric shells coupled with internal fluid

Walid Larbi, Jean-François Deü, Roger Ohayon

► To cite this version:

Walid Larbi, Jean-François Deü, Roger Ohayon. Vibration of axisymmetric composite piezoelectric shells coupled with internal fluid. International Journal for Numerical Methods in Engineering, 2007, 71, pp.1412 - 1435. 10.1002/nme.1987. hal-01422438

HAL Id: hal-01422438

https://hal.science/hal-01422438

Submitted on 25 Dec 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Vibration of axisymmetric composite piezoelectric shells coupled with internal fluid

W. Larbi, J.-F. Deu and R. Ohayon

*Structural Mechanics and Coupled Systems Laboratory, Conservatoire National des Arts et Métiers, 292 rue Saint-Martin, 75141 Paris Cedex 03, France*

This paper presents the theoretical and finite element formulations of piezoelectric composite shells of revolution filled with compressible fluid. The originality of this work lies (i) in the development of a variational formulation for the fully coupled fluid/piezoelectric structure system, and (ii) in the finite element implementation of an inexpensive and accurate axisymmetric adaptive laminated conical shell element. Various modal results are presented in order to validate and illustrate the efficiency of the proposed fluid-structure finite element formulation.

**KEY WORDS:** axisymmetric conical shells; piezoelectric laminates; fluid-structure interaction; free vibrations; finite element method

## 1. INTRODUCTION

In the context of noise reduction techniques by active treatments, this work presents the theoretical formulation and the finite element implementation of the specific problem of piezoelectric composite shells of revolution filled with compressible fluid. The originality of this work lies (i) in the modal formulation of the fully coupled fluid/piezoelectric structure system, and (ii) in the development of an inexpensive and accurate adaptive axisymmetric laminated conical shell element.

For the resolution of structural-acoustic vibration problems, two classes of methods can be distinguished according to whether they are based on analytical or numerical discrete approaches. Analytical methods are classically restricted to simple geometries and boundary conditions. For
example, the vibration of cylindrical thin shells filled with a compressible fluid has been studied using a wave propagation approach in [1] for elastic structures, or using a state-space method in [2] for piezoelectric structures. Three-dimensional exact solutions have recently been proposed for the free vibration analysis of arbitrary thick simply supported piezoelectric hollow cylinders filled with compressible fluid [3, 4]. Methods based on discretization of the structural-acoustic governing equations, such as the finite element method, are not restricted to specific boundary conditions but require a higher computational effort. For example, a mixed formulation based on a piezoelectric shell finite element for the structure and a three-dimensional boundary element for the fluid, is developed in [5] to calculate the coupled response of smart structural-acoustic systems. For problems involving axisymmetric geometries, semi-analytic approaches are frequently used because they are computationally more efficient and accurate than fully discrete methods approaches. These methods, which employed a discretization in Fourier series for the circumferential direction, have recently been developed for vibrations analysis of axisymmetric piezoelectric shells [6–8], as well as axisymmetric laminated composite cylinders filled with fluid [9–11].

In this paper, a semi-analytical approach is used for adaptive conical laminated shells filled with compressible fluid. The piezoelectric shell model is based on the Kirchhoff–Love theory and is supposed to be polarized in thickness direction. The associated finite element model combines an equivalent single-layer approach for the mechanical behaviour with a layerwise representation of the electric potential in the thickness direction. In this work, we consider that the fluid is inviscid, compressible and barotropic, gravity effect being neglected. Extension to vibration damping of more complex systems involving surface gravity effects in the inviscid incompressible case [12] or internal gravity effects in the inviscid compressible case [13] will be the purpose of further investigations. Considering the above assumptions, the fluctuating fluid forces are determined by means of the linearized Euler equations and reduced to the Helmholtz equation for the pressure. The dynamic formulation of the coupled system is derived from a variational principle involving shell displacements and rotation, electric potential and fluid pressure. Concerning the finite element discretization, a one-dimensional element with two nodes is used for the adaptive conical shell structure and a quadrangular element with four nodes for the fluid. Moreover, an appropriate interface element is developed to ensure the interaction between the fluid and structure. For the shell element, each node has (i) four mechanical degrees-of-freedom to describe axial, radial and circumferential displacements, as well as the rotation about the circumferential axis, and (ii) one electrical degree-of-freedom related to the sensor voltage per piezoelectric layer. For the fluid element, the pressure is the unique nodal unknown. It should be mentioned that the variables associated with these two elements are expanded in Fourier series in the circumferential direction. Finally, some numerical examples are presented in order to validate and demonstrate the effectiveness of the proposed formulation.

## 2. FINITE ELEMENT FORMULATION OF THE PIEZOELECTRIC STRUCTURE COUPLED WITH AN INTERNAL FLUID

Let us consider the linear vibrations of an elastic/piezoelectric structure completely filled with a homogeneous, inviscid and compressible fluid, neglecting gravity effects. We establish in this section the variational formulation of the coupled problem and the corresponding matrix equations.
<figure><img src="image_2.png" /><figcaption>Figure 1. Fluid/piezoelectric structure-coupled system.</figcaption></figure>

## 2.1. Governing equations of the spectral-coupled problem

We consider a piezoelectric structure occupying the domain $\Omega_S$ at the equilibrium. The structure is clamped on a part $\Gamma_u$ and traction free on the complementary part $\Gamma_\sigma$ of its external boundary. The interior fluid domain is denoted by $\Omega_F$ and the fluid-structure interface by $\Sigma$ (see Figure 1). For the considered spectral problem, electric boundary conditions are defined by null electric potential on $\Gamma_\phi$ and null surface density of electric charge on the remaining part $\Gamma_D$. Thus, the total structure boundary is $\partial\Omega_S = \Gamma_u \cup \Gamma_\sigma \cup \Sigma = \Gamma_\phi \cup \Gamma_D$ with $\Gamma_u \cap \Gamma_\sigma \cap \Sigma = \Gamma_\phi \cap \Gamma_D = \emptyset$.

The linearized deformation tensor is denoted by $\varepsilon$ and the corresponding stress tensor by $\sigma$. Moreover, $\mathbf{D}$ denotes the electric displacement and $\mathbf{E}$ the electric field. $\rho_S$ is the mass density of the structure and $\mathbf{n}^S$ is the unit normal external to $\Omega_S$.

Since the compressible fluid is assumed to be inviscid, instead of describing its motion by a fluid displacement vector field, which requires an appropriate discretization of the fluid irrotationality constraint [14], we use the pressure scalar field $p$. Let us consider $c_F$ as the constant speed of sound in the fluid, $\rho_F$ the mass density of the fluid, and $\mathbf{n}$ the unit normal external to $\Omega_F$.

The local equations describing the spectral fluid-piezoelectric structure-coupled problem are

$$ \mathrm{div} \, \boldsymbol{\sigma} + \omega^2 \rho_S \mathbf{u} = \mathbf{0} \quad \text{in } \Omega_S \qquad (1a) $$

$$ \boldsymbol{\sigma} \mathbf{n}^S = \mathbf{0} \quad \text{on } \Gamma_\sigma \qquad (1b) $$

$$ \mathbf{u} = \mathbf{0} \quad \text{on } \Gamma_u \qquad (1c) $$

$$ \boldsymbol{\sigma} \mathbf{n}^S = p \mathbf{n} \quad \text{on } \Sigma \qquad (1d) $$

$$ \mathrm{div} \, \mathbf{D} = \mathbf{0} \quad \text{in } \Omega_S \qquad (2a) $$

$$ \mathbf{D} \cdot \mathbf{n}^S = \mathbf{0} \quad \text{on } \Gamma_D \qquad (2b) $$

$$ \phi = \mathbf{0} \quad \text{on } \Gamma_\phi \qquad (2c) $$
$$ \Delta p + \frac{\omega^2}{c_F^2} p = 0 \qquad \text{in } \Omega_F \qquad (3a) $$

$$ \frac{\partial p}{\partial n} = \omega^2 \rho_F \mathbf{u} \cdot \mathbf{n} \quad \text{on } \Sigma \qquad (3b) $$

where $\omega$ is the angular frequency.

Equation (1a) corresponds to the elastodynamic equation in the absence of body force; Equations (1b) and (1c) are the classical mechanical boundary conditions; Equation (1d) results from the action of pressure forces exerted by the fluid on the structure; Equation (2a) corresponds to the electric charge equation for a dielectric medium; Equations (2b) and (2c) are the previously described electric boundary conditions; Equation (3a) is the Helmholtz equation; and Equation (3b) is the contact condition for the fluid on $\Sigma$.

The stress tensor $\sigma$ and electric displacement vector $D$ are related to the linear strain tensor $\epsilon$ and electric field vector $E$ through the converse and direct linear piezoelectric constitutive equations

$$ \sigma = c\epsilon - e^T E \qquad (4) $$

$$ D = e\epsilon + dE \qquad (5) $$

where $c$ is the elastic stiffness matrix, $d$ is the dielectric permittivity matrix, and $e$ is the piezoelectric coupling matrix.

Moreover, we have the following gradient relations between the linearized strain tensor $\epsilon$ and the displacement $\mathbf{u}$, and between the electric field $E$ and the electric potential $\phi$:

$$ \epsilon = \frac{1}{2}(\nabla \mathbf{u} + \nabla^T \mathbf{u}) \qquad (6) $$

$$ E = -\nabla\phi \qquad (7) $$

For a detailed derivation of those classical equations, we refer the reader, for example, to [15] for piezoelectric aspects and to [16] for fluid-structure aspects.

## 2.2. Variational formulation in terms of ($\mathbf{u}$, $\phi$, $p$)

The local equations of Section 2.1 are expressed in terms of the chosen unknown fields of the piezoelectric structure/fluid boundary value problem, i.e. the structural mechanical displacement $\mathbf{u}$, the electric potential in the structure $\phi$, and the fluid pressure $p$.

In order to obtain the variational formulation associated with the local equations of the coupled fluid-piezoelectric structure system given in Equations (1)-(3), the test-function method is applied. We proceed in three steps, successively considering the equations relating to the structure (subject to fluid pressure actions), the electric charge equation for a dielectric medium, and the equations relating to the fluid (subject to a wall displacement). We recall that the chosen unknown fields are the structure displacement $\mathbf{u}$, the electric potential $\phi$ and the fluid pressure $p$.

First, we introduce the space $C_u^*$ of sufficiently regular functions $\mathbf{u}$ defined in $\Omega_S$ verifying $\mathbf{u}=\mathbf{0}$ on $\Gamma_u$. Multiplying Equation (1a) by any function $\delta\mathbf{u} \in C_u^*$, then applying Green's formula, and
finally taking Equations (1b) and (1d) into account, leads to

$$ \int_{\Omega_S} \sigma(\mathbf{u}, \phi) : \varepsilon(\delta\mathbf{u}) \, dv - \int_{\Sigma} p \mathbf{n} \cdot \delta\mathbf{u} \, ds - \omega^2 \rho_S \int_{\Omega_S} \mathbf{u} \cdot \delta\mathbf{u} \, dv = 0 \quad \forall \delta\mathbf{u} \in C_u^* \quad (8) $$

Secondly, we consider the space $C_\phi^*$ of sufficiently regular functions $\phi$ in $\Omega_S$ verifying $\phi=0$ on $\Gamma_\phi$. Multiplying Equation (2a) by any function $\delta\phi \in C_\phi^*$ and integrating over $\Omega_S$, we have

$$ \int_{\Omega_S} \operatorname{div} \mathbf{D}(\mathbf{u}, \phi) \delta\phi \, dv = 0 \quad \forall \delta\phi \in C_\phi^* \quad (9) $$

Finally, we consider the space $C_p$ of sufficiently regular functions $p$ defined in $\Omega_F$. Multiplying Equation (3a) by $\delta p \in C_p$, applying Green's formula, and taking Equation (3b) into account, we obtain

$$ \frac{1}{\rho_F} \int_{\Omega_F} \nabla p \cdot \nabla \delta p \, dv - \frac{\omega^2}{\rho_F c_F^2} \int_{\Omega_F} p \delta p \, dv - \omega^2 \int_{\Sigma} \mathbf{u} \cdot \mathbf{n} \delta p \, ds = 0 \quad \forall \delta p \in C_p \quad (10) $$

Using the piezoelectric constitutive Equations (4) and (5), and taking the electric boundary condition (2b) into account, Equations (8) and (9) can be rewritten in the following forms:

$$ \int_{\Omega_S} \mathbf{c}\varepsilon(\mathbf{u}) : \varepsilon(\delta\mathbf{u}) \, dv - \int_{\Omega_S} \mathbf{e}^T \mathbf{E}(\phi) : \varepsilon(\delta\mathbf{u}) \, dv - \int_{\Sigma} p \mathbf{n} \cdot \delta\mathbf{u} \, ds - \omega^2 \rho_S \int_{\Omega_S} \mathbf{u} \cdot \delta\mathbf{u} \, dv = 0 \quad (11) $$

and

$$ \int_{\Omega_S} \mathbf{e}\varepsilon(\mathbf{u}) \cdot \mathbf{E}(\delta\phi) \, dv + \int_{\Omega_S} \mathbf{d}\mathbf{E}(\phi) \cdot \mathbf{E}(\delta\phi) \, dv = 0 \quad (12) $$

Thus, the variational formulation of the coupled fluid–piezoelectric structure spectral problem consists in finding $\omega \in \mathbb{R}^+$ and $(\mathbf{u}, \phi, p) \in (C_u^*, C_\phi^*, C_p)$ satisfying Equations (10)–(12).

## 2.3. Finite element formulation

Let us introduce $\mathbf{U}$, $\mathbf{\Phi}$ and $\mathbf{P}$ corresponding to the vectors of nodal values of $\mathbf{u}$, $\phi$ and $p$, respectively, and the submatrices corresponding to the various bilinear forms involved in Equations (10)–(12) defined by

$$ \int_{\Omega_S} \mathbf{c}\varepsilon(\mathbf{u}) : \varepsilon(\delta\mathbf{u}) \, dv \Rightarrow \delta\mathbf{U}^T \mathbf{K}_u \mathbf{U} \quad (13a) $$

$$ \rho_S \int_{\Omega_S} \mathbf{u} \cdot \delta\mathbf{u} \, dv \Rightarrow \delta\mathbf{U}^T \mathbf{M}_u \mathbf{U} \quad (13b) $$

$$ \int_{\Omega_S} \mathbf{e}^T \mathbf{E}(\phi) : \varepsilon(\delta\mathbf{u}) \, dv \Rightarrow \delta\mathbf{U}^T \mathbf{C}_{u\phi} \mathbf{\Phi} \quad (13c) $$

$$ \int_{\Omega_S} \mathbf{e}\varepsilon(\mathbf{u}) \cdot \mathbf{E}(\delta\phi) \, dv \Rightarrow \delta\mathbf{\Phi}^T \mathbf{C}_{u\phi}^T \mathbf{U} \quad (13d) $$
$$ \int_{\Sigma} p \mathbf{n} \cdot \delta \mathbf{u} \, d\mathbf{s} \Rightarrow \delta \mathbf{U}^{\mathrm{T}} \mathbf{C}_{\mathrm{up}} \mathbf{P} \qquad (13e) $$

$$ \frac{1}{\rho_F} \int_{\Omega_F} \nabla p \cdot \nabla \delta p \, dv \Rightarrow \delta \mathbf{P}^{\mathrm{T}} \mathbf{K}_p \mathbf{P} \qquad (13f) $$

$$ \frac{1}{\rho_F c_F^2} \int_{\Omega_F} p \delta p \, dv \Rightarrow \delta \mathbf{P}^{\mathrm{T}} \mathbf{M}_p \mathbf{P} \qquad (13g) $$

$$ \int_{\Omega_S} \mathbf{d}\mathbf{E}(\phi) \cdot \mathbf{E}(\delta\phi) \, d\mathbf{v} \Rightarrow \delta\mathbf{\Phi}^{\mathrm{T}}\mathbf{K}_{\phi}\mathbf{\Phi} \qquad (13h) $$

The variational equations (10)–(12) for the fluid–piezoelectric structure-coupled problem can be written in discretized form

$$ \begin{bmatrix} \mathbf{K}_u & -\mathbf{C}_{u\phi} & -\mathbf{C}_{up} \\ \mathbf{C}_{u\phi}^{\mathrm{T}} & \mathbf{K}_{\phi} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{K}_p \end{bmatrix} \begin{bmatrix} \mathbf{U} \\ \mathbf{\Phi} \\ \mathbf{P} \end{bmatrix} = \omega^2 \begin{bmatrix} \mathbf{M}_u & \mathbf{0} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{0} \\ \mathbf{C}_{up}^{\mathrm{T}} & \mathbf{0} & \mathbf{M}_p \end{bmatrix} \begin{bmatrix} \mathbf{U} \\ \mathbf{\Phi} \\ \mathbf{P} \end{bmatrix} \qquad (14) $$

### Remarks

* The matrices involved in Equation (14) are unsymmetric. A symmetric formulation with mass coupling can be established by introducing the fluid displacement potential as additional variable [16]. For the sake of brevity, this procedure is not detailed here and will be the subject of a future work.

* Using the second line of Equation (14), the degrees-of-freedom associated with the electric potential can be expressed in terms of structure displacements as

$$ \mathbf{\Phi} = - \mathbf{K}_{\phi}^{-1} \mathbf{C}_{u\phi}^{\mathrm{T}} \mathbf{U} \qquad (15) $$

Thus, after substitution of $\mathbf{\Phi}$ into Equation (14), we get the following spectral problem in terms of $\mathbf{U}$ and $\mathbf{P}$:

$$ \begin{bmatrix} \mathbf{K}_u + \mathbf{C}_{u\phi} \mathbf{K}_{\phi}^{-1} \mathbf{C}_{u\phi}^{\mathrm{T}} & -\mathbf{C}_{up} \\ \mathbf{0} & \mathbf{K}_p \end{bmatrix} \begin{bmatrix} \mathbf{U} \\ \mathbf{P} \end{bmatrix} = \omega^2 \begin{bmatrix} \mathbf{M}_u & \mathbf{0} \\ \mathbf{C}_{up}^{\mathrm{T}} & \mathbf{M}_p \end{bmatrix} \begin{bmatrix} \mathbf{U} \\ \mathbf{P} \end{bmatrix} \qquad (16) $$

It should be noted that the term $\mathbf{C}_{u\phi} \mathbf{K}_{\phi}^{-1} \mathbf{C}_{u\phi}^{\mathrm{T}}$ represents an ‘added-stiffness’ contribution due to the electromechanical coupling [17].

* A standard unsymmetric system, associated to the fluid–elastic structure coupled problem, can be obtained from Equation (16) by choosing a zero elastic-electric coupling matrix $\mathbf{C}_{u\phi}$ (i.e. when piezoelectric coupling constants are zero)

$$ \begin{bmatrix} \mathbf{K}_u & -\mathbf{C}_{up} \\ \mathbf{0} & \mathbf{K}_p \end{bmatrix} \begin{bmatrix} \mathbf{U} \\ \mathbf{P} \end{bmatrix} = \omega^2 \begin{bmatrix} \mathbf{M} & \mathbf{0} \\ \mathbf{C}_{up}^{\mathrm{T}} & \mathbf{M}_p \end{bmatrix} \begin{bmatrix} \mathbf{U} \\ \mathbf{P} \end{bmatrix} \qquad (17) $$
### 3. APPLICATION TO AXISYMMETRIC COMPOSITE CONICAL SHELLS FILLED WITH COMPRESSIBLE FLUID

In this section, the proposed finite element formulation is applied to the vibration analysis of axisymmetric conical elastic/piezoelectric shells coupled with an internal compressible fluid.

#### 3.1. Elastic laminated conical shell

Consider a hollow truncated conical shell of length $L$, slope angle $α$, having an arbitrary constant thickness $h$ and made of $N$ orthotropic layers (see Figures 2 and 3). The particular geometry of the structure is used to reduce the dimensions of the problem through a semi-analytical procedure whereby variables are expanded as Fourier series in the circumferential direction. The orthogonality property of harmonic functions decouples the contributions of the components of the series and allows a separate analysis for each circumferential harmonic.

<figure><img src="image_3.png" /><figcaption>Figure 2. Conical shell: geometry and notations.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3. Cross-sectional view of the laminated composite conical shell.</figcaption></figure>
### 3.1.1. Strain-displacement relations.
Since a classical shell element based on Kirchhoff-Love hypotheses is used in this work, the components of the displacement vector can be written as

$$u_x = u(x, \theta) + z\beta(x, \theta) \qquad (18a)$$

$$u_\theta = v(x, \theta) + z\gamma(x, \theta) \qquad (18b)$$

$$u_z = w(x, \theta) \qquad (18c)$$

$u$, $v$ and $w$ being the axial, circumferential and radial displacement components at the middle surface of the shell. $\beta$ and $\gamma$ are the rotations about the $\theta$ and $x$ axes, respectively, given by

$$\beta = -\frac{\partial w}{\partial x} \qquad (19)$$

$$\gamma = \frac{\cos \alpha}{R} v - \frac{1}{R} \frac{\partial w}{\partial \theta} \qquad (20)$$

where $R$ denotes the local radius of the shell written by

$$R = x \sin \alpha \qquad (21)$$

With these notations, the state of strain at any point $x$, $\theta$ (on the shell) and $z$ (distance from the neutral surface), is equal to the sum of the middle surface strains (membrane) and the strains due to the changes of curvature (flexural). Denoting the middle surface extensional strains in the lengthwise and circumferential directions by $\epsilon_{xx}$ and $\epsilon_{\theta\theta}$, the middle surface changes in curvature in the lengthwise and circumferential directions by $\kappa_{xx}$ and $\kappa_{\theta\theta}$, and the middle surface in-plane shear strain and change in twist by $\epsilon_{x0}$ and $\kappa_{x0}$, the state of strain at any point $(x, \theta, z)$ within the shell can be written as

$$\epsilon_{xx} = \epsilon_{xx} + z\kappa_{xx} \qquad (22a)$$

$$\epsilon_{\theta\theta} = \epsilon_{\theta\theta} + z\kappa_{\theta\theta} \qquad (22b)$$

$$2\epsilon_{x0} = \epsilon_{x0} + z\kappa_{x0} \qquad (22c)$$

where the factor 2 in the left-hand side of Equation (22c) comes from the definition of the strain tensor (Equation (6)) and the middle surface extensional strains are defined by

$$\epsilon_{xx} = \frac{\partial u}{\partial x} \qquad (23a)$$

$$\epsilon_{\theta\theta} = \frac{\sin \alpha}{R} u + \frac{1}{R} \frac{\partial v}{\partial \theta} + \frac{\cos \alpha}{R} w \qquad (23b)$$

$$\epsilon_{x0} = \frac{1}{R} \frac{\partial u}{\partial \theta} - \frac{\sin \alpha}{R} v + \frac{\partial v}{\partial x} \qquad (23c)$$
and the middle surface curvatures are expressed by

$$ \kappa_{xx} = -\frac{\partial^2 w}{\partial x^2} \qquad (24a) $$

$$ \kappa_{00} = \frac{\cos \alpha \partial v}{R^2 \partial \theta} - \frac{1}{R^2} \frac{\partial^2 w}{\partial \theta^2} - \frac{\sin \alpha \partial w}{R \partial x} \qquad (24b) $$

$$ \kappa_{x\theta} = -\frac{2 \cos \alpha \sin \alpha}{R^2} v + \frac{\cos \alpha \partial v}{R \partial x} + \frac{2 \sin \alpha \partial w}{R^2 \partial \theta} - \frac{2}{R} \frac{\partial^2 w}{\partial x \partial \theta} \qquad (24c) $$

Equations (23) and (24) correspond to strain-displacement relations for thin conical shells in Love theory [18].

### 3.1.2. Fourier series expansion.
The geometric and material rotational symmetry of the structure allows a semi-analytical formulation for the shell finite element. Thus, the displacement field is expended in Fourier series in the circumferential direction according to the following expressions:

$$ u = u_0(x) + \sum_{n=1}^{\infty} u_n(x) \cos n\theta \qquad (25a) $$

$$ v = v_0(x) + \sum_{n=1}^{\infty} v_n(x) \sin n\theta \qquad (25b) $$

$$ w = w_0(x) + \sum_{n=1}^{\infty} w_n(x) \cos n\theta \qquad (25c) $$

where *n* identifies the harmonic component.

In fact, due to the axisymmetry, the modes for $n \ge 1$ appear by pair, the two members of a pair being denoted as companion modes. Those companion modes have exactly the same frequency and their shapes have the same geometry, the only difference being in their angular position. Moreover, displacement components appear grouped in $(u, \partial v/\partial \theta, w)$ or in $(\partial u/\partial \theta, v, \partial w/\partial \theta)$ in Equations (23) and (24). As a consequence, the displacement components of one companion mode are of the form $(\cos n\theta, \sin n\theta, \cos n\theta)$ or $(\sin n\theta, \cos n\theta, \sin n\theta)$. To obtain both companion modes, one would have to introduce both sine and cosine functions in each displacement component $(u, v, w)$. However, for sake of brevity, only one of the companion modes is considered in Equations (25) and consequently in our simulations.

Using Equations (23)–(25), middle surface extensional strains and curvatures can be written as

$$ \epsilon_{xx} = \epsilon_{xx}^{0}(x) + \sum_{n=1}^{\infty} \epsilon_{xx}^{n}(x) \cos n\theta, \quad \kappa_{xx} = \kappa_{xx}^{0}(x) + \sum_{n=1}^{\infty} \kappa_{xx}^{n}(x) \cos n\theta $$

$$ \epsilon_{00} = \epsilon_{00}^{0}(x) + \sum_{n=1}^{\infty} \epsilon_{00}^{n}(x) \cos n\theta, \quad \kappa_{00} = \kappa_{00}^{0}(x) + \sum_{n=1}^{\infty} \kappa_{00}^{n}(x) \cos n\theta $$

$$ \epsilon_{x\theta} = \epsilon_{x\theta}^{0}(x) + \sum_{n=1}^{\infty} \epsilon_{x\theta}^{n}(x) \sin n\theta, \quad \kappa_{x\theta} = \kappa_{x\theta}^{0}(x) + \sum_{n=1}^{\infty} \kappa_{x\theta}^{n}(x) \sin n\theta $$
where the strain components can be rewritten, for any circumferential harmonic $n$, in the following matrix form:

$$ \mathbf{e}_n = \mathbf{D}_m \begin{bmatrix} u_n \\ v_n \\ w_n \end{bmatrix}, \quad \mathbf{k}_n = \mathbf{D}_b \begin{bmatrix} u_n \\ v_n \\ w_n \end{bmatrix} \qquad (26) $$

with $\mathbf{e}_n = [\epsilon_{xx}^n \ \epsilon_{00}^n \ \epsilon_{x0}^n]^T$, $\mathbf{k}_n = [\kappa_{xx}^n \ \kappa_{00}^n \ \kappa_{x0}^n]^T$ and the gradient operators

$$ \mathbf{D}_m = \begin{bmatrix} \frac{\partial}{\partial x} & 0 & 0 \\ \frac{S}{R} & \frac{n}{R} & \frac{C}{R} \\ -\frac{n}{R} & -\frac{S}{R} + \frac{\partial}{\partial x} & 0 \end{bmatrix}, \quad \mathbf{D}_b = \begin{bmatrix} 0 & 0 & -\frac{\partial^2}{\partial x^2} \\ 0 & \frac{nC}{R^2} & \frac{n^2}{R^2} - \frac{S}{R} \frac{\partial}{\partial x} \\ 0 & -2\frac{CS}{R^2} + \frac{C}{R} \frac{\partial}{\partial x} & -2\frac{nS}{R^2} + 2\frac{n}{R} \frac{\partial}{\partial x} \end{bmatrix} $$

where $C = \cos \alpha$, $S = \sin \alpha$ and the subscripts m and b state for membrane and bending.

### 3.1.3. Degrees-of-freedom and shape functions.
The displacement components associated with each harmonic $n$ are predicted by discretizing the shell with one-dimensional two nodes axisymmetric elements. Each node has four degrees-of-freedom to describe axial, circumferential and radial displacements, as well as the rotation about the circumferential axis. One element can thus be described in terms of the following nodal degrees-of-freedom vector:

$$ \mathbf{U}_n^e = [u_{n1} \ v_{n1} \ w_{n1} \ \beta_{n1} \ | \ u_{n2} \ v_{n2} \ w_{n2} \ \beta_{n2}]^T \qquad (27) $$

where superscript $e$ denotes an elementary quantity, and subscripts 1 and 2 correspond to nodes bounding the element.

The longitudinal and circumferential displacements are assumed to vary linearly along the axial co-ordinate $x$, while the radial displacement is described by Hermite cubic polynomial in order to maintain C1 continuity due to the application of Kirchhoff-Love hypothesis. Thus, the elementary middle surface displacements and rotations of the shell are given in terms of the nodal degrees-of-freedom by

$$ \begin{bmatrix} u_n^e \\ v_n^e \\ w_n^e \end{bmatrix} = \mathbf{N}_d(x) \mathbf{U}_n^e, \quad \begin{bmatrix} \beta_n^e \\ \gamma_n^e \\ 0 \end{bmatrix} = \begin{bmatrix} -\frac{\partial w_n^e}{\partial x} \\ \frac{C}{R} v_n^e + \frac{n}{R} w_n^e \\ 0 \end{bmatrix} = \mathbf{N}_r(x) \mathbf{U}_n^e \qquad (28) $$

where the interpolation matrices are defined by

$$ \mathbf{N}_d = \begin{bmatrix} N_1^1 & 0 & 0 & 0 & N_2^1 & 0 & 0 & 0 \\ 0 & N_1^1 & 0 & 0 & 0 & N_2^1 & 0 & 0 \\ 0 & 0 & N_1^c & N_2^c & 0 & 0 & N_3^c & N_4^c \end{bmatrix} $$
and

$$ \mathbf{N}_r = \begin{bmatrix} 0 & 0 & -\frac{\partial N_1^c}{\partial x} & -\frac{\partial N_2^c}{\partial x} \\ 0 & \frac{C}{R} N_1^1 & \frac{n}{R} N_1^c & \frac{n}{R} N_2^c \\ 0 & 0 & 0 & 0 \end{bmatrix} \begin{bmatrix} 0 & 0 & -\frac{\partial N_3^c}{\partial x} & -\frac{\partial N_4^c}{\partial x} \\ 0 & \frac{C}{R} N_2^1 & \frac{n}{R} N_3^c & \frac{n}{R} N_4^c \\ 0 & 0 & 0 & 0 \end{bmatrix} $$

with the following linear $N_i^1$ ($i=1,2$) and cubic $N_i^c$ ($i=1...4$) shape functions given in the reference element

$$ N_1^1 = \frac{1}{2}(1 - \xi), \quad N_2^1 = \frac{1}{2}(1 + \xi) $$

$$ N_1^c = \frac{1}{4}(2 - 3\xi + \xi^3), \quad N_2^c = -\frac{L_e}{8}(1 - \xi - \xi^2 + \xi^3) $$

$$ N_3^c = \frac{1}{4}(2 + 3\xi - \xi^3), \quad N_4^c = -\frac{L_e}{8}(-1 - \xi + \xi^2 + \xi^3) $$

in which $\xi$ is the elementary co-ordinate that varies from $\xi = -1$ at node 1 to $\xi = 1$ at node 2 and $L_e$ is the length of the element.

Moreover, elementary extensional strains and curvatures can be expressed by

$$ \mathbf{e}_n^e = \mathbf{B}_m(x) \mathbf{U}_n^e, \quad \mathbf{k}_n^e = \mathbf{B}_b(x) \mathbf{U}_n^e \qquad (29) $$

with the following discretized gradient operators:

$$ \mathbf{B}_m = \mathbf{D}_m \mathbf{N}_d, \quad \mathbf{B}_b = \mathbf{D}_b \mathbf{N}_d $$

Note that the third line of zeros in the $\mathbf{N}_r$ expression is used in order to define the products $\mathbf{N}_d^T \mathbf{N}_r$ and $\mathbf{N}_r^T \mathbf{N}_d$ in the mass matrix defined below.

### 3.1.4. Elementary mass and stiffness matrices.
The interpolations of displacements and strains presented in the previous section are used to express the elementary mass and stiffness matrices of the composite laminated shell.

The mass matrix is evaluated without neglecting rotational inertia. Using relation (13b) and combining Equations (18), (25) and (28), the elementary mass matrix is defined, for any harmonic $n$, by

$$ \mathbf{M}_u^e = \pi(1 + \delta_{0n}) \int_{-1}^{1} [I_0 \mathbf{N}_d^T \mathbf{N}_d + I_1 (\mathbf{N}_d^T \mathbf{N}_r + \mathbf{N}_r^T \mathbf{N}_d) + I_2 \mathbf{N}_r^T \mathbf{N}_r] R(\xi) J_S d\xi $$

where $J_S$ is the Jacobian determinant which is simply given by $J_S = L_e/2$, $\delta$ is the Kronecker symbol and the zero, first and second mass moments of inertia are given by

$$ I_0 = \sum_{k=1}^{N} \rho_S^k (z_k - z_{k-1}), \quad I_1 = \frac{1}{2} \sum_{k=1}^{N} \rho_S^k (z_k^2 - z_{k-1}^2), \quad I_2 = \frac{1}{3} \sum_{k=1}^{N} \rho_S^k (z_k^3 - z_{k-1}^3) $$
In the same way, using definition (13a) and combining Equations (22), (26) and (29), the elementary elastic stiffness matrix is given, for any harmonic n, by

$$ \mathbf{K}_u^e = \pi(1 + \delta_{0n}) \int_{-1}^{1} (\mathbf{B}_m^T \mathbf{A} \mathbf{B}_m + \mathbf{B}_m^T \mathbf{B} \mathbf{B}_b + \mathbf{B}_b^T \mathbf{B} \mathbf{B}_m + \mathbf{B}_b^T \mathbf{D} \mathbf{B}_b) R(\xi) J_S d\xi $$

The matrices A, D and B are extensional, bending and extensional-bending coupling stiffness of the laminated composite (see for example [19]). The components of these 3 × 3 symmetric matrices are given by

$$ A_{ij} = \sum_{k=1}^{N} \overline{Q}_{ij}^{k} (z_k - z_{k-1}), \quad B_{ij} = \frac{1}{2} \sum_{k=1}^{N} \overline{Q}_{ij}^{k} (z_k^2 - z_{k-1}^2), \quad D_{ij} = \frac{1}{3} \sum_{k=1}^{N} \overline{Q}_{ij}^{k} (z_k^3 - z_{k-1}^3) $$

where i, j = 1, 2, 6 and $\overline{Q}_{ij}^k$ represent the reduced material stiffness constants for each layer k in the global co-ordinate system deduced from the assumption of zero normal stress in the thickness direction.

## 3.2. Finite element modelling of the fluid domain

The aim of this section is to introduce the finite element discretization of the fluid domain considering the axisymmetry of the problem (see Figure 4). As for the structure displacements, fluid pressure can be expressed in terms of Fourier series by

$$ p(r, \theta, Z) = p_0(r, Z) + \sum_{n=1}^{\infty} p_n(r, Z) \cos n\theta \quad (30) $$

where $p_0$ and $p_n$ are, respectively, the axisymmetric and the nth harmonic component of the pressure.

The considered fluid domain is discretized into quadrilateral elements. Each element is bounded by four nodes and the nodal pressure is considered as the only unknown variable.

For any harmonic n, the elementary pressure $p_n^e(r, Z)$ can be expressed in terms of the vector of nodal pressures $\mathbf{P}_n^e$

$$ p_n^e(r, Z) = \mathbf{N}_p(r, Z) \mathbf{P}_n^e \quad (31) $$

where the interpolation matrix is defined by

$$ \mathbf{N}_p = [N_1 \ N_2 \ N_3 \ N_4] $$

<figure><img src="image_5.png" /><figcaption>Figure 4. Fluid domain: geometry and notations.</figcaption></figure>
with the following linear shape functions, given in the reference element

$$N_1 = \frac{1}{4}(1 - \eta)(1 - \zeta), \quad N_2 = \frac{1}{4}(1 + \eta)(1 - \zeta), \quad N_3 = \frac{1}{4}(1 + \eta)(1 + \zeta), \quad N_4 = \frac{1}{4}(1 - \eta)(1 + \zeta)$$

From definition (13g), and using Equations (30) and (31), the elementary mass matrix of the fluid can be expressed, for each harmonic $n$, by

$$\mathbf{M}_p^e = \frac{\pi(1 + \delta_{0n})}{\rho_F c_F^2} \int_{-1}^{1} \int_{-1}^{1} \mathbf{N}_p^T \mathbf{N}_p r(\eta, \zeta) J_F d\eta \, d\zeta \quad (32)$$

where $J_F$ is the Jacobian determinant of the transformation from the real to the reference element.

In order to compute the elementary stiffness matrix of the fluid, we have to evaluate the pressure gradient, for each harmonic $n$, using

$$\nabla p_n^e(r, Z) = \mathbf{B}_p(r, Z) \mathbf{P}_n^e \quad (33)$$

with the following discretized gradient operators:

$$\mathbf{B}_p(r, Z) = \left[ \frac{\partial}{\partial r} \frac{n}{r} \frac{\partial}{\partial Z} \right]^T \mathbf{N}_p(r, Z) \quad (34)$$

From Equations (13f), (30) and (33), and after transforming the discretized gradient operator in the reference co-ordinate system, the elementary stiffness matrix of the fluid can be expressed, for each harmonic $n$, by

$$\mathbf{K}_p^e = \frac{\pi(1 + \delta_{0n})}{\rho_F} \int_{-1}^{1} \int_{-1}^{1} \mathbf{B}_p^T \mathbf{B}_p r(\eta, \zeta) J_F d\eta \, d\zeta \quad (35)$$

### 3.3. Elementary fluid-structure coupling matrix

Using the appropriate interface conditions at the common boundary between the fluid and the structure, the coupled fluid-structure matrix is derived. In first approach, a compatible mesh is considered at the interface so that the axisymmetric interface element can be easily defined. The interface element is bounded by two nodes, each having two degrees-of-freedom to describe, for any harmonic $n$, the normal structure displacement $w_n^e$ and the fluid pressure $p_n^e$.

The elementary normal displacement $w_n^e$ is discretized by cubic shape functions and can be written, in terms of the nodal structure degrees-of-freedom, by

$$w_n^e = [0 \ 0 \ N_1^c \ N_2^c \ 0 \ 0 \ N_3^c \ N_4^c] \mathbf{U}_n^e = \mathbf{N}_w \mathbf{U}_n^e \quad (36)$$

The elementary pressure $p_n^e$ is discretized by linear shape functions and can be written, in terms of fluid pressure in nodes located at the interface, by

$$p_n^e = [N_1^c \ N_2^c]^T \begin{bmatrix} p_{n1} \\ p_{n2} \end{bmatrix} = \mathbf{N}_{pi} \begin{bmatrix} p_{n1} \\ p_{n2} \end{bmatrix} \quad (37)$$

Thus, the elementary 8 × 2 coupling matrix can be written as

$$\mathbf{C}_{up}^e = \pi(1 + \delta_{0n}) \int_{-1}^{1} \mathbf{N}_w^T \mathbf{N}_{pi} R(\xi) J_S d\xi \quad (38)$$
### 3.4. Extension to piezoelectric conical shell

In this section, we present an extension of the laminated composite conical shell to the case of piezoelectric layers with radial polarization. The chosen multilayer model combines an equivalent single-layer assumption for the mechanical displacement and a layerwise representation of the transverse electric potential (the electric potential is assumed to vary linearly in the thickness of each piezoelectric layer). The advantages of this mixed laminate theory are linked to its effectiveness to model thin composite shells and to capture the through-thickness electric heterogeneity induced by the piezoelectric layers [20, 21].

From the previous electric potential assumption and neglecting the in-plane components, the electric field is defined, for each piezoelectric layer $k$, by its transverse component

$$E^k = -V^k/h_k \quad (39)$$

where $V^k$ is the electrical potential difference between top and bottom surfaces ($V^k = \phi^k - \phi^{k-1}$) and $h_k$ is the layer thickness ($h_k = z_k - z_{k-1}$).

Instead of describing the electric state by the electric potentials $\phi^k$, we use in the finite element discretization the electric potential differences $V^k$. Moreover, as for the displacement and pressure fields, the potential difference is developed in Fourier series by

$$V^k(x, \theta) = V_0^k(x) + \sum_{n=1}^{\infty} V_n^k(x) \cos n\theta \quad (40)$$

Thus, the degrees-of-freedom of the axisymmetric shell element described in the previous section are augmented by the electric potential difference (voltage) of each layer. For a circumferential wave number $n$, the electric field of one multilayer piezoelectric shell element can be written in the following form:

$$\mathbf{E}_n^e = - \mathbf{B}_\phi \mathbf{\Phi}_n^e \quad (41)$$

with $\mathbf{E}_n^e = [E_n^1 \dots E_n^k \dots E_n^N]^T$, $\mathbf{\Phi}_n^e = [V_n^1 \dots V_n^k \dots V_n^N]^T$, and

$$\mathbf{B}_{\phi} = \begin{bmatrix} \frac{1}{h_1} & 0 & \cdots & 0 \\ 0 & \ddots & 0 & \vdots \\ \vdots & 0 & \frac{1}{h_k} & 0 \\ \vdots & \vdots & 0 & \ddots & 0 \\ 0 & \cdots & \cdots & 0 & \frac{1}{h_N} \end{bmatrix} \quad (42)$$

where $N$ is the total number of layers.

With these considerations, the elementary electromechanical coupling stiffness matrix, given by Equation (13c), can be written, for each circumferential wave number $n$, by

$$\mathbf{C}_{u\phi}^{e} = \pi(1 + \delta_{0n}) \int_{-1}^{1} (\mathbf{B}_{\mathrm{m}}^{\mathrm{T}} \mathbf{F} \mathbf{B}_{\phi} + \mathbf{B}_{\mathrm{b}}^{\mathrm{T}} \mathbf{G} \mathbf{B}_{\phi}) R(\xi) J_{\mathrm{s}} \, d\xi \quad (43)$$
F and G being the membrane-electric and bending-electric coupling matrices, defined by

$$ \mathbf{F} = \begin{bmatrix} h_1 \bar{e}_{31}^1 & \cdots & h_k \bar{e}_{31}^k & \cdots & h_N \bar{e}_{31}^N \\ h_1 \bar{e}_{32}^1 & \cdots & h_k \bar{e}_{32}^k & \cdots & h_N \bar{e}_{32}^N \\ \vdots & \ddots & \vdots & \ddots & \vdots \\ h_1 \bar{e}_{36}^1 & \cdots & h_k \bar{e}_{36}^k & \cdots & h_N \bar{e}_{36}^N \end{bmatrix} \quad (44) $$

and

$$ \mathbf{G} = \begin{bmatrix} \frac{z_1^2 - z_0^2}{2} \bar{e}_{31}^1 & \cdots & \frac{z_k^2 - z_{k-1}^2}{2} \bar{e}_{31}^k & \cdots & \frac{z_N^2 - z_{N-1}^2}{2} \bar{e}_{31}^N \\ \frac{z_1^2 - z_0^2}{2} \bar{e}_{32}^1 & \cdots & \frac{z_k^2 - z_{k-1}^2}{2} \bar{e}_{32}^k & \cdots & \frac{z_N^2 - z_{N-1}^2}{2} \bar{e}_{32}^N \\ \frac{z_1^2 - z_0^2}{2} \bar{e}_{36}^1 & \cdots & \frac{z_k^2 - z_{k-1}^2}{2} \bar{e}_{36}^k & \cdots & \frac{z_N^2 - z_{N-1}^2}{2} \bar{e}_{36}^N \end{bmatrix} \quad (45) $$

where $\bar{e}_{3j}^k$ ($j=1, 2, 6$) are the reduced piezoelectric constants for each layer $k$ in the global co-ordinate system deduced from the assumption of zero normal stress in the thickness direction.

The previously defined matrices are derived from the global constitutive equations of the laminate, which connect, for any harmonic $n$, the resultant in-plane forces $N_n$ and bending moments $M_n$, to the middle surface extensional strains $e_n$ and curvatures $\kappa_n$, and electric field $E_n$

$$ \begin{bmatrix} \mathbf{N}_n \\ \mathbf{M}_n \end{bmatrix} = \begin{bmatrix} \mathbf{A} & \mathbf{B} \\ \mathbf{B} & \mathbf{D} \end{bmatrix} \begin{bmatrix} \mathbf{e}_n \\ \mathbf{\kappa}_n \end{bmatrix} - \begin{bmatrix} \mathbf{F} \\ \mathbf{G} \end{bmatrix} \mathbf{E}_n \quad (46) $$

Finally, from Equation (13h), the electric matrix $\mathbf{K}_\phi^e$ is given, for each harmonic $n$, by

$$ \mathbf{K}_\phi^e = \pi(1 + \delta_{0n}) \int_{-1}^{1} \mathbf{B}_\phi^T \mathbf{H} \mathbf{B}_\phi R(\xi) J_S d\xi \quad (47) $$

where $\mathbf{H}$ is defined by

$$ \mathbf{H} = \begin{bmatrix} h_1 \bar{d}_{33}^1 & 0 & \cdots & 0 \\ 0 & \ddots & 0 & \vdots \\ \vdots & 0 & h_k \bar{d}_{33}^k & 0 \\ \vdots & \vdots & 0 & \ddots \\ 0 & \cdots & 0 & h_N \bar{d}_{33}^N \end{bmatrix} \quad (48) $$

with $\bar{d}_{33}^k$ the dielectric permittivity in the thickness direction.

The matrix $\mathbf{H}$ is derived by integrating the transverse electric displacement $D^k$ over the thickness of each layer such that

$$ \mathbf{Q}_n = - \mathbf{F}^T \mathbf{e}_n - \mathbf{G}^T \mathbf{\kappa}_n - \mathbf{H} \mathbf{E}_n \quad (49) $$
where the electric displacement resultant vector $Q_n$ is defined by

$$Q_n = \left[ \int_{z_0}^{z_1} D^1 dz \cdots \int_{z_{k-1}}^{z_k} D^k dz \cdots \int_{z_{N-1}}^{z_N} D^N dz \right]^T \quad (50)$$

and where $D^k$ is given for each layer in the global co-ordinate system by

$$D^k = \bar{e}_{31}^k \varepsilon_{xx} + \bar{e}_{32}^k \varepsilon_{00} + 2 \bar{e}_{36}^k \varepsilon_{x0} + \bar{d}_{33}^k E^k \quad (51)$$

## 4. NUMERICAL RESULTS

In this last section, numerical results, obtained with a Matlab program developed by the authors, are proposed in order to validate and analyse the previously described formulation.

Section 4.1 presents the performances of the composite shell finite element for elastic or piezoelectric structures without internal fluid. Section 4.2 concerns elastic and piezoelectric cylindrical shells coupled with an internal compressible fluid.

### 4.1. Free vibration of shells without internal fluid

#### 4.1.1. Isotropic conical shell.
Consider an isotropic single-layer conical shell of length $L = 7.5$ m, slope angle $\alpha = \pi/6$, thickness $h = 0.05$ m, and radius at its large edge $R_2 = 5$ m. The isotropic material used in this first example is aluminium, with Young's modulus $E = 69000$ MPa, Poisson ratio $\nu = 0.3$, and mass density $\rho = 2700$ kg m$^{-3}$. Moreover, the shell is discretized with 50 axisymmetric elements.

Table I presents a comparison of the non-dimensional frequency parameter $\lambda = \omega R_2 \sqrt{\rho(1-\nu^2)/E}$ for various set of boundary conditions at both the base and top of the cone: free-free (F-F), simply supported-simply supported (SS-SS) $\nu = w = 0$, and clamped-clamped (C-C) $u = v = w = \beta = 0$. Moreover, the results are compared with those given by Shu [22] using a generalized differential quadrature technique, and those obtained by the finite element code Nastran using quadrangular shell elements. In this example, the circumferential wave number $n$ ranges from 1 to 9. An excellent agreement can be observed for all modes. From this comparison, it can be concluded

<table><caption>Table I. Frequency parameter $\lambda = \omega R_2 \sqrt{\rho(1-\nu^2)/E}$ of an isotropic conical shell for $n = 1-9$.</caption><thead><tr><th rowspan="2">n</th><th colspan="2">F-F</th><th colspan="3">SS-SS</th><th colspan="3">C-C</th></tr><tr><th>Nastran</th><th>Present</th><th>Nastran</th><th>Present</th><th>Shu [22]</th><th>Nastran</th><th>Present</th><th>Shu [22]</th></tr></thead><tbody><tr><td>1</td><td>—</td><td>—</td><td>0.679</td><td>0.679</td><td>0.695</td><td>0.752</td><td>0.753</td><td>0.752</td></tr><tr><td>2</td><td>0.014</td><td>0.015</td><td>0.358</td><td>0.359</td><td>0.359</td><td>0.466</td><td>0.466</td><td>0.466</td></tr><tr><td>3</td><td>0.037</td><td>0.037</td><td>0.210</td><td>0.211</td><td>0.212</td><td>0.315</td><td>0.316</td><td>0.315</td></tr><tr><td>4</td><td>0.065</td><td>0.066</td><td>0.190</td><td>0.190</td><td>0.181</td><td>0.260</td><td>0.260</td><td>0.261</td></tr><tr><td>5</td><td>0.096</td><td>0.096</td><td>0.216</td><td>0.216</td><td>0.216</td><td>0.264</td><td>0.265</td><td>0.265</td></tr><tr><td>6</td><td>0.132</td><td>0.132</td><td>0.252</td><td>0.253</td><td>0.254</td><td>0.293</td><td>0.294</td><td>0.293</td></tr><tr><td>7</td><td>0.173</td><td>0.174</td><td>0.298</td><td>0.298</td><td>0.298</td><td>0.334</td><td>0.335</td><td>0.335</td></tr><tr><td>8</td><td>0.221</td><td>0.222</td><td>0.351</td><td>0.352</td><td>0.353</td><td>0.385</td><td>0.387</td><td>0.387</td></tr><tr><td>9</td><td>0.275</td><td>0.276</td><td>0.413</td><td>0.415</td><td>0.414</td><td>0.446</td><td>0.448</td><td>0.448</td></tr></tbody></table>
<table><caption>Table II. Mechanical and geometrical data of an elastic laminated cylinder [23].</caption><thead><tr><th>Material properties</th><th>Geometric data</th></tr></thead><tbody><tr><td>E<sub>1</sub> = 19 GPa, E<sub>1</sub>/E<sub>2</sub> = 2.5</td><td>Inner layer thickness = h/3</td></tr><tr><td>G<sub>12</sub> = 4.1 GPa</td><td>Middle layer thickness = h/3</td></tr><tr><td>v<sub>12</sub> = 0.26</td><td>Outer layer thickness = h/3</td></tr><tr><td>ρ = 1643 kg m<sup>-3</sup></td><td>h/R = 0.002</td></tr></tbody></table>

<table><caption>Table III. Frequency parameter $\lambda = \omega\sqrt{\rho R^2/E_2}$ of a composite cylinder.</caption><thead><tr><th rowspan="2">n</th><th colspan="3">L/R = 1</th><th colspan="3">L/R = 10</th></tr><tr><th>Present</th><th>Lam and Loy [23]</th><th>Zhang [24]</th><th>Present</th><th>Lam and Loy [23]</th><th>Zhang [24]</th></tr></thead><tbody><tr><td>1</td><td>1.061284</td><td>1.061284</td><td>1.061283</td><td>0.083910</td><td>0.083908</td><td>0.083908</td></tr><tr><td>2</td><td>0.804055</td><td>0.804054</td><td>0.804052</td><td>0.030011</td><td>0.030009</td><td>0.030008</td></tr><tr><td>3</td><td>0.598332</td><td>0.598331</td><td>0.598328</td><td>0.015194</td><td>0.015193</td><td>0.015191</td></tr><tr><td>4</td><td>0.450145</td><td>0.450144</td><td>0.450140</td><td>0.012177</td><td>0.012176</td><td>0.012174</td></tr><tr><td>5</td><td>0.345254</td><td>0.345253</td><td>0.345248</td><td>0.015232</td><td>0.015231</td><td>0.015230</td></tr><tr><td>6</td><td>0.270755</td><td>0.270754</td><td>0.270747</td><td>0.021179</td><td>0.021179</td><td>0.021178</td></tr></tbody></table>

that the present axisymmetric two-node shell element is numerically accurate and computationally efficient for thin isotropic conical shells.

### 4.1.2. Laminated composite cylinder.
This example gives a comparison of the frequency parameter $\lambda = \omega R\sqrt{\rho E_2}$ for a three-layer cross-ply [0, $\pi/2$, 0] cylindrical shell with simply supported (SS-SS) boundary conditions. The mechanical and geometrical shell parameters are given in Table II. Moreover, 100 finite elements are used for the discretization of the cylinder.

The results presented in Table III show the excellent performance of our finite element model compared to the analytical solution proposed by Lam and Loy [23] and the wave propagation approach used by Zhang [24]. It can be observed that for shells with small length-to-radius ratios, as in the case $L/R = 1$, the minimum frequency occurs for the larger $n$. Nevertheless, for $L/R = 10$, the natural frequencies decrease when $n$ varies from 1 to 4, and increase for $n$ larger than 4.

### 4.1.3. Laminated composite conical shell.
In order to validate our finite element model for composite non-cylindrical structures, the natural frequencies of a three-layer [0, $\pi/2$, 0] conical shell are compared with those given by Nastran with quadrangular composite shell elements. The geometrical data, boundary conditions and spatial discretization are the same as those given for the isotropic conical shell example (see Section 4.1.1). Here, the thickness of each layer is $h/3$ and a graphite-epoxy material is used (see Table IV).

As in the previous case, a very good agreement between the frequencies can be observed in Table V, thus validating our composite conical shell element. It is also observed that the influence of boundary conditions on the variation of frequency with circumferential wave number $n$ is significant. For the F-F conical shell, the frequencies increase as the circumferential wave number
<table><caption>Table IV. Properties of graphite-epoxy and PZT-5H materials.</caption><thead><tr><th>Properties</th><th>Graphite-epoxy</th><th>PZT-5H</th></tr></thead><tbody><tr><td>c<sub>11</sub> (GPa)</td><td>183.443</td><td>126</td></tr><tr><td>c<sub>12</sub> (GPa)</td><td>4.363</td><td>79.5</td></tr><tr><td>c<sub>13</sub> (GPa)</td><td>4.363</td><td>84.1</td></tr><tr><td>c<sub>22</sub> (GPa)</td><td>11.662</td><td>126</td></tr><tr><td>c<sub>23</sub> (GPa)</td><td>3.918</td><td>84.1</td></tr><tr><td>c<sub>33</sub> (GPa)</td><td>11.662</td><td>117</td></tr><tr><td>c<sub>44</sub> (GPa)</td><td>2.87</td><td>23</td></tr><tr><td>c<sub>55</sub> (GPa)</td><td>7.17</td><td>23</td></tr><tr><td>c<sub>66</sub> (GPa)</td><td>7.17</td><td>23.3</td></tr><tr><td>e<sub>15</sub> (C m<sup>-2</sup>)</td><td>0</td><td>17</td></tr><tr><td>e<sub>24</sub> (C m<sup>-2</sup>)</td><td>0</td><td>17</td></tr><tr><td>e<sub>31</sub> (C m<sup>-2</sup>)</td><td>0</td><td>-6.5</td></tr><tr><td>e<sub>32</sub> (C m<sup>-2</sup>)</td><td>0</td><td>-6.5</td></tr><tr><td>e<sub>33</sub> (C m<sup>-2</sup>)</td><td>0</td><td>23.3</td></tr><tr><td>d<sub>11</sub> (10<sup>-10</sup> F m<sup>-1</sup>)</td><td>153</td><td>150.3</td></tr><tr><td>d<sub>22</sub> (10<sup>-10</sup> F m<sup>-1</sup>)</td><td>153</td><td>150.3</td></tr><tr><td>d<sub>33</sub> (10<sup>-10</sup> F m<sup>-1</sup>)</td><td>153</td><td>130</td></tr><tr><td>ρ (kg m<sup>-3</sup>)</td><td>1590</td><td>7500</td></tr></tbody></table>

<table><caption>Table V. Frequencies (Hz) of a laminated [0, π/2, 0] conical shell.</caption><thead><tr><th rowspan="2">n</th><th colspan="2">F-F</th><th colspan="2">SS-SS</th><th colspan="2">C-C</th></tr><tr><th>Nastran</th><th>Present</th><th>Nastran</th><th>Present</th><th>Nastran</th><th>Present</th></tr></thead><tbody><tr><td>1</td><td>—</td><td>—</td><td>93.531</td><td>93.562</td><td>95.536</td><td>95.636</td></tr><tr><td>2</td><td>1.485</td><td>1.538</td><td>53.452</td><td>53.456</td><td>59.639</td><td>59.786</td></tr><tr><td>3</td><td>3.936</td><td>3.970</td><td>35.277</td><td>35.291</td><td>42.951</td><td>43.090</td></tr><tr><td>4</td><td>7.175</td><td>7.199</td><td>28.422</td><td>28.451</td><td>35.612</td><td>35.721</td></tr><tr><td>5</td><td>10.870</td><td>10.892</td><td>28.519</td><td>28.560</td><td>34.345</td><td>34.438</td></tr><tr><td>6</td><td>14.852</td><td>14.880</td><td>31.610</td><td>31.667</td><td>36.758</td><td>36.880</td></tr><tr><td>7</td><td>19.292</td><td>19.332</td><td>35.995</td><td>36.085</td><td>41.151</td><td>41.337</td></tr><tr><td>8</td><td>24.325</td><td>24.386</td><td>41.428</td><td>41.570</td><td>46.883</td><td>47.156</td></tr><tr><td>9</td><td>29.988</td><td>30.079</td><td>47.822</td><td>48.036</td><td>53.701</td><td>54.089</td></tr><tr><td>10</td><td>36.289</td><td>36.421</td><td>55.094</td><td>55.405</td><td>61.453</td><td>61.985</td></tr></tbody></table>

n increases. For the SS-SS and C-C conical shells, the frequencies decrease for circumferential wave number from n = 1 to 4 and increase from n = 5 to 9.

For illustration purpose, the mode shapes are presented in Figure 5. These three-dimensional deformations have been constructed from the one-dimensional element and using the Fourier series expansion of the displacement field in the circumferential direction.

### 4.1.4. Piezoelectric cylindrical shell.
In order to validate the formulation for piezoelectric shells, we present in this example a comparative study for a simply supported (SS-SS) piezoelectric
<figure><img src="image_6.png" /><figcaption>Figure 5. The five first mode shapes of composite conical shell: (a) free-free; (b) simply supported-simply supported; and (c) clamped-clamped.</figcaption></figure>

<table><caption>Table VI. Frequencies (Hz) of a piezoelectric cylindrical shell.</caption><thead><tr><th rowspan="2">n</th><th rowspan="2">m</th><th colspan="3">Short circuited</th><th colspan="3">Open circuited</th></tr><tr><th>Present</th><th>Exact</th><th>Error %</th><th>Present</th><th>Exact</th><th>Error %</th></tr></thead><tbody><tr><td rowspan="4">1</td><td>1</td><td>87.806</td><td>88.231</td><td>0.482</td><td>93.370</td><td>93.790</td><td>0.448</td></tr><tr><td>2</td><td>216.127</td><td>216.260</td><td>0.062</td><td>228.162</td><td>228.227</td><td>0.028</td></tr><tr><td>0</td><td>278.711</td><td>275.956</td><td>-0.998</td><td>278.711</td><td>275.956</td><td>-0.998</td></tr><tr><td>3</td><td>312.564</td><td>311.647</td><td>-0.294</td><td>331.697</td><td>330.649</td><td>-0.317</td></tr><tr><td>4</td><td>368.995</td><td>366.968</td><td>-0.553</td><td>383.703</td><td>392.134</td><td>2.150</td></tr><tr><td rowspan="5">2</td><td>1</td><td>36.705</td><td>37.316</td><td>1.639</td><td>39.489</td><td>40.109</td><td>1.546</td></tr><tr><td>2</td><td>113.540</td><td>114.210</td><td>0.587</td><td>121.979</td><td>122.677</td><td>0.569</td></tr><tr><td>3</td><td>193.357</td><td>193.791</td><td>0.224</td><td>207.472</td><td>207.916</td><td>0.213</td></tr><tr><td>4</td><td>259.388</td><td>259.346</td><td>-0.016</td><td>278.508</td><td>278.434</td><td>-0.027</td></tr><tr><td>5</td><td>308.942</td><td>308.464</td><td>-0.155</td><td>332.121</td><td>331.551</td><td>-0.172</td></tr></tbody></table>

cylindrical shell. The geometrical properties are $L = 5$ m, $R = 1$ m and $h = 0.02$ m. Moreover, the piezoelectric material is the PZT-5H (see Table IV).

Table VI presents the frequencies calculated by our method and those given by an exact three-dimensional solution proposed by Deü and Larbi [4]. This exact solution is based on a mixed state-space approach previously developed for the free vibration analysis of laminated piezoelectric plates actuated by transverse shear mechanisms [25]. In this table, the modes are classified, for two values of circumferential harmonic ($n$), in order of increasing frequency. To specify the mode type, the corresponding axial wave number ($m$) obtained with the exact solution is also given.
<table><caption>Table VII. Frequencies of a cylindrical shell coupled with a compressible fluid.</caption><thead><tr><th>m, n</th><th>Rigid cavity</th><th>Empty shell</th><th>Empty shell [1]</th><th>Coupled pb</th><th>Coupled pb [1]</th></tr></thead><tbody><tr><td>1,2</td><td>651.04</td><td>12.17</td><td>12.17</td><td>4.77</td><td>4.93</td></tr><tr><td>1,3</td><td>822.10</td><td>19.62</td><td>19.61</td><td>8.53</td><td>8.94</td></tr><tr><td>2,3</td><td>953.94</td><td>23.45</td><td>23.28</td><td>10.24</td><td>10.64</td></tr><tr><td>2,2</td><td>701.41</td><td>27.72</td><td>28.06</td><td>10.93</td><td>11.48</td></tr><tr><td>3,3</td><td>955.90</td><td>32.46</td><td>31.98</td><td>14.25</td><td>14.66</td></tr><tr><td>1,4</td><td>934.42</td><td>36.44</td><td>36.47</td><td>17.21</td><td>18.26</td></tr><tr><td>2,4</td><td>1192.97</td><td>37.42</td><td>37.37</td><td>17.72</td><td>18.73</td></tr><tr><td>3,4</td><td>1194.38</td><td>40.09</td><td>39.78</td><td>19.08</td><td>19.96</td></tr></tbody></table>

As it can be observed from this table, there is a good agreement between finite element (with 100 conical shell elements) and exact solutions for different electric boundary conditions corresponding to short-circuited ($\phi = 0$) or open-circuited ($D = 0$) piezoelectric materials. As expected, the natural frequencies are higher in the open-circuit case than in the closed-circuit one. It is noted that the error committed by the finite element approximation is lower than 3%. This comparison enables us to validate the piezoelectric aspect of the shell. It is important to note that this validation is restricted to the first modes, i.e. when a linear piezoelectric potential along the thickness of the piezoelectric layer is sufficient.

## 4.2. Free vibration of cylindrical shells filled with a compressible fluid

### 4.2.1. Elastic case: comparison with a wave propagation approach.
In order to validate the proposed formulation in a fluid-structure case, we consider an elastic isotropic cylinder filled with compressible fluid. The shell is clamped at both ends and has the following properties: length $L = 20$ m, radius $R = 1$ m, thickness $h = 0.01$ m, Young's modulus $E = 210$ GPa, Poisson ratio $\nu = 0.3$, and mass density $\rho_S = 7850$ kg m$^{-3}$. The considered fluid is a liquid with a mass density $\rho_F = 1000$ kg m$^{-3}$ and a speed of sound $c_F = 1500$ m s$^{-1}$. Moreover, 40 conical shell elements and 200 quadrangular axisymmetric fluid elements are used in the simulations.

Table VII gives the frequencies of (i) the fluid in rigid cavity, (ii) the empty shell, and (iii) the fluid-structure-coupled system. The structural frequencies in the coupled and uncoupled cases are compared to those obtained with a wave propagation method [1]. This comparison enables us to check the validity of the fluid-structure proposed formulation. It can be observed that the coupled frequencies correspond to half of the empty shell frequencies. This is due to the added mass effect of the fluid which is, in this particular case of a liquid, quasi-incompressible for the first low-frequency modes. In order to illustrate the results in the coupled case, the mode shapes, in terms of fluid pressure and structure displacement, are presented in Figure 6.

### 4.2.2. Elastic case: comparison with an exact state-space solution.
In our second fluid-structure example, we consider a simply supported isotropic cylindrical shell containing a compressible fluid. The shell has the following geometrical and material properties: length $L = 20$ m, radius $R = 1$ m, thickness $h = 0.002$ m, Young's modulus $E = 2.07788$ GPa, Poisson ratio $\nu = 0.317756$ and mass density $\rho_S = 8166$ kg m$^{-3}$. The considered fluid and the finite element discretization are the same as in the preceding example.
<figure><img src="image_7.png" /><figcaption>Figure 6. Fluid-structure-coupled modes: fluid pressure and shell total displacement.</figcaption></figure>

<table><caption>Table VIII. Frequencies (Hz) of an empty cylindrical shell for *n* varying from 1 to 10.</caption><thead><tr><th>n</th><th>FE</th><th>Exact [4]</th><th>Loy et al. [26]</th></tr></thead><tbody><tr><td>1</td><td>13.587</td><td>13.561</td><td>13.548</td></tr><tr><td>2</td><td>4.637</td><td>4.595</td><td>4.592</td></tr><tr><td>3</td><td>4.288</td><td>4.257</td><td>4.263</td></tr><tr><td>4</td><td>7.234</td><td>7.211</td><td>7.225</td></tr><tr><td>5</td><td>11.546</td><td>11.519</td><td>11.542</td></tr><tr><td>6</td><td>16.899</td><td>16.863</td><td>16.897</td></tr><tr><td>7</td><td>23.245</td><td>23.197</td><td>23.244</td></tr><tr><td>8</td><td>30.574</td><td>30.511</td><td>30.573</td></tr><tr><td>9</td><td>38.882</td><td>38.802</td><td>38.881</td></tr><tr><td>10</td><td>48.168</td><td>48.069</td><td>48.168</td></tr></tbody></table>

In the case where the fluid is not taken into account, Table VIII shows a comparison, for *n* varying from 1 to 10, between our finite element results, an exact state-space solution based on the three-dimensional theory of elasticity [4], and the results given by Loy et al. [26] using a Love shell theory and a Ritz method. A very good agreement between the three approaches can be observed.
<table><caption>Table IX. Frequencies (Hz) of an elastic cylindrical shell filled with a compressible fluid: comparison with an exact state-space solution [4].</caption><thead><tr><th rowspan="2">m</th><th colspan="2">n = 1</th><th colspan="2">n = 2</th><th colspan="2">n = 3</th><th colspan="2">n = 4</th></tr><tr><th>FE</th><th>Exact</th><th>FE</th><th>Exact</th><th>FE</th><th>Exact</th><th>FE</th><th>Exact</th></tr></thead><tbody><tr><td>1</td><td>2.408</td><td>2.441</td><td>0.899</td><td>0.912</td><td>0.947</td><td>0.968</td><td>1.779</td><td>1.839</td></tr><tr><td>2</td><td>9.073</td><td>9.179</td><td>3.409</td><td>3.453</td><td>2.040</td><td>2.068</td><td>2.142</td><td>2.193</td></tr><tr><td>3</td><td>18.633</td><td>18.817</td><td>7.484</td><td>7.559</td><td>4.243</td><td>4.280</td><td>3.254</td><td>3.290</td></tr><tr><td>4</td><td>29.547</td><td>29.849</td><td>12.902</td><td>12.987</td><td>7.357</td><td>7.391</td><td>5.144</td><td>5.155</td></tr><tr><td>5</td><td>40.536</td><td>41.147</td><td>19.405</td><td>19.470</td><td>11.283</td><td>11.283</td><td>7.701</td><td>7.666</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>Figure 7. Variation of frequency with: (a) circumferential wave number *n* and (b) axial wave number *m*.</figcaption></figure>

<table><caption>Table X. Frequencies (Hz) of a piezoelectric cylindrical shell filled with a compressible fluid: comparison with an exact solution [4].</caption><thead><tr><th rowspan="2">n</th><th rowspan="2">m</th><th colspan="3">Short circuited</th><th colspan="3">Open circuited</th></tr><tr><th>Present</th><th>Exact</th><th>Error %</th><th>Present</th><th>Exact</th><th>Error %</th></tr></thead><tbody><tr><td rowspan="5">1</td><td>1</td><td>45.445</td><td>44.760</td><td>-1.531</td><td>49.491</td><td>48.271</td><td>-2.528</td></tr><tr><td>2</td><td>108.920</td><td>107.390</td><td>-1.424</td><td>115.311</td><td>115.614</td><td>0.263</td></tr><tr><td>3</td><td>154.565</td><td>154.131</td><td>-0.281</td><td>166.237</td><td>166.065</td><td>-0.103</td></tr><tr><td>4</td><td>191.016</td><td>189.327</td><td>-0.892</td><td>209.635</td><td>204.211</td><td>-2.656</td></tr><tr><td>5</td><td>220.018</td><td>217.190</td><td>-1.302</td><td>236.253</td><td>234.469</td><td>-0.761</td></tr><tr><td rowspan="5">2</td><td>1</td><td>19.637</td><td>19.916</td><td>1.404</td><td>21.001</td><td>21.464</td><td>2.157</td></tr><tr><td>2</td><td>62.497</td><td>62.597</td><td>0.159</td><td>66.928</td><td>67.623</td><td>1.027</td></tr><tr><td>3</td><td>108.109</td><td>108.575</td><td>0.430</td><td>116.416</td><td>117.248</td><td>0.710</td></tr><tr><td>4</td><td>148.610</td><td>149.022</td><td>0.276</td><td>161.251</td><td>160.912</td><td>-0.211</td></tr><tr><td>5</td><td>182.892</td><td>182.817</td><td>-0.041</td><td>197.330</td><td>197.400</td><td>0.036</td></tr></tbody></table>

Table IX presents the natural frequencies for the fluid-structure-coupled problem. As can be seen in this table, the coupled frequencies obtained with our finite element approach agree very well with those obtained from the exact solution [27] for different values of circumferential ($n = 1–4$) and axial ($m = 1–5$) harmonics. It is also observed in Figure 7(a) that for all axial
wave number *m*, the frequencies decrease as the circumferential wave number *n* increases except for *m* = 1. Moreover, Figure 7(b) shows, for all circumferential wave number *n*, that the frequencies increase as the axial wave number *m* increases.

4.2.3. *Piezoelectric case: comparison with a state-space exact solution.* In this last example, the piezoelectric cylinder studied in Section 4.1.4 is filled with a compressible fluid characterized by $\rho_F = 1000 \text{ kg m}^{-3}$ and $c_F = 1500 \text{ m s}^{-1}$. Table X presents the eigenfrequencies of the coupled system computed by our finite element method and those given by the exact state-space approach [4]. The results, given for the two first circumferential harmonics, show a good agreement between the two methods. However, as previously mentioned, the difference between exact and finite element solutions can increase for higher modes, in particular due to the linearity assumption of the electric potential.

5. CONCLUSIONS

This paper concerns variational formulation and the finite element implementation of the fluid/
piezoelectric structure vibration problems. The coupled system consists of a piezoelectric structure
(described by its displacement field and its electric potential) containing a compressible fluid
(described by its pressure field). With this description, the variational formulation of the coupled
problem is established and the corresponding matrix equations are presented. This formulation
is then applied to the vibration of axisymmetric conical elastic/piezoelectric composite shells
coupled with an internal fluid. The piezoelectric shell is based on Kirchhoff-Love theory and
is supposed to be polarized in thickness direction. The associated finite element formulation
combines an equivalent single-layer approach for the mechanical behaviour with a layerwise
representation of the electric potential in the thickness direction. The fluid is modelled by the
Helmholtz equation in terms of pressure. The particular geometry of the structure is used to
reduce the dimensions of the problem through a semi-analytical procedure whereby variables
are expanded as Fourier series. Following these considerations, the shell is discretized by a one-
dimensional element with two nodes, and the fluid by a quadrilateral axisymmetric element with
four nodes. Moreover, an appropriate interface element is developed to ensure the interaction
between fluid and structure. Finally, numerical examples are presented validating the development
of the composite axisymmetric shell element, the piezoelectric aspect and the fluid-structure finite
element implementation. In this paper, we have used a non-symmetric fluid/piezoelectric structure
formulation. The symmetrization can be carried out using procedures as described in [16] and will
be the purpose of a further publication.

ACKNOWLEDGEMENTS

This research is partially carried out within the 6th Framework Programme of the European Commission (InMAR project) under Grant No. NMP2-CT-2003-501084.

REFERENCES

1. Zhang XM, Liu GR, Lam KY. Coupled vibration analysis of fluid-filled cylindrical shells using the wave propagation approach. *Applied Acoustics* 2001; **62**(3):229-243.

2. Ding HJ, Chen WQ, Guo YM, Yang QD. Free vibration of piezoelectric cylindrical shells filled with compressible fluid. *International Journal of Solids and Structures* 1997; **34**(16):2025-2034.
3. Chen WQ, Bian ZG, Lv CF, Ding HJ. 3D free vibration analysis of a functionally graded piezoelectric hollow cylinder filled with compressible fluid. *International Journal of Solids and Structures* 2004; **41**(3):947-964.

4. Deü J-F, Larbi W. A state space method for free-vibration analysis of radially polarized laminated piezoelectric cylinder filled with fluid. In *Proceedings of the Eighth International Conference on Computational Structures Technology*, CST2006, Las Palmas de Gran Canaria, Spain, 12-15 September 2006.

5. Kaljevic I, Saravanos DA. Steady-state response of acoustic cavities bounded by piezoelectric composite shell structures. *Journal of Sound and Vibration* 1997; **204**(3):459-479.

6. Hussein M, Heyliger PR. Discrete layer analysis of axisymmetric vibrations of laminated piezoelectric cylinders. *Journal of Sound and Vibration* 1996; **192**(3):995-1013.

7. Benjeddou A, Gorge V, Ohayon R. Use of piezoelectric shear response in adaptive sandwich shells of revolution—part 1: theoretical formulation. *Journal of Intelligent Material Systems and Structures* 2001; **12**(4):235-245.

8. Pinto Correia IF, Mota Soares CM, Mota Soares CA, Herskovits J. Active control of axisymmetric shells with piezoelectric layers: a mixed laminated theory with a high order displacement field. *Computers and Structures* 2002; **80**(27-30):2265-2275.

9. Xi ZC, Yam LH, Leung TP. Free vibration of a laminated composite circular cylindrical shell partially filled with fluid. *Composites Part B: Engineering* 1997; **28**(4):359-374.

10. Zhang YL, Reese JM, Gorman DG. A comparative study of axisymmetric finite elements for the vibration of thin cylindrical shells conveying fluid. *International Journal for Numerical Methods in Engineering* 2002; **54**(1):89-110.

11. Kochupillai J, Ganesan N, Padmanabhan C. Model reduction for parametric instability analysis in shells conveying fluid. *Journal of Sound and Vibration* 2003; **262**(3):633-649.

12. Schotté J-S, Ohayon R. Incompressible hydroelastic vibrations, finite element modelling of the elastogravity operator. *Computers and Structures* 2005; **83**(2-3):209-219.

13. Andrianarison O, Ohayon R. Compressibility and gravity effects in internal fluid-structure vibrations: basic equations and appropriate variational formulations. *Computer Methods in Applied Mechanics and Engineering* 2006; **195**(17-18):1958-1972.

14. Bermúdez A, Rodríguez R. Finite element computation of the vibration modes of a fluid-solid system. *Computer Methods in Applied Mechanics and Engineering* 1994; **119**(3-4):355-370.

15. Tzou HS. *Piezoelectric Shells: Distributed Sensing and Control of Continua*. Kluwer Academic Publishers: Dordrecht, Boston, London, 1993.

16. Morand HJ-P, Ohayon R. *Fluid-Structure Interaction*. Wiley: New York, 1995.

17. Galucio AC, Deü J-F, Ohayon R. A fractional derivative viscoelastic model for hybrid active-passive damping treatments in time domain—application to sandwich beams. *Journal of Intelligent Material Systems and Structures* 2005; **16**(1):33-45.

18. Soedel W. *Vibrations of Shells and Plates* (3rd edn). Marcell Dekker: New York, 2004.

19. Reddy JN. *Mechanics of Laminated Composite Plates: Theory and Analysis*. CRC Press: Boca Raton, FL, 1997.

20. Saravanos DA. Mixed laminate theory and finite element for smart piezoelectric composite shell structures. *AIAA Journal* 1997; **35**(8):1327-1333.

21. Lammering R, Mesecke-Rischmann S. Multi-field variational formulations and related finite elements for piezoelectric shells. *Smart Materials and Structures* 2003; **12**(6):904-913.

22. Shu C. Free vibration analysis of composite laminated conical shells by generalized differential quadrature. *Journal of Sound and Vibration* 1996; **194**(4):587-604.

23. Lam KY, Loy CT. Analysis of rotating laminated cylindrical shells by different thin shell theories. *Journal of Sound and Vibration* 1995; **186**(1):23-35.

24. Zhang XM. Vibration analysis of cross-ply laminated composite cylindrical shells using the wave propagation approach. *Applied Acoustics* 2001; **62**(11):1221-1228.

25. Deü J-F, Benjeddou A. Free-vibration analysis of laminated plates with embedded shear-mode piezoceramic layers. *International Journal of Solids and Structures* 2005; **42**(7):2059-2088.

26. Loy CT, Lam KY, Reddy JN. Vibration of functionally graded cylindrical shells. *International Journal of Mechanical Sciences* 1999; **41**(3):309-324.

27. Larbi W, Deü J-F. Solution tridimensionnelle pour l'analyse des vibrations d'une structure cylindrique avec fluide interne. In *Proceedings of the First International Congress Design and Modelling of Mechanical Systems, CMSM'2005*, Hammamet, Tunisia, 23-25 March 2005.