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

HAL
open science

# Wave Transmission and Reflection Analysis Based on the Three-dimensional Second Strain Gradient Theory

Bo Yang, Mohamed Ichchou, Christophe Droz, Abdelmalek Zine

► To cite this version:

Bo Yang, Mohamed Ichchou, Christophe Droz, Abdelmalek Zine. Wave Transmission and Reflection Analysis Based on the Three-dimensional Second Strain Gradient Theory. WMVC 2022: 10th International Conference on Wave Mechanics and Vibrations, Jul 2022, Lisbon, Portugal. 10.1007/978-3-031-15758-5_78. http://hal-03811882

HAL Id: hal-03811882

https://hal.science/hal-03811882

Submitted on 12 Oct 2022

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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Wave Transmission and Reflection Analysis Based on the
Three-dimensional Second Strain Gradient Theory

Bo Yang¹, Mohamed Ichchou¹*, Christophe Droz², and Abdelmalek Zine³

¹ LTDS - CNRS UMR 5513, Vibroacoustics & Complex Media Research Group, École Centrale de Lyon,
France

boyang-chn@hotmail.com, mohamed.ichchou@ec-lyon.fr

² Univ. Gustave Eiffel, Inria, COSYS/SII, I4S Team, France

christophe.droz@inria.fr

³ Institute Camille Jordan - CNRS UMR 5208, École Centrale de Lyon, France

abdel-malek.zine@ec-lyon.fr

**Abstract** The scattering of guided waves through a coupling region is a crucial information when studying waveguides. In this paper, the second strain gradient theory (SSG) is used to describe wave transmission and reflection in a three-dimensional micro-sized medium. First, the constitutive relation of 3D SSG model is derived while six quintic Hermite polynomial shape functions are used for the displacement field. Then Hamilton's principle is used for the weak formulation of the unit-cell's stiffness matrix finite element stiffness, mass matrices and force vector. Eventually the wave diffusion (i.e. including reflection and transmission coefficients) are computed and discussed for various coupling conditions.

**Keywords** Second strain gradient theory · Wave finite element method · Wave transmission · Wave reflection.

# 1 Introduction

The dynamical properties of guided waves such as wave transmission and reflection have been widely studied over the past decade especially in the field of acoustics, earthquake and electromagnetic. Initially, studies focused on the guided waves interaction at interfaces between different macro-medias. But for the micro-sized structure with size effects, the micro-particles such as atoms with associated energy on the free surface of the structure has a significant influence on the structure's behavior. This energy related to surface atoms is called surface free energy which produces surface tension. The surface tension can not be ignored due to the very large ratio between the surface and the volume of structure. On the other hand, long-range or non-local interaction between micro-particles has also an indispensable effect on the micro-sized structure's dynamical behaviors [1,2,20]. The wave propagation and diffusion in micro-medias can no longer be reasonably predicted by Classical Theory (CT) of continuum mechanics[3].

Therefore, the non-classical continuum theories of elasticity that can interpret the properties of micro-sized structures have been proposed. Generally, these non-classical theories can be categorized into non-local elasticity theory [4], micro-continuum theory [5,6], surface elasticity theory [7] and strain gradient family [8]. The strain gradient family is composed of the couple stress theory, the first and second strain gradient theories and the modified couple stress theory. Mindlin established one of the strain gradient family called First Strain Gradient (SG) theory [9] in which the constitutive relations is composed of strain and the first gradient of strain. The atomic structure with the nearest and next nearest interactions between different particles is used to describe the SG theory in the framework of lattice spring model, but only in noncentro-symmetric materials [10]. In order to explore the properties of centro-symmetric materials, the Second Strain Gradient (SSG) theory [11] was put forward, which offers a reasonable description of the strain and surface tension properties on the micro-structure's surface by introducing the high-order parameters. The constitutive relations in SSG theory is a function of strain, first gradient of strain and second gradient of strain. The connections between the SSG theory and lattice spring model with the nearest, next nearest and next-next nearest neighbor interactions for 1D structures can be confirmed through the Fourier series transform [3,22].

On the other hand, in order to study the wave diffusion in complex structures, the numerical
methods such as Spectral Finite Mlement (SFE) method [12] and Semi-Analytical Finite Element (SAFE) [13] can be used. In the past decade, the Wave Finte Element Method (WFEM) [14,15,23] has attracted many works. The advantage of WFEM is the convenient application in engineering field. Since it can be developed from the Finte Element Method (FEM) packages which allows the current element library and grid generation procedures to be applied for the modelling of different waveguide structures. In addition, WFEM can reduce a global periodic structure into a single substructure or unit cell based on the periodic structures theory [16,24]. The resulting stiffness and mass matrices are post processed to offer the dynamic stiffness matrix. The dynamical properties of the periodic structure can be reflected through the spectral analysis of the unit cell [17,21]. The main objective of this work is, firstly, to calculate the multi mode propagation in a 3D periodic waveguide by SSG theory, and, secondly, to confirm the wave diffusion under a complex coupling condition.

The article's structure is the following: in section 2, the constitutive relations of 3D micro-sized model are introduced in the SSG theory framework and the weak formulations including element stiffness, mass matrices and force vector are calculated. Afterwards, in section 3, free wave propagation characteristics are expressed by solving eigenvalue problems in the direct WFEM framework, diffusion matrix for a complex coupling condition are confirmed. In section 4, wave diffusion is introduced. Ultimately, some useful conclusions are presented in section 5.

