# Wave Finite Element Formulation of the Acoustic Transmission Through Complex Infinite Plates

Jean-Loup Christen, Mohamed Ichchou, Abdelmalek Zine, Bernard Troclet

► To cite this version:

Jean-Loup Christen, Mohamed Ichchou, Abdelmalek Zine, Bernard Troclet. Wave Finite Element Formulation of the Acoustic Transmission Through Complex Infinite Plates. Acta Acustica united with Acustica, 2016, 102 (6), pp.984-991. 10.3813/AAA.919013. hal-01716013

HAL Id: hal-01716013

https://hal.science/hal-01716013

Submitted on 3 Nov 2020

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
# Wave Finite Element Formulation of the Acoustic Transmission Through Complex Infinite Plates

Jean-Loup Christen<fnref n="1" />), Mohamed Ichchou<fnref n="1" />), Abdelmalek Zine<fnref n="2" />), Bernard Troclet<fnref n="3" />,<fnref n="4" />)

## Summary

A finite element-based derivation of the transmission loss (TL) of anisotropic layered infinite plates is presented in this paper. The wave-finite element method (WFE) is used to represent the plate with a finite element model of a single unit cell. The incident acoustic field is a kown plane wave, and the reflected and transmitted pressures are supposed to be plane waves with unknown amplitudes and phases. The periodicity conditions on the unit cell allow to find a simple matrix equation linking the amplitudes of the transmitted and reflected fields as a function of the incident one. This approach is validated for several cases against classical analytical models for thin plates and sandwich constructions, where the results agree perfectly for a reasonable mesh size. The method is then used to study the effect of stacking order in a laminated composite plate. The main interest of the method is the use of finite elements, which enables a relative easy modelling since most packages readily include different formulations, compared to analytical models, where different formulations have to be implemented for every kind of material.

## 1. Introduction

Sound transmission modelling is a very important topic in the industry, yet it is often still difficult to calculate for complex geometries or materials. Several methods have been proposed for the computation of transmission loss through infinite plate-like structures. Analytic models based on the constitutive equations for infinite thin plates [1], sandwich constructions [2, 3, 4], and even double plates with periodical stiffeners [5] have been proposed in the past. The Transfer Matrix Method (TMM) [6, 7] is a general framework used to compute acoustic transmission features of infinite plane isotropic layerings comprising elastic solids, poroelastic materials and air gaps. It relies on an analytical formulation of the problem in each layer, thus making it necessary to derive new equations for every kind of material that can be encountered (isotropic or anisotropic, elastic or poroelastic).

Finite Element (FE) modelling has the practical advantage over analytical models that most constitutive laws and couplings between materials are already implemented in commercial codes. However, meshing a large and simple structure such as a plate can lead to high computational times. The Wave-Finite Element (WFE) method is a response to this issue in the case of large periodic structures.

It couples finite elements and the periodic structure theory (also called Bloch's theory), that has been successfully used in the past for the analysis of free wave propagation in mono- and bi-dimensional waveguides [8]. The WFE was also used to compute vibroacoustic indicators such as group velocities and modal densities for use in Statistical Energy Analysis (SEA) [9, 10]. Extensions of this method towards the calculation of forced responses of waveguides have also been proposed [11], without considering fluid-structure coupling. An application of 1D WFE formulation was used by Serra et al. [12] for finite multilayers comprising poroelastic materials, but the transfer approach necessitated meshing the whole surface of the plate. An assessment of wave methods can be found in [13], comparing 3 methods, namely TMM, the WFE+SEA procedure of [9] and a Rayleigh-Ritz procedure.

Since the periodic structure theory deals primarily with infinite structures, the WFE can be used to compute the transmission loss of infinite plane structures excited by a superposition of plane waves directly, without having to resort to SEA. The derivation of this application is proposed in this paper, making use of the WFE framework for a plate with homogeneous properties in its plane. This enables to model an arbitrary layering with the finite element model of a single unit cell, which is just 1 element wide in each of the plane's dimensions, as in [14] and [15].

<footnote n="1">1) LTDS, École Centrale de Lyon, 36 avenue Guy de Collongue, 69134 Écully, France.
jean-loup.christen@notk.org</footnote>

<footnote n="2">2) ICJ, École Centrale de Lyon</footnote>

<footnote n="3">3) Airbus Defence and Space, 66 Route de Verneuil, 78133 Les Mureaux Cedex, France</footnote>

<footnote n="4">4) LMT, ENS Cachan, Université Paris-Saclay, 61 Avenue du Président Wilson, 94230 Cachan, France</footnote>
This paper is structured as follows: the WFE model under forced plane wave is presented in Section 2, with the derivation of the formulation of the TL and absorption coefficient. The method is validated on several cases in Section 3: an isotropic thin plate, an orthotropic plate, a sandwich construction with a soft viscoelastic core and a stiff sandwich plate with honeycomb core. A discussion on the element size closes this section. Finally, Section 4 presents an application case in which the proposed method is used to study the effect of the order of the layers in a composite plate.

