Modeling Bloch Waves in Prestressed Phononic Crystal Plates

1. Introduction

The propagation of elastic waves in periodic structures is governed by the elastomorphic and material parameters of its unit cell (Kushwaha and Halevi, 1994; Sainidou et al., 2005). Properly conceiving their design in terms of size, shape, and arrangement, as well as choosing their density and elastic properties demonstrated great potential for attaining exceptional dynamic behavior, such as frequency bandgaps (Kushwaha and Halevi, 1994; Mártinez-Sala et al., 1995; Liu et al., 2000; Miniaci et al., 2018a), negative refraction (Morvan et al., 2010; Zhu et al., 2014; Zhu and Semperlotti, 2016), topological protection (Mousavi et al., 2015; Süsstrunk and Huber, 2015; Pal et al., 2016; Miniaci et al., 2018b), etc. This opened up new perspectives in many fields, ranging from microelectromechanical systems to nondestructive evaluation (Pennec et al., 2010; Craster and Guenneau, 2012; Deymier, 2013), including but not limited to wave filters, waveguiding and beam splitting, sensing devices, wave splitters, vibration shielding, subwavelength imaging (Sukhovich et al., 2009; Li et al., 2015; Trainiti et al., 2015; Colombi et al., 2016; Miniaci et al., 2016, 2017; Su et al., 2016; Miniaci et al., 2019).

However, most of the configurations proposed so far are limited by the fact that their unconventional dynamic properties are accessible to specific frequencies only, once the unit cell design is conceived, being its geometrical and mechanical properties constant with respect to time. To overcome this limitation, periodic systems with adaptive elastic properties have received great attention in the recent years. Indeed, introducing an additional degree of freedom, enables the reversible tailoring of the dispersive diagram of the structure, considerably enhancing their functionality and technological applications. To achieve this goal, several approaches have been proposed and explored so far. They include active systems, such as piezoelectric materials (Casadei et al., 2010; Bergamini et al., 2014; Kherraz et al., 2016), temperature (Jim et al., 2009; Airoldi et al., 2011; Cheng et al., 2011; Wu et al., 2018) and magneto-based techniques (Robillard et al., 2009; Matar et al., 2012; Guo and Wei, 2016; Zhang et al., 2017), actuated polymers via electromagnetic waves (Walker et al., 2014, 2017) as well as the application of an external mechanical load. Among them, the latter has shown to be particularly promising to achieve tunable dispersive properties. Bigoni et al. (2008) formulated a theoretical model for an orthotropic, prestressed (compressible) elastic layer vibrating on an elastic half space assuming long-wave asymptotics for the solution. It was found that the influence of the prestress over bandgaps and passbands, so that it could be exploited as a tuning parameter to shift the dispersion curves of the system. A modeling tool for the prediction of controlled bandgap structures responding to flexural vibrations was provided by Gei et al. (2009), also relaxing the hypothesis of perfect periodicity, allowing thus to consider quasi-periodic structures (Gei, 2010). Wang and Bertoldi (2012) showed how mechanical deformations enabled the tuning of the phononic bandgaps in 3D periodic elastomeric structures. The ability of the elastomers to undergo small as well as large strain deformations guaranteed the reversibility and repeatability of the process. Both the linear and nonlinear regimes of elastic deformation were explored, including different geometrical topology for triggering mechanical instability. Galich et al. (2016) showed how local deformations experienced by a composite structure can bring to local stiffening influencing the instability-induced interfaces on elastic wave propagation in finitely deformed layered materials. Finally, Norris and Parnell (2012) tried to exploit the effects of the prestress to theoretically show the feasibility of a hyperelastic cloaking via transformation elasticity.

In this paper, the effects of an applied mechanical load on the dispersive diagram for phononic crystal plates is investigated in the case of elastic deformations, so to have a complete reversibility of the phenomena. The paper is organized as follows: section 2 provides the description of an Updated Lagrangian scheme, in which a representative unit cell of the phononic crystal plate is studied in its static and dynamic deformed configurations. In section 3, the procedure to extract the band diagrams and in particular the energy velocity using the equations described in section 2 is detailed. Finally, section 4 investigates the effect of a reversible compressive and tensile state of prestress on an epoxy plate with cylindrical inclusions.

2. Floquet-Bloch Analysis of a Prestressed Phononic Crystal

2.1. Static Analysis

In what follows, the equations of static and dynamic equilibrium are derived with reference to the generic infinite phononic crystal (PC) shown in Figure 1, where C₀ indicates the static undeformed configuration, C the static configuration resulting from the application of an external load, and C′ the dynamic configuration undergoing a harmonic motion. The so called Updated Lagrangian (UL) scheme is employed to analyze the PC in the three different equilibrium configurations. The UL approach consists of two steps: in the first step, C₀ is used as the reference configuration to calculate the displacement and stress fields relative to the configuration C, while in the second step, C is used as the new reference configuration to evaluate the dynamic equilibrium of the unit cell in the configuration C′.

FIGURE 1

Figure 1. Schematic representation of the Updated Lagrangian (UL) scheme for the infinite phononic crystal (PC) plate.

Due to the periodicity of the system under consideration, in the undeformed configuration C₀ any scalar, vector or tensor function ϕ satisfies the condition ϕ(x₀ + r₀m) = ϕ(x₀), where m ∈ ℤ² and r₀ = [r₀₁, r₀₂] is the matrix of lattice vectors. This allows to restrict the computational domain to a unit cell of domain Ω₀ and boundary ∂Ω₀ (Collet et al., 2011; Mazzotti et al., 2017). Similarly, the unit cell in the configuration C satisfies the periodic condition ϕ(x + rm) = ϕ(x), with r = [r₁, r₂]. For this cell, the mapping between the coordinates ${\text{x}}_0={\left\{x_{01},x_{02},x_{03}\right\}}^{\text{T}}$ of a material point in C₀ and the coordinates $\text{x}={\left\{x_1,x_2,x_3\right\}}^{\text{T}}$ of the same point in C is established by the deformation gradient F(x₀) = ∇₀x = I + ∇₀u₀, where ∇₀ denotes the gradient operator defined with respect to C₀, u₀ is the pseudo-static displacement resulting from the application of external, Ω₀-periodic, volume (f_V) and surface (t_S) loads. Following Zhang and Parnell (2017), the deformation gradient can be decomposed as F(x₀) = F_L(x₀)F_P(x₀), where F_L indicates the affine deformation gradient of the lattice points such that r = F_Lr₀, while F_P denotes a periodic non-affine deformation.