## 2 A brief of second strain gradient elasticity

### 2.1 3D constitutive relations calculation

The strain energy density $\bar{U}$ composed of strain $\varepsilon=sym(\nabla W)$, first gradient of strain $\xi=\nabla\nabla W$ and second gradient of strain $\zeta=\nabla\nabla\nabla W$ in the SSG theory framework was put forward by Mindlin [11], as below:

$$ \bar{U} = \frac{1}{2} \lambda \varepsilon_{ii} \varepsilon_{jj} + \mu \varepsilon_{ij} \varepsilon_{ij} + a_1 \xi_{ijj} \xi_{ikk} + a_2 \xi_{iik} \xi_{kjj} + a_3 \xi_{iik} \xi_{jjk} + a_4 \xi_{ijk} \xi_{ijk} + a_5 \xi_{ijk} \xi_{jki} + b_1 \zeta_{ijjj} \zeta_{kkll} \\ + b_2 \zeta_{ijkl} \zeta_{ijll} + b_3 \zeta_{ijkl} \zeta_{jkll} + b_4 \zeta_{ijkl} \zeta_{llkj} + b_5 \zeta_{ijkl} \zeta_{lljk} + b_6 \zeta_{ijkl} \zeta_{ijkl} + b_7 \zeta_{ijkl} \zeta_{jkli} + c_1 \varepsilon_{ii} \zeta_{jjkk} \\ + c_2 \varepsilon_{ij} \zeta_{ijkl} + c_3 \varepsilon_{ij} \zeta_{kkij}, \quad (1) $$

where $W$ is the displacement vector, symbol $\nabla$ means the gradient operator, $\lambda$ and $\mu$ represent the Lamé parameters which are related to the Young's modulus $E$, the Poisson's ratio $\nu$ and the shear modulus $G$ [3]. $a_i$, $b_i$ and $c_i$ denote the higher order parameters [10] in SSG theory. The higher order parameters for Al and Cu are shown in Table 1. Based on the 3D elasticity theory, the vector of

<table><caption>**Table 1.** Higher order material parameters $a_i$ (eV/Å), $b_i$ (eV·Å), $c_i$ (eV/Å).</caption><thead><tr><th>Material</th><th>a<sub>1</sub></th><th>a<sub>2</sub></th><th>a<sub>3</sub></th><th>a<sub>4</sub></th><th>a<sub>5</sub></th><th>b<sub>1</sub></th><th>b<sub>2</sub></th><th>b<sub>3</sub></th><th>b<sub>4</sub></th><th>b<sub>5</sub></th><th>b<sub>6</sub></th><th>b<sub>7</sub></th><th>c<sub>1</sub></th><th>c<sub>2</sub></th><th>c<sub>3</sub></th></tr></thead><tbody><tr><td>Al</td><td>0.140</td><td>0.002</td><td>-0.008</td><td>0.096</td><td>0.258</td><td>0.792</td><td>0.064</td><td>-0.194</td><td>-0.001</td><td>0.001</td><td>16.156</td><td>48.529</td><td>0.504</td><td>0.357</td><td>0.178</td></tr><tr><td>Cu</td><td>0.183</td><td>0.010</td><td>0.001</td><td>0.072</td><td>0.189</td><td>0.661</td><td>0.066</td><td>0.206</td><td>0.002</td><td>0.002</td><td>12.625</td><td>37.940</td><td>0.845</td><td>0.573</td><td>0.347</td></tr></tbody></table>

displacement field defined in the Cartesian coordinate system $(x, y, z)$ is given as:

$$ \mathbf{W}(x, y, z, t) = (u(x, y, z, t), v(x, y, z, t), w(x, y, z, t))^{\mathrm{T}}, \quad (2) $$

where $u$, $v$ and $w$ are the displacements along $x$, $y$ and $z$ direction.

The relations between strains and displacement components can be defined by introducing the vectors of first, second and third order derivatives of displacement components:

$$ \varepsilon = \mathbf{Q}_1 \mathbf{W}, \quad \xi = \mathbf{Q}_2 \mathbf{W}, \quad \zeta = \mathbf{Q}_3 \mathbf{W}, \quad (3) $$

where

$$ \mathbf{Q}_1 = \begin{bmatrix} \partial_1 & 0 & 0 & 0 & \partial_3 & 0 \\ 0 & \partial_2 & 0 & \partial_3 & 0 & \partial_2 \\ 0 & 0 & \partial_3 & \partial_2 & \partial_1 & \partial_1 \end{bmatrix}^{\mathrm{T}}, \quad \mathbf{Q}_2 = \begin{bmatrix} \mathbf{t}_1 & 0 & 0 \\ 0 & \mathbf{t}_1 & 0 \\ 0 & 0 & \mathbf{t}_1 \end{bmatrix} \otimes (\partial_{11} \partial_{22} \partial_{33} 2\partial_{12} 2\partial_{13} 2\partial_{23})^{\mathrm{T}}, \quad \mathbf{Q}_3 = \begin{bmatrix} \mathbf{t}_2 & 0 & 0 \\ 0 & \mathbf{t}_2 & 0 \\ 0 & 0 & \mathbf{t}_2 \end{bmatrix} \otimes $$
($\partial_{111} \partial_{222} \partial_{333} 3\partial_{112} 3\partial_{113} 3\partial_{221} 3\partial_{223} 3\partial_{331} 3\partial_{332} 6\partial_{123}$)ᵀ, symbol ⊗ stands for the Kronecker product, **t**₁ with size 6 × 1 and **t**₂ with size 10 × 1 are the matrices whose element value is 1.