## 2. Forced response to an incident plane wave

We consider a multilayered plate with one surface lying in the plane $z = 0$. Each layer is considered homogeneous, while this is of course not the case for the assembly. As a homogeneous structure is periodic for any period, the unit cell is a parallelepiped. It is modelled with standard (8 nodes) brick elements. The FE discretisation has only a single element in each of the x and y directions, and as many as needed in the z direction.

### 2.1. Equilibrium formulation

Let **D** denote the dynamic stiffness matrix of the unit cell. This matrix is computed from the finite element stiffness and mass matrices **K**, **M** with the relationship

$$ \mathbf{D} = \mathbf{K} - \omega^2 \mathbf{M}. \qquad (1) $$

Viscoelastic damping can be taken into account through a complex part of the stiffness matrix **K**.

The equilibrium equation of the unit cell then writes

$$ \mathbf{D}\mathbf{u} = \mathbf{f} + \mathbf{e}, \qquad (2) $$

where **f** represents the forces due to the interaction with the neighbouring cells, and **e** the forces imposed on the plate by the acoustic pressure fields on each side. As can be seen in Figure 1, the displacement vector **u** can be partitioned as

$$ \mathbf{u} = \begin{pmatrix} \mathbf{u}_1 \\ \mathbf{u}_2 \\ \mathbf{u}_3 \\ \mathbf{u}_4 \end{pmatrix}. \qquad (3) $$

We now search for solutions as plane waves with the same wave characteristics as the projection of the incident wave on the plate, hence

$$ \begin{aligned} \mathbf{u}_2 &= \lambda_x \mathbf{u}_1, \\ \mathbf{u}_3 &= \lambda_y \mathbf{u}_1, \\ \mathbf{u}_4 &= \lambda_x \lambda_y \mathbf{u}_1, \end{aligned} $$

<figure><img src="image_2.png" /><figcaption>Figure 1. FE model of a unit cell and notation for the applied forces.</figcaption></figure>

where $\lambda_x = \exp(ik_x d_x)$ and $\lambda_y = \exp(ik_y d_y)$. The reduced variable $\mathbf{u}_1$ can be used with the transformation

$$ \mathbf{u} = \begin{pmatrix} \mathbf{I} \\ \lambda_x \mathbf{I} \\ \lambda_y \mathbf{I} \\ \lambda_x \lambda_y \mathbf{I} \end{pmatrix} \mathbf{u}_1 = \mathbf{T} \mathbf{u}_1, \qquad (4) $$

where **I** is the identity matrix.

The imposed pressure is condensed to forces exerted on the nodes, and can be partitioned in the same way,

$$ \mathbf{e} = \mathbf{T} \mathbf{e}_1. \qquad (5) $$

The internal forces exerted on each vertical segment are derived from the stress field applied to the normal to the two faces related to that segment. The linearity of the stress tensor with respect to displacement field **u** and the previous periodicity relations lead to

$$ \mathbf{f}_1 + \lambda_x^{-1} \mathbf{f}_2 + \lambda_y^{-1} \mathbf{f}_3 + \lambda_x^{-1} \lambda_y^{-1} \mathbf{f}_4 = 0, \qquad (6) $$

or $\tilde{\mathbf{T}}\mathbf{f} = 0$, where $\tilde{\mathbf{T}} = (\mathbf{I}, \lambda_x^{-1}\mathbf{I}, \lambda_y^{-1}\mathbf{I}, \lambda_x^{-1}\lambda_y^{-1}\mathbf{I})$. Equation (2) can then be premultiplied by $\tilde{\mathbf{T}}$ and rewritten as

$$ \tilde{\mathbf{T}}\mathbf{D}\mathbf{T}\mathbf{u}_1 = \tilde{\mathbf{T}}\mathbf{T}\mathbf{e}_1 = 4\mathbf{e}_1. \qquad (7) $$

Let's denote $\mathbf{A} = \frac{1}{4}\tilde{\mathbf{T}}\mathbf{D}\mathbf{T}$. We finally obtain the reduced equation,

$$ \mathbf{A}\mathbf{u}_1 = \mathbf{e}_1, \qquad (8) $$

linking the forces and the displacements on one single segment instead of the initial four.

### 2.2. Acoustic coupling

We are concerned with an oblique incident plane wave impinging the plate on one side, with a known amplitude $p_I$. The interaction between this incident wave and the plate create a reflected wave on the same side, and a transmitted wave on the other one, whose respective amplitudes $p_R$ and $p_T$ we want to calculate. The sound field on the
incident side is then the superposition of the incident and reflected waves,