The specific material density of the unit cell in C₀ is ρ₀(x₀), while, assuming a hyperelastic material behavior described by the Murnaghan's model (Murnaghan et al., 1937; Pau and Lanza di Scalea, 2015; Dubuc et al., 2017, 2018), the tensor of tangential elastic moduli with respect to C₀ is expressed by ${\text{D}}_0=4{\partial }^2\text{Ψ}/(\partial \text{C}\partial \text{C})$, where the elastic energy density Ψ has the form

$\begin{array}{l}\text{Ψ}=\frac{1}{2}(\lambda +2\mu )I_1^2(\text{E})-2\mu I_2(\text{E})+\frac{1}{3}(l+2m)I_1^3(\text{E})\\ -2m{I}_1(\text{E})I_2(\text{E})+n{I}_3(\text{E}),& (1)\end{array}$

in which λ and μ denote the first and second Lamé parameters, respectively, (l, m, n) the third order Murnaghan parameters, C = F^T F the right Cauchy-Green deformation tensor, $\text{E}=\frac{1}{2}({\text{F}}^{\text{T}}\text{F}-\text{I})$ the Green-Lagrange strain tensor and I₁(E), I₂(E), and I₃(E) its first, second and third invariants, respectively.

Finding the coordinates of the deformed configuration C and the associated stress fields requires the solution of a variational functional of the form (Bonet and Wood, 2008)

$\begin{array}{ll}{\displaystyle {\int }_{\Omega }}\text{S}(\text{x}):\delta \text{E}(\text{x})\text{d}\Omega ={\displaystyle {\int }_{\Omega }}{\text{f}}_{V0}(\text{x})·\delta {\text{u}}_0\text{d}\Omega +{\displaystyle {\int }_{\partial \Omega }}{\text{t}}_{S0}(\text{x})·\delta {\text{u}}_0\text{d}(\partial \Omega ),& (2)\end{array}$

subjected to the Dirichlet boundary conditions

$\begin{array}{ll}{\text{u}}_0(\text{x}+\text{r}\text{n}(\text{x}))={\text{u}}_0(\text{x})\text{ on }\partial \Omega ,& (3)\end{array}$

in which δu₀ is an arbitrary admissible kinematic variation of u₀(x₀), S = D₀ : E the second Piola-Kirchhoff stress tensor and n(x₀) the outward pointing normal at x₀ ∈ ∂Ω₀. The application of a standard Galerkin approach to Equation (2) results in the generalized system of equations

$\begin{array}{ll}\left[{\Gamma }_0^{\text{T}}\text{K}({\text{Q}}_0){\Gamma }_0\right]{\text{Q}}_0(\text{X})={\text{P}}_0(\text{X}),& (4)\end{array}$

where K(Q₀) is a static stiffness matrix, P₀ represents the global vector of nodal forces, Q₀ denotes a global vector of independent nodal displacements and Γ₀ is a mapping operator resulting from Equation (3) and realizing the condition U₀ = Γ₀Q₀, in which U₀ indicates the full vector of nodal displacements. In this work, the solution of Equation (4) is carried out using the weak form module implemented in Comsol Multiphysics 5.3 (Comsol, 2017).

2.2. Geometry and Mechanical Properties Updating

Once the set of independent static displacements Q₀ is obtained, the reference configuration is updated from C₀ to C by calculating the corresponding nodal coordinates x = x₀ + Γ₀Q₀(x₀). The updated material properties in C are given by $\rho ={\rho }_0{(det \text{F})}^{-1}$ and $D_{ijkl}={(\mathrm{\text{det }}\text{F})}^{-1}{F}_{iI}{F}_{jJ}{F}_{kK}{F}_{lL}{(D_0)}_{IJKL}$, while the Cauchy stress tensor is obtained from the relation σ = (det F)⁻¹ FSF^T. The geometry of the unit cell in the configuration C is then re-meshed and used as the basis for the linear dynamic analysis presented in the next section.

2.3. Dynamic Analysis Using the Floquet-Bloch Decomposition

Assuming C as the new reference configuration, the position vector for the unit cell in the dynamic deformed configuration C′ is given by x′ = x + u(x, t), where u(x, t) = u(x)exp(−iωt) is a time-harmonic perturbation superimposed on C, being t the time and ω the angular frequency. In the rest of the paper, the time-dependence exp(−iωt) is dropped for conciseness. According to the small-on-large displacement hypotheses usually assumed in acousto-elasticity (Mazzotti et al., 2012; Pau and Lanza di Scalea, 2015; Shim et al., 2015), only small perturbations u(x, t) are considered. In this case, the approximations x′ ≈ x and C′ ≈ C hold and any displacement-dependent vector and tensor field can be obtained in a linearized incremental form. Following this procedure, the stress tensor is linearized in the direction of u(x, t) by applying a first order Taylor series expansion of σ(x) about x, resulting in (Mazzotti et al., 2012)

$\begin{array}{ll}{{D}}_{\text{u}}\sigma (\text{x})={\sigma }_0(\text{x})\otimes \nabla \text{u}(\text{x})+\text{D}(\text{x}):\text{e}(\text{x}),& (5)\end{array}$

where $\text{e}(\text{x})=\frac{1}{2}\left[\nabla \text{u}(\text{x})+{(\nabla \text{u}(\text{x}))}^{\text{T}}\right]$ denotes the linearized Green-Lagrange strain tensor.

From the application of the Floquet-Bloch theorem, any small harmonic perturbation u(x) can be expressed in the form (Collet et al., 2011)

$\begin{array}{ll}\text{u}(\text{x})=\tilde{\text{u}}(\text{x})\mathrm{\text{exp}}(\mathrm{\text{i}}\text{k}·\text{x}),& (6)\end{array}$