Then, the constitutive relations for 3D model by SSG theory can be defined as:

$$ \tau_1 = \mathcal{L}\varepsilon + \mathcal{C}\zeta, \quad \tau_2 = \mathcal{A}\xi, \quad \tau_3 = \mathcal{B}\zeta + \mathcal{C}^T\varepsilon, \qquad (4) $$

in which $\mathcal{L}$ is matrix including classical parameters, $\mathcal{A}$, $\mathcal{B}$ and $\mathcal{C}$ are matrices including higher order parameters [18]. Finally, the strain energy density for SSG theory can be rewritten as the matrix form:

$$ \bar{u} = \frac{1}{2} (\varepsilon^T \mathcal{L} \varepsilon + \xi^T \mathcal{A} \xi + \zeta^T \mathcal{B} \zeta + 2 \zeta^T \mathcal{C}^T \varepsilon). \qquad (5) $$

Here should be noted that, Eq. 5 is the basic form of building a 3D model using partial differential equations (PDE) weak form in some commercial numerical simulation software (e.g., COMSOL).

## 2.2 Finite element discretization

After obtaining the strain energy density, the next step is to calculate the weak form including stiffness and mass matrices and force vector. The definition of node degree of freedoms (DOFs) for 1D and 3D Hermite elements. Firstly, in order to ensure the continuity of higher derivatives up to the second order between 1D elements, the six-term polynomial function is considered to interpolate the scalar field $W_1=u(x, t)$ inside a 1D element, as follows:

$$ W_1 = \left[ 1 \ x \ x^2 \ x^3 \ x^4 \ x^5 \right] \left[ s_0 \ s_1 \ s_2 \ s_3 \ s_4 \ s_5 \right]^T = \mathbf{x}\mathbf{s}. \qquad (6) $$

The evaluation of the nodal displacement vector $\bar{\mathbf{w}}_1^{(e)}$ of 1D element gives:

$$ \bar{\mathbf{w}}_1^{(e)} = (\mathbf{d}_1 \ \mathbf{d}_2 \ \mathbf{d}_3 \ \mathbf{d}_4 \ \mathbf{d}_5 \ \mathbf{d}_6)^T \mathbf{s} = \mathbf{d}\mathbf{s}, \qquad (7) $$

where $\mathbf{d}_1=(1 - l_e \ l_e^2 - l_e^3 \ l_e^4 - l_e^5)$, $\mathbf{d}_2=(0 \ 1 - 2l_e \ 3l_e^2 - 4l_e^3 \ 5l_e^4)$, $\mathbf{d}_3=(0 \ 0 \ 2 - 6l_e \ 12l_e^2 - 20l_e^3)$, $\mathbf{d}_4=(1 \ l_e^2 \ l_e^3 \ l_e^4 \ l_e^5)$, $\mathbf{d}_5=(0 \ 1 \ 2l_e \ 3l_e^2 \ 4l_e^3 \ 5l_e^4)$, $\mathbf{d}_6=(0 \ 0 \ 2 \ 6l_e \ 12l_e^2 \ 20l_e^3)$. Then, submitting Eq.7 into Eq.6, the displacement vector within the 1D element can be derived by employing the six quin-tic Hermite polynomial shape function and nodal displacement vector, as follows:

$$ W_1 = \mathbf{x}\mathbf{d}^{-1}\bar{\mathbf{w}}_1^{(e)} = \mathbf{S}(x)\bar{\mathbf{w}}_1^{(e)}, \qquad (8) $$

in which the shape function $\mathbf{S}(x)$ is written as:

$$ \mathbf{S}(x) = (S_1^0(x), S_1^1(x), S_1^2(x), S_2^0(x), S_2^1(x), S_2^2(x)), \qquad (9) $$

where $S_1^0(x)=\frac{5x^3}{8l_e^3}-\frac{15x}{16l_e}-\frac{3x^5}{16l_e^5}+\frac{1}{2}$, $S_1^1(x)=\frac{5l_e}{16}-\frac{7x}{16}-\frac{3x^2}{8l_e}+\frac{5x^3}{8l_e^2}+\frac{x^4}{16l_e^3}-\frac{3x^5}{16l_e^4}$, $S_1^2(x)=\frac{l_e^2}{16}-\frac{l_e x}{16}-\frac{x^2}{8}+\frac{x^3}{8l_e}+\frac{x^4}{16l_e^2}-\frac{x^5}{16l_e^3}$, $S_2^0(x)=\frac{15x}{16l_e}-\frac{5x^3}{8l_e^3}+\frac{3x^5}{16l_e^5}+\frac{1}{2}$, $S_2^1(x)=\frac{3x^2}{8l_e}-\frac{7x}{16}-\frac{5l_e}{16}+\frac{5x^3}{8l_e^2}-\frac{x^4}{16l_e^3}-\frac{3x^5}{16l_e^4}$, $S_2^2(x)=\frac{l_e x}{16}+\frac{l_e^2}{16}-\frac{x^2}{8}-\frac{x^3}{8l_e}+\frac{x^4}{16l_e^2}+\frac{x^5}{16l_e^3}$. $\mathbf{S}(y) = \mathbf{S}(x)|_{x=y} = [N_1^0(y), S_1^1(y), S_1^2(y), S_2^0(y), S_2^1(y), S_2^2(y)]$, $\mathbf{S}(z) = \mathbf{S}(x)|_{x=z} = [S_1^0(z), S_1^1(z), S_1^2(z), S_2^0(z), S_2^1(z), S_2^2(z)]$. The shape function of hexahedral element:

$$ \mathcal{S}(x, y, z) = (\mathbf{S}_1(x, y, z) \otimes \mathbf{e}_1, \mathbf{S}_2(x, y, z) \otimes \mathbf{e}_2, \mathbf{S}_3(x, y, z) \otimes \mathbf{e}_3)^T. \qquad (10) $$

The element in $\mathbf{S}_p(x, y, z)$ and $\mathbf{e}_p$ ($p = 1, 2, 3$) are defined as:

$$ S_p^{i(j,k,l)}(x, y, z) = S_{i'}^j(x)S_{i''}^k(y)S_{i'''}^l(z), \quad \mathbf{e}_p = (\epsilon_{p1} \ \epsilon_{p2} \ \epsilon_{p3}), \qquad (11) $$

where $i = 1, ..., 8$. $j, k, l = 0, 1, 2$. $S_p^{i(j,k,l)}(x, y, z)$ is associated with the DOFs $\partial^{j+k+l}u_1^i/(\partial x^j \partial y^k \partial z^l)$ of node $i$ of the hexahedron element. $i', i'', i''' = 1, 2$ relate to the node number in the corresponding 1D element and they take values of 1 or 2 if the coordinate value of node $i$ is $-l_e$ or $l_e$. The displacement vector $\mathbf{W}(x, y, z)$ within 3D element can be expressed as:

$$ \mathbf{W}(x, y, z, t) = \mathcal{S}(x, y, z)\mathbf{w}^{(e)}(t), \qquad (12) $$
where $\mathbf{w}^{(e)}=[(\mathbf{w}_1^{(e)})^T, (\mathbf{w}_2^{(e)})^T, (\mathbf{w}_3^{(e)})^T]^T$, $(\mathbf{w}_p^{(e)})^T=[(\mathbf{w}_p^{1(e)})^T, (\mathbf{w}_p^{2(e)})^T, ..., (\mathbf{w}_p^{8(e)})^T]^T$ ($p=1,2,3$). Then, integrating the strain energy density over its volume to obtain the strain potential energy $\mathcal{U}$ as:

$$ \mathcal{U} = \int_V \bar{\mathcal{U}} dV = \frac{1}{2} \mathbf{w}^{(e)T} \int_V (\mathbb{S}^T \mathbf{Q}_1^T \mathcal{L} \mathbf{Q}_1 \mathbb{S} + \mathbb{S}^T \mathbf{Q}_2^T \mathcal{A} \mathbf{Q}_2 \mathbb{S} + \mathbb{S}^T \mathbf{Q}_3^T \mathcal{B} \mathbf{Q}_3 \mathbb{S} + 2 \mathbb{S}^T \mathbf{Q}_3^T \mathcal{C} \mathbf{Q}_1 \mathbb{S}) dV \mathbf{w}^{(e)}. \quad (13) $$

On the other hand, the beam kinetic energy is expressed as:

$$ \mathcal{T} = \frac{1}{2} \int_V \left(\frac{\partial \mathbf{W}}{\partial t}\right)^T \rho \left(\frac{\partial \mathbf{W}}{\partial t}\right) dV = \frac{1}{2} \left(\frac{\partial \mathbf{w}^{(e)}}{\partial t}\right)^T \int_V (\mathbb{S}^T \rho \mathbb{S}) dV \left(\frac{\partial \mathbf{w}^{(e)}}{\partial t}\right), \quad (14) $$

where $\rho$ denotes the linear mass density, $l_1$ and $l_2$ are higher-order length-scale parameters. Meanwhile, the work done $\delta W$ by external force can be expressed as:

$$ \delta W = \int_V \delta \mathbf{W}^T \mathbf{f}_V dV + \int_S \delta \mathbf{W}^T \mathbf{f}_S dS = \delta \mathbf{w}^{(e)T} \left( \int_V (\mathbb{S}^T \mathbf{f}_V) dV + \int_S (\mathbb{S}^T \mathbf{f}_S) dS \right), \quad (15) $$

where $\mathbf{f}_V$ is the volume force, $\mathbf{f}_S$ means the face force. The 3D element stiffness, mass matrices and force vector can be confirmed as:

$$ \mathbf{K}^{(e)} = \int_V (\mathbb{S}^T \mathbf{Q}_1^T \mathcal{L} \mathbf{Q}_1 \mathbb{S} + \mathbb{S}^T \mathbf{Q}_2^T \mathcal{A} \mathbf{Q}_2 \mathbb{S} + \mathbb{S}^T \mathbf{Q}_3^T \mathcal{B} \mathbf{Q}_3 \mathbb{S} + 2 \mathbb{S}^T \mathbf{Q}_3^T \mathcal{C} \mathbf{Q}_1 \mathbb{S}) dV, \quad (16) $$