$$ 
\begin{aligned}
P_I &= p(x, y, z)|_{z>0} \\
&= p_I \exp(-ik_{z,I}z) + p_R \exp(+ik_{z,I}z),
\end{aligned} 
\quad (9) $$

where the common factor $\exp(i(\omega t - k_x x - k_y y))$ has been omitted for legibility. The wavenumber components satisfy the relationship $k_x^2 + k_y^2 + (k_{z,I})^2 = k^2 = (\omega/c_I)^2$ and $c_I$ is the speed of sound in the fluid on the incident side. In the same way, the sound field on the transmission side contains only the transmitted wave propagating in the same direction as the incident wave,

$$ P_T = p(x, y, z)|_{z<0} = p_T \exp(-ik_{z,T}z), \quad (10) $$

with the relationship $k_x^2 + k_y^2 + (k_{z,T})^2 = (\omega/c_T)^2$, where $c_T$ is the speed of sound in the fluid on the transmission side. Potential phase difference between the pressure fields are accounted for through the fact that the amplitudes $p_R$ and $p_T$ may be complex. The wavenumbers $k_x$ and $k_y$ are conserved across the plate, so only the $k_z$ component may vary with the nature of the fluid. We will consider in the following that the fluid is the same on both sides, thus $\rho_I = \rho_{hO} = \rho_0$ and $c_I = c_T = c_0$.

In this case, the incidence angle $\theta$ is the angle between the wave vector and the normal of the plate, while the azimuthal angle $\phi$ gives its orientation in the plane. We therefore have $k_x = k \sin\theta \cos\phi$, $k_y = k \sin\theta \sin\phi$ and $k_z = k \cos\theta$.

The load imposed on the plate can be written from these two pressure fields lumped on the nodes of the finite element model. As the pressure force is exerted along the normal to the structure, the only non-zero terms in the vector $\mathbf{e}_1$ of forces imposed on segment 1 will be those relative to the DOFs in z-direction on the incident side, $e_I$, and on the transmission side $e_T$. For an elastic material, these quantities are scalar. Let $\mathbf{u}_O$ and $\mathbf{e}_O$ denote the displacement and force vectors on all other degrees of freedom in the segment. The force exerted on the nodes of segment 1 is then, as shown in Figure 1

$$ \mathbf{e}_1 = \begin{pmatrix} e_I \\ \mathbf{e}_O \\ e_T \end{pmatrix} = \begin{pmatrix} S(p_I + p_R) \\ \mathbf{0} \\ Sp_T \end{pmatrix}, \quad (11) $$

where $S$ is the free surface of the element, which is identical on both sides. All other vectors and matrices can be written following the same decomposition, allowing to rewrite (8) as

$$ \begin{pmatrix} a_{II} & A_{IO} & a_{IT} \\ A_{OI} & A_{OO} & A_{OT} \\ a_{TI} & A_{TO} & a_{TT} \end{pmatrix} \begin{pmatrix} u_I \\ \mathbf{u}_O \\ u_T \end{pmatrix} = \begin{pmatrix} S(p_I + p_R) \\ \mathbf{0} \\ Sp_T \end{pmatrix}. \quad (12) $$

The second line of (12) leads to

$$ \mathbf{u}_O = -\mathbf{A}_{OO}^{-1} (\mathbf{A}_{OI} \mathbf{u}_I + \mathbf{A}_{OT} \mathbf{u}_T), \quad (13) $$

which allows to condense this equation into

$$ \begin{pmatrix} b_{II} & b_{IT} \\ b_{TI} & b_{TT} \end{pmatrix} \begin{pmatrix} u_I \\ u_T \end{pmatrix} = \begin{pmatrix} S(p_I + p_R) \\ Sp_T \end{pmatrix}. \quad (14) $$

The fluid-structure interaction is characterized by the continuity of normal particle velocity at the interface. This writes for the incident side,

$$ \rho_0 \omega^2 u_I = \frac{\partial P_I}{\partial z} = -ik_z (p_I - p_R), \quad (15) $$

and for the transmission side:

$$ \rho_0 \omega^2 u_T = \frac{\partial P_T}{\partial z} = -ik_z p_T. \quad (16) $$

Introducing acoustic admittance $Y_0 = \cos(\theta)/(i\omega\rho_0 c_0)$, we obtain

$$ u_I = Y_0 (p_I - p_R), $$

$$ u_T = Y_0 p_T. $$

These expressions for $u_I$ and $u_T$ can be re-injected into (14), leading to two scalar equations linking the unknowns $p_R$ and $p_T$ and the incident pressure $p_I$,