in which $\tilde{\text{u}}(\text{x})$ is a Ω-periodic displacement amplitude and k ∈ Λ is the Bloch wavenumber vector, being Λ the reciprocal unit cell defined in C by the reciprocal lattice vector basis g_j, j = 1, 2, satisfying the condition r_i · g_j = 2πδ_ij, where δ_ij is the Kronecker delta. These are related to the reciprocal lattice vector basis g_0j in C₀ by ${\text{g}}_j={\text{F}}_L^{-\text{T}}{\text{g}}_{0j}$, j = 1, 2 (Zhang and Parnell, 2017). In the reciprocal lattice domain, the Bloch wavenumber vector is described in terms of its orientation angle ϑ as k(ϑ) = kΦ(ϑ), where ϑ is defined with respect to g₁, k = ‖k(ϑ)‖₂ and Φ(ϑ) = {cosϑ, sinϑ, 0}^T. The relation between the orientation vector in C₀ and C writes $\Phi ={\text{F}}_L^{-\text{T}}{\Phi }_0$.

From Equation (6), and by defining the k-shifted gradient of a generic Ω-periodic vector field $\tilde{\varphi }(\text{x})\in {ℂ}^3$ as

$\begin{array}{ll}{\nabla }_{\text{k}}\tilde{\varphi }(\text{x},\vartheta )={\nabla }_{\text{x}}\tilde{\varphi }(\text{x})+\text{i}k\tilde{\varphi }(\text{x})\otimes \Phi (\vartheta ),& (7)\end{array}$

the solution of the elastodynamic problem for free vibrations of the unit cell in C subjected to an initial stress σ₀ can be obtained from the variational statement

$\begin{array}{ll}\begin{array}{c}-{\omega }^2{\displaystyle {\int }_{\Omega }}\rho (\text{x})\delta {\tilde{\text{u}}}^{*}(\text{x})·\tilde{\text{u}}(\text{x})\text{d}\Omega \hfill \\ +{\displaystyle {\int }_{\Omega }}\delta {\tilde{\text{e}}}_{\text{k}}^{*}(\text{x},\vartheta ):\text{D}(\text{x}):{\tilde{\text{e}}}_{\text{k}}(\text{x},\vartheta )\text{d}\Omega ,\hfill \\ +{\displaystyle {\int }_{\Omega }}{\sigma }_0(\text{x}):\left[{({\nabla }_{\text{k}}\delta {\tilde{\text{u}}}^{*}(\text{x}))}^{\text{T}}{\nabla }_{\text{k}}\tilde{\text{u}}(\text{x})\right]\text{d}\Omega =0,\hfill \end{array}\hfill & (8)\hfill \end{array}$

subjected to the Dirichlet boundary condition

$\begin{array}{ll}\tilde{\text{u}}(\text{x}+\text{r}\text{n}(\text{x}))=\tilde{\text{u}}(\text{x})\text{ on}\partial \Omega ,& (9)\end{array}$

in which (·)* stands for the conjugate of a complex vector or tensor field, ${\tilde{\text{e}}}_{\text{k}}(\text{x},\vartheta )=\frac{1}{2}\left[{\nabla }_{\text{k}}\tilde{\text{u}}(\text{x})+{({\nabla }_{\text{k}}\tilde{\text{u}}(\text{x}))}^{\text{T}}\right]$ follows from Equation (7) and n(x) represents the outpointing surface normal at x ∈ ∂Ω. After the application of a Galerkin discretization scheme, Equations (8) and (9) lead to the following generalized linear eigenvalue problem

$\begin{array}{ll}\left\{{\Gamma }^{\text{T}}\left[k^2{\text{K}}_3(\vartheta )+\text{i}k({\text{K}}_2(\vartheta )-{\text{K}}_2^{\text{T}}(\vartheta ))+{\text{K}}_1-{\omega }^2\text{M}\right]\Gamma \right\}\tilde{\text{Q}}(\omega ,\vartheta )=0,& (10)\end{array}$

which forms the basis of the band structure analysis for the prestressed PC. In Equation (10), Γ is a mapping operator implementing the Dirichlet boundary condition in Equation (9) such that $\tilde{\text{U}}(\vartheta )=\Gamma \tilde{\text{Q}}(\vartheta )$ is verified, where $\tilde{\text{U}}(\vartheta )$ is the global vector of nodal displacement amplitudes and $\widetilde{\text{Q}(\vartheta )}$ a subvector of $\tilde{\text{U}}(\vartheta )$ collecting only its independent components. The mass operator M and the stiffness operators K₁, K₂, and K₃ are given by

$\begin{array}{llll}\text{M}& =& {\displaystyle \bigcup _e}{\displaystyle {\int }_{{\Omega }^{(e)}}}{\text{N}}^{\text{T}}(\text{x})\rho (\text{x})\text{N}(\text{x})\text{d}\Omega ,& (11)\end{array}$

$\begin{array}{llll}{\text{K}}_1& =& {\displaystyle \bigcup _e}{\displaystyle {\int }_{{\Omega }^{(e)}}}{\text{N}}^{\text{T}}(\text{x})\left[{\text{B}}^{\text{T}}\text{D}(\text{x})\text{B}+{\text{B}}_0^{\text{T}}{\Sigma }_0(\text{x}){\text{B}}_0\right]\text{N}(\text{x})\text{d}\Omega & (12)\end{array}$

$\begin{array}{llll}{\text{K}}_2(\vartheta )& =& {\displaystyle \bigcup _e}{\displaystyle {\int }_{{\Omega }^{(e)}}}{\text{N}}^{\text{T}}(\text{x})\left[{\text{B}}^{\text{T}}\text{D}(\text{x})\text{H}(\vartheta )+{\text{B}}_0^{\text{T}}{\Sigma }_0(\text{x}){\text{H}}_0(\vartheta )\right]\text{N}(\text{x})\text{d}\Omega ,& (13)\end{array}$

$\begin{array}{llll}{\text{K}}_3(\vartheta )& =& {\displaystyle \bigcup _e}{\displaystyle {\int }_{{\Omega }^{(e)}}}{\text{N}}^{\text{T}}(\text{x})\left[{\text{H}}^{\text{T}}(\vartheta )\text{D}(\text{x})\text{H}(\vartheta )+{\text{H}}_0^{\text{T}}(\vartheta ){\Sigma }_0(\text{x}){\text{H}}_0(\vartheta )\right]\text{N}(\text{x})\text{d}\Omega ,& (14)\end{array}$