$$ \mathbf{M}^{(e)} = \int_V (\mathbb{S}^T \rho \mathbb{S}) dV, \quad \mathbf{F}^{(e)} = \int_V (\mathbb{S}^T \mathbf{f}_V) dV + \int_S (\mathbb{S}^T \mathbf{f}_S) dS. $$

### 3 Diffusion analysis through wave finite element method

<figure><img src="image_3.png" /><figcaption>Fig. 1. Two coupled waveguides through a coupling element ($\bar{x}$ is the local coordinate).</figcaption></figure>

This section is concerned with the characterization of coupling conditions between two semi-infinite periodic waveguides which are connected through an elastic coupling element at surfaces 1 and 2 as shown in Fig.1. The reflection coefficients ($\mathbf{R}$) and transmission coefficients ($\mathbf{T}$) are confirmed through a complex coupling condition: the wave modes calculation in waveguides is based on the CT but SSG theory in coupling element (CT-SSG). The propagation constants $\mathbf{\Lambda}$ and eigenvectors $\mathbf{\Psi}_u$ in a unit cell can be solved by direct Bloch formulation [17,3,15] as:

$$ [\mathbf{D}_{RL}(\omega)\mathbf{\Lambda}^{-1} + (\mathbf{D}_{RR}(\omega) + (\mathbf{D}_{LL}(\omega)) + \mathbf{D}_{LR}(\omega)\mathbf{\Lambda})] \mathbf{\Psi}_u = 0, \quad (17) $$

where $\mathbf{\Lambda}=\text{diag}\{\lambda_j\}_{j=1,...,2p}$, $\mathbf{\Psi}_u = \{\phi_j\}_{j=1,...,2p}$ which can be divided into $\mathbf{\Psi}_u^+ = \{\phi_j^+\}_{j=1,...,p}$ and $\mathbf{\Psi}_u^- = \{\phi_j^-\}_{j=p+1,...,2p}$, in which $p=m$ for SSG, $p=n$ for CT. The waves propagate to positive when $|1/\lambda_j|<1$. The waves propagate to negative when $|\lambda_j|>1$. Here, $\lambda_j$ take the form $\lambda_j=\exp(-i\kappa_j L(q))$, subscript $q=c$ for coupling element, $q=k, k+1$ for unit cells. It should be noted that the spectral analysis method of the coupling element and unit cell $k+1$ is the same as that of unit cell $k$. The displacement field of coupling element ($c$) and unit cells ($k, k+1$) can be represented by the superposition of the eigenmodes:

$$ \hat{\mathbf{w}}^{(q)\pm} = \mathbf{\Psi}_{u}^{(q)\pm} \mathbf{v}^{(q)\pm} \mathbf{P}^{(q)\pm}, \quad (18) $$
where $\mathbf{P}^{(q)\pm}$ is the amplitudes of wave modes. $\mathbf{v}^{(q)\pm}=\text{diag}\{\exp(\mp i\kappa_j^{(q)}\bar{x})\}$ ($0 \le \bar{x} \le L_{(q)}$), symbol $q=c$ for coupling element, $q=k, k+1$ for unit cells. On the other hand, the force components of coupling element $(c)$ and unit cells $(k, k+1)$ can be expressed as:

$$ \hat{\mathbf{F}}^{(q)\pm} = \mathbf{\Psi}_F^{(q)\pm} \mathbf{v}^{(q)\pm} \mathbf{P}^{(q)\pm}, \quad (19) $$

in which $\mathbf{\Psi}_F^{(q)\pm}=\mathbf{D}_{LL}^{(q)}\mathbf{\Psi}_u^{(q)\pm}+\mathbf{D}_{LR}^{(q)}\mathbf{\Psi}_u^{(q)\pm}\text{diag}\{\exp(\mp i\kappa_j^{(q)}L_{(q)})\}$, symbol $q=c$ for coupling element, $q=k, k+1$ for unit cells. The eigensolutions will be the same when waveguides 1 and 2 have the same cross-section and material. The state vector on the right side of surface 1, represented by $R_1$:

$$ \begin{pmatrix} \hat{\mathbf{w}}_{R_1}^{(c)} \\ \hat{\mathbf{F}}_{R_1}^{(c)} \end{pmatrix} = \begin{bmatrix} \mathbf{\Psi}_u^{(c)+} \mathbf{P}^{(c)+} + \mathbf{\Psi}_u^{(c)-} \mathbf{P}^{(c)-} \\ \mathbf{\Psi}_F^{(c)+} \mathbf{P}^{(c)+} + \mathbf{\Psi}_F^{(c)-} \mathbf{P}^{(c)-} \end{bmatrix}. \quad (20) $$

In addition, the state vector on the left side of surface 2, represented by $L_2$, is written as:

$$ \begin{pmatrix} \hat{\mathbf{w}}_{L_2}^{(c)} \\ \hat{\mathbf{F}}_{L_2}^{(c)} \end{pmatrix} = \begin{bmatrix} \mathbf{\Psi}_u^{(c)+} \mathbf{v}^{(c)+}|_{\bar{x}=L_c} \mathbf{P}^{(c)+} + \mathbf{\Psi}_u^{(c)-} \mathbf{v}^{(c)-}|_{\bar{x}=L_c} \mathbf{P}^{(c)-} \\ \mathbf{\Psi}_F^{(c)+} \mathbf{v}^{(c)+}|_{\bar{x}=L_c} \mathbf{P}^{(c)+} + \mathbf{\Psi}_F^{(c)-} \mathbf{v}^{(c)-}|_{\bar{x}=L_c} \mathbf{P}^{(c)-} \end{bmatrix}. \quad (21) $$