$$ \begin{pmatrix} b_{II} + \frac{S}{Y_0} & -b_{IT} \\ b_{TI} & -b_{TT} - \frac{S}{Y_0} \end{pmatrix} \begin{pmatrix} p_R \\ p_T \end{pmatrix} = p_I \begin{pmatrix} b_{II} - \frac{S}{Y_0} \\ b_{TI} \end{pmatrix}. \quad (17) $$

Solving this equation gives the acoustic transparency $\tau = |p_T/p_I|^2$ and the absorption coefficient $\alpha = 1 - |p^R/p^I|^2$. The transmission loss for a plane wave is then finally

$$ \mathrm{TL} = -10 \log_{10} \tau. \quad (18) $$

This leads also to the diffuse field transmission loss by integrating over all possible incidences between minimum and maximum bounds $\theta \in [\theta_{min}, \theta_{max}]$ and directions $\phi \in [0, 2\pi]$. It is recommended to avoid the full range $[0, \frac{\pi}{2}]$ for $\theta$ because numerical errors occur close to grazing incidences $\theta = \pi/2$. However, $\theta_{min} = 0$ causes no issue. The diffuse field transparency is then

$$ \tau_d = \frac{\int_0^{2\pi} \int_0^{\theta_{max}} \tau(\omega\theta, \phi) \sin\theta \cos\theta d\theta d\phi}{\int_0^{2\pi} \int_0^{\theta_{max}} \sin\theta \cos\theta d\theta d\phi}, \quad (19) $$

and the diffuse field transmission loss

$$ \mathrm{TL}_d(\omega) = -10 \log_{10} \tau_d(\omega). \quad (20) $$

### 3. Validation results

In order to assess the validity of the method, a comparison with various analytical models is performed in this section. A discussion on the size of the elements to be used is proposed in Section 3.5.
<figure><img src="image_3.png" /><figcaption>Figure 2. Transmission loss of an infinite steel plate for an oblique plane wave with $\theta = \pi/4$ o: analytic model (plate), +:TMM, -: current method, $d_x = 0.1$ mm.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3. Diffuse field TL of an infinite steel plate +: analytical model (thin plate), -: current method, $d_x = 0.1$ mm.</figcaption></figure>

### 3.1. Isotropic plate

The first case is an isotropic aluminium plate, under two different excitations: a 45° plane wave and a diffuse field. The material parameters are the Young modulus $E = 70$ GPa, the Poisson ratio $\nu = 0.33$, the density $\rho = 2700$ kg.m$^{-3}$. The loss factor is $\eta = 0.5\%$, which leads to a complex Young's modulus $E = E(1 + i\eta)$. The unit cell consists of 10 elements through the thickness, with an overall thickness of 2 mm and lengths $d_x = d_y = 0.1$ mm.

Two kinds of analytical models are used for comparison. The first one is based on the Kirchhoff-Love plate assumption, while the other is the Transfer Matrix Method (TMM) described in [7]. All three model agree very well with each other, with a difference of less than 1 dB over most of the frequency range, except in a tiny region around the coincidence frequency, which is slightly underestimated in the thin plate model. Because the TMM and the current model agree in this region, this difference must be due to the fact that the Kirchhoff-Love assumptions neglect shear, while the other two use a 3D elasticity theory.

The diffuse field integration has been performed over the range $[0, 0.999\pi/2]$ to avoid numerical issues with grazing incidence. The range was discretised with 1000 points. The three models are found in this case to agree up to 1 dB over the whole frequency range, the small difference at coincidence disappearing because of the averaging over 1000 incidence angles.

### 3.2. Orthotropic plate

The diffuse field transmission loss of an orthotropic plate is compared against an analytical model [16]. The material parameters are given in Table I, where $E$ stands for Young's modulus, $\nu$ for Poisson's ratios and $G$ for shear modulus. The TL results are presented in Figure 4 for frequencies ranging between 100 Hz and 20 kHz. The mesh of the unit cell has 10 elements in the $z$ direction and lengths $d_x = d_y = 0.1$ mm. Here again, the analytical and current models agree well with each other, with an error

<table><caption>Table I. Parameters of the orthotropic plate.</caption><tr><td>E<sub>x</sub><br>(GPa)</td><td>E<sub>y</sub> = E<sub>z</sub><br>(GPa)</td><td>v<sub>xy</sub> = v<sub>xz</sub><br>(-)</td><td>v<sub>yz</sub><br>(-)</td><td>G<sub>xy</sub> = G<sub>xz</sub><br>(GPa)</td></tr><tr><td>224</td><td>6.9</td><td>0.25</td><td>0.3</td><td>56.58</td></tr><tr><td>G<sub>yz</sub><br>(GPa)</td><td>ρ<br>(kg.m<sup>-3</sup>)</td><td>η<br>(-)</td><td>h<br>(mm)</td><td></td></tr><tr><td>1.38</td><td>1578</td><td>0.5%</td><td>50</td><td></td></tr></table>