where Ω^(e) denotes the domain of the e-th finite element of the mesh, ⋃_e(·) stands for the standard direct stiffness assembling procedure, N(x) is a matrix of shape functions for the e-th element, Σ₀(x) is a block-diagonal matrix of the form

$\begin{array}{ll}{\Sigma }_0(\text{x})=\left[\begin{array}{ccc}{\sigma }_0(\text{x})& 0& 0\\ 0& {\sigma }_0(\text{x})& 0\\ 0& 0& {\sigma }_0(\text{x})\end{array}\right],& (15)\end{array}$

while the different compatibility operators are expressed as

$\begin{array}{lll}\text{B}& ={\displaystyle \sum _{i=1}^3}\frac{\partial }{\partial x_i}{\text{L}}_i,\text{H}(\vartheta )={\displaystyle \sum _{i=1}^3}\left({\text{z}}_i^{\text{T}}\Phi (\vartheta )\right){\text{L}}_i,& (16)\end{array}$

$\begin{array}{lll}{\text{B}}_0& ={\displaystyle \sum _{i=1}^3}\frac{\partial }{\partial x_i}{\text{L}}_{0i},{\text{H}}_0(\vartheta )={\displaystyle \sum _{i=1}^3}\left({\text{z}}_i^{\text{T}}\Phi (\vartheta )\right){\text{L}}_{0i},& (17)\end{array}$

being z_i a unit vector identifying the i-th coordinate in the Cartesian frame of reference. The operators L_i and L_0i are given in the Appendix. It should be noted that, since the mass and stiffness matrices in Equations (11)–(14) are evaluated with respect to the deformed configuration C, a new mesh for the deformed geometry of the unit cell needs to be generated after the static analysis described in section 2.1 has been completed and before the eigenvalue analysis is carried out.

3. Band Structure Analysis and Energy Velocity Extraction

The homogeneous problem in Equation (10) can be solved in the wavenumber k(ϑ, ω) and corresponding Floquet eigenvectors $\tilde{\text{Q}}(\vartheta ,\omega )$ for any fixed direction ϑ and circular frequency ω ∈ ℝ, ω > 0, from which the band structures of the propagation and attenuation constants of a specific mode m are obtained by taking the real and imaginary components of the corresponding wavenumber k_m(ϑ, ω), respectively. For the particular case of a lossless structure (Im(D) = 0) immersed in vacuum, the modes supported by the crystal can be either propagative (Im(k) = 0) or evanescent (Re(k) = 0), the latter belonging to the so called deaf frequency range (bandgap).

In addition to the Bloch wavenumber k(ϑ, ω), from the computed set of solutions (k_m(ϑ, ω), Q_m(ϑ, ω)) it is possible to extract the band structure of the energy velocity, which corresponds to the velocity of propagation of packets of waves having close central frequency (Brillouin, 1953; Mazzotti et al., 2012). The energy velocity of a specific Bloch mode can be found as the ratio between the time-averaged energy flux in the direction ϑ over one period T = 2π/ω and the time-averaged total mechanical energy over the same period, i.e.

$\begin{array}{ll}{c}_e(\vartheta ,\omega )=\frac{\langle {I}(\vartheta ,\omega )\rangle }{\langle {K}(\vartheta ,\omega )\rangle +\langle {W}(\vartheta ,\omega )\rangle +\langle {{W}}_{{\sigma }_0}(\vartheta ,\omega )\rangle },& (18)\end{array}$

where $\langle \varphi \rangle =\frac{\omega }{2\pi }{\int }_t^{t+2\pi /\omega }\varphi dt$ denotes the time-averaging operation, ${I}$(ϑ, ω) denotes the total energy flux along the orientation ϑ of the Bloch wavenumber k, ${K}$(ϑ, ω) indicates the total kinetic energy while ${W}$(ϑ, ω) and ${W}$_σ0(ϑ, ω) represent the total stored elastic energy related to the harmonic motion and prestress, respectively. The time-averaged expressions for these quantities are given by

$\begin{array}{llll}\langle {I}(\vartheta ,\omega )\rangle & =& \frac{\omega }{2}\text{Im}\left\{{\displaystyle {\int }_{\Omega }}\left[({{D}}_{\text{u}}\tilde{\sigma }(\text{x},\vartheta ,\omega ){\tilde{\text{u}}}^{*}(\text{x},\vartheta ,\omega ))·\Phi (\vartheta )\right]\mathrm{\text{d}}\Omega \right\},& (19)\end{array}$

$\begin{array}{llll}\langle {K}(\vartheta ,\omega )\rangle & =& \frac{{\omega }^2}{4}\text{Re}\left\{{\displaystyle {\int }_{\Omega }}{\tilde{\text{u}}}^{*}(\text{x},\vartheta ,\omega )·\rho (\text{x})\tilde{\text{u}}(\text{x},\vartheta ,\omega )\mathrm{\text{d}}\Omega \right\},& (20)\end{array}$

$\begin{array}{llll}\langle {W}(\vartheta ,\omega )\rangle & =& \frac{1}{4}\text{Re}\left\{{\displaystyle {\int }_{\Omega }}\left[{\tilde{\text{e}}}_{\text{k}}^{*}(\text{x},\vartheta ,\omega ):\text{D}(\text{x}):{\tilde{\text{e}}}_{\text{k}}(\text{x},\vartheta ,\omega )\right]\mathrm{\text{d}}\Omega \right\},& (21)\end{array}$

$\begin{array}{llll}\langle {{W}}_{{\sigma }_0}(\vartheta ,\omega )\rangle & =& \frac{1}{4}\text{Re}\left\{{\displaystyle {\int }_{\Omega }}\left[{\left({({\nabla }_{\text{k}}\tilde{\text{u}}(\text{x},\vartheta ,\omega ))}^{\text{T}}{\nabla }_{\text{k}}\tilde{\text{u}}(\text{x},\vartheta ,\omega )\right)}^{*}:{\sigma }_0(\text{x})\right]\mathrm{\text{d}}\Omega \right\}.& (22)\end{array}$

Equation (18) can be evaluated at any given solution k_m(ϑ, ω) by means of a Gauss quadrature scheme over the finite element mesh of the unit cell in the configuration C, in which the displacement, strain and stress fields can be post-processed from Q_m(ϑ, ω) using nodal interpolations.

4. Numerical Applications: Epoxy Plate With Cylindrical Steel Inclusions

The numerical method presented in the previous section is here applied to study a PC plate made of steel cylinders embedded in an epoxy matrix. In the undeformed configuration, the plate is 2.5 mm thick, while the cylindrical inclusions have a radius of 3.0 mm and are arranged in a square lattice of 10.0 mm side length, as shown in Figure 2A. The material properties for the steel and the epoxy are reported in Table 1.

FIGURE 2

Figure 2. Unit cell of steel cylinders embedded in an epoxy matrix: undeformed geometry (A) and deformed geometry corresponding to isotropic traction (B) and compression (C). The deformation scale is amplified of 20 times.

TABLE 1

Table 1. Material constants for the steel and the epoxy.

To study how an initial state of stress affects the passbands and bandgaps of the PC plate, two different deformed configuration are considered. In the first configuration, the plate is subjected to a state of compressive stress representative of a normal displacement u₀(x₀) · n₀(x₀) = −0.025 mm applied at each point x₀ belonging to a lateral face of the unit cell. In the second case, a generic tensile state is applied by assuming a normal displacement u₀(x₀) · n₀(x₀) = 0.025 mm for each point x₀ belonging to a lateral face of the unit cell. Figures 2B,C report the von Mises stress consequent to the compression and tension prestress condition, respectively. It is worth noticing that the maximum magnitude of the applied prestretch has been limited to the typical strength at yield of both epoxy and steel in order to ensure the elastic behavior of the materials. This choice enables a full reversibility of the undeformed configuration, and thus of the dispersive behavior of the PC plate, once the load is removed, which is of more practical interest with respect to the case of a permanently induced deformation.

The band diagrams in terms of wavenumber k and energy velocity c_e versus frequency are shown in Figure 3. In these diagrams, the blue and red dots denote the dispersion curves relative to the unit cell subjected to the compressive stress and tensile stress shown in Figures 2B,C, respectively. It should be noted that, since the deformation due to prestress applied to each face of the cell is isotropic in the xy-plane in both cases, the deformation gradient F and its affine component F_L are diagonal and, as a consequence, the orientation of the reciprocal lattice vectors g_j in the deformed configurations does not change with respect to that in the deformed configuration. This implies that the orientation angle θ of the Bloch wavevector also remains unchanged between the undeformed and deformed configurations. The results presented in Figure 3 have been obtained for an angle θ = 0 measured with respect to g₀₁ and g₁. In classical Bloch analysis, this angle corresponds to the Γ − X direction of the irreducible Brillouin zone.

FIGURE 3

Figure 3. Band diagrams of the wavenumber k (Left) and energy velocity c_e (Right) versus frequency for the PC plate under isotropic compression and isotropic traction reported in Figures 2B,C, respectively.

From the dispersion analysis of the undeformed unit cell, it can be found that, at frequencies below 110 kHz, only one bandgap is present along the direction ϑ = 0, which is located in the [68.3–88.6] kHz frequency range (see Figure A1 in Appendix B). This result, which has not been reported in the Figure 3 for the sake of clarity, can be readily obtained by carrying out the analysis of sections 2.3 and 3 on the undeformed geometry and by setting σ₀ = 0 in Equation (15). However, when the cell is subjected to a compressive stress as shown in Figure 2B, the bandgap experiences a shift towards higher frequencies, being it located, in this case, in the [70.7 − 87.3] kHz frequency range. On the other hand, when the cell is subjected to traction, the lower bound of the bandgap is observed at 68.3 kHz while its upper bound at 84.3 kHz, which corresponds to a frequency downshift of approximately 3.5% of the bandgap with respect to the case of the PC plate under compression. These results suggest that, even in the elastic regime, a deformation due to prestress can lead to significant changes in the passband and bandgap behavior.

Furthermore, looking at the behavior of the energy velocity for the two deformed configurations, it is possible to infer that, similarly to the case of homogeneous plates (Pau and Lanza di Scalea, 2015), a generic state of compression leads to an increase of the energy velocity, whereas a state of traction leads, in general, to its decrease. This situation can be readily deduced from the inspection of the energy velocity band diagrams in Figure 3, from which it can be noted that the largest variations take place in proximity of the cutoff frequencies. Moreover, in the case of prestressed homogeneous plates (Dubuc et al., 2017), characteristic spikes in the energy velocity dispersion curves appear due to mode coupling. This behavior can also be observed for the third fundamental mode at around 36 kHz when the PC plate is subjected to a compression.

Finally, in the case of the PC plate under compression, it is possible to observe the presence of a cutoff wavenumber for the first fundamental mode at about 108 rad/m and null frequency. This value can be derived by solving the eigenvalue problem reported in Equation (10) in k(ϑ) by imposing ω = 0. In this case, the eigenvalue problem corresponds to that of the buckling load and, consequently, the lowest eigenvalue (the cutoff wavenumber) corresponds to a bifurcation point that indicates the onset of buckling for the infinite PC plate. It should be noted that the specific value of the cutoff wavenumber depends on the magnitude and distribution of the compressive stress in the PC.

5. Conclusions

An Updated Lagrangian computational scheme has been presented for the calculation of the band diagrams of phononic crystal plates subjected to a generic state of prestress. The scheme involves the solution of a geometrically nonlinear static problem for a representative unit cell and the application of the Floquet-Bloch theorem to the linearized equations of acousto-elasticity on the deformed configuration of the unit cell undergoing the static load. For the case of an epoxy plate with cylindrical steel inclusions, it has been demonstrated that the existence of a prestress state of compression or tension can lead to significant changes in the passbands and bandgaps of a phononic crystal plate, even in the case of prestress inducing full reversible deformation state, which is of more practical interest with respect to the case of a permanently induced deformation. It has been observed that when the cell is subjected to a compressive stress, the bandgap experiences a shift toward higher frequencies, whereas when the cell is subjected to traction, we observe a frequency downshift of approximately the 3.5% of the bandgap with respect to the case of the phononic plate under compression. Similarly to the case of homogeneous plates, a generic state of compression leads to an increase of the energy velocity, whereas a state of traction tends to lower its value. The largest variations have been observed in proximity of the cutoff frequencies. Characteristic spikes in the energy velocity dispersion curves appear due to mode coupling for the third fundamental mode at around 36 kHz when the phononic plate is subjected to a compression.