Combining Eq.20 and Eq.21, the relation between state vector on left and right side of coupling element will be conformed:

$$ (\hat{\mathbf{w}}_{L_2}^{(c)}, \hat{\mathbf{F}}_{L_2}^{(c)})^T = \mathbf{X}^{(c)} (\hat{\mathbf{w}}_{R_1}^{(c)}, \hat{\mathbf{F}}_{R_1}^{(c)})^T, \quad (22) $$

with

$$ \mathbf{X}^{(c)} = \begin{bmatrix} \mathbf{\Psi}_u^{(c)+} \mathbf{v}^{(c)+}|_{\bar{x}=L_c}, \mathbf{\Psi}_u^{(c)-} \mathbf{v}^{(c)-}|_{\bar{x}=L_c} \\ \mathbf{\Psi}_F^{(c)+} \mathbf{v}^{(c)+}|_{\bar{x}=L_c}, \mathbf{\Psi}_F^{(c)-} \mathbf{v}^{(c)-}|_{\bar{x}=L_c} \end{bmatrix} \begin{bmatrix} \mathbf{\Psi}_u^{(c)+} \mathbf{\Psi}_u^{(c)-} \\ \mathbf{\Psi}_F^{(c)+} \mathbf{\Psi}_F^{(c)-} \end{bmatrix}^{-1}. \quad (23) $$

Next, assuming that the incident waves come from the infinity of waveguide 1 and there is no reflection from the end of waveguide 2. The initial boundary of waveguide 1 is also non-reflecting. The state vector on the left side of surface 1, represented by $L_1$, is expressed as:

$$ \begin{pmatrix} \hat{\mathbf{w}}_{L_1}^{(k)} \\ \hat{\mathbf{F}}_{L_1}^{(k)} \end{pmatrix} = \begin{bmatrix} \mathbf{\Psi}_u^{(k)+} \mathbf{v}^{(k)+}|_{\bar{x}=L_k} \mathbf{P}^{(k)+} + \mathbf{\Psi}_u^{(k)-} \mathbf{v}^{(k)-}|_{\bar{x}=L_k} \mathbf{P}^{(k)-} \\ \mathbf{\Psi}_F^{(k)+} \mathbf{v}^{(k)+}|_{\bar{x}=L_k} \mathbf{P}^{(k)+} + \mathbf{\Psi}_F^{(k)-} \mathbf{v}^{(k)-}|_{\bar{x}=L_k} \mathbf{P}^{(k)-} \end{bmatrix}. \quad (24) $$

The state vector on the right side of surface 2, represented by $R_2$, is written as:

$$ \begin{pmatrix} \hat{\mathbf{w}}_{R_2}^{(k+1)} \\ \hat{\mathbf{F}}_{R_2}^{(k+1)} \end{pmatrix} = \begin{bmatrix} \mathbf{\Psi}_u^{(k+1)+} \mathbf{v}^{(k+1)+}|_{\bar{x}=0} \mathbf{P}^{(k+1)+} \\ \mathbf{\Psi}_F^{(k+1)+} \mathbf{v}^{(k+1)+}|_{\bar{x}=0} \mathbf{P}^{(k+1)+} \end{bmatrix}. \quad (25) $$

Here should be noted that the size of state vector for coupling element is $2m \times 1$, but $2n \times 1$ for unit cells $k$ and $k+1$. The higher order parts in state vector for coupling element is $2(m-n) \times 1$. In order to ensure the continuity on surfaces 1 and 2, defining new state vectors including higher order parts $(2(m-n) \times 1)$ for unit cells $k$ and $k+1$:

$$ \begin{gathered} (\hat{\mathbf{w}}_{L_1}^{*(k)}, \hat{\mathbf{F}}_{L_1}^{*(k)})^T = (\hat{\mathbf{w}}_{L_1}^{(k)}, \hat{\mathbf{w}}_{L_1}^{'(k)}, \hat{\mathbf{F}}_{L_1}^{(k)}, \hat{\mathbf{F}}_{L_1}^{'(k)})^T, \\ (\hat{\mathbf{w}}_{R_2}^{*(k+1)}, \hat{\mathbf{F}}_{R_2}^{*(k+1)})^T = (\hat{\mathbf{w}}_{R_2}^{(k+1)}, \hat{\mathbf{w}}_{R_2}^{'(k+1)}, \hat{\mathbf{F}}_{R_2}^{(k+1)}, \hat{\mathbf{F}}_{R_2}^{'(k+1)})^T, \end{gathered} \quad (26) $$

where $\hat{\mathbf{w}}_{L_1}^{(k)}, \hat{\mathbf{F}}_{L_1}^{(k)}$ and $\hat{\mathbf{w}}_{L_1}^{(k+1)}, \hat{\mathbf{F}}_{L_1}^{(k+1)}$ are unknown higher order displacements and forces vectors for unit cell $k$ and $k+1$ respectively. The continuity on surfaces 1 and 2 is:

$$ (\hat{\mathbf{w}}_{R_1}^{(c)}, \hat{\mathbf{F}}_{R_1}^{(c)})^T = (\hat{\mathbf{w}}_{L_1}^{*(k)}, \hat{\mathbf{F}}_{L_1}^{*(k)})^T, \quad (\hat{\mathbf{w}}_{L_2}^{(c)}, \hat{\mathbf{F}}_{L_2}^{(c)})^T = (\hat{\mathbf{w}}_{R_2}^{*(k+1)}, \hat{\mathbf{F}}_{R_2}^{*(k+1)})^T. \quad (27) $$

Combining Eq.22, Eq.26 and Eq.27, assume that higher order parts for $k$ and $k+1$ are 0. Define $\mathbf{P}^{(k)+}=\mathbf{I}$, $\mathbf{P}^{(k)-}=\mathbf{R}$, $\mathbf{P}^{(k+1)+}=\mathbf{T}$, the $\mathbf{R}$ and $\mathbf{T}$ can be conformed as:

$$ \mathbf{G} = \mathbf{X}^{(c)}\mathbf{H}, \quad (28) $$

where $\mathbf{G}=[\mathbf{\Psi}_u^{(k+1)+}\mathbf{v}^{(k+1)+}|_{\bar{x}=0}\mathbf{T}, \hat{\mathbf{w}}_{R_2}^{(k+1)}, \mathbf{\Psi}_F^{(k+1)+}\mathbf{v}^{(k+1)+}|_{\bar{x}=0}\mathbf{T}, \mathbf{0}]^T$, $\mathbf{H}=[\mathbf{\Psi}_u^{(k)+}\mathbf{v}^{(k)+}|_{\bar{x}=L_k}\mathbf{I} + \mathbf{\Psi}_u^{(k)-}\mathbf{v}^{(k)-}|_{\bar{x}=L_k}\mathbf{R}, \hat{\mathbf{w}}_{L_1}^{(k)}, \mathbf{\Psi}_F^{(k)+}\mathbf{v}^{(k)+}|_{\bar{x}=L_k}\mathbf{I} + \mathbf{\Psi}_F^{(k)-}\mathbf{v}^{(k)-}|_{\bar{x}=L_k}\mathbf{R}, \mathbf{0}]^T$.
## 4 Numerical applications

<figure><img src="image_4.png" /><figcaption>Fig. 2. Finite element model of two waveguides coupled by a coupling element.</figcaption></figure>

In this part, the WFEM is applied to analyze the multi-mode diffusion. The unit cells $k$ and $k+1$ with $L_{k/k+1} = 50a_0$, $L_y = 300a_0$ and $L_z = 300a_0$ ($a_0$ is the lattice parameter) as presented in Fig.2. The coupling element with $L_c = 100a_0$ as shown in Fig.2. Materials Aluminum (Al) and Copper (Cu) are used here. The Young's modulus $E$ is 70 GPa for Al and 110 GPa for Cu, linear mass density $\rho$ is $2.7 \text{ g/cm}^3$ for Al and $8.96 \text{ g/cm}^3$ for Cu. The damping loss factor $\eta=1e^{-4}$. Unit cells $k$ and $k+1$ are meshed into 16 3D elements, coupling element is meshed 64 3D elements.

In order to study the wave diffusion under a complex coupling condition, as shown in Fig.2, first

<figure><img src="image_5.png" /><figcaption>Fig. 3. Absolute values of diffusion coefficients (the materials of coupling element, unit cell $k$ and $k+1$ are Al).</figcaption></figure>

of all, defining the materials of coupling element, unit cell $k$ and $k+1$ are Al, the diffusion model of coupling element is built by SSG theory with higher-order parameters, the diffusion models of unit cell $k$ and $k+1$ are built by CT. The **R** and **T** coefficients can be calculated from Eq. 28 including bending, tension and torsion modes. As shown in Fig.3, the black lines denote the WFEM results. The value
of **R** representing the non-classical part of reflection is no longer 0, and the value of **T** representing the non-classical part of transmission is no longer 1. The influence of non-local interactions caused by higher-order parameters can be reflected by this model. On the other hand, an analytical method [19], shown by the red lines, is used to valid the WFEM results for bending and tension. As we can see, for the reflection coefficient, the result obtained by WFEM is very close to the one by the analytical method at low frequency, but the results are different at high frequencies. For the transmission coefficient, the results by WFEM matches the results by analytical method well.

<figure><img src="image_6.png" /><figcaption>**Fig. 4.** Absolute values of diffusion coefficients (the material of coupling element is Cu, unit cell *k* and *k*+1 are Al).</figcaption></figure>

In addition, the joint influence of classical parameters (i. g., Young's modulus, Poisson's ratio and mass density) and higher-order parameters of material on diffusion is also a very meaningful study. Defining the material of coupling element is Cu, unit cell *k* and *k*+1 are Al, the diffusion model of coupling element is built by SSG theory with higher-order parameters, the diffusion models of unit cell *k* and *k*+1 are built by CT. The **R** and **T** coefficients including bending, tension and torsion modes can be illustrated by Eq. 28 as well. As shown in Fig.4, the black lines denote the WFEM results and the red lines represent the results from an analytical method [19]. The diffusion is different from the case presented in Fig.3, this shows that the impedance mismatch is not only due to the non-local interactions caused by higher order parameters in the SSG theory model but also the local interactions caused by classical parameters.

**5 Conclusions**

In this paper, SSG theory is used for the multi-mode diffusion analysis within the WFEM framework. The diffusion is confirmed through two different cases. For the first case, the value of **R** representing the non-classical part of reflection is no longer 0, and the value of **T** representing the non-classical part of transmission is no longer 1. The influence of non-local interactions caused by higher-order parameters can be reflected by this model. The second case shows that the impedance mismatch is not only due to the non-local interactions caused by higher order parameters in the SSG theory model but also the local interactions caused by classical parameters such as Young's modulus, Poisson's ratio and mass density.
Acknowledgements

This work is supported by LabEx CeLyA (Centre Lyonnais d'Acoustique, ANR-10-LABX-0060) of Université de Lyon. The research of B. Yang is funded by the China Scholarship Council (CSC).

References

1. Miller, R.E., Shenoy, V.B.: Size-dependent elastic properties of nano-sized structural elements. Nanotechnology 11(3), 139 (2000)

2. Lim, C.W., He, L.H.: Size-dependent nonlinear response of thin elastic films with nano-scale thickness. International Journal of Solids and Structures 46(11), 15–26 (2004)

3. Yang, B., Droz, C., Zine, A.M., Ichchou, M.N.: Dynamic analysis of second strain gradient elasticity through a wave finite element approach. Composite Structures (2020). https://doi.org/10.1016/j.compstruct.2020.11342

4. Kroner, E.: Elasticity theory of materials with long range cohesive forces. International Journal of Solids and Structures 3, 731–742 (1967)

5. Eringen, A.: Simple microfluids. International Journal of Engineering Science 2, 205–217 (1964)

6. Eringen, A.: Linear theory of micropolar elasticity. Journal of Applied Mathematics and Mechanics 15, 909–923 (1966)

7. Gurtin, M.E., Murdoch, A.I.: A continuum theory of elastic material surfaces. Archive for rational mechanics and analysis 57(4), 291–323 (1975)

8. Lazar, M., Maugin, G.A., Aifantis, E.C.: Dislocation in second strain gradient elasticity. International Journal of Solids and Structures pp. 1787–1817 (2006)

9. Mindlin, R.D.: Micro-structure in linear elasticity. Archive for Rational Mechanics and Analysis 16, 51–78 (1964)

10. Shodja, H.M., Ahmadpoor, F., A., T.: Calculation of the additional constants for fcc materials in second strain gradient elasticity: behavior of a nano-size bernoulli-euler beam with surface effects. Applied Mechanics 72(2), 021008 (2010)

11. Mindlin, R.D.: Second gradient of strain and surface tension in linear elasticity. International Journal of Solids and Structures pp. 147–438 (1965)

12. Mahapatra, D., Gopalakrishnan, S.: A spectral finite element for analysis of wave propagation in uniform composite tubes. Journal of Sound and Vibration 268, 429–463 (2003)

13. innveden, S., Fraggstedt, M.: Waveguide finite elements for curved structures. Journal of Sound and Vibration 312, 644–671 (2008)

14. Mencik, J.M., Ichchou, M.N.: Wave finite elements in guided elastodynamics with internal fluid. International Journal of Solids and Structures 44, 2148–2167 (2007)

15. Droz, C., Lainé, J.P., Ichchou, M.N., Inquiété, G.: A reduced formulation for the free-wave propagation analysis in composite structures. Composite Structures 113, 134–144 (2014)

16. Mead, D.: A general theory of harmonic wave propagation in linear periodic systems with multiple coupling. Journal of Sound and Vibration 27, 235–260 (1973)

17. Duhamel, D., Mace, B.R., Brennan, M.: Finite element analysis of the vibrations of waveguides and periodic structures. Journal of Sound and Vibration 294, 205–220 (2006)

18. Torabi, J., Ansari, R., Bazdid-Vahdati, M., Darvizeh, M.: Second strain gradient finite element analysis of vibratory nanostructures based on the three-dimensional elasticity theory. Iranian Journal of Science and Technology, Transactions of Mechanical Engineering 44(3), 631–645 (2020)

19. Zhu, G., Zine, A., Droz, C., Ichchou, M.: Wave transmission and reflection analysis through complex media based on the second strain gradient theory. European Journal of Mechanics/A Solids (2021). https://doi.org/10.1016/j.euromechsol.2021.104326

20. Rosi, G., Placidi, L., Auffray, N.: On the validity range of strain-gradient elasticity: a mixed static-dynamic identification procedure. European Journal of Mechanics-A/Solids 69, 179–191 (2018)

21. Ahsani, S., Boukadia, R., Droz, C., Claeys, C., Deckers, E., Desmet, W.: Diffusion based homogenization method for 1d wave propagation. Mechanical Systems and Signal Processing (2020). https://doi.org/10.1016/j.ymssp.2019.106515

22. Khakalo, S., Niiranen, J.: Form ii of mindlin's second strain gradient theory of elasticity with a simplification: For materials and structures from nano- to macro-scales. European Journal of Mechanics / A Solids 71, 292–319 (2018)

23. Zhao, J. Q., Zeng, P., Pan, B.: Improved hermite finite element smoothing method for full-field strain measurement over arbitrary region of interest in digital image correlation. Optics and Lasers in Engineering 50, 1662–1671 (2012)

24. Singh, R., Droz, C., Ichchou, M., Franco, Bareille, F., De Rosa, O., S.: Stochastic wave finite element quadratic formulation for periodic media: 1d and 2d. Mechanical Systems and Signal Processing 136, 106431 (2020)