<figure><img src="image_5.png" /><figcaption>Figure 4. Diffuse field TL of an orthotropic plate + analytic, -: current method. The difference is less than 1 dB on the whole range (100 Hz, 20 kHz).</figcaption></figure>

smaller than 1 dB over the whole range. There are two critical frequencies due to the very different Young's moduli in the $x$ and $y$ directions. An expression for these frequencies can classically be derived from the analytical model, and is given by the formula

$$f_{\text{crit},x} = \frac{c_0^2}{2\pi} \sqrt{\frac{12(1 - v_{xy}v_{yx})\rho}{E_x h^2}}, \quad (21)$$

where $v_{xy} = (E_y/E_x)v_{xy}$. The other critical frequency $f_{\text{crit},y}$ is obtained by replacing $E_x$ by $E_y$ in (21). In this
case, the frequencies obtained by these formulas agree well with the dips observed in the TL curve.

### 3.3. Stiff sandwich plate

The 6th order constitutive equation for the bending of a sandwich beam was given by Mead and Markus [17] and later extended to a plate by Narayanan and Shanbhag [3]. The main assumptions are that the predominant deformation of the core is due to shear, while that of the skins is due to pure bending. This implies that the Young modulus of the core is high compared to the out-of-plane shear modulus. The model has five main parameters, namely skin bending stiffness $D_t$, overall bending stiffness $B$, damping $\eta$, surface mass $m$ and shear parameter $g$. The constitutive equation reads

$$D_t(1 + i\eta)\nabla^6 w - g(D_t + B)(1 + i\eta)\nabla^4 w + m\frac{\partial^2 w}{\partial t^2} - \frac{mB}{N}\frac{\partial^2}{\partial t^2}\nabla^2 w = (\nabla^2 - g)(p_L - p_R), \quad (22)$$

where $w = v/i\omega$ is the normal displacement of the plate. In the considered frame where a forced wave is imposed on the plate with a wavenumber $k = \frac{\omega}{c} \sin\theta$, the spatial derivative operator $\nabla$ can be replaced by $-ik$. This leads to the following expression of the impedance:

$$Z(\omega, \theta) = \frac{D_t k^6 + g(D_t + B)k^4 - m\omega^2 (k^2 + g(1 - v_s^2))}{i\omega(k^2 + g)} \quad (23)$$

For sandwiches made of isotropic materials and identical skins, the skin bending stiffness is $D_t = E_s h_s^3 / [6(1 - v_s^2)]$ the overall bending stiffness is $B = E_s h_c^2 h_s (1 + h_s / h_c)^2 / 2$ and the shear parameter is $g = 2G_c / (E_s h_s h_c)$, where $G_c = G_{c,xz} = G_{c,yz}$ is the core transverse shear modulus, $E_s$ and $v_s = v_{s,xy}$ are the in-plane Young's modulus and Poisson ratio of the skins, and $h_s$ (resp. $h_c$) is the thickness of one skin (resp. the core). In-plane core shear modulus and all Young's moduli of the core are neglected in the analytical model.

The two models are compared between 1 kHz and 20 kHz on a sandwich plate with isotropic skins and orthotropic core impinged by a 45° plane wave. The parameters of the sandwich construction are given in Table II, and the TL results shown in Figure 6. With a mesh using 3 elements for each skin and 45 for the core, and again $d_x = d_y = 0.1$ mm, the two models differ also by less than 1 dB in the whole range. Since the analytical model takes shear into account, the agreement is also excellent at coincidence.

### 3.4. Sandwich plate with soft core

The same calculation was made with a thick sandwich plate made of a rather soft isotropic core and thick skins. The material properties are described in Table III. The transmission loss curve is shown on Figure 6. Again this correlates well with the analytical TMM, with a maximum error of 1 dB over the considered range. Three dips in TL

<table><caption>Table II. Parameters of the stiff sandwich plate.</caption><thead><tr><th></th><th><i>E</i><sub>x</sub> (GPa)</th><th><i>E</i><sub>z</sub> (GPa)</th><th><i>ρ</i> (kg.m<sup>-3</sup>)</th><th><i>ν</i></th></tr></thead><tbody><tr><td>Skin</td><td>70</td><td>70</td><td>2700</td><td>0.33</td></tr><tr><td>Core</td><td>0.08</td><td>8</td><td>8</td><td>0.1</td></tr></tbody><thead><tr><th></th><th><i>G</i><sub>xy</sub> (MPa)</th><th><i>G</i><sub>xz</sub> (MPa)</th><th><i>G</i><sub>yz</sub> (MPa)</th><th><i>h</i> (mm)</th></tr></thead><tbody><tr><td>Skin</td><td>26310</td><td>26310</td><td>26310</td><td>2</td></tr><tr><td>Core</td><td>0.04</td><td>4</td><td>4</td><td>30</td></tr></tbody></table>

<table><caption>Table III. Parameters of the soft sandwich plate.</caption><thead><tr><th></th><th><i>E</i> (MPa)</th><th><i>ν</i></th><th><i>ρ</i> (kg.m<sup>-3</sup>)</th><th><i>h</i> (mm)</th></tr></thead><tbody><tr><td>Skin</td><td>8300</td><td>0.15</td><td>629.9</td><td>6.35</td></tr><tr><td>Core</td><td>8.3</td><td>0</td><td>16</td><td>38.1</td></tr></tbody></table>

<figure><img src="image_6.png" /><figcaption>Figure 5. TL of a sandwich plate impinged by a 45° plane wave. +: Analytical model, — current method.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 6. Transmission loss for the sandwich plate impinged by a 45° plane wave o TMM (analytic) — current method.</figcaption></figure>

can be observed between 1 kHz and 6 kHz, corresponding to coincidence frequencies. The first two drops around 1730 Hz and 5470 Hz are due to coincidences with a wave having a symmetric motion of the skin and compression
<figure><img src="image_8.png" /><figcaption>Figure 7. Dispersion curves of the out-of-plane waves of the soft sandwich structure and trace wavenumber of the incident acoustic plane wave at 45°.</figcaption></figure>

in the core, while the third one around 5620 Hz is due to coincidence with the bending wave. The dispersion curves for these two waves are presented in Figure 7.

The displacement field can be calculated back from (13). They are displayed at several frequencies on Figure 8. The first one is related to the symmetric motion of the skins (Figure 8a), while the second is related to a global bending motion of the sandwich, or antisymmetric motion of the skins (Figure 8b). In both cases, it can be observed that a shearing motion of the skins is involved, due to their rather high thickness. In the antisymmetric case, the cross-section remain straight and vertical, indicating that only shear is at stake, with practically no bending motion.

The field represented on Figure 8c is at a frequency between the two coincidences, showing a decoupling of the two skins. On the latter example, it may be noted that the amplitude of the z-displacement is higher on the transmission side than on the incident side, yet the amplitude of the transmitted wave is of course smaller than that of the incident wave. This is due to the fact that the pressure field on the incident side is result of interference between incident and reflected waves, while the transmitted field consists of a single plane wave.

### 3.5. Discussion on element size

In the validation cases presented in this section, the mesh parameters are chosen so as to ensure convergence with the TMM model, which is based on 3D isotropic elastic theory and therefore exact in this case. A coarser mesh leads to errors in the coincidence frequency region. Of all the waves considered in acoustic transmission, the acoustic wave in the fluid has the shortest wavelength, at the highest frequency (20 kHz) considered in all these cases. Since it is imposed on the structure at grazing incidence, it is relevant to compare the mesh size to it. The element sizes in x and y directions must therefore satisfy about 100 elements per wavelength, which is much stricter than the usual 10 elements per wavelength used for finite elements.

<figure><img src="image_9.png" /><figcaption>Figure 8. Displacement fields at several frequencies. The color indicates the absolute value of the z-displacement, in meters, for an incident amplitude of 1Pa. The incident wave comes from the top. (a) Symmetric coincidence at 5472 Hz. (b) Antisymmetric coincidence at 5623 Hz. (c) Displacement field at 5548 Hz.</figcaption></figure>

Since there is only one element in the in-plane directions, this has no impact on the number of degrees of freedom in the mesh. However, the results are also found to be quite sensitive to the aspect ratio of the elements, which means that the element size in the z-direction must be about the same as that in the plane. For thick layers, a small $d_x$ may thus lead to a high number of elements through the thickness. On the other hand, it also sets a maximum size for $d_x$ and $d_y$, since there must be at least one element per layer. Above 10 kHz, this latter limit is usually less restrictive than the “100 element per wavelength” criterion for plates with thickness in the order of 1cm or less.

### 4. Influence of stacking on the transmission loss

The current approach can be used to study the influence of the stacking sequence of a multilayer plate on the TL. We consider a 1.6 mm thick composite plate made of 24 layers of the same transverse isotropic material described in Table IV in order to study the effect of the stacking sequence on the TL. Three different stacking sequences have been studied, which are described in Table V. All of these have the same proportion of layers with fibres oriented towards each direction, at -45, 0, 45 or 90 degrees of the global
<table><caption>Table IV. Material parameters of a single layer. For a transverse isotropic material, we have $E_y = E_z$, $v_{xz} = v_{xy}$ and $G_{xz} = G_{xy}$.</caption><thead><tr><th>E<sub>x</sub><br>(GPa)</th><th>E<sub>z</sub><br>(GPa)</th><th>v<sub>xy</sub><br>(-)</th><th>v<sub>yz</sub><br>(-)</th><th>G<sub>xy</sub><br>(GPa)</th></tr></thead><tbody><tr><td>175.0</td><td>6.90</td><td>0.2542</td><td>0.4689</td><td>4.18</td></tr><tr><th>G<sub>yz</sub><br>(GPa)</th><th>ρ<br>(kg.m<sup>-3</sup>)</th><th>η<br>(%)</th><th>h<br>(μm)</th><td></td></tr><tr><td>2.35</td><td>1520</td><td>0.1</td><td>67</td><td></td></tr></tbody></table>

<table><caption>Table V. The three stacking configurations.</caption><tr><td>Initial</td><td>[45<sub>3</sub>; -45<sub>3</sub>; 0<sub>3</sub>; 90<sub>3</sub>]<sub>s</sub></td></tr><tr><td>Sorted</td><td>[-45<sub>6</sub>; 0<sub>6</sub>; 45<sub>6</sub>; 90<sub>6</sub>]</td></tr><tr><td>Random</td><td>[0; 90; -45; 90; 45; 0<sub>3</sub>; -45; 90; -45; 45; 90; -45; 0; 90; -45; 45<sub>2</sub>; 90; -45; 45; 0; 45]</td></tr></table>

coordinate system, orientation 0 being aligned with the x axis. The layers are supposed to be perfectly bonded to each other, so the FE model ensures continuity of displacements at the interfaces. The three configuration are a classical symmetric quasi-isotropic stacking (“initial” configuration in Table V), a random permutation of it (“random” case), and a grouping of the same orientation in ascending order (“sorted” case).

The TL was computed for these three sequences with the method described above, using 2 elements per layer and a mesh size of 0.1 mm in the x and y directions. The results are presented in Figure 9.

Concerning the low frequency behaviour, it can be seen that the stacking sequence has practically no influence for frequencies below 2 kHz, and little up to 4 kHz. This is due to the fact that below this frequency, the plate behaves according to the mass law, and the actual stiffness has no effect on TL.

In high frequencies (above 8 kHz), the random and initial stacking sequences have the same asymptotic behaviour. This is due to the dominant effect of damping on the high-frequency TL, above the highest critical frequency. These effects have been studied in [18]. The sorted case exhibits the same asymptotic behaviour, which is only visible for much higher frequencies, above 20 kHz and thus not represented on Figure 9.

Between these, all configurations present a coincidence region, characterized by a lower TL. Because of the presence of different layer orientations, the coincidence zone of all plates is strictly contained between the two extreme critical frequencies of a plate where all layers are identically oriented. These bounds are given by (21), which yields for this material $f_{crit,x} = 3770$ Hz and $f_{crit,y} = 19$ kHz. The random and initial configurations exhibit only one dip, indicating that they behave closely to an isotropic plate in this respect. This is due to a good repartition of the layer orientations throughout the thickness. On the other hand, the sorted stacking exhibits a wide coincidence region with two critical frequencies, a behaviour that is similar to that of an orthotropic plate.

<figure><img src="image_10.png" /><figcaption>Figure 9. Diffuse field transmission loss for 3 different stacking sequences.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Figure 10. Wavenumbers as a function of direction angle $\phi$ at 7 kHz for the 3 stacking sequences.</figcaption></figure>

The coincidence phenomenon is illustrated in Figure 10, which shows the free propagating wavenumbers as a function of direction $\phi$ at a frequency of 7 kHz, which is in the coincidence region for all 3 cases. These results were obtained with a WFE calculation using the same mesh. The acoustic trace wavenumber $k_a = 2\pi f/c_0$ is represented for a grazing incidence, and does not depend on $\phi$. The random stacking k-space has a shape close to a circle, and lies mostly below $k_a$. The initial stacking has a more distorted k-space, partly above and partly below the acoustic wavenumber, but the deviation to the circle remains small, hence the narrow coincidence region. Finally, the sorted case has a very stretched k-space, close to an ellipse with major axis oriented around $-\pi/6$ radians, which also crosses the $k_a$ circle. This behaviour is related to the
strongly inhomogeneous nature of this stacking case: for
a wave with direction $\phi = 0$, the stiffest layer (oriented
along the fibers) is on the outer part of the plate, while a
wave with $\phi = \pi/4$ “sees” the stiffest layer close to the
central plane.

**5. Conclusion**

A new method for computing the transmission loss of
infinite plates made of arbitrarily layered materials has
been presented in this paper. The originality of the ap-
proach is that it couples the periodic structure theory to
finite elements as in the WFE method, while taking into
account the coupling with an acoustic fluid analytically.
The method is validated against classical analytic Transfer
Matrix Method computation for isotropic, orthotropic and
sandwich constructions. The main interest of the method
is that it natively takes into account complex behaviour of
the structure, such as symmetric motion of the skins with
respect to the neutral plane. Besides, the use of finite ele-
ments makes it easy to retrieve the displacement field in-
side the structure. The computation time is much higher
than that of analytical methods, but remains low, typically
up to about 10ms per frequency for a plane wave. A possi-
ble use of this method is the validation of analytical model
for multilayers.

An application case is proposed to study the effect of
the stacking sequence of a layered composite plate on the
diffuse field TL. It is shown that the stacking sequence has
an effect on the width of the coincidence range, and little
influence outside it, the low frequencies being governed by
mass law and the higher by damping effects. The most ex-
treme case of a sorted stacking shows a wide coincidence
zone limited by two critical frequencies, while the random
permutation exhibits a dip in a narrower zone around the
lowest critical frequency (7 kHz in this case). In applicat
ions where noise control is important at frequencies cor-
responding to the coincidence range, stacking sequence
should therefore be considered for acoustic criteria too.

Acknowledgements

The authors gratefully acknowledge Airbus Defence and
Space for their financial support. This work was partly
conducted under the SMART-NEST project funded by the
European Commission.

References

[1] F. J. Fahy, P. Gardonio: Sound and structural vibration.
Academic Press, 2007.

[2] G. Kurtze, B. G. Watters: New wall design for high trans-
mission loss or high damping. The Journal of the Acousti-
cal Society of America **31** (1959) 739–748.

[3] S. Narayanan, R. L. Shanbhag: Sound transmission through a damped sandwich panel. Journal of Sound and Vibration **80** (1982) 315–327.

[4] S. Kumar, L. Feng, U. Orrenius: Predicting the sound trans-
mission loss of honeycomb panels using the wave propaga-
tion approach. Acta Acustica unite with Acustica **97** (2011)
869–876.

[5] G.-F. Lin, J. M. Garrelick: Sound transmission through pe-
riodically framed parallel plates. The Journal of the Acous-
tical Society of America **61** (1977) 1014–1018.

[6] B. Brouard, D. Lafarge, J.-F. Allard: A general method of
modelling sound propagation in layered media. Journal of
Sound and Vibration **183** (1995) 129–142.

[7] J.-F. Allard, N. Atalla: Modelling multilayered systems with porous materials using the transfer matrix method. – In: Propagation of sound in porous media. John Wiley & Sons, Ltd., 2009, Kap. 11.

[8] C. Droz, C. Zhou, M. Ichchou, J.-P. Lainé: A hybrid wave-
mode formulation for the vibro-acoustic analysis of 2D
periodic structures. Journal of Sound and Vibration **363**
(2016) 285–302.

[9] V. Cotoni, R. Langley, P. Shorter: A statistical energy anal-
ysis subsystem formulation using finite element and peri-
odic structure theory. Journal of Sound and Vibration **318**
(2008) 1077–1108.

[10] D. Chronopoulos, B. Troclet, M. Ichchou, J.-P. Lainé: A unified approach for the broadband vibroacoustic response of composite shells. Composites: Part B **43** (2012) 1837–1846.

[11] J. M. Renno, B. R. Mace: Calculating the forced response of two-dimensional homogeneous media using the wave and finite element method. Journal of Sound and Vibration **330** (2011) 5913–5927.

[12] Q. Serra, M. N. Ichchou, J.-F. Deü: On the use of trans-
fer approaches to predict the vibroacoustic response of
poroelastic media. Journal of Computational Acoustics **23**
(2015) 1550020.

[13] U. Orrenius, H. Liu, A. Wareing, S. Finnveden, V. Cotoni:
Wave modelling in predictive vibro-acoustics: Applications
to rail vehicles and aircraft. Wave Motion **51** (2014) 635–
649.

[14] E. Manconi: Modelling wave propagation in two-dimen-
sional structures using a wave/finite element technique.
Dissertation. University of Parma, 2008.

[15] B. R. Mace, E. Manconi: Modelling wave propagation in
two-dimensional structures using finite element analysis.
Journal of Sound and Vibration **318** (2008) 884–902.

[16] J.-L. Guyader, C. Lesueur: Transmission of reverberant
sound through orthotropic viscoelastic multilayered plates.
Journal of Sound and Vibration **70** (1980) 319–332.

[17] D. J. Mead, S. Markus: The forced vibration of a three-
layer damped sandwich beam with arbitrary boundary con-
ditions. Journal of Sound and Vibration **10** (1969) 163–
175.

[18] J.-L. Christen, M. Ichchou, B. Troclet, O. Bareille, M.
Ouisse: Global sensitivity analysis of analytical vibroa-
coustic transmission models. Journal of Sound and Vibra-
tion **368** (2016) 121 – 134.