The results presented herein provide insights in the behavior of phononic crystal plates with tunable dispersive properties, and suggest new leverages for wave manipulation valuable in many application fields such as wave filters, waveguiding and beam splitting, sensing devices, and vibration shielding.

Author Contributions

MMa and MMi conceived the idea and wrote the first draft of the paper. MMa performed numerical calculations. All the authors discussed the manuscript.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The handling editor declared a past co-authorship with MMa.

Acknowledgments

MMi has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement N. 754364.

References

Airoldi, L., and Ruzzene, M. (2011). Design of tunable acoustic metamaterials through periodic arrays of resonant shunted piezos. New J. Phys. 13:113010. doi: 10.1088/1367-2630/13/11/113010

CrossRef Full Text | Google Scholar

Bergamini, A., Delpero, T., De Simoni, L., Di Lillo, L., Ruzzene, M., and Ermanni, P. (2014). Phononic crystal with adaptive connectivity. Adv. Mater. 26, 1343–1347. doi: 10.1002/adma.201305280

PubMed Abstract | CrossRef Full Text | Google Scholar

Bigoni, D., Gei, M., and Movchan, A. (2008). Dynamics of a prestressed stiff layer on an elastic half space: filtering and band gap characteristics of periodic structural models derived from long-wave asymptotics. J. Mech. Phys. Solids 56, 2494–2520. doi: 10.1016/j.jmps.2008.02.007

CrossRef Full Text | Google Scholar

Bonet, J., and Wood, R. D. (2008). Nonlinear Continuum Mechanics for Finite Element Analysis. New York, NY: Cambridge University Press.

Google Scholar

Brillouin, L. (1953). Wave Propagation in Periodic Structures: Electric Filters and Crystal Lattices, 2nd Edn. New York, NY: Dover Publications (1953).

Google Scholar

Casadei, F., Ruzzene, M., Dozio, L., and Cunefare, K. A. (2010). Broadband vibration control through periodic arrays of resonant shunts: experimental investigation on plates. Smart Mater. Struct. 19:015002. doi: 10.1088/0964-1726/19/1/015002

CrossRef Full Text | Google Scholar

Cattani, C., and Rushchitsky, J. (2007). Wavelet and Wave Analysis as Applied to Materials with Micro or Nanostructure. Singapore: World Scientific.

Cheng, Y., Liu, X. J., and Wu, D. J. (2011). Temperature effects on the band gaps of lamb waves in a one-dimensional phononic-crystal plate (l). J. Acoust. Soc. Am. 129, 1157–1160. doi: 10.1121/1.3543970

PubMed Abstract | CrossRef Full Text | Google Scholar

Collet, M., Ouisse, M., Ruzzene, M., and Ichchou, M. (2011). Floquet-bloch decomposition for the computation of dispersion of two-dimensional periodic, damped mechanical systems. Int. J. Solids Struct. 48, 2837–2848. doi: 10.1016/j.ijsolstr.2011.06.002

CrossRef Full Text | Google Scholar

Colombi, A., Colquitt, D., Roux, P., Guenneau, S., and Craster, R. V. (2016). A seismic metamaterial: the resonant metawedge. Sci. Rep. 6:27717. doi: 10.1038/srep27717

PubMed Abstract | CrossRef Full Text | Google Scholar

Comsol (2017). Comsol Multiphysics 5.3. Available online at: https://www.comsol.com/

Craster, R. V., and Guenneau, S. (2012). Acoustic Metamaterials: Negative Refraction, Imaging, Lensing and Cloaking. Springer Series in Materials Science. London: Springer (2012).

Deymier, P. A. (2013). Acoustic Metamaterials and Phononic Crystals, Springer Series in Solid-State Sciences. Berlin; Heidelberg: Springer (2013).

Dubuc, B., Ebrahimkhanlou, A., and Salamone, S. (2017). The effect of applied stress on the phase and group velocity of guided waves in anisotropic plates. J. Acoust. Soc. Am. 142, 3553–3563. doi: 10.1121/1.5016969

PubMed Abstract | CrossRef Full Text | Google Scholar

Dubuc, B., Ebrahimkhanlou, A., and Salamone, S. (2018). Higher order longitudinal guided wave modes in axially stressed seven-wire strands. Ultrasonics 84, 382–391. doi: 10.1016/j.ultras.2017.12.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Galich, P., Slesarenko, V., and Rudykh, S. (2016). Elastic wave propagation in soft microstructured composites undergoing finite deformations. PAMM 16, 627–628. doi: 10.1002/pamm.201610302

CrossRef Full Text | Google Scholar

Gei, M. (2010). Wave propagation in quasiperiodic structures: stop/pass band distribution and prestress effects. Int. J. Solids Struct. 47, 3067–3075. doi: 10.1016/j.ijsolstr.2010.07.008

CrossRef Full Text | Google Scholar

Gei, M., Movchan, A., and Bigoni, D. (2009). Band-gap shift and defect-induced annihilation in prestressed elastic structures. J. Appl. Phys. 105:063507. doi: 10.1063/1.3093694

CrossRef Full Text | Google Scholar

Guo, X., and Wei, P. (2016). Dispersion relations of elastic waves in one-dimensional piezoelectric/piezomagnetic phononic crystal with initial stresses. Ultrasonics 66, 72–85. doi: 10.1016/j.ultras.2015.11.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Jim, K. L., Leung, C. W., Lau, S. T., Choy, S. H., and Chan, H. L. W. (2009). Thermal tuning of phononic bandstructure in ferroelectric ceramic/epoxy phononic crystal. Appl. Phys. Lett. 94:193501. doi: 10.1063/1.3136752

CrossRef Full Text | Google Scholar

Kherraz, N., Haumesser, L., Levassort, F., Benard, P., and Morvan, B. (2016). Controlling bragg gaps induced by electric boundary conditions in phononic piezoelectric plates. Appl. Phys. Lett. 108:093503. doi: 10.1063/1.4943138

CrossRef Full Text | Google Scholar

Kushwaha, M. S., and Halevi, P. (1994). Band gap engineering in periodic elastic composites. Appl. Phys. Lett. 64, 1085–1087. doi: 10.1063/1.110940

CrossRef Full Text | Google Scholar

Li, J., Wu, F., Zhong, H., Yao, Y., and Zhang, X. (2015). Acoustic beam splitting in two-dimensional phononic crystals using self-collimation effect. J. Appl. Phys. 118:144903. doi: 10.1063/1.4932138

CrossRef Full Text | Google Scholar

Liu, Z., Zhang, X., Mao, Y., Zhu, Y. Y., Yang, Z., Chan, C. T., et al. (2000). Locally resonant sonic materials. Science 289, 1734–1736. doi: 10.1126/science.289.5485.1734

PubMed Abstract | CrossRef Full Text | Google Scholar

Mártinez-Sala, R., Sancho, J., Sánchez, J. V., Gomez, V., Llinares, J., and Meseguer, F. (1995). Sound attenuation by sculpture. Nature 378:241. doi: 10.1038/378241a0

CrossRef Full Text | Google Scholar

Matar, O. B., Robillard, J. F., Vasseur, J. O., Hladky-Hennion, A.-C., Deymier, P. A., Pernod, P., et al. (2012). Band gap tunability of magneto-elastic phononic crystal. J. Appl. Phys. 111:054901. doi: 10.1063/1.3687928

CrossRef Full Text | Google Scholar

Mazzotti, M., Marzani, A., Bartoli, I., and Viola, E. (2012). Guided waves dispersion analysis for prestressed viscoelastic waveguides by means of the safe method. Int. J. Solids Struct. 49, 2359–2372. doi: 10.1016/j.ijsolstr.2012.04.041

CrossRef Full Text | Google Scholar

Mazzotti, M., Miniaci, M., and Bartoli, I. (2017). Band structure analysis of leaky bloch waves in 2d phononic crystal plates. Ultrasonics 74, 140–143. doi: 10.1016/j.ultras.2016.10.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Miniaci, M., Gliozzi, A. S., Morvan, B., Krushynska, A., Bosia, F., Scalerandi, M., et al. (2017). Proof of concept for an ultrasensitive technique to detect and localize sources of elastic nonlinearity using phononic crystals. Phys. Rev. Lett. 118:214301. doi: 10.1103/PhysRevLett.118.214301

PubMed Abstract | CrossRef Full Text | Google Scholar

Miniaci, M., Krushynska, A., Bosia, F., and Pugno, N. M. (2016). Large scale mechanical metamaterials as seismic shields. N. J. Phys. 18:083041. doi: 10.1088/1367-2630/18/8/083041

CrossRef Full Text | Google Scholar

Miniaci, M., Mazzotti, M., Radzieski, M., Kherraz, N., Kudela, P., Ostachowicz, W., et al. (2018a). Experimental observation of a large low-frequency band gap in a polymer waveguide. Front. Mater. 5:8. doi: 10.3389/fmats.2018.00008

CrossRef Full Text | Google Scholar

Miniaci, M., Pal, R. K., Manna, R., and Ruzzene, M. (2019). Valley based splitting of topologically protected helical waves in elastic plates. arXiv 1811.04814.

Google Scholar

Miniaci, M., Pal, R. K., Morvan, B., and Ruzzene, M. (2018b). Experimental observation of topologically protected helical edge modes in patterned elastic plates. Phys. Rev. X 8:031074. doi: 10.1103/PhysRevX.8.031074

CrossRef Full Text | Google Scholar

Morvan, B., Tinel, A., Hladky-Hennion, A.-C., Vasseur, J., and Dubus, B. (2010). Experimental demonstration of the negative refraction of a transverse elastic wave in a two-dimensional solid phononic crystal. Appl. Phys. Lett. 96:101905. doi: 10.1063/1.3302456

CrossRef Full Text | Google Scholar

Mousavi, S. H., Khanikaev, A. B., and Wang, Z. (2015). Topologically protected elastic waves in phononic metamaterials. Nat. Comm. 6:8682. doi: 10.1038/ncomms9682

PubMed Abstract | CrossRef Full Text | Google Scholar

Murnaghan, F. D. (1937). Finite deformations of an elastic solid. Amer. J. Math. 59, 235–260. doi: 10.2307/2371405

CrossRef Full Text | Google Scholar

Norris, A. N., and Parnell, W. J. (2012). Hyperelastic cloaking theory: transformation elasticity with pre-stressed solids. Proc. R. Soc. Lond. A 468, 2881–2903. doi: 10.1098/rspa.2012.0123

CrossRef Full Text | Google Scholar

Pal, R. K., Schaeffer, M., and Ruzzene, M. (2016). Helical edge states and topological phase transitions in phononic systems using bi-layered lattices. J. Appl. Phys. 119:084305. doi: 10.1063/1.4942357

CrossRef Full Text | Google Scholar

Pau, A., and Lanza di Scalea, F. (2015). Nonlinear guided wave propagation in prestressed plates. J. Acoust. Soc. Am. 137, 1529–1540. doi: 10.1121/1.4908237

PubMed Abstract | CrossRef Full Text | Google Scholar

Pennec, Y., Vasseur, J. O., Djafari-Rouhani, B., Dobrzyśki, L., and Deymier, P. A. (2010). Two-dimensional phononic crystals: examples and applications. Surf. Sci. Rep. 65, 229–291. doi: 10.1016/j.surfrep.2010.08.002

CrossRef Full Text | Google Scholar

Robillard, J.-F., Matar, O. B., Vasseur, J. O., Deymier, P. A., Stippinger, M., Hladky-Hennion, A.-C., et al. (2009). Tunable magnetoelastic phononic crystals. Appl. Phys. Lett. 95:124104. doi: 10.1063/1.3236537

CrossRef Full Text | Google Scholar

Sainidou, R., Stefanou, N., Psarobas, I., and Modinos, A. (2005). A layer-multiple-scattering method for phononic crystals and heterostructures of such. Comput. Phys. Commun. 166, 197–240. doi: 10.1016/j.cpc.2004.11.004

CrossRef Full Text | Google Scholar

Shim, J., Wang, P., and Bertoldi, K. (2015). Harnessing instability-induced pattern transformation to design tunable phononic crystals. Int. J. Solids Struct. 58, 52–61. doi: 10.1016/j.ijsolstr.2014.12.018

CrossRef Full Text | Google Scholar

Su, X., and Norris, A. N. (2016). Focusing, refraction, and asymmetric transmission of elastic waves in solid metamaterials with aligned parallel gaps. J. Acoust. Soc. Am. 139, 3386–3394. doi: 10.1121/1.4950770

PubMed Abstract | CrossRef Full Text | Google Scholar

Sukhovich, A., Merheb, B., Muralidharan, K., Vasseur, J. O., Pennec, Y., Deymier, P. A., et al. (2009). Experimental and theoretical evidence for subwavelength imaging in phononic crystals. Phys. Rev. Lett. 102:154301. doi: 10.1103/PhysRevLett.102.154301

PubMed Abstract | CrossRef Full Text | Google Scholar

Süsstrunk, R., and Huber, S. D. (2015). Observation of phononic helical edge states in a mechanical topological insulator. Science 349, 47–50. doi: 10.1126/science.aab0239

PubMed Abstract | CrossRef Full Text | Google Scholar

Trainiti, G., Rimoli, J. J., and Ruzzene, M. (2015). Wave propagation in periodically undulated beams and plates. Int. J. Solids Struct. 75–76, 260–276. doi: 10.1016/j.ijsolstr.2015.08.019

CrossRef Full Text | Google Scholar

Walker, E., Reyes, D., Rojas, M. M., Krokhin, A., Wang, Z., and Neogi, A. (2014). Tunable ultrasonic phononic crystal controlled by infrared radiation. Appl. Phys. Lett. 105:143503. doi: 10.1063/1.4894489

CrossRef Full Text | Google Scholar

Walker, E. L., Wang, Z., and Neogi, A. (2017). Radio-frequency actuated polymer-based phononic meta-materials for control of ultrasonic waves. NPG Asia Mater. 9:e350. doi: 10.1038/am.2016.209

CrossRef Full Text | Google Scholar

Wang, L., and Bertoldi, K. (2012). Mechanically tunable phononic band gaps in three-dimensional periodic elastomeric structures. Int. J. Solids Struct. 49, 2881–2885. doi: 10.1016/j.ijsolstr.2012.05.008

CrossRef Full Text | Google Scholar

Wu, Y., Yu, K., Yang, L., Zhao, R., Shi, X., and Tian, K. (2018). Effect of thermal stresses on frequency band structures of elastic metamaterial plates. J. Sound Vibrat. 413, 101–119. doi: 10.1016/j.jsv.2017.10.014

CrossRef Full Text | Google Scholar

Zhang, P., and Parnell, W. J. (2017). Soft phononic crystals with deformation-independent band gaps. Proc. R. Soc. A Math. Phys. Eng. Sci. 473:20160865. doi: 10.1098/rspa.2016.0865

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, S., Shi, Y., and Gao, Y. (2017). Tunability of band structures in a two-dimensional magnetostrictive phononic crystal plate with stress and magnetic loadings. Phys. Lett. A 381, 1055–1066. doi: 10.1016/j.physleta.2017.01.044

CrossRef Full Text | Google Scholar

Zhu, H., and Semperlotti, F. (2016). Anomalous refraction of acoustic guided waves in solids with geometrically tapered metasurfaces. Phys. Rev. Lett. 117:034302. doi: 10.1103/PhysRevLett.117.034302

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhu, Q., Burtin, C., and Binetruy, C. (2014). Acoustoelastic effect in polyamide 6: linear and nonlinear behaviour. Polym. Testing 40, 178–186. doi: 10.1016/j.polymertesting.2014.09.007

CrossRef Full Text | Google Scholar

Zhu, R., Liu, X. N., Hu, G. K., Sun, C. T., and Huang, G. L. (2014). Negative refraction of elastic waves at the deep-subwavelength scale in a single-phase metamaterial. Nat. Commun. 5, 241–241. doi: 10.1038/ncomms6510

PubMed Abstract | CrossRef Full Text | Google Scholar

Appendix A

The operators L_i and L_0i are expressed as

$\begin{array}{lll}{\text{L}}_1& =\left[\begin{array}{ccc}1& 0& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 0& 1\\ 0& 1& 0\end{array}\right],{\text{L}}_2=\left[\begin{array}{ccc}0& 0& 0\\ 0& 1& 0\\ 0& 0& 0\\ 0& 0& 1\\ 0& 0& 0\\ 1& 0& 0\end{array}\right],{\text{L}}_3=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0& 0& 1\\ 0& 1& 0\\ 1& 0& 0\\ 0& 0& 0\end{array}\right],& (23)\end{array}$

$\begin{array}{lll}{\text{L}}_{01}& =\left[\begin{array}{ccc}1& 0& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 1& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 0& 1\\ 0& 0& 0\\ 0& 0& 0\end{array}\right],{\text{L}}_{02}=\left[\begin{array}{ccc}0& 0& 0\\ 1& 0& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 1& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 0& 1\\ 0& 0& 0\end{array}\right],{\text{L}}_{03}=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 1& 0& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 1& 0\\ 0& 0& 0\\ 0& 0& 0\\ 0& 0& 1\end{array}\right].& (24)\end{array}$

Appendix B

The band diagram in terms of wavenumber k and energy velocity c_e for the undeformed unit cell of Figure 2 at ϑ = 0 is reported in Figure A1.

FIGURE A1

Figure A1. Band diagrams of the wavenumber k (Left) and energy velocity c_e (Right) vs. frequency for the PC plate in the undeformed configuration reported in Figure 2A.