Modeling a labyrinthine acoustic metamaterial through an inertia-augmented relaxed micromorphic approach

Abstract

We present an inertia-augmented relaxed micromorphic model that enriches the relaxed micromorphic model previously introduced by the authors via a term Curl $\stackrel{·}{P}$ in the kinetic energy density. This enriched model allows us to obtain a good overall fitting of the dispersion curves while introducing the new possibility of describing modes with negative group velocity that are known to trigger negative refraction effects. The inertia-augmented model also allows for more freedom on the values of the asymptotes corresponding to the cut-offs. In the previous version of the relaxed micromorphic model, the asymptote of one curve (pressure or shear) is always bounded by the cut-off of the following curve of the same type. This constraint does not hold anymore in the enhanced version of the model. While the obtained curves’ fitting is of good quality overall, a perfect quantitative agreement must still be reached for very small wavelengths that are close to the size of the unit cell.

1. Introduction

Metamaterials are materials whose mechanical properties go beyond those of classical materials thanks to their heterogeneous microstructure. They can show unusual static/dynamic responses such as negative Poisson’s ratio [1], twist or bend in response to being pushed or pulled [2,3], band gaps [4–7], cloaking [8,9], focusing [10,11], channeling [12,13], negative refraction [14–16], etc. The working frequency of each metamaterial strongly depends on the characteristic size and the geometry of the underlying unit cell, as well as on the choice of the base material. In this paper, we present a labyrinthine metamaterial that, thanks to the use of a polymer-based material and an optimized distribution of mass inside the unit cell (see Figure 1), gives rise to a wide acoustic band gap with characteristic unit cell’s size of the order of centimeters. Conceptually, similar geometries have already been explored in literature [17–22].

Figure 1.

Geometry of the unit cell. Left: details of one unit cell (rotated by ${45}^{°}$ ) showing the tetragonal symmetry. For $a=20\mathrm{mm}$ as the size of the unit size we consider later, both the bars and holes have a thickness of 0.4mm each. Right: a $3\times 3$ section of the metamaterial made up of this unit cell (red dashed square).

The direct finite element modeling of structures build up with this labyrinthine metamaterial is unfeasible due to the extremely tight meshing that would be needed to correctly cover the narrow strips of material inside each unit cell. It is thus apparent the need for a homogenized model to use this type of very promising metamaterials in actual engineering designs. Various homogenization techniques have been developed with the purpose of providing rigorous predictions of the macroscopic metamaterial’s mechanical response when the properties of the base materials and their spatial distribution are known. These homogenization approaches have been shown to be useful in describing the overall behavior of metamaterials in the static and quasi-static regimes [23–35] as well as, more recently, in the dynamic regime [36–47]. However, these models are often unsuited to deal with finite-size metamaterials, because they are based on upscaling techniques valid for unbounded media. Because of that finite-size metamaterials’ structures are mostly investigated via finite element simulations which are performed using directly the microstructured material [48]. The downside of this approach is that the computational cost quickly becomes unsustainable (especially for unit cells as the one presented in this paper), although the propagation patterns obtained are very accurate. This heavily limits the possibility of exploring large-scale or very convoluted geometric meta-structures.

To overcome this problem and open up the possibility of designing complex meta-structures using the metamaterial presented in this paper as a basic building block, we propose to use an inertia-augmented relaxed micromorphic model. This model is based on the relaxed micromorphic model that we previously established [49–53] and has been augmented with a new inertia term accounting for coupled space-time derivatives of the micro-distortion tensor. The relaxed micromorphic model has extensively proven its efficacy in describing the broadband behavior of many infinite and finitesize metamaterials [52–56] and is extended in this paper so as to be able to account for negative group velocity which was not the case before. We will show that the proposed model is able to describe well the labyrinthine metamaterial’s response for a large range of frequencies (going beyond the first band gap) and wave numbers (approaching the size of the unit cell) and for all directions of propagation with a limited number of frequency- and scale-independent constitutive parameters. The new inertia-augmented term will be shown to trigger modes with negative group velocities that are known to be associated with negative refraction phenomena. The results presented in this paper will allow us to shortly present new designs of finite-size labyrinthine metamaterials’ structures that can control elastic energy in the acoustic regime for eventual subsequent re-use.

1.1. A polyethylene-based metamaterial for acoustic control

In this section, we present a new unit cell’s design that gives rise to a metamaterial for acoustic control. This unit cell is designed to achieve a band gap at relatively low frequencies $(600-2000\mathrm{Hz})$ so that application for acoustic control can be targeted. The unit cell considered is made out of polyethylene, cf. Table 1. Compared to aluminium or titanium, which we used for the metamaterials studied in literature [54–56], polyethylene gives rise to lower wave speeds, thus allowing band-gap phenomena to appear at lower frequencies.

Table 1.

Size, bulk material constants and wave speed for polyethylene unit cell where ${\overline{\kappa }}_{\mathrm{B}}$ is the plane strain bulk modulus. ^a

$a$	$\rho$	${\overline{\kappa }}_{\mathrm{B}}$	${\mu }_{\mathrm{B}}$	$c_{\mathrm{p},\mathrm{B}}$	$c_{\mathrm{s},\mathrm{B}}$
(mm)	$(\mathrm{kg}/{\mathrm{m}}^2)$	(MPa)	(MPa)	(m/s)	(m/s)
20	900	3160	262	1950	540

a

The relations between the wave speed and the elastic constants are $c_{\mathrm{p}}=\sqrt{\frac{{\overline{\kappa }}_{\mathrm{B}}+{\mu }_{\mathrm{B}}}{\rho }}=\sqrt{\frac{{\lambda }_{\mathrm{B}}+2{\mu }_{\mathrm{B}}}{\rho }}$ and $c_{\mathrm{s}}=\sqrt{\frac{{\mu }_{\mathrm{B}}}{\rho }}$ , while the plane strain bulk modulus is ${\overline{\kappa }}_{\mathrm{B}}={\lambda }_{\mathrm{B}}+{\mu }_{\mathrm{B}}$ .

A further lowering of the band gap is obtained through the adoption of a labyrinth-type geometry, cf. Figure 1. This structure presents a tetragonal symmetry and thus features a reduced number of parameters with respect to a fully anisotropic system. The circular center of the unit cell is connected by thin bars allowing the heavier center to move easily, thus giving rise to local resonance phenomena of relatively low frequencies while additionally providing a very soft macro-material behavior.

2. Relaxed micromorphic modeling of finite-size metamaterials

We briefly recall the weak and strong form of the governing equations of both classical Cauchy and relaxed micromorphic continua. The Lagrangian ${\mathcal{L}}_{\mathrm{c}}$ for the classical Cauchy model is

$${\mathcal{L}}_{\mathrm{c}}(\nabla u,\stackrel{·}{u})=\frac{1}{2}\rho \parallel \stackrel{·}{u}{\parallel }^2-\frac{1}{2}〈\mathbb{C}\mathrm{sym}\nabla u,\mathrm{sym}\nabla u〉,$$

(1)

where $u$ is the displacement field, $〈·,·〉:{\mathbb{R}}^{n\times n}\times {\mathbb{R}}^{n\times n}\to \mathbb{R}$ is the scalar product, $\rho$ is the apparent mass density, and $\mathbb{C}$ is the classical fourth-order elasticity tensor. The Lagrangian ${\mathcal{L}}_{\mathrm{m}}$ for the relaxed micromorphic model enhanced with the micro-inertia term $\mathrm{Curl}\stackrel{·}{P}$ is [52–55] ¹

$$\begin{array}{c}{\mathcal{L}}_{\mathrm{m}}(\stackrel{·}{u},\nabla \stackrel{·}{u},\stackrel{·}{P},\mathrm{Curl}\stackrel{·}{P},\nabla u,P,\mathrm{Curl}P)=\frac{1}{2}\rho \parallel \stackrel{·}{u}{\parallel }^2+\frac{1}{2}〈{\mathbb{J}}_{\mathrm{m}}\mathrm{sym}\stackrel{·}{P},\mathrm{sym}\stackrel{·}{P}〉+〈{\mathbb{J}}_{\mathrm{c}}\mathrm{skew}\stackrel{·}{P},\mathrm{skew}\stackrel{·}{P}〉\hfill \\ +\frac{1}{2}〈{\mathbb{T}}_{\mathrm{e}}\mathrm{sym}\nabla \stackrel{·}{u},\mathrm{sym}\nabla \stackrel{·}{u}〉+\frac{1}{2}〈{\mathbb{T}}_{\mathrm{c}}\mathrm{skew}\nabla \stackrel{·}{u},\mathrm{skew}\nabla \stackrel{·}{u}〉\hfill \\ +\frac{1}{2}〈{\mathbb{M}}_{\mathrm{s}}\mathrm{sym}\mathrm{Curl}\stackrel{·}{P},\mathrm{sym}\mathrm{Curl}\stackrel{·}{P}〉+\frac{1}{2}〈{\mathbb{M}}_{\mathrm{a}}\mathrm{skew}\mathrm{Curl}\stackrel{·}{P},\mathrm{skew}\mathrm{Curl}\stackrel{·}{P}〉\hfill \\ -\frac{1}{2}〈{\mathbb{C}}_{\mathrm{e}}\mathrm{sym}(\nabla u-P),\mathrm{sym}(\nabla u-P)〉\hfill \\ -\frac{1}{2}〈{\mathbb{C}}_{\mathrm{c}}\mathrm{skew}(\nabla u-P),\mathrm{skew}(\nabla u-P)〉-\frac{1}{2}〈{\mathbb{C}}_{\mathrm{micro}}\mathrm{sym}P,\mathrm{sym}P〉\hfill \\ -\frac{1}{2}〈{\mathbb{L}}_{\mathrm{s}}\mathrm{sym}\mathrm{Curl}P,\mathrm{sym}\mathrm{Curl}P〉-\frac{1}{2}〈{\mathbb{L}}_{\mathrm{a}}\mathrm{skew}\mathrm{Curl}P,\mathrm{skew}\mathrm{Curl}P〉,\hfill \end{array}$$

(2)

where $u\in {\mathbb{R}}^3$ is the macroscopic displacement field; $P\in {\mathbb{R}}^{3\times 3}$ is the non-symmetric micro-distortion tensor; $\rho$ is the macroscopic apparent density; ${\mathbb{J}}_{\mathrm{m}}$ , ${\mathbb{J}}_{\mathrm{c}}$ , ${\mathbb{T}}_{\mathrm{e}}$ , ${\mathbb{T}}_{\mathrm{c}}$ , ${\mathbb{M}}_{\mathrm{s}}$ , ${\mathbb{M}}_{\mathrm{a}}$ are fourth-order micro-inertia tensors; and ${\mathbb{C}}_{\mathrm{e}}$ , ${\mathbb{C}}_{\mathrm{c}}$ , ${\mathbb{C}}_{\mathrm{micro}}$ , ${\mathbb{L}}_{\mathrm{s}}$ , and ${\mathbb{L}}_{\mathrm{a}}$ are fourth-order elastic tensors (for more details see Appendix 1). ² The action functional ${\mathcal{A}}_i$ of the considered continuum can be defined based on the Lagrangian function ${\mathcal{L}}_i$ and the external work ${\mathcal{U}}_i$ as

$${\mathcal{A}}_{\mathrm{c}}=\underset{\mathrm{\Omega }\times [0,T]}{\int }{\mathcal{L}}_{\mathrm{c}}\mathrm{d}x\mathrm{d}t-\underset{0}{\overset{T}{\int }}{\mathcal{U}}_{\mathrm{c}}^{\mathrm{ext}}\mathrm{d}t,\mathrm{and}{\mathcal{A}}_{\mathrm{m}}=\underset{\mathrm{\Omega }\times [0,T]}{\int }{\mathcal{L}}_{\mathrm{m}}\mathrm{d}x\mathrm{d}t-\underset{0}{\overset{T}{\int }}{\mathcal{U}}_{\mathrm{m}}^{\mathrm{ext}}\mathrm{d}t,$$

(3)

for classical Cauchy and relaxed micromorphic media, respectively. In equation (3), $\mathrm{\Omega }\subset {\mathbb{R}}^3$ is the domain of the considered continuum in its reference configuration and $[0,T]$ is a time interval during which the deformation of the continuum is observed. In the case of a conservative system, the virtual work of internal actions ${\mathcal{W}}_i^{\mathrm{int}}=\delta {\mathcal{L}}_i$ can be defined as the first variation of the Lagrangian function, while the virtual work of external actions ${\mathcal{W}}_i^{\mathrm{ext}}=\delta {\mathcal{U}}_i^{\mathrm{ext}}$ is given as the first variation of the external work. In formulas, we have

$$\delta {\mathcal{A}}_{\mathrm{c}}=\underset{0}{\overset{T}{\int }}{\mathcal{W}}_{\mathrm{c}}^{\mathrm{int}}\mathrm{d}t-\underset{0}{\overset{T}{\int }}{\mathcal{W}}_{\mathrm{c}}^{\mathrm{ext}}\mathrm{d}t,\mathrm{and}\delta {\mathcal{A}}_{\mathrm{m}}=\underset{0}{\overset{T}{\int }}{\mathcal{W}}_{\mathrm{m}}^{\mathrm{int}}\mathrm{d}t-\underset{0}{\overset{T}{\int }}{\mathcal{W}}_{\mathrm{m}}^{\mathrm{ext}}\mathrm{d}t,$$

(4)

for classical Cauchy and relaxed micromorphic media, respectively. In equation (4), the variation operator $\delta$ indicates that the variation must be taken with respect to the unknown kinematics fields ( $u$ for Cauchy media and $(u,P)$ for the relaxed micromorphic model). Following classical variational calculus, the strong form of the bulk equations of motion and the Neumann boundary conditions for the Cauchy and the relaxed micromorphic model can be obtained via a least action principle stating that the first variation of the action functional must be vanishing. In absence of external body loads, the application of a least-action principle to Cauchy and relaxed micromorphic models gives the following equilibrium equation: [49, 53, 54]

$$\rho \stackrel{··}{u}=\mathrm{Div}\sigma ,\sigma :=\mathbb{C}\mathrm{sym}\nabla u,$$

(5)

for the classical Cauchy model, and

$$\rho \stackrel{··}{u}-=\mathrm{Div}\hat{\sigma }=\mathrm{Div}\tilde{\sigma },\overline{\sigma }=\tilde{\sigma }-s-\mathrm{Curl}m-\mathrm{Curl}\hat{m},$$

(6)

for the relaxed micromorphic model, where we set

$$\begin{array}{c}\tilde{\sigma }:={\mathbb{C}}_{\mathrm{e}}\mathrm{sym}(\nabla u-P)+{\mathbb{C}}_{\mathrm{c}}\mathrm{skew}(\nabla u-P),\hat{\sigma }:={\mathbb{T}}_{\mathrm{e}}\mathrm{sym}\nabla \stackrel{··}{u}+{\mathbb{T}}_{\mathrm{c}}\mathrm{skew}\nabla \stackrel{··}{u},\hfill \\ \overline{\sigma }:={\mathbb{J}}_{\mathrm{m}}\mathrm{sym}\stackrel{··}{P}+{\mathbb{J}}_{\mathrm{c}}\mathrm{skew}\stackrel{··}{P},s:={\mathbb{C}}_{\mathrm{micro}}\mathrm{sym}P,\hfill \\ m:={\mathbb{L}}_{\mathrm{s}}\mathrm{sym}\mathrm{Curl}P+{\mathbb{L}}_{\mathrm{a}}\mathrm{skew}\mathrm{Curl}P,\hat{m}:={\mathbb{M}}_{\mathrm{s}}\mathrm{sym}\mathrm{Curl}\stackrel{··}{P}+{\mathbb{M}}_{\mathrm{a}}\mathrm{skew}\mathrm{Curl}\stackrel{··}{P}.\hfill \end{array}$$

(7)

The Neumann boundary condition for the classical Cauchy model are

$$t_{\mathrm{c}}:=\sigma n=t_{\mathrm{c}}^{\mathrm{ext}},$$

(8)

with $t_{\mathrm{c}}^{\mathrm{ext}}$ as the externally applied traction vector and $n$ as the outward-pointing normal to the boundary, and for the relaxed micromorphic model

$$t_{\mathrm{m}}:=(\tilde{\sigma }+\hat{\sigma })n=t_{\mathrm{m}}^{\mathrm{ext}}\mathrm{and}\tau :=(m+\hat{m})\times n={\tau }^{\mathrm{ext}},$$

(9)

where $t_{\mathrm{m}}^{\mathrm{ext}}$ is the generalized traction vector, ${\tau }^{\mathrm{ext}}$ is the double traction second-order tensor, and the cross product × is understood row-wise. In absence of curvature terms $({\mathbb{L}}_{\mathrm{s}}={\mathbb{L}}_{\mathrm{a}}={\mathbb{M}}_{\mathrm{s}}={\mathbb{M}}_{\mathrm{a}}=0)$ , the boundary condition (9)₂ involving $\tau$ must not be assigned on the boundary.

As it is well known, Dirichlet-type boundary conditions can also be alternatively considered in the form $u=u^{\mathrm{ext}}$ for Cauchy media, while $u=u^{\mathrm{ext}}$ and $P\times n={\varphi }^{\mathrm{ext}}$ for the relaxed micromorphic model where $n$ is the outward-pointing normal vector to the surface and ${\varphi }^{\mathrm{ext}}$ is an assigned second-order tensor.

2.1. Tetragonal symmetry / shape of elastic tensors (in Voigt notation)

In the following equation (10), we report the elastic tensors expressed in Voigt notation for the tetragonal class of symmetry ³ , where only the parameters involved under the plane strain hypothesis are explicitly presented. Thereby, the symbol ☆ indicates that the specific entry do not intervene under the plane strain hypothesis. ⁴ The class of symmetry has been chosen accordingly with the symmetry of the unit cell presented in Figure 1, under the assumption that the same class of symmetry applies both at the micro- and the macro-scale.

$$\begin{array}{c}{\mathbb{C}}_{\mathrm{e}}=\left(\begin{array}{ccccc}{\kappa }_{\mathrm{e}}+{\mu }_{\mathrm{e}}\hfill & {\kappa }_{\mathrm{e}}-{\mu }_{\mathrm{e}}\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{{\kappa }_{\mathrm{e}}-{\mu }_{\mathrm{e}}}& {\kappa }_{\mathrm{e}}+{\mu }_{\mathrm{e}}\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & ⋮\hfill & \ddots \hfill & \hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & \hfill & {\mu }_{\mathrm{e}}^{*}\hfill \end{array}\right),{\mathbb{C}}_{\mathrm{micro}}=\left(\begin{array}{ccccc}{\kappa }_{\mathrm{m}}+{\mu }_{\mathrm{m}}\hfill & {\kappa }_{\mathrm{m}}-{\mu }_{\mathrm{m}}\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{{\kappa }_{\mathrm{m}}-{\mu }_{\mathrm{m}}}& {\kappa }_{\mathrm{m}}+{\mu }_{\mathrm{m}}\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & ⋮\hfill & \ddots \hfill & \hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & \hfill & {\mu }_{\mathrm{m}}^{*}\hfill \end{array}\right),\hfill \\ {\mathbb{J}}_{\mathrm{m}}=\rho L_{\mathrm{c}}^2\left(\begin{array}{ccccc}{\kappa }_{\gamma }+{\gamma }_1\hfill & {\kappa }_{\gamma }-{\gamma }_1\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{{\kappa }_{\gamma }-{\gamma }_1}& {\kappa }_{\gamma }+{\gamma }_1\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & ⋮\hfill & \ddots \hfill & \hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & \hfill & {\gamma }_1^{*}\hfill \end{array}\right),{\mathbb{T}}_{\mathrm{e}}=\rho L_{\mathrm{c}}^2\left(\begin{array}{ccccc}{\overline{\kappa }}_{\gamma }+{\overline{\gamma }}_1\hfill & {\overline{\kappa }}_{\gamma }-{\overline{\gamma }}_1\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{{\overline{\kappa }}_{\gamma }-{\overline{\gamma }}_1}& {\overline{\kappa }}_{\gamma }+{\overline{\gamma }}_1\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & ⋮\hfill & \ddots \hfill & \hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & \hfill & {\overline{\gamma }}_1^{*}\hfill \end{array}\right),\hfill \\ {\mathbb{L}}_{\mathrm{s}}=L_{\mathrm{c}}^2\left(\begin{array}{ccccc}★\hfill & ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & ⋮\hfill & \begin{array}{cc}{\alpha }_1\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& {\alpha }_1\hfill \end{array}\hfill & \begin{array}{c}0\hfill \\ \multicolumn{1}{c}{0}\end{array}\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & \begin{array}{cc}0\hfill & 0\hfill \end{array}\hfill & ★\hfill \end{array}\right),{\mathbb{M}}_{\mathrm{s}}=\rho L_{\mathrm{c}}^4\left(\begin{array}{ccccc}★\hfill & ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{★}& ★\hfill & ★\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & ⋮\hfill & \begin{array}{cc}{\beta }_1\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& {\beta }_1\hfill \end{array}\hfill & \begin{array}{c}0\hfill \\ \multicolumn{1}{c}{0}\end{array}\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & \begin{array}{cc}0\hfill & 0\hfill \end{array}\hfill & ★\hfill \end{array}\right),\hfill \\ {\mathbb{C}}_{\mathrm{c}}=\left(\begin{array}{ccc}★\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& ★\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 4{\mu }_{\mathrm{c}}\hfill \end{array}\right),{\mathbb{J}}_{\mathrm{c}}=\rho L_{\mathrm{c}}^2\left(\begin{array}{ccc}★\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& ★\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 4{\gamma }_2\hfill \end{array}\right),{\mathbb{T}}_{\mathrm{c}}=\rho L_{\mathrm{c}}^2\left(\begin{array}{ccc}★\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& ★\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 4{\overline{\gamma }}_2\hfill \end{array}\right),\hfill \\ {\mathbb{L}}_{\mathrm{a}}=L_{\mathrm{c}}^2\left(\begin{array}{ccc}4{\alpha }_2\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 4{\alpha }_2\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & ★\hfill \end{array}\right),{\mathbb{M}}_{\mathrm{a}}=\rho L_{\mathrm{c}}^4\left(\begin{array}{ccc}4{\beta }_2\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 4{\beta }_2\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & ★\hfill \end{array}\right).\hfill \end{array}$$

(10)

It is emphasized that the matrix representation of ${\mathbb{L}}_{\mathrm{s}}$ , ${\mathbb{L}}_{\mathrm{a}}$ , ${\mathbb{M}}_{\mathrm{s}}$ , and ${\mathbb{M}}_{\mathrm{a}}$ differs from the others since the non-zero terms in $(\mathrm{Curl}P{)}_{\mathrm{ij}}$ and in $(\mathrm{Curl}\stackrel{·}{P}{)}_{\mathrm{ij}}$ are the “out-of-plane curvatures” related to the indexes $i=1,2$ and $j=3$ . To better explain this fact, we consider the shape of $\mathrm{Curl}m$ given by

$$\begin{array}{c}P=\left(\begin{array}{ccc}•\hfill & •\hfill & 0\hfill \\ \multicolumn{1}{c}{•}& •\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill \end{array}\right),\mathrm{Curl}P=\left(\begin{array}{ccc}0\hfill & 0\hfill & •\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & •\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill \end{array}\right),\hfill \\ m=({\mathbb{L}}_{\mathrm{s}}+{\mathbb{L}}_{\mathrm{a}})\mathrm{Curl}P=\left(\begin{array}{ccc}0\hfill & 0\hfill & •\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & •\hfill \\ \multicolumn{1}{c}{•}& •\hfill & 0\hfill \end{array}\right),\mathrm{Curl}m=\left(\begin{array}{ccc}•\hfill & •\hfill & 0\hfill \\ \multicolumn{1}{c}{•}& •\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & •\hfill \end{array}\right).\hfill \end{array}$$

(11)

Note that the structure of $\mathrm{Curl}\hat{m}=\mathrm{Curl}(({\mathbb{M}}_{\mathrm{s}}+{\mathbb{M}}_{\mathrm{a}})\mathrm{Curl}\stackrel{··}{P})$ is identical.

As it can be seen from equation (11)₄, even if the sym/skew decomposition is enforced in the constitutive law (we recall that ${\mathbb{L}}_{\mathrm{s}}$ , ${\mathbb{L}}_{\mathrm{a}}$ , ${\mathbb{M}}_{\mathrm{s}}$ , and ${\mathbb{M}}_{\mathrm{a}}$ have minor symmetries, see Appendix 1), it is not guaranteed that the term $\mathrm{Curl}m$ appearing in the equilibrium equations (6) will retain only in-plane components for the most general constitutive tensor belonging to the tetragonal symmetry class. In order to avoid this out-of-plane contribution, we must impose one of the following conditions:

$${\alpha }_2={\alpha }_1\wedge {\beta }_2={\beta }_1\mathrm{or}{P}_{22}=0.$$

(12)

Since the condition (12)₂ is not desirable because it prevents the presence of pressure waves, the condition (12)₁ will be adopted. This choice, while reducing the number of parameters in the model, does not reduce the number of the independent ones for $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ as shown in Appendix 1.

Furthermore, we can always compute the parameters ⁵ of the meso-scale depending on the micro-parameters and a new set of macro-parameters [51,54,57]:

$${\mathbb{C}}_{\mathrm{e}}={\mathbb{C}}_{\mathrm{micro}}({\mathbb{C}}_{\mathrm{micro}}-{\mathbb{C}}_{\mathrm{macro}}{)}^{-1}{\mathbb{C}}_{\mathrm{macro}},$$

(13)

where ${\mathbb{C}}_{\mathrm{macro}}$ has the same structure of ${\mathbb{C}}_{\mathrm{micro}}$ and ${\mathbb{C}}_{\mathrm{e}}$ . For the tetragonal class of symmetry in 2D, relation (13) particularizes to

$${\mu }_{\mathrm{e}}=\frac{{\mu }_{\mathrm{m}}{\mu }_{\mathrm{M}}}{{\mu }_{\mathrm{m}}-{\mu }_{\mathrm{M}}},{\kappa }_{\mathrm{e}}=\frac{{\kappa }_{\mathrm{m}}{\kappa }_{\mathrm{M}}}{{\kappa }_{\mathrm{m}}-{\kappa }_{\mathrm{M}}},{\mu }_{\mathrm{e}}^{*}=\frac{{\mu }_{\mathrm{m}}^{*}{\mu }_{\mathrm{M}}^{*}}{{\mu }_{\mathrm{m}}^{*}-{\mu }_{\mathrm{M}}^{*}},$$

(14)

In the limit of infinitesimal small unit cells or rather of an indefinitely large body, these macro-parameters are obtained as the homogenization to a classical Cauchy material [51,57] and can directly be obtained by the slope of the acoustic curves at the origin, cf. section 5.2. Throughout this paper, all material parameters introduced will be positive to guarantee the positive definiteness of their corresponding tensors.

3. Dispersion curves

We assume a plane strain ⁶ time-harmonic ansatz for the displacement $u$ and the micro-distortion tensor $P$

$$v_j={\mathrm{\Psi }}_j{e}^{i(k_1{x}_1+k_2{x}_2-\omega t)},$$

(15)

where $v_j$ represents the generic component of $u$ or $P$ , ${\mathrm{\Psi }}_j$ is a scalar amplitude, $(k_1,k_2{)}^T=k(\sin \varphi ,\cos \varphi {)}^T$ are the wave vector components with $\varphi$ as the angle giving the direction of propagation, $k$ the wave vector length, and $\omega$ is the frequency.

Substituting the ansatz (15) in the equilibrium equations (6), we obtain the homogeneous algebraic linear system

$$A\mathrm{\Psi }=0$$

(16)

where $A=A(\omega ,k,\varphi )\in {\mathbb{C}}^{6\times 6}$ is the acoustic tensor which depends on the frequency $\omega$ , the wave vector length $k$ , the angle of propagation $\varphi$ , all the constitutive parameters in equation (10), and $\mathrm{\Psi }\in {\mathbb{R}}^6$ is the vector of amplitudes. ⁷ The non-trivial solutions of the system (16) are obtained when $A$ is singular, i.e. when $detA=0$ , which provides relations between $k=k(\omega ,\varphi )$ (or $\omega =\omega (k,\varphi )$ ), the so-called dispersion relations.

The acoustic tensor can be written as follows:

$$A=\left(\begin{array}{cc}{B}_1\hfill & C_1\hfill \\ \multicolumn{1}{c}{C_2}& B_2\hfill \end{array}\right),$$

(17)

where $(B_1,B_2,C_1,C_2)\in {\mathbb{C}}^{3\times 3}$ .

Since $detA=det(B_1)det(B_2-C_2{B}_1^{-1}{C}_1)$ [58], it is clear that the determinant of $A$ becomes the product of the two independent factors $det{B}_1$ and $det{B}_2$ if either $C_1=0$ or $C_2=0$ . If the reference system is chosen to be aligned with the direction of the wave vector (i.e. $\varphi =\theta$ , with $\theta$ the angle of rotation of the reference system), the sub-matrices $C_1$ and $C_2$ expressions are as follows:

$$C_1=\frac{\sin (4\theta )}{2}\left(\begin{array}{ccc}{k}^2(q_{\mathrm{e}}-L_{\mathrm{c}}^2\rho {\overline{q}}_{\gamma }{\omega }^2)\hfill & \mathrm{ik}{q}_{\mathrm{e}}\hfill & \mathrm{ik}{q}_{\mathrm{e}}\hfill \\ \multicolumn{1}{c}{-\mathrm{ik}{q}_{\mathrm{e}}}& -\rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2+q_{\mathrm{e}}+q_{\mathrm{m}}\hfill & -\rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2+q_{\mathrm{e}}+q_{\mathrm{m}}\hfill \\ \multicolumn{1}{c}{\mathrm{ik}{q}_{\mathrm{e}}}& \rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2-q_{\mathrm{e}}-q_{\mathrm{m}}\hfill & \rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2-q_{\mathrm{e}}-q_{\mathrm{m}}\hfill \end{array}\right),$$

(18)

$$C_2=\frac{\sin (4\theta )}{2}\left(\begin{array}{ccc}{k}^2(q_{\mathrm{e}}-L_{\mathrm{c}}^2\rho {\overline{q}}_{\gamma }{\omega }^2)\hfill & \mathrm{ik}{q}_{\mathrm{e}}\hfill & -\mathrm{ik}{q}_{\mathrm{e}}\hfill \\ \multicolumn{1}{c}{-\mathrm{ik}{q}_{\mathrm{e}}}& -\rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2+q_{\mathrm{e}}+q_{\mathrm{m}}\hfill & \rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2-q_{\mathrm{e}}-q_{\mathrm{m}}\hfill \\ \multicolumn{1}{c}{-\mathrm{ik}{q}_{\mathrm{e}}}& -\rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2+q_{\mathrm{e}}+q_{\mathrm{m}}\hfill & \rho L_{\mathrm{c}}^2{q}_{\gamma }{\omega }^2-q_{\mathrm{e}}-q_{\mathrm{m}}\hfill \end{array}\right),$$

where $q_{\mathrm{e}}={\mu }_{\mathrm{e}}-{\mu }_{\mathrm{e}}^{*},q_{\mathrm{m}}={\mu }_{\mathrm{m}}-{\mu }_{\mathrm{m}}^{*},q_{\gamma }={\gamma }_1-{\gamma }_1^{*},\mathrm{and}{\overline{q}}_{\gamma }={\overline{\gamma }}_1-{\overline{\gamma }}_1^{*}$ . It is highlighted that, since the constitutive law for the curvature terms $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ depends on just one parameter, they are isotropic, and this makes the matrices $C_1$ and $C_2$ independent with respect ${\alpha }_i$ and ${\beta }_i$ , while the expressions of the latters in $B_1$ and $B_2$ are not affected by the rotation of the reference system.

From equation (18) it is possible to deduce that the condition for which $C_1=C_2=0$ (one of the two would be already enough) is $\theta =(\pi n)/4$ with $n\in \mathbb{N}$ . This means that, when the reference system is aligned with the direction of wave propagation, and both are aligned with a symmetry axes of the material, the determinant $detA=det(B_1)det(B_2)$ is the product of two independent factors.

These two independent factors can be associated with pure-pressure waves $det(B_1)$ and pure-shear waves $det(B_2)$ .

In particular, this allows us to reduce the order of the dispersion polynomial by

$$detA=p(k^2,k^4,k^6,k^8,{\omega }^2,{\omega }^4,{\omega }^6,{\omega }^8,{\omega }^{10},{\omega }^{12})=p_1(k^2,k^4,{\omega }^2,{\omega }^4,{\omega }^6)p_2(k^2,k^4,{\omega }^2,{\omega }^4,{\omega }^6),$$

(19)

where $p_1$ and $p_2$ are of the form

$$p_i=c_0{k}^2+c_1{k}^4-(c_2+c_3{k}^2+c_4{k}^4){\omega }^2+(c_5+c_6{k}^2+c_7{k}^4){\omega }^4-(c_8+c_9{k}^2+c_{10}{k}^4){\omega }^6,$$

(20)

easing the calculation of the roots of these polynomials significantly, i.e. the expressions of the dispersion curves. It is noted that we can always consider $k(\omega )$ , instead of $\omega (k)$ , for shorter analytical expression as we must only solve a quadratic equation (in $k^2$ ) and not a third-order polynomial (in ${\omega }^2$ ) but use $\omega (k)$ for easier plotting and its natural split in three distinct expressions for each dispersion curve.

4. New considerations on the relaxed micromorphic parameters

In this section, we draw some useful considerations about the consistency of the relaxed micromorphic model with respect to a change of unit cell’s size and of the material properties of the base material. The model’s consistency is checked against a standard Bloch–Floquet analysis of the wave propagation performed using the unit cell described in section 1.1 with built in periodic Bloch–Floquet boundary conditions from Comsol Multiphysics®.

The following two connections between the properties of the unit cell and the behavior of the dispersion curves can be drawn:

• The dispersion curves scale proportionally in $\omega$ with respect to the speed of the wave of the bulk material composing the unit cell;

• The dispersion curves scale inversely in both $\omega$ and $k$ with respect to the size of the unit cell.

Both results are useful to avoid repeating the time-consuming fitting procedure when changing the size of the cell and the base material’s properties while keeping the unit cell’s geometry unchanged.

4.1. Consistency of the relaxed micromorphic model with respect to a change in the unit cell’s bulk material properties

The dispersive properties of a microstructured isotropic Cauchy material depend exclusively and linearly on the wave speeds of the bulk material once its geometry is fixed. Indeed, as it is well known, the dispersion relatives can be always written as $\omega =k{c}_i$ , $i=\{p,s\}$ , where $c_{\mathrm{p}}$ and $c_{\mathrm{s}}$ are the pressure and shear wave speeds, respectively, defined as follows:

$$c_{\mathrm{p}}:=\sqrt{\frac{{\kappa }_{\mathrm{M}}+{\mu }_{\mathrm{M}}}{\rho }},c_{\mathrm{s}}:=\sqrt{\frac{{\mu }_{\mathrm{M}}}{\rho }},k=\frac{\omega }{c_i},i=\{p,s\}.$$

(21)

This implies that by scaling the elastic coefficients by a constant $a$ and the density by another constant $b$ , we have

$${\tilde{c}}_{\mathrm{p}}:=c_{\mathrm{p}}\sqrt{\frac{a}{b}},{\tilde{c}}_{\mathrm{s}}:=c_{\mathrm{s}}\sqrt{\frac{a}{b}},k=\sqrt{\frac{b}{a}}\frac{\omega }{c_i},i=\{p,s\}.$$

(22)

Therefore, the response of an effective model should also change accordingly. Thus, we observe that by multiplying all the relaxed micromorphic elastic coefficients $({\mathbb{C}}_{\mathrm{e}},{\mathbb{C}}_{\mathrm{c}},{\mathbb{C}}_{\mathrm{micro}},{\mathbb{L}}_{\mathrm{s}},{\mathbb{L}}_{\mathrm{a}})$ by a constant $a$ and the apparent density $\rho$ by another constant $b$ , we can rewrite equation (20) as ⁸

$$\begin{array}{c}{a}^3(c_0{k}^2+c_1{k}^4)-a^{2}b(c_2+c_3{k}^2+c_4{k}^4){\omega }^2+a{b}^2(c_5+c_6{k}^2+c_7{k}^4){\omega }^4\hfill \\ -b^3(c_8+c_9{k}^2+c_{10}{k}^4){\omega }^6=0.\hfill \end{array}$$

(23)

By collecting $a^3$ we obtain

$$\begin{array}{c}{a}^3[(c_0{k}^2+c_1{k}^4)-(c_2+c_3{k}^2+c_4{k}^4){\left(\omega \sqrt{\frac{b}{a}}\right)}^2+(c_5+c_6{k}^2+c_7{k}^4){\left(\omega \sqrt{\frac{b}{a}}\right)}^4\hfill \\ -(c_8+c_9{k}^2+c_{10}{k}^4){\left(\omega \sqrt{\frac{b}{a}}\right)}^6]=0,\hfill \end{array}$$

(24)

from which we can introduce a scaled frequency $\tilde{\omega }=\omega \sqrt{\frac{b}{a}}$ . Increasing the stiffness (a) or decreasing the density (b) of the base material will cause an overall shifting of the dispersion curves toward lower frequencies.

This is consistent with what is observed looking at the dispersion properties of a microstructured material’s unit cell obtained via a Bloch–Floquet analysis. In the case for which $a=b$ , the roots of equation (24) do not change at all. Thanks to this identification, we can now easily change the material constituting the unit cell (without changing the geometry) by scaling the material parameters accordingly without repeating the whole fitting process. In particular, all the cut-offs and asymptotes will be scaled by a quantity $\sqrt{\frac{a}{b}}$ .

We explicitly remark again that scaling the macroscopic apparent density $\rho$ of the unit cell by a factor $b>0$ will change the frequency $\omega$ by the factor $\frac{1}{\sqrt{b}}$ , i.e. the frequency is inversely proportional to the square root of the density of the unit cell. The wavenumber $k$ is invariant under changing the macroscopic apparent density since the periodicity of the unit cell remains unaltered.

4.2. Consistency of the relaxed micromorphic model with respect to a change in the unit cell’s size

While keeping the geometry and the material unaltered, the dispersion properties of a microstructured isotropic Cauchy material are inversely proportional to the size of its unit cell, meaning that halving the size of the unit cell will double the frequency response for each value of the length $k$ of the wave vector, which also changes with the same inverse proportionality since it represents the spatial periodicity of the structure. This can be easily retrieved by performing standard Bloch–Floquet analysis.

In order to obtain this behavior with the relaxed micromorphic model, we must scale the elastic curvature tensors ( ${\mathbb{L}}_{\mathrm{s}},{\mathbb{L}}_{\mathrm{a}}$ ) and all the micro-inertia tensors ( ${\mathbb{J}}_{\mathrm{m}},{\mathbb{J}}_{\mathrm{c}}$ , ${\mathbb{T}}_{\mathrm{e}},{\mathbb{T}}_{\mathrm{c}}$ ) by the square of the size of the unit cell, and the micro-inertia curvature tensors ( ${\mathbb{M}}_{\mathrm{s}},{\mathbb{M}}_{\mathrm{a}}$ ) by the fourth power of the size of the unit cell.

To prove this, we consider a change in the size of the unit cell by some arbitrary factor $t>0$ , which requires the scaling of the characteristic length $L_{\mathrm{c}}$ (which is now considered to be equal to the size of the unit cell) by the same factor $t$ . Assuming that all the other material parameters used remain constant, in equation (20) we substitute $L_{\mathrm{c}}\to t{L}_{\mathrm{c}}$ (see the coefficients in Appendices 2–4)

$$\begin{array}{c}{c}_0{k}^2+t^2{c}_1{k}^4-c_2{\omega }^2+t^2{c}_3{k}^2{\omega }^2+t^4{c}_4{k}^4{\omega }^2+t^2{c}_5{\omega }^4\hfill \\ +t^4{c}_6{k}^2{\omega }^4+t^6{c}_7{k}^4{\omega }^4-t^4{c}_8{\omega }^6+t^6{c}_9{k}^2{\omega }^6+t^8{c}_{10}{k}^4{\omega }^6=0.\hfill \end{array}$$

(25)

We can now collect $\frac{1}{t^2}$ in equation (25) arriving at

$$\begin{array}{c}\frac{1}{t^2}[c_0{(\mathrm{tk})}^2+c_1{(\mathrm{tk})}^4-c_2{(t\omega )}^2+c_3{(\mathrm{tk})}^2{(t\omega )}^2+c_4{(\mathrm{tk})}^4{(t\omega )}^2+c_5{(t\omega )}^4\hfill \\ +c_6{(\mathrm{tk})}^2{(t\omega )}^4+c_7{(\mathrm{tk})}^4{(t\omega )}^4-c_8{(t\omega )}^6+c_9{(\mathrm{tk})}^2{(t\omega )}^6+c_{10}{(\mathrm{tk})}^4{(t\omega )}^6]=0.\hfill \end{array}$$

(26)

It is clear to see from the comparison between equations (25) and (26) that their roots are simply linearly scaled with respect the factor $t$ . This result enables the choice of having $L_{\mathrm{c}}$ equal to the size of the unit cell and thus allows us to change the latter in the microstructured material considered without needing to repeat the whole fitting procedure.

In particular, scaling the size $L_{\mathrm{c}}$ of the unit cell by the factor $t>0$ will change the frequency $\omega$ and wavenumber $k$ by the factor $\frac{1}{t}$ , i.e. the frequency and wavenumber are reverse proportional to the size of the unit cell. This simple observation allows to perform the fitting procedure for the relaxed micromorphic model only once for each geometry of the unit cell: changing the size of the unit cell will result in an automatic fitting when suitably rescaling $\omega$ and $k$ where none of the material parameters (except $L_{\mathrm{c}}$ ) must be changed. Note that the slopes at the origin of the dispersion curves do not change when changing the size of the unit cell while keeping the geometry fixed.

4.3. Relaxed micromorphic cut-offs

The cut-offs of the dispersion curves play an important role in fitting the material parameters of the relaxed micromorphic model [51, 54, 57]. For the convenience of the reader, we show the calculations of the analytic expressions again. In the case $k=0$ , the dispersion relation (20) simplifies into

$$\begin{array}{c}-c_2{\omega }^2+c_5{\omega }^4-c_8{\omega }^6=0\iff {\omega }^2({\omega }^4-\frac{c_5}{c_8}{\omega }^2+\frac{c_2}{c_8})=0\hfill \\ \iff {\omega }_1^2=0,{\omega }_{2,3}^2=\frac{c_5\pm \sqrt{c_5^2-4c_2{c}_8}}{2c_8}.\hfill \end{array}$$

(27)

The coefficients $c_2,c_5,c_8$ depend on the elastic parameters ${\mu }_{\mathrm{m}},{\kappa }_{\mathrm{m}},{\mu }_{\mathrm{m}}^{*},{\mu }_{\mathrm{c}},{\mu }_{\mathrm{e}},{\kappa }_{\mathrm{e}},{\mu }_{\mathrm{e}}^{*}$ ; the micro-inertia parameters ${\kappa }_{\gamma },{\gamma }_1,{\gamma }_1^{*},{\gamma }_2$ ; the macroscopic apparent density $\rho$ ; and characteristic length $L_{\mathrm{c}}$ but are independent of the parameters for $\mathrm{Curl}P$ , $\nabla u$ , and $\mathrm{Curl}\stackrel{·}{P}$ , cf. Appendix 2.

Equations (27) can be simplified as in Table 2 with

$${\omega }_{\mathrm{r}}=\sqrt{\frac{{\mu }_{\mathrm{c}}}{\rho L_{\mathrm{c}}^2{\gamma }_2}},{\omega }_{\mathrm{s}}=\sqrt{\frac{{\mu }_{\mathrm{e}}+{\mu }_{\mathrm{m}}}{\rho L_{\mathrm{c}}^2{\gamma }_1}},{\omega }_{\mathrm{ss}}=\sqrt{\frac{{\mu }_{\mathrm{e}}^{*}+{\mu }_{\mathrm{m}}^{*}}{\rho L_{\mathrm{c}}^2{\gamma }_1^{*}}},{\omega }_{\mathrm{p}}=\sqrt{\frac{{\kappa }_{\mathrm{e}}+{\kappa }_{\mathrm{m}}}{\rho L_{\mathrm{c}}^2{\kappa }_{\gamma }}}.$$

(28)

Table 2.

Cut-off expressions for the pressure waves (left) and for the shear waves (right).

	$0^{°}$	${45}^{°}$		$0^{°}$	${45}^{°}$
Shear	${\omega }_2={\omega }_{\mathrm{ss}}$	${\omega }_2={\omega }_{\mathrm{s}}$	Pressure	${\omega }_2={\omega }_{\mathrm{s}}$	${\omega }_2={\omega }_{\mathrm{ss}}$
	${\omega }_3={\omega }_{\mathrm{r}}$	${\omega }_3={\omega }_{\mathrm{r}}$		${\omega }_3={\omega }_{\mathrm{p}}$	${\omega }_3={\omega }_{\mathrm{p}}$

We recognize that the expressions ${\omega }_{\mathrm{s}}$ and ${\omega }_{\mathrm{ss}}$ change from pressure to shear and shear to pressure, respectively, when going from $0^{°}$ to ${45}^{°}$ of incidence. Since the dispersion curves of the unit cell have two cut-offs that coincide, we chose them to be ${\omega }_{\mathrm{s}}$ and ${\omega }_{\mathrm{ss}}$ . Therefore, we introduce the following relation:

$${\omega }_{\mathrm{s}}={\omega }_{\mathrm{ss}}\iff \frac{{\mu }_{\mathrm{e}}+{\mu }_{\mathrm{m}}}{{\gamma }_1}=\frac{{\mu }_{\mathrm{e}}^{*}+{\mu }_{\mathrm{m}}^{*}}{{\gamma }_1^{*}}.$$

(29)

The values of these cut-offs have been fixed according to Comsol Multiphysics® simulations as

$${\omega }_{\mathrm{r}}=554.61\mathrm{Hz},{\omega }_{\mathrm{s}}={\omega }_{\mathrm{ss}}=2011.83\mathrm{Hz},{\omega }_{\mathrm{p}}=2048.19\mathrm{Hz},$$

(30)

and the values of last points from Comsol Multiphysics® are used to fix the asymptotes, cf. Table 3.

Table 3.

Numerical values of the asymptotes via Bloch–Floquet analysis using Comsol Multiphysics®.

	$0^{°}$	${45}^{°}$		$0^{°}$	${45}^{°}$
Shear	340.42 Hz	606.31 Hz	Pressure	581.95 Hz	689.20 Hz
	637.79 Hz	689.21 Hz		2040.81 Hz	1883.97 Hz
	2078.39 Hz	2375.78 Hz		2156.21 Hz	2376.45 Hz

5. Fitting of the relaxed micromorphic parameters: the particular case of vanishing curvature (without Curl $P$ and Curl $\stackrel{·}{P}$ )

In the numerical applications considered in this work, we start without considering the tensors ${\mathbb{L}}_{\mathrm{s}},{\mathbb{L}}_{\mathrm{c}}$ and ${\mathbb{M}}_{\mathrm{s}},{\mathbb{M}}_{\mathrm{c}}$ , i.e. we neglect the effect of $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ on the dynamic regime. This fundamentally changes the shape of the analytical expression (20) resulting in a new polynomial $p(k^2,{\omega }^2,{\omega }^4,{\omega }^6)$ with reduced order of $k$ . Assuming that the wavenumber is always positive, we only have one single expression $k(\omega )$ describing all three dispersion curves:

$$c_0{k}^2-(c_2+c_3^{*}{k}^2){\omega }^2+(c_5+c_6^{*}{k}^2){\omega }^4-(c_8+c_9^{*}{k}^2){\omega }^6=0,k=\omega \sqrt{\frac{c_2-c_5{\omega }^2+c_8{\omega }^4}{c_0-c_3^{*}{\omega }^2+c_6^{*}{\omega }^4-c_9^{*}{\omega }^6}}.$$

(31)

The coefficients $c_0,c_2,c_3^{*},c_5,c_6^{*},c_8,c_9^{*}$ used ⁹ have a different expression at $0^{°}$ and ${45}^{°}$ angle of incidence as well as for pressure and shear waves and are included in Appendix 2. The structure of $k(\omega )$ in equation (31)₂ can be also reported as a non-linear classic dispersion relation $k=\frac{\omega }{c_{\mathrm{p}/\mathrm{s}}(\omega )}$ with a frequency dependent group velocity $c_{\mathrm{p}/\mathrm{s}}$ .

5.1. Asymptotes

Instead of fitting dispersion curves pointwise, we focus on using the analytical expression of the cut-offs ( $k=0$ ) and of the asymptotes ( $k\to \infty$ ). The explicit expression of the cut-offs is already discussed in section 4.3. We use a similar approach to calculate the asymptotes as well by considering the limit $k\to \infty$ where only the terms with the highest order of $k$ are important. Thus, we arrive at

$$c_0-c_3^{*}{\omega }^2+c_6^{*}{\omega }^4-c_9^{*}{\omega }^6=0\iff {\omega }^6-\frac{c_6^{*}}{c_9^{*}}{\omega }^4+\frac{c_3^{*}}{c_9^{*}}{\omega }^2-\frac{c_0}{c_9}=0.$$

(32)

In contrast to the analytical expression of the cut-offs, we are not able to simplify the asymptotes’ expression in a feasible way. This is mainly due to the fact that we must solve a third-order polynomial while we only had two non-zero cut-offs each before.

Since the dispersion curves of the unit cell obtained via Bloch–Floquet analysis are by nature periodic, the limit for $k\to \infty$ is per se meaningless when considering the Bloch–Floquet approach. The value for $k=\frac{1}{2L_{\mathrm{c}}}$ (with $L_{\mathrm{c}}$ the size of the unit cell) is the periodicity limit. On the other hand, this limit has, of course, meaning for a continuum model like the relaxed micromorphic model. To reconcile these two limits in the fitting procedure, we will impose that the limit $k\to \infty$ for our continuum model will coincide with the periodicity limit of the Bloch–Floquet curves $k=\frac{1}{2L_{\mathrm{c}}}$ . This strategy allows us to preserve the width of the band gap, cf. Figure 2.

Figure 2.

The fitting of the dispersion curves (black line) includes information about the slope at zero and the asymptotes (blue lines), while the remaining shape comes automatically. The Bloch–Floquet analysis used in Comsol Multiphysics® can only evaluate up to the periodicity limit $k=\frac{1}{2L_{\mathrm{c}}}$ (black dots) which causes a slight gap between the dispersion curves at $k=\frac{1}{2L_{\mathrm{c}}}$ and their corresponding asymptotes (red part).

5.2. Fitting

We start the fitting with the macroscopic apparent density $\rho$ and values of the macro parameters ${\kappa }_{\mathrm{M}},{\mu }_{\mathrm{M}},{\mu }_{\mathrm{M}}^{*}$ , i.e. the material constants necessary for the classical homogenization of an infinite large micromorphic material. For an anisotropic Cauchy material, the speed of the acoustic waves is

$$\begin{array}{c}{c}_{\mathrm{p}}=\sqrt{\frac{{\kappa }_{\mathrm{M}}+{\mu }_{\mathrm{M}}}{\rho }},c_{\mathrm{s}}=\sqrt{\frac{{\mu }_{\mathrm{M}}^{*}}{\rho }},\hfill \\ {\overline{c}}_{\mathrm{p}}=\sqrt{\frac{{\kappa }_{\mathrm{M}}+{\mu }_{\mathrm{M}}^{*}}{\rho }},{\overline{c}}_{\mathrm{s}}=\sqrt{\frac{{\mu }_{\mathrm{M}}}{\rho }},\hfill \end{array}$$

(33)

where $c_p,c_s$ are the speed of pressure and shear wave, respectively, for $0^{°}$ of incidence, while ${\overline{c}}_p,{\overline{c}}_s$ describe an incidence angle of ${45}^{°}$ . For the tetragonal class of symmetry we choose, it holds

$$c_{\mathrm{p}}^2+c_{\mathrm{s}}^2=\frac{{\kappa }_{\mathrm{M}}+{\mu }_{\mathrm{M}}+{\mu }_{\mathrm{M}}^{*}}{\rho }={\overline{c}}_{\mathrm{p}}^2+{\overline{c}}_{\mathrm{s}}^2,$$

(34)

reducing the system of equations to just three independent quantities. For the relaxed micromorphic model, we fit these macro parameters by the slope of the corresponding acoustic dispersion curves at $k=0$

$${\mu }_{\mathrm{M}}={\overline{c}}_{\mathrm{s}}^2\rho ,{\mu }_{\mathrm{M}}^{*}=c_{\mathrm{s}}^2\rho ,{\kappa }_{\mathrm{M}}=(c_{\mathrm{p}}^2-{\overline{c}}_{\mathrm{s}}^2)\rho =({\overline{c}}_{\mathrm{p}}^2-c_{\mathrm{s}}^2)\rho .$$

(35)

The remaining unknown density $\rho$ is given by the material and geometry of the cell we used (see Figure 1) and is directly computed as $\rho =\overline{\rho }\frac{A_{\mathrm{total}}-A_{\mathrm{voids}}}{A_{\mathrm{total}}}$ with $\overline{\rho }$ as the density of polyethylene and $A_{\mathrm{total}}$ , $A_{\mathrm{voids}}$ being the area of the whole cell and its voids, respectively. We list the numerical values in Table 4.

Table 4.

Numerical values of the macroscopic apparent density $\rho$ and the elastic macro-parameters ${\kappa }_{\mathrm{M}},{\mu }_{\mathrm{M}},{\mu }_{\mathrm{M}}^{*}$ for our unit cell when considering polyethylene as base material.

$L_{\mathrm{c}}$	$\rho$	${\kappa }_{\mathrm{M}}$	${\mu }_{\mathrm{M}}$	${\mu }_{\mathrm{M}}^{*}$
(mm)	(kg/m 3)	(kPa)	(kPa)	(kPa)
20	361.2	352.2	200	99.36

As a second step, we use the analytical expressions for the cut-offs (28) and calculate

$$\begin{array}{c}{\kappa }_{\gamma }=\frac{{\kappa }_{\mathrm{e}}+{\kappa }_{\mathrm{m}}}{\rho L_{\mathrm{c}}^2{\omega }_p^2},{\gamma }_1=\frac{{\mu }_{\mathrm{e}}+{\mu }_{\mathrm{m}}}{\rho L_{\mathrm{c}}^2{\omega }_s^2},\hfill \\ {\gamma }_1^{*}=\frac{{\mu }_{\mathrm{e}}^{*}+{\mu }_{\mathrm{m}}^{*}}{\rho L_{\mathrm{c}}^2{\omega }_{\mathrm{ss}}^2},{\gamma }_2=\frac{{\mu }_{\mathrm{c}}}{\rho L_{\mathrm{c}}^2{\omega }_r^2}.\hfill \end{array}$$

(36)

Thus, we can reduce the system of independent variables by the four inertia parameters ${\kappa }_{\gamma },{\gamma }_1,{\gamma }_1^{*},{\gamma }_2$ . Note again that we implied ${\omega }_{\mathrm{s}}={\omega }_{\mathrm{s}s}$ , see section 4.3.

As the third step, we use the analytic expression of the asymptotes for $0^{°}$ and ${45}^{°}$ resulting in 12 expressions in total which is more than the number of independent parameters remaining. Thus, we fit the last eight remaining parameters, namely the micro parameters ${\kappa }_{\mathrm{m}},{\mu }_{\mathrm{m}},{\mu }_{\mathrm{m}}^{*},{\mu }_{\mathrm{c}}$ and inertia parameters ${\overline{\kappa }}_{\gamma },{\overline{\gamma }}_1,{\overline{\gamma }}_1^{*},{\overline{\gamma }}_2$ , numerically, by minimizing the square error of the analytical expressions (32) and their corresponding numerical values from Comsol Multiphysics® (for $k=\frac{1}{2L_{\mathrm{c}}}$ for $0^{°}$ and $k=\frac{1}{\sqrt{2}{L}_{\mathrm{c}}}$ for ${45}^{°}$ ). Thereby, we utilize the fact that each group of asymptotes only depend on four independent parameters to speed up the calculations, as summarized in Table 5.

Table 5.

Dependence of the asymptotes of the dispersion curves on the free material parameters as function of the direction of propagation ( $0^{°}$ / ${45}^{°}$ ) and type of wave (shear/pressure).

	Pressure	Shear
$0^{°}$	${\kappa }_{\mathrm{m}},{\overline{\kappa }}_{\gamma },{\mu }_{\mathrm{m}},{\overline{\gamma }}_1$	${\mu }_{\mathrm{c}},{\overline{\gamma }}_2,{\mu }_{\mathrm{m}}^{*},{\overline{\gamma }}_1^{*}$
${45}^{°}$	${\kappa }_{\mathrm{m}},{\overline{\kappa }}_{\gamma },{\mu }_{\mathrm{m}}^{*},{\overline{\gamma }}_1^{*}$	${\mu }_{\mathrm{c}},{\overline{\gamma }}_2,{\mu }_{\mathrm{m}},{\overline{\gamma }}_1$

We calculate these values with Mathematica using the NMinimize-algorithm with the inbuilt method RandomSearch running in a loop for multiple times. Here, we must start with reasonable initial values for all parameters involved. Note that it is always possible to imply ${\mathbb{C}}_{\mathrm{micro}}=2{\mathbb{C}}_{\mathrm{macro}}$ as a conservative first guess implying

$${\mathbb{C}}_{\mathrm{e}}={\mathbb{C}}_{\mathrm{micro}}({\mathbb{C}}_{\mathrm{micro}}-{\mathbb{C}}_{\mathrm{macro}}{)}^{-1}{\mathbb{C}}_{\mathrm{macro}}=2{\mathbb{C}}_{\mathrm{macro}}(2{\mathbb{C}}_{\mathrm{macro}}-{\mathbb{C}}_{\mathrm{macro}}{)}^{-1}{\mathbb{C}}_{\mathrm{macro}}=2{\mathbb{C}}_{\mathrm{macro}},$$

(37)

i.e. simplifying the starting expressions ¹⁰ by ${\mathbb{C}}_{\mathrm{e}}={\mathbb{C}}_{\mathrm{micro}}$ .

We list the numerical values of all parameters obtained via the fitting procedure presented in this section of the relaxed micromorphic model in Table 6.

Table 6.

Obtained numerical values of the relaxed micromorphic model without curvature fitted for the metamaterial whose unit cell is given in Figure 1.

${\kappa }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}^{*}$	${\mu }_{\mathrm{c}}$	${\kappa }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}^{*}$
(kPa)	(kPa)	(kPa)	(kPa)	(kPa)	(kPa)	(kPa)
546.0	286.9	175.04	15.24	992.8	664.4	229.9

${\overline{\kappa }}_{\gamma }$	${\overline{\gamma }}_1$	${\overline{\gamma }}_1^{*}$	${\overline{\gamma }}_2$	${\kappa }_{\gamma }$	${\gamma }_1$	${\gamma }_1^{*}$	${\gamma }_2$
6.928	0.496	0.523	2.092	2.541	0.407	1.628	0.692

5.3. Discussion

The fitting shown in Figure 3 behaves well for all frequencies $\omega$ and wavenumber $k$ for $0^{°}$ of incidence but looses some precision for an incidence angle of ${45}^{°}$ especially for higher values of $k$ . This calls for a further generalization of the relaxed micromorphic model which will be object of following papers. In any case, the achieved overall precision already allows us to explore the dispersive metamaterial’s characteristics at a satisfactory level.

Figure 3.

Dispersion curves $\omega (k)$ for $0^{°}$ (left) and ${45}^{°}$ (right) with pressure curves colored in yellow and shear in blue. The dots are the points computed with Comsol Multiphysics® while the smooth curves show the analytical expression of the dispersion curves for the relaxed micromorphic model without curvature, i.e. for ${\alpha }_1={\beta }_1=0$ . The values of the curve’s horizontal asymptotes are also shown with dashed lines.

The absence of higher-order terms ( $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ ) caused the reduction to a single expression $k(\omega )$ describing all three dispersion curves in one, cf. equation (31). Thus, for every frequency $\omega$ , there is exactly one wavenumber $k$ which may be imaginary if the term inside the root is negative. For the plots here, we only show $k(\omega )$ where the expression is real-valued and ignore imaginary $k(\omega )$ which arise in the band gap and for higher frequencies. Moreover, we cannot have two distinct wavenumbers with the same frequency which implies that all curves are monotonic. For every group of dispersion curves, e.g. the three pressure waves for ${45}^{°}$ incidence, each individual curve is bounded by the others. Starting with the acoustic curves, their asymptote must be below the cut-off of the lower optic curve of the same type (pressure or shear), while the asymptote of the lower optic curves is bounded from above by the cut-off of the highest optic curve. In particular, self-intersection between two pressure or two shear curves is not possible with this simplified version of the relaxed micromorphic model.

On the other hand, we observe that for the numerical values from Comsol Multiphysics® the asymptote of the acoustic shear wave at ${45}^{°}$ should be slightly higher than the cut-off of the lower optic curve with $581.95$ and $554.61$ , respectively. In addition, assuming that all micro parameters ${\kappa }_{\mathrm{m}},{\mu }_{\mathrm{m}},{\mu }_{\mathrm{m}}^{*}$ (see Table 6) are larger than their corresponding macro counterparts ${\kappa }_{\mathrm{M}},{\mu }_{\mathrm{M}},{\mu }_{\mathrm{M}}^{*}$ , we did not manage to generate decreasing dispersion curves as observed for the lower optic pressure wave for an angle of incidence of ${45}^{°}$ . This effect can instead be achieved when considering mixed space-time derivatives on the micro-distortion tensor $P$ , cf. section 7.

When the relative positions of the curves allow to fit the cut-offs and the asymptotes of each curve separately (e.g. for an incidence angle of $0^{°}$ shown here), the simplified version of the relaxed micromorphic model used does indeed show very good results. We want to emphasize that we only used the limit cases $k=0$ (cut-offs) and $k\to \infty$ (asymptotes) for the fitting procedure but have an appreciable approximation for all values of $k,\omega$ .

6. Fitting of the relaxed micromorphic parameters with curvature (with $\mathbf{Curl}\mathit{P}$ )

To gain more freedom on the dispersion curves shape, we will now consider the addition of $\mathrm{Curl}P$ resulting in a higher-order polynomial in $k$ which enables up to two distinct values of $k$ for every $\omega \in$ . At first, we still neglect $\mathrm{Curl}\stackrel{·}{P}$ resulting in the slightly reduced characteristic polynomial

$$c_0{k}^2+c_1{k}^4-(c_2+c_3{k}^2+c_4^{*}{k}^4){\omega }^2+(c_5+c_6^{*}{k}^2+c_7^{*}{k}^4){\omega }^4-(c_8+c_9^{*}{k}^2){\omega }^6=0.$$

(38)

The characteristic polynomial is of second order regarding $k^2$ which, assuming $k>0$ , results in two distinct roots

$$k_{1,2}=\frac{-c_0+c_3{\omega }^2-c_6{\omega }^4-c_9^{*}{\omega }^6\pm \sqrt{-4(c_1+c_4^{*}{\omega }^2+c_7{\omega }^4)(-c_2+c_5{\omega }^4+c_8{\omega }^6)+{(c_0-c_3{\omega }^2+c_6^{*}{\omega }^4+c_9{\omega }^6)}^2}}{2(c_1+c_4^{*}{\omega }^2+c_7^{*}{\omega }^4)}$$

(39)

describing the dispersion relation. The coefficients $c_0,c_1,c_2,c_3,c_4^{*},c_5,c_6^{*},c_7^{*},c_8,c_9^{*}$ used ¹¹ for the expression above change from $0^{°}$ to ${45}^{°}$ as well as for considering pressure and shear waves are included in the Appendix 3.

6.1. Asymptotes

Because the cut-offs are independent of the coefficients with higher order of $k$ , they do not change with the addition of $\mathrm{Curl}P$ . Instead, the expressions of the asymptotes hugely differ compared to the expression without $\mathrm{Curl}P$ discussed before. We only include the terms with the highest order of $k$ available and compute

$$\begin{array}{c}{c}_1-c_4^{*}{\omega }^2+c_7^{*}{\omega }^4=0\iff {\omega }^4-\frac{c_4^{*}}{c_7^{*}}{\omega }^2+\frac{c_1}{c_7^{*}}=0\hfill \\ \iff {\omega }_{1,2}^2=\frac{c_4^{*}\pm \sqrt{{(c_4^{*})}^2-4c_1{c}_7^{*}}}{2c_7^{*}}.\hfill \end{array}$$

(40)

Surprisingly, the asymptotes with $\mathrm{Curl}P$ are significantly simpler because we must only solve a second-order polynomial instead of a third-order polynomial needed for the cut-offs and the asymptotes without $\mathrm{Curl}P$ . We now only have four distinct horizontal asymptotes (two shear and two pressure) in contrast to six before, which means that we must allow that the two curves (one shear and one pressure) will tend to infinity for high values of $k$ , and our choice falls on the two highest optic curves. The same reasoning about the use of the asymptote in section 5.1 is applied here besides for the two highest optic curves that do not have a horizontal asymptote.

6.2. Fitting

The fitting procedure starts exactly as described in section 5.2. As a first step, we fit the macro parameters ${\kappa }_{\mathrm{M}},{\mu }_{\mathrm{M}},{\mu }_{\mathrm{M}}^{*}$ using the slopes of the dispersion curves for $k=0$ which results in the same values as before. Following this, we express the micro-inertia ${\kappa }_{\gamma },{\gamma }_1,{\gamma }_1^{*},{\gamma }_2$ as a function of the numerical values of the cut-offs and the remaining independent material parameters. Overall, we arrive at four unknown micro parameters ${\kappa }_{\mathrm{m}},{\mu }_{\mathrm{m}},{\mu }_{\mathrm{m}}^{*},{\mu }_{\mathrm{c}}$ ; four unknown inertia parameters ${\overline{\kappa }}_{\gamma },{\overline{\gamma }}_1,{\overline{\gamma }}_1^{*},{\overline{\gamma }}_2$ ; and the new elastic parameter ${\alpha }_1$ belonging to $\mathrm{Curl}P$ .

Although the expressions of the asymptotes are different from the ones without the $\mathrm{Curl}P$ , we still have the same split between the parameters, resulting in four independent parameters for every group of asymptotes, cf. Table 7.

Table 7.

Dependence of the asymptotes of the dispersion curves on the free material parameters as function of the direction of propagation $(0^{°}/{45}^{°})$ and type of wave (shear/pressure).

	Pressure	Shear
$0^{°}$	${\kappa }_{\mathrm{m}},{\overline{\kappa }}_{\gamma },{\mu }_{\mathrm{m}},{\overline{\gamma }}_1$	${\mu }_{\mathrm{c}},{\overline{\gamma }}_2,{\mu }_{\mathrm{m}}^{*},{\overline{\gamma }}_1^{*}$
${45}^{°}$	${\kappa }_{\mathrm{m}},{\overline{\kappa }}_{\gamma },{\mu }_{\mathrm{m}}^{*},{\overline{\gamma }}_1^{*}$	${\mu }_{\mathrm{c}},{\overline{\gamma }}_2,{\mu }_{\mathrm{m}},{\overline{\gamma }}_1$

In particular, all eight expressions are independent on ${\alpha }_1$ and thus $\mathrm{Curl}P$ . Hence, we have no analytic expressions to assign a value to this parameter. Note that the dispersion curves itself (in particular, their shape) still depend on ${\alpha }_1$ but the cut-offs and asymptotes remain independent of this parameter. Because we miss the higher optic curves asymptotes, we introduce four new expressions as a replacement, by considering $k=\frac{1}{2L_{\mathrm{c}}}$ at $0^{°}$ and $\frac{1}{\sqrt{2}{L}_{\mathrm{c}}}$ at ${45}^{°}$ (instead of $k\to \infty$ ) for the corresponding curves, i.e. the highest value for $k$ for which we still have numerical values using Bloch–Floquet analysis. Note again that $L_{\mathrm{c}}$ is always the size of the unit cell. These “pseudo-asymptotes” depend on the same parameters as their corresponding acoustic and lower optic curve asymptotes (cf. Table 7) and the additional independent parameter ${\alpha }_1$ . The analytical expression is again too large to be included here but depend on all the coefficients $c_0,\dots ,c_9^{*}$ of the dispersion polynomial (38), cf. Appendix 3.

We list the numerical values of all parameters used for the fitting of the micromorphic model in Table 8.

Table 8.

Numerical values of the relaxed micromorphic model with $\mathrm{Curl}P$ fitted for the metamaterial whose unit cell is given in Figure 1.

${\kappa }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}^{*}$	${\mu }_{\mathrm{c}}$	${\kappa }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}^{*}$	${\alpha }_1$
(MPa)	(MPa)	(MPa)	(MPa)	(MPa)	(MPa)	(MPa)	(MPa)
1.674	20.04	8.748	15.33	4465	202.4	100.5	14.34
${\overline{\kappa }}_{\gamma }$	${\overline{\gamma }}_1$	${\overline{\gamma }}_1^{*}$	${\overline{\gamma }}_2$	${\kappa }_{\gamma }$	${\gamma }_1$	${\gamma }_1^{*}$	${\gamma }_2$
6.932	5.82	404.4	0.672	3.496	34.61	345	15.13

6.3. Discussion

The fitting shown in Figure 4 including the curvature $\mathrm{Curl}P$ is worse compared to the one without, cf. Figure 3 and Table 8. This is mainly because we lost the higher optic curves asymptotes. Additionally, the shape of the curves which comes automatically by fitting all material parameters using only the cut-offs and asymptotes does not match properly with the data computed numerically with Comsol Multiphysics®. Most importantly, the fitting with $\mathrm{Curl}P$ shows the additional challenge of asymptotes that are approached at a very high wavenumber $k>>\frac{1}{\sqrt{2}{L}_{\mathrm{c}}}$ resulting in a poor fit for values of $k$ between zero and the size of the unit cell, even if the slopes of the acoustic curves close to zero are well fitted. In particular, the acoustic shear wave is notably too slow at approaching its limit. Increasing ${\alpha }_1$ substantially helps for lower frequencies but causes a much higher “pseudo-asymptote”, i.e. a worse fit for higher frequencies. To show this, we did a second fit only using the asymptotes of the acoustic and lower optic curves with ${\alpha }_1$ as an independent parameter set by hand, cf. Figure 5.

Figure 4.

Dispersion curves $\omega (k)$ for $0^{°}$ (left) and ${45}^{°}$ (right) with pressure curves colored in yellow and shear in blue. The dots are the points computed with Comsol Multiphysics® while the smooth curves show the analytical expression of the dispersion curves for the relaxed micromorphic model for ${\beta }_1=0$ . The values of the curve’s horizontal asymptotes are also shown with dashed lines.

Figure 5.

Dispersion curves $\omega (k)$ for $0^{°}$ (left) and ${45}^{°}$ (right) with pressure curves colored in yellow and shear in blue. The dots are the points computed with Comsol Multiphysics® while the smooth curves show the analytical expression of the dispersion curves for the relaxed micromorphic model for ${\beta }_1=0$ . The values of the curve’s horizontal asymptotes are also shown with dashed lines.

Most other values remain at the same magnitude but are slightly higher, cf. Table 9. In future works, we will consider an enhanced relaxed micromorphic model to better describe these effects.

Table 9.

Second set of alternative numerical values of the relaxed micromorphic model with $\mathrm{Curl}P$ fitted for the metamaterial whose unit cell is given in Figure 1.

${\kappa }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}^{*}$	${\mu }_{\mathrm{c}}$	${\kappa }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}^{*}$	${\alpha }_1$
(MPa)	(MPa)	(MPa)	(MPa)	(MPa)	(MPa)	(MPa)	(MPa)
1.418	20.04	8.936	15.61	468.4	202.4	100.5	240.0
${\overline{\kappa }}_{\gamma }$	${\overline{\gamma }}_1$	${\overline{\gamma }}_1^{*}$	${\overline{\gamma }}_2$	${\kappa }_{\gamma }$	${\gamma }_1$	${\gamma }_1^{*}$	${\gamma }_2$
3.616	9.552	408.0	4.288	3.116	34.61	351.2	15.45

7. Fitting of the relaxed micromorphic parameters with enhanced kinetic energy (with $\mathrm{Curl}\stackrel{·}{P}$ )

We now include $\mathrm{Curl}\stackrel{·}{P}$ which reintroduces the third asymptote by considering the full dispersion polynomial

$$c_0{k}^2+c_1{k}^4-(c_2+c_3{k}^2+c_4{k}^4){\omega }^2+(c_5+c_6{k}^2+c_7{k}^4){\omega }^4-(c_8+c_9{k}^2+c_{10}{k}^4){\omega }^6=0,$$

(41)

with coefficients $c_0,\dots ,c_{10}$ depending on all the material parameters described before and ${\beta }_1$ belonging to $\mathrm{Curl}\stackrel{·}{P}$ , cf. Appendix 4.

7.1. Asymptotes

Again, the cut-offs are independent on the coefficients with higher order of $k$ and thus they do not change with respect to the two previous cases. For the asymptotes we only consider the terms with the highest order of $k$ available and compute

$$c_1-c_4{\omega }^2+c_7{\omega }^4-c_{10}{\omega }^6=0\iff {\omega }^6-\frac{c_7}{c_{10}}{\omega }^4+\frac{c_4}{c_{10}}{\omega }^2-\frac{c_1}{c_{10}}=0.$$

(42)

We have again three asymptotes (the roots of a third-order polynomial) which, in general, causes the analytical expressions to be impractical rather quickly. However, in this case, it is possible to find one root by hand

$$c_1-c_4{\omega }^2+c_7{\omega }^4-c_{10}{\omega }^6=0\iff (1-\rho L_{\mathrm{c}}^2\frac{{\beta }_1}{{\alpha }_1}{\omega }^2)(c_1-c_4^{*}{\omega }^2+c_7^{*}{\omega }^4)$$

(43)

with $c_4^{*}$ and $c_7^{*}$ from the dispersion relation (38) without $\mathrm{Curl}\stackrel{·}{P}$ which holds because

$$\rho L_{\mathrm{c}}^2\frac{{\beta }_1}{{\alpha }_1}{c}_1=c_4-c_4^{*}\mathrm{and}\rho L_{\mathrm{c}}^2\frac{{\beta }_1}{{\alpha }_1}{c}_4^{*}=c_7-c_7^{*}\mathrm{and}\rho L_{\mathrm{c}}^2\frac{{\beta }_1}{{\alpha }_1}{c}_7^{*}=c_{10}$$

(44)

for all combinations of shear/pressure and $0^{°}$ / ${45}^{°}$ angle of incidence. This is remarkable because it allows us to use the same analytical expressions for the acoustic and lower optic asymptotes from the calculations without $\mathrm{Curl}\stackrel{·}{P}$ of section 6, while the general expressions of the dispersion curves differ because of the addition of new terms. Second, the third asymptote

$${\omega }_3=\sqrt{\frac{{\alpha }_1}{\rho L_{\mathrm{c}}^2{\beta }_1}}$$

(45)

is identical for shear and pressure and invariant under change of angle of incidence. The same reasoning about the use of the asymptote in section 5.1 is applied here.

7.2. Fitting

The fitting starts as described in the sections before, using cut-offs and the slope of acoustic waves for $k=0$ to reduce the number of independent parameters. Here we arrive at the same four unknown micro parameters ${\kappa }_{\mathrm{m}},{\mu }_{\mathrm{m}},{\mu }_{\mathrm{m}}^{*},{\mu }_{\mathrm{c}}$ and four unknown inertia parameters ${\overline{\kappa }}_{\gamma },{\overline{\gamma }}_1,{\overline{\gamma }}_1^{*},{\overline{\gamma }}_2$ with the additional elastic parameter ${\alpha }_1,{\beta }_1$ belonging to $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ , respectively. For the latter, we choose the acoustic pressure and lower optic shear curves as the one with the superimposed asymptote. The former are determined numerically using the remaining asymptotes (eight in total), see also Table 10.

Table 10.

Dependence of the asymptotes of the dispersion curves on the free material parameters as function of the direction of propagation $(0^{°}/{45}^{°})$ and type of wave (shear/pressure).

	Pressure	Shear
$0^{°}$	${\kappa }_{\mathrm{m}},{\overline{\kappa }}_{\gamma },{\mu }_{\mathrm{m}},{\overline{\gamma }}_1$	${\mu }_{\mathrm{c}},{\overline{\gamma }}_2,{\mu }_{\mathrm{m}}^{*},{\overline{\gamma }}_1^{*}$
${45}^{°}$	${\kappa }_{\mathrm{m}},{\overline{\kappa }}_{\gamma },{\mu }_{\mathrm{m}}^{*},{\overline{\gamma }}_1^{*}$	${\mu }_{\mathrm{c}},{\overline{\gamma }}_2,{\mu }_{\mathrm{m}},{\overline{\gamma }}_1$
Superimposed	${\alpha }_1,{\beta }_1$	${\alpha }_1,{\beta }_1$

However, because the analytical expression (45) of these four highest asymptotes is identical, while, in general, the numerical values from Comsol Multiphysics® can be distinct, we will arrive at the average of these four asymptotes. Moreover, we can only fix the ratio $\frac{{\alpha }_1}{{\beta }_1}$ which leaves us with one last free parameter, i.e. independent of the values of all cut-offs and asymptotes, which we fit by hand.

We list the numerical values of all parameters used in the micromorphic model in Table 11.

Table 11.

Numerical values for the relaxed micromorphic model with enhanced kinetic and potential energy terms $(\mathrm{Curl}P,\mathrm{Curl}\stackrel{·}{P})$ for the metamaterial whose unit cell is given in Figure 1.

${\kappa }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}$	${\mu }_{\mathrm{m}}^{*}$	${\mu }_{\mathrm{c}}$	${\kappa }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}$	${\mu }_{\mathrm{e}}^{*}$	${\alpha }_1$
(MPa)	(kPa)	(kPa)	(kPa)	(kPa)	(kPa)	(kPa)	(Pa)
35.22	228.4	27.67	110.7	355.8	1635	973.6	597.4

${\overline{\kappa }}_{\gamma }$	${\overline{\gamma }}_1$	${\overline{\gamma }}_1^{*}$	${\overline{\gamma }}_2$	${\kappa }_{\gamma }$	${\gamma }_1$	${\gamma }_1^{*}$	${\gamma }_2$	${\beta }_1$
1.168	2.004	2.276	0.820	58.68	3.188	1.184	1.856	1.345·10−2

7.3. Discussion

The introduction of $\mathrm{Curl}\stackrel{·}{P}$ increases the quality of fitting again (see Figure 6) with the recovery of the higher optic curve asymptote, while the remaining free parameter ${\alpha }_1,{\beta }_1$ can improve the shape of the curves. In particular, we can fit a dispersion curve with negative group velocity perfectly, the lower optic pressure curve for ${45}^{°}$ , for the first time.

Figure 6.

Dispersion curves $\omega (k)$ for $0^{°}$ (left) and ${45}^{°}$ (right) with pressure curves colored in yellow and shear in blue. The dots are the points computed with Comsol Multiphysics® while the smooth curves show the analytical expression of the dispersion curves for the relaxed micromorphic model. The values of the curve’s horizontal asymptotes are also shown with dashed lines.

The fitting of the other curves is slightly worse compared to the calculations without curvature terms, and it is mainly due to the degenerated shape of $\mathrm{Curl}\stackrel{·}{P}$ and $\mathrm{Curl}P$ in a plane problem, resulting in an isotropic behavior with only a single independent parameter. This causes four asymptotes to coincide.

We will not discuss in a fourth section the fitting with an enhanced kinetic energy with $\mathrm{Curl}\stackrel{·}{P}$ but without the curvature $\mathrm{Curl}P$ , since it leads to a degenerate dispersion polynomial $p(k,\omega )$ similar to the one described in section 6 resulting in only four configurable horizontal asymptotes for shear and pressure each. However, while we lost the higher optic asymptotes by adding $\mathrm{Curl}P$ without $\mathrm{Curl}\stackrel{·}{P}$ , the case with $\mathrm{Curl}\stackrel{·}{P}$ and without curvature $\mathrm{Curl}P$ instead forces the limits of the one curve to zero as it can be seen as a limit case of expression (45), i.e. ${\alpha }_1\to 0$ .

8. Summary of the obtained results

Comparing the results of sections 5–7, the first approach which does not include $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ shows the best agreement to the numerical values from Comsol Multiphysics® overall. At the same time, it must be noted that this simplified model without any curvature terms is only capable of describing monotonic dispersion curves with three disjointed domains for the three pressure waves and three disjointed domains for the three shear waves (crossing between pressure and shear waves is allowed, while it is not allowed between curves of the same type). The latter property can be easily deduced from equation (31)₂ which guarantees that for each value of $\omega$ there can only be one value of $k$ .

In order to remove this constraint, the dispersion relation $p(k,\omega )=0$ must contain higher order terms of the wavenumber $k$ to allow the overlapping of the domains of curves of the same kind. To this aim, we started considering the full relaxed micromorphic model (with $\mathrm{Curl}P$ ) and then augmented it with a new inertia term ( $\mathrm{Curl}\stackrel{·}{P}$ ). The comparison of the analytic expression of the asymptotes with or without $\mathrm{Curl}\stackrel{·}{P}$ shown in equation (43) suggests that, when a new term is added to the elastic energy density, it is always better to include the corresponding dynamic part as well. In particular, considering $\mathrm{Curl}P$ without its counterpart $\mathrm{Curl}\stackrel{·}{P}$ or vice versa causes missing terms in the dispersion polynomial $p(k,\omega )=0$ resulting in fewer horizontal asymptotes.

While the model without curvature terms gives the best quantitative agreement (see Figure 3), the augmented relaxed micromorphic model still gives a good agreement and opens up the important possibility of decreasing modes with negative group velocity, cf. Figure 6. Note that $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ degenerate for planar problems: for both terms, it just remains a single independent constitutive parameter ( ${\alpha }_1$ and ${\beta }_1$ ) restricting the class of symmetry for their constitutive tensors to the isotropic one.

We want to emphasize that the main focus of this work is not the result of the fitting of the three different approaches per se, but the semi-analytical fitting algorithm itself and the underlying consistency of the relaxed micromorphic model with respect to the material properties and some of the geometrical characteristic of the metamaterial that it represents. Using the complex but analytically defined expressions of the asymptotes, we can find a numerical fit of all material parameters by only giving the numerical values computed with Comsol Multiphysics® and the apparent mass density $\rho$ of the unit cell. Note that we only use the cut-offs $k=0$ and asymptotes $k\to \infty$ for calculating the material parameters while the shape of the curves for intermediate values of $k$ comes automatically.

The routine is completely written with Mathematica allowing us to use symbolic calculations. The essential part of the fitting procedure uses the inbuilt algorithm NMinimize (with the Method RandomSearch) to minimize the mean square error of the asymptotes between the relaxed micromorphic model and the numerical values of the finite element approach in Comsol Multiphysics®. Therefore, in general, if a local minimum is found, it is not guaranteed that it corresponds to a global optimum as well.

In order to better understand the nature of the minimization problem, we visualize the impact of each parameter thanks to Figure 7 for the case without $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ of section 5.2: given a set of material constants values, we move one parameter at the time and plot the impact on the error for the asymptotes.

Figure 7.

Mean squared error between the analytic expressions from the relaxed micromorphic model and their corresponding numerical values from Comsol Multiphysics® while changing the scaling coefficient $t$ for different material parameters. Left: the static parameters ${\kappa }_{\mathrm{m}},{\mu }_{\mathrm{m}},{\mu }_{\mathrm{m}}^{*},{\mu }_{\mathrm{c}}$ are colored in orange, blue, green, and purple, while the dashed red line shows the impact of scaling ${\mathbb{C}}_{\mathrm{micro}}$ as a whole. Right: the dynamic parameters ${\overline{\kappa }}_{\gamma },{\overline{\gamma }}_1,{\overline{\gamma }}_1^{*},{\overline{\gamma }}_2$ are colored in orange, blue, green, and purple, while the dashed red line shows the impact of scaling ${\mathbb{T}}_{\mathrm{e}}$ as a whole.

It is crucial to state here that the expression of the coefficients shown in Appendices 2–4 are a priori constrained by the relations described in section 5.2 which reduces the number of independent parameters. This guarantees that the values of the cut-offs will not change if the micro parameters change because the corresponding micro-inertia parameters ${\kappa }_{\gamma },{\gamma }_1,{\gamma }_1^{*},{\gamma }_2$ are automatically scaled as well.

When considering the impact of the different independent parameters on the fitting error of the asymptotes, visible dissimilarity arises. While most coefficients show a simple behavior with one distinct minimum, the impact on the dispersion curves of the parameters ${\mu }_{\mathrm{m}}^{*}$ and ${\overline{\gamma }}_1^{*}$ is very small. In addition to the numerical difficulty caused by a vanishing gradient, it suggests finding additional quantities to constrain these material parameters, e.g. use a static test to fix ${\mu }_{\mathrm{m}}^{*}$ beforehand. Overall, the minimum problem behaves remarkably well given its simple shape, which is in accordance with the fact that the fitting procedure always converges to very similar values of the parameters regardless of their initial guess.

9. Conclusion and perspectives

We presented an inertia-augmented relaxed micromorphic model that enhances the model proposed in [52–56] via the addition of a term $\mathrm{Curl}\stackrel{·}{P}$ in the kinetic energy. We then used this model to describe the dynamical behavior of a labyrinthine metamaterial and compared its predictability when setting all curvature terms to zero, when setting only $\mathrm{Curl}\stackrel{·}{P}$ to zero and when considering the full model (including $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ ). The model’s efficiency is tested by comparing the obtained dispersion curves to the ones issued via a standard Bloch–Floquet analysis. We find that the model without curvature terms gives the best average behavior (when considering the whole set of frequencies and wave numbers). However, this reduced model shows two main drawbacks: (1) it cannot give rise to curves with negative group velocity and (2) the horizontal asymptote of each curve is bounded by the cut-off of the following one. This implies a limitation of the fitting quality at higher wave numbers for some directions of propagation. To remove these constraints, the full model (including $\mathrm{Curl}P$ and $\mathrm{Curl}\stackrel{·}{P}$ ) can be used. The fitting that we obtained with the full model remains of good quality (although precision is sometimes lowered point-wise) and we could perfectly describe (for the first time in a micromorphic framework) a higher-frequency mode with negative group velocity. The constraint concerning the asymptotes’ boundedness is also removed in the full model, but a perfect fitting at higher wavenumbers cannot still be achieved due to a lack of extra microscopic degrees of freedom that are needed (at least for some directions) for a perfect fitting at wavelengths approaching the size of the unit cell.

We also tested the model’s performances by considering only $\mathrm{Curl}P$ and not $\mathrm{Curl}\stackrel{·}{P}$ . The fitting quality is visibly worsening, thus suggesting that, when introducing a term in the strain energy, its dynamical counterpart in the kinetic energy should be always considered as well.

In addition to the fitting comparison presented in the present paper, one main result that we present here is the fitting procedure itself that has been automatized to a big extent by asking only the cut-offs and asymptotes to be imposed a priori. This has been done by imposing the exact value of the cut-offs and minimizing the asymptotes’ mean square error compared to the exact numerical values issued via Bloch–Floquet analysis. The rest of the curves’ fitting follows directly.

The routine for the fitting procedure is written with Mathematica allowing us to use symbolic calculations. The essential part of the routine uses the built-in algorithm NMinimize (with the Method RandomSearch) to minimize the asymptote’s mean square error. We finally checked whether the found local minimum is also a global minimum. We found that all the elastic parameters achieve a global minimum in the computed configuration, but ${\mu }_{\mathrm{m}}^{*}$ and ${\overline{\gamma }}_1^{*}$ seem to have little effect on the overall fitting after a certain threshold. This result seems to suggest that the values of these parameters should be eventually fixed beforehand with the help of extra independent static and/or dynamic tests.

Based on the findings of this paper, we will briefly present some insight that will give directions to the follow-up research:

• Further enhance the relaxed micromorphic model via the addition of extra microscopic degrees of freedom to increase its precision at very small wavelengths (approaching the unit cell’s size);

• Design complex large-scale meta-structures that control elastic energy using the new labyrinthine metamaterial as a basic building block. This design would not be otherwise possible due to the huge number of degrees of freedom resulting from the meshing of all the tiny elements contained in the labyrinthine unit cells;

• Study negative refraction phenomena in meta-structures including the new labyrinthine metamaterial as a basic building block;

• Design complex structures for wave control simultaneously including the different metamaterials that were characterized via the relaxed micromorphic model until now.