An asymptotic hyperbolic-elliptic model for flexural-seismic metasurfaces

TL;DR

This paper develops an asymptotic hyperbolic-elliptic model for elastic metasurfaces with flexural resonators, simplifying dispersion relations and revealing how junction conditions influence wave propagation and stop bands.

Contribution

It introduces a new asymptotic scalar hyperbolic model for flexural resonators on elastic surfaces, extending previous compressional models and highlighting the impact of junction conditions.

Findings

The model provides explicit dispersion relations with closed-form solutions.

Flexural resonators significantly alter wave stop bands.

Junction conditions dramatically affect wave propagation features.

Abstract

We consider a periodic array of resonators, formed from Euler-Bernoulli beams, attached to the surface of an elastic half-space. Earlier studies of such systems have concentrated on compressional resonators. In this paper we consider the effect of the flexural motion of the resonators, adapting a recently established asymptotic methodology that leads to an explicit scalar hyperbolic equation governing the propagation of Rayleigh-like waves. Compared with classical approaches, the asymptotic model yields a significantly simpler dispersion relation, with closed form solutions, shown to be accurate for surface wave-speeds close to that of the Rayleigh wave. Special attention is devoted to the effect of various junction conditions joining the beams to the elastic half-space which arise from considering flexural motion and are not present for the case of purely compressional resonators. Such effects are shown to provide significant and interesting features and, in particular, the choice of junction conditions dramatically changes the distribution and sizes of stop bands. Given that flexural vibrations in thin beams are excited more readily than compressional modes and the ability to model elastic surface waves using the scalar wave equation (i.e. waves on a membrane), the paper provides new pathways toward novel experimental set-ups for elastic metasurfaces.

Tables

Table 1. Table 1 : The numerical system parameters for the beam resonator and half-space used to produce the dispersion curves shown in Figs. 3 , 4 , 5 and 6 .

Symbol	Definition	Value
$l$	Lattice spacing	2 m
$\rho$	Half-space density	13000 kg m-3
$\mu$	Half-space shear modulus	325 MPa
$\lambda$	Half-space first Lamè parameter	702 MPa
$L$	Resonator length	14 m
$h$	Resonator diameter	0.3m
${\rho }_R$	Resonator density	450 kg m-3
$E_R$	Resonator stiffness	1.70 GPa

Equations

$$u_1=\frac{\partial \varphi }{\partial x_1}-\frac{\partial \psi }{\partial x_3},u_3=\frac{\partial \varphi }{\partial x_3}+\frac{\partial \psi }{\partial x_1}.$$

$$\frac{{\partial }^2\varphi }{\partial x_1^2}+\frac{{\partial }^2\varphi }{\partial x_3^2}-\frac{1}{c_1^2}\frac{{\partial }^2\varphi }{\partial t^2}=0,\frac{{\partial }^2\psi }{\partial x_1^2}+\frac{{\partial }^2\psi }{\partial x_3^2}-\frac{1}{c_2^2}\frac{{\partial }^2\psi }{\partial t^2}=0,$$

$$E_R{I}_R\frac{{\partial }^4\mathbf{U}}{\partial x_3^4}=-M_R\frac{{\partial }^2\mathbf{U}}{\partial t^2},i=1,2,$$

$$E_R=\frac{E}{1-{\nu }^2}$$

$$I_R=\frac{\pi h^4}{64},$$

$$\frac{{\partial }^3\mathbf{U}}{\partial x_3^3}=0,\frac{{\partial }^2\mathbf{U}}{\partial x_3^2}=0.$$

$${\sigma }_{3i}=0,i=1,2,3.$$

$$\begin{array}{rl}{\sigma }_{31}& =2\mu {\varphi }_{,13}+\mu {\psi }_{,11}-\mu {\psi }_{,33},\\ {\sigma }_{33}& =(\lambda +2\mu ){\varphi }_{,33}+\lambda {\varphi }_{,11}+2\mu {\psi }_{,13}.\end{array}$$

$$\varphi =A_{\varphi }{e}^{i(k{x}_1-\omega t)-k\alpha x_3},\psi =A_{\psi }{e}^{i(k{x}_1-\omega t)-k\beta x_3},$$

$$\alpha =\sqrt{1-{\left(\frac{\omega }{c_{1}k}\right)}^2},\beta =\sqrt{1-{\left(\frac{\omega }{c_{2}k}\right)}^2}.$$

$$4\alpha \beta -(1+{\beta }^2{)}^2=0,$$

$$\alpha =\sqrt{1-{\left(\frac{c_R}{c_1}\right)}^2}\triangleq {\alpha }_R,\beta =\sqrt{1-{\left(\frac{c_R}{c_2}\right)}^2}\triangleq {\beta }_R.$$

$${\sigma }_{31}=\frac{H}{l^2},{\sigma }_{31}=\frac{V}{l^2}$$

$$ϵ=\left|1-\frac{\omega }{c_{R}k}\right|\ll 1,$$

$$\begin{array}{rl}\frac{{\partial }^2\varphi }{\partial x_3^2}+{\alpha }_R^2\frac{{\partial }^2\varphi }{\partial x_1^2}& =0,\\ \frac{{\partial }^2\psi }{\partial x_3^2}+{\beta }_R^2\frac{{\partial }^2\psi }{\partial x_1^2}& =0,\end{array}$$

$$\varphi =A_{\varphi }{e}^{i(k{x}_1-\omega t)-k{\alpha }_R{x}_3},\psi =A_{\psi }{e}^{i(k{x}_1-\omega t)-k{\beta }_R{x}_3}.$$

$${\psi }_{,11}-\frac{1}{c_R^2}{\psi }_{,tt}$$

$$B=\frac{{\beta }_R}{{\alpha }_R}(1-{\alpha }_R^2)+\frac{{\alpha }_R}{{\beta }_R}(1-{\beta }_R^2)-(1-{\beta }_R^4),$$

$${\varphi }_{,11}-\frac{1}{c_R^2}{\varphi }_{,tt}$$

$${\psi }_{,1}=-\frac{2}{1+{\beta }_R^2}{\varphi }_{,3},$$

$$\Delta \varphi -\frac{1}{c_R^2}{\varphi }_{,tt}$$

$${\varphi }_{,i}=\frac{2}{1+{\beta }_R^2}{\psi }_{i,3},{\varphi }_{,3}=-\frac{1+{\beta }_R^2}{2}\left({\psi }_{1,1}+{\psi }_{2,2}\right),i=1,2.$$

$$\begin{array}{rl}{U}_1(t,0)& =u_1{|}_{x_3=0},\\ U_{1,33}(t,0)& =0,\end{array}$$

$$H(x_1,x_2,t)$$

$$H(x_1,x_2,t)=E_R\frac{\pi h^4}{64}{K}^3\frac{\cos (KL)\cosh (KL)-1}{\cosh (KL)\sin (KL)-\cos (KL)\sinh (KL)}{u}_1{|}_{x_3=0},$$

$$K={\left(\frac{64M_R}{E_R\pi h^4}\right)}^{\frac{1}{4}}\sqrt{\omega }.$$

$${\sigma }_{31}$$

$$-2i\mu k^2\alpha \varphi -\mu k^2(1+{\beta }^2)\psi$$

$$\left((\lambda +2\mu )k^2{\alpha }^2-\lambda k^2\right)\varphi -2i\mu k^2\beta \psi$$

$$\left((\lambda +2\mu )k^2{\alpha }^2-\lambda k^2\right)\left(-\mu k^2\left(1+{\beta }^2\right)-k{K}^3\hat{H}\beta \right)=\left(2\mu k^2\alpha +k{K}^3\hat{H}\right)\left(-2\mu k^2\beta \right),$$

$${\psi }_{,11}-\frac{1}{c_R^2}{\psi }_{,tt}$$

Full text

1 Introduction

In recent years there has been an explosion of interest, both scholarly and popular, in the use of so-called locally resonant metamaterials in the control of wave propagation. The concept has its origins in optics [1], but has since found application across a broad range of fields [2], including acoustics [3, 4, 5, 6], elasticity [7, 8, 9], flexural waves in thin plates [10, 11, 12], and seismic waves [13, 14]. Despite the term metamaterial first appearing at the turn of the millennium [15], the development of this relatively new field continues apace. The essential concept is that, using arrays of local resonators, metamaterials can be designed to have properties that do not occur in materials found in nature; the canonical example is negative refraction [16].

At their inception, metamaterials were designed for the control of waves in bulk media. This idea can be extended to surface waves, creating so-called metasurfaces [17]; these meta surfaces have the substantial advantage over metamaterials in that they are relatively easy to manufacture [18]. The majority of work on metasurfaces has focused on electromagnetism [18] and, for example, spoof-surface-plasmons [19]. At present, however, the use of metasurfaces has begun to gain traction in other fields, such as acoustics [20, 21] and, very recently, elastodynamics [22, 23] and fluid-solid interactions [24]. In comparison to, say, electromagnetism and acoustics, the treatment of metasurfaces in the framework of elastodynamics is far more challenging and offers a greater array of unique phenomena.

Despite the singular challenges associated with mechanical metamaterials, they have a substantial advantage in that the range of scales over which they can be utilized is large compared with Electromagnetic metamaterials: from the order of millimeters [25] for cloaking of small-scale structures to several meters [13] for seismic metamaterials. The second notable feature is the existence of two bulk waves (shear and pressure), which propagate at two different speeds; this phenomena is unique to elastic waves and has no analogue in electromagnetics or acoustics. The existence of two distinct classes of wave adds both challenges and novel features. Elastic media also support surface waves, again travelling at two different speeds, that are localised to the surface of the medium, but do not decay in the directions parallel to the surface. This is in contrast to electromagnetism where surface plasmons [26] decay in all directions. The closest analogue of “true” surface waves in electromagnetism are so-called spoof-surface-plasmons [19], which are supported on structured surfaces.

Recently, there have been several analytical [22, 23], numerical [14, 27], and experimental [28] studies on the use of sub-wavelength resonators for the control of surface waves on elastic substrates. Notable features of such resonant arrays include filtering, focusing, rainbow trapping, and mode-conversion of surface waves to bulk waves [29].

Interestingly, a recent numerical and experimental study [30] has shown that the hybridisation of the two zeroth-order modes in thin plates can be controlled by using symmetric and asymmetric distributions of resonators. The experimental studies [28, 30] are incumbent on prior completion of numerical or analytical studies. The numerical approaches are computationally intensive and require extensive GPU computing facilities [14], while complete analytical studies present their own significant challenges [23].

In the present paper we develop a novel asymptotic formulation, initially introduced in [31], to analyse the interaction of surface waves with a resonant array resting on an elastic half-space. This asymptotic approach is used to analyse the dynamic effects near the surface of semi-infinite elastic bodies and consists of an explicit hyperbolic equation governing the Rayleigh wave propagation as well as a pseudo-static elliptic equation governing decay into the interior. This approach has proven a powerful tool for a variety of near-surface elastodynamic problems, including 3-dimensional moving load problems [32] and the effect of adding pre-stress to the elastic half-space [33]. For more general works, the reader is referred  [34, 35, 36]. Recently this hyperbolic-elliptic model was applied to a resonant array of rods attached to a half-space  [37]; these rods were assumed to undergo only longitudinal vibrations as was the case in previous analytical [23] and numerical [22] studies. However recent experimental data obtained for a dense forest of pine trees, acting as subwavelength resonators, suggests that the flexural resonances as well as the junction conditions [38] may play an important role; yet, little attention has thus far been given to the junction condition joining the resonators to the substrate.

Here we extend this approach to include the effects of flexural motion in the thin rods and investigate the effect of these flexural resonances as well as the junction conditions coupling the resonators to the elastic half-space. In order to better see the effect of the flexural motion, and because in the frequency ranges considered flexural waves dominate, the longitudinal motion of the resonators is neglected. These new interactions bring new features, including the ability to control the width of the deep sub-wavelength stop bands. This is an important feature as, usually, sub-wavelength resonant arrays are only capable of producing very narrow band gaps [7] and, as a result, it is challenging to make use of these band gaps in the design of filters and other wave control devices. While a non-dimensionalisation of this problem would be useful to consider, this manuscript is following from an existing body of work with strong practical motivations; in the literature of such problems it is more typical to see the results presented in dimensional form.

As is usual for the treatment of such low amplitude bending waves, the effect of the weight of the beams on vertical force into the half-space and the bending motion of the resonators is considered to be negligible compared to the bending forces.

The present paper is be organised as follows. The problem formulation is given in Section 2 with the formulation for Rayleigh waves in subsection 2.2. The sections which follow then each detail the treatment of a different junction condition between the resonators and the half-space. Section 3 considers the resonators as simply supported, with the resonators able to freely rotate but conserve displacements and horizontal stresses. In Section 4 the resonators are considered as ideally attached at the base to rails which can move freely but matches the bending moment and the gradient at the surface. For the junction considered in Section 5 we consider the resonators as ideally attached to the half-space and so all variables between the resonator and half-space are fully matched. In each section, the solution obtained from the asymptotic formulation is compared with an implicit evaluation of the full unimodal solution as used in the previous treatment of similar systems [23].

2 Statement of the problem

Before considering the effects of the various interface conditions, we begin with a statement of the problem. We consider the propagation of surface waves on an infinite elastic half-space upon which rests an array of thin Euler-Bernoulli beams, as illustrated in Figure 1. These surface waves propagate with constant amplitude along the free surface of the half-space and decay, exponentially, into the bulk. We begin by formulating the three-dimensional problem.

2.1 An Array of Beams embedded into an Elastic Half-Space

We consider the linearly elastic motion of a three-dimensional half-space ${-\infty<x_{1}<\infty}$ , ${-\infty<x_{2}<\infty}$, ${0\leq x_{3}<\infty}$. For the usual displacement vector $\boldsymbol{u}=[u_{1},u_{2},u_{3}]^{\mathrm{T}}$ we consider only motion in the $x_{1}$ and $x_{3}$ directions and introduce two displacement potentials $\phi$ and $\psi$ such that

[TABLE]

Then for linear elasticity, where $\rho$ is the volume mass density and where $\lambda$ and $\mu$ are Lamé’s first and second parameters respectively, the equations of motion can be expressed as

[TABLE]

where $t$ is time and $c_{1}=\sqrt{(\lambda+2\mu)/\rho}$ and $c_{2}=\sqrt{\mu/\rho}$ are the longitudinal and shear wave speed respectively.

The flexural motion of the thin elastic beams that pattern the surface of the elastic half-space and form the resonant array obey the Euler-Bernoulli beam equation

[TABLE]

where $\boldsymbol{U}(t,x_{3})=[U_{1},U_{2}]^{\mathrm{T}}$ is the displacement vector parallel to the surface of the half-space ${x_{3}=0}$, for simplicity $E_{R}$ is introduced to represent the resonator stiffness,

[TABLE]

where $E$ and $\nu$ respectively are the Young Modulus and Poisson ratio of the resonators, $M_{R}$ is the mass per unit length of the resonators, and $I_{R}$ is the area moment of inertia of the beams, which for resonators with cylindrical cross-section is given by,

[TABLE]

where $h$ is the diameter of a beam. The ends of the beams at $x_{3}=-L$ are assumed to be free such that the moments and transverse forces vanish

[TABLE]

For the moment, we will leave the interface condition at the base of each resonator ($x_{3}=0$), connecting it to the half-space, undefined. These junction conditions will be defined and studied in subsequent sections. In general, the interface between one-dimensional and three-dimensional bodies results in a singularly perturbed problem which requires careful consideration of the boundary layer; see, for example, [39].

For sufficiently thin resonators the resistance to bending motion is, typically, much smaller than the resistance to compression and therefore, in the linear regime, the flexural and compressional deformations are decoupled at leading order (see, for example, [40]). Moreover [23] suggest that, for the case of resonant arrays on thin plates, the effects of flexural and compressional resonances can be considered independently. Resonant arrays composed of purely compressional resonators atop a half-space have previously been considered in [22, 23, 37]. However, a detailed analysis of the influence of flexural resonances on the propagation of surface waves is lacking. Therefore, in the present paper, we focus on the flexural motion of the resonators and consider pure bending only such that vertical motion on the surface of the half-space does not couple to compressional deformation of the resonators. It is however emphasised that vertical tractions may arise at the surface of the half-space as a result of rotations in the flexural resonators [40] and these are taken into account.

2.2 Asymptotic Formulation for Rayleigh Waves

Considering only waves propagating in the positive $x_{1}$ direction leads to a simplification of the governing equations of elasticity and subsequently a simplification to the asymptotic model, for which the full 3-dimensional model is given in [32].

Using the displacement potentials (1) and governing equations (2) from before, we introduce the usual stress tensor $\boldsymbol{\sigma}$, where in the case of a free surface, stress free junction conditions are imposed yielding,

[TABLE]

In terms of the displacement potentials the stresses are given by

[TABLE]

where, as usual, a subscript of the variable denotes a derivative. We introduce the following Ansätze corresponding to time-harmonic waves of radian frequency $\omega$

[TABLE]

where $A_{\phi}$, $A_{\psi}$, $k$ and are constants and $\alpha$ and $\beta$ are given by

[TABLE]

Substitution of these Ansätze into the free surface junction conditions yields the condition for the existence of a non-trivial solution [23]

[TABLE]

from which an implicit expression for the Rayleigh wave speed, denoted by $c_{R}$, can be obtained. The dispersion equation for such waves is linear, $\omega=c_{R}k$, and the Rayleigh decay constants are subsequently given by

[TABLE]

In the parameter regime of interest, where the characteristic wavelength of the travelling wave is sufficiently large compared with both the cross-section and separation of the resonators, then, as in previous treatments of such problems [37, 41], the surface stress associated with the resonant array may be distributed in the form

[TABLE]

where $H$ and $V$ are the amplitudes of the horizontal and vertical forces respectively, and $l$ is the separation between the resonators.

Each type of junction condition treated in this paper results from the matching of variables between the resonators and the surface of the half-space at $x_{3}=0$. The same approach will be used in the treatment of the moments arising from the resonators. However, we defer this discussion until the need arises in §4.

We now develop the asymptotic model, initially introduced by Kaplunov et al. for a near-Rayleigh surface waves (see [32] and references therein for further details). This asymptotic model is a leading order perturbation of the equations of motion (2) around the Rayleigh wave eigensolution above (11). A perturbation in the form of a small surface stress is imposed which causes a small deviation from the Rayleigh speed. This expansion is valid in the regime where

[TABLE]

and, in particular, the solution is known to be valid only if $\frac{\omega}{k}\sim c_{R}$. From inspection of the surface conditions this is equivalent to when both of the surface stresses are small. The equations for the bulk then reduce to

[TABLE]

which then lead to the travelling wave solutions (9) becoming

[TABLE]

Depending on the applied surface stress, the model is formulated in terms of either the longitudinal or shear wave potential. If $\sigma_{33}=0$, then it is convenient to work with the shear wave potential and the junction condition along the surface $x_{3}=0$ is

[TABLE]

where

[TABLE]

If, on the other hand, $\sigma_{31}=0$ then the junction condition at $x_{3}=0$ is written in terms of the pressure potential such that

[TABLE]

Moreover, on the surface $z=0$, the displacement potentials are related by

[TABLE]

allowing for the full solution to be obtained. Not only is this form of the surface equations much simpler than the full unimodal form, but for a known stress they will each only require the solving of a single one-dimensional hyperbolic problem along the surface with the known Rayleigh wave speed.

It is straight-forward to recover the three-dimensional field from the two-dimensional asymptotic model. In particular, it is necessary to introduce a second shear potential $\psi_{2}$, corresponding to the stress $\sigma_{32}$, and then simply replace any second order derivatives with respect to $x_{1}$ with the two dimensional Laplace operator, $\Delta$. For instance, the vertical surface relation (19) in the full three-dimensional form becomes

[TABLE]

The final step in recovering the full three-dimensional formulation is to replace the surface potential relation (20) by its three-dimensional equivalents

[TABLE]

The problem under consideration concerns the effect of a resonant array which, as described earlier, can be treated as a distributed stress on a travelling Rayleigh wave; this makes the asymptotic model introduced here an ideal tool to investigate the problem. For completeness, however, the solution to the full problem (cf. §22.1) will also be computed numerically and used as a comparison to the asymptotic solution as verification. We shall now proceed to consider the effect of various interface conditions on the propagation of surface waves over the half-space.

3 Simply Supported Beams

We begin by considering the resonators as simply supported at the surface of the half-space. In this system, there is no need to conserve bending angle and there is no overall bending moment at the coupling. Instead only the horizontal displacement and the transverse force require matching, as illustrated in Figure 2.

We emphasise that, as discussed in §2.1 on p.2.1, we consider the influence of the resonators under pure bending only and, since there are no rotations at the junction point, the vertical component of traction does not couple into the beam. Note also that, since the beam rotations are not conserved in the half-plane, any rigid body rotations of the beams are also neglected. The junction conditions at the end of the beam atop the half-space are then

[TABLE]

and the stress at the end of the beam is,

[TABLE]

where $U_{1}$ is the horizontal displacement of the the beam and $H(x_{1},y,t)$ is the horizontal force at the surface of the half-space caused by the resonators.

By solving the beam equation (3) with these junction conditions and the free end conditions (6) in the usual way we obtain the relation for horizontal force,

[TABLE]

where

[TABLE]

This force is applied as a point force in the regular array. If the size of the array and wavelength of the wave are assumed to be large in comparison with the distance between the point forces then we can homogenise these point forces into a continuous surface stress using the distribution (13). This stress is then given by,

[TABLE]

where $l$ is the distance between two adjacent resonators.

We first consider the full unimodal formulation for the system. From (8), the surface stresses can be expressed in terms of displacement potentials such that

[TABLE]

From (1), we have that $U_{1}(0)=k(i\phi+\beta\psi)$ and, by substituting from (25), we obtain the full unimodal dispersion equation

[TABLE]

where $\hat{H}=\sigma_{31}K^{-3}/U_{1}(0)$. It is clear that in this form this expression will be difficult to manipulate. If the asymptotic model produces a close approximation which can be used to investigate the behaviour of the dispersion relation, then this will enable us to understand this system on a more fundamental level than the full unimodal solution can.

For the asymptotic model, we use equation (17) along the surface,

[TABLE]

which, when combined with (1) and (20), yields the asymptotic dispersion relation

[TABLE]

The asymptotic dispersion relation can therefore be given explicitly, in terms of the roots of the quadratic equation,

[TABLE]

This asymptotic solution can also be obtained directly from the full dispersion relation obtained earlier by expanding around the Rayleigh solution. In particular, the Rayleigh condition (11) corresponds to the roots of the function

[TABLE]

where $r=\omega^{2}c_{2}^{-2}k^{-2}$. For $0<r\ll 1$, the leading order term of (34) is

[TABLE]

Therefore, in the case of simply supported beams and small $r$, the full unimodal dispersion relation, (30) can be expanded in the form

[TABLE]

whence

[TABLE]

which is identical to the asymptotic dispersion relation (33) given that $\beta=\beta_{R}$.

It is clear that the asymptotic dispersion relation (33) is much simpler than the full unimodal equation (30) obtained previously. Importantly, the asymptotic dispersion relation can be solved explicitly for $k$ as a function of $\omega$ and, therefore, manipulated and interpreted with greater ease. We compare the asymptotic and full unimodal dispersion relations in Figure 3 using the parameter values in Table 1.

The homogenisation scheme from (13) assumes that the wavelength is on a larger scale than the separation of the resonators, or equivalently $kl\ll 2\pi$. Subtituting $l$ from the the system parameters in Table 1, this condition is $k\ll\pi\mathrm{m}^{-1}$ for each of the illustrative figures, Figures 3, 4, 5 and 6. In addition, while the asymptotic method gives a solution which is valid for all $k$, we note that in the full unimodal system $\beta$ will become purely imaginary when $\omega/k>c_{2}$ and so no propagating solution can exist in this region of the plot.

Figure 3 shows good agreement between the dispersion curves for the asymptotic and full unimodal approaches; particularly near the Rayleigh sound-line and for lower frequencies, as expected. In addition, the asymptotes created by the flexural resonances of the beams extends the range of validity of the asymptotic approach beyond the regime of small $k$ and small stresses where one would expect the method to remain accurate. The flexural resonances associated with these asymptotes create fine deep sub-wavelength band gaps which have been shown to have a range of physical applications (see, for example, [23, 22, 28, 10]).

4 Beams on a Rail

Suppose now that the resonators are held at $z=0$ to a rail parallel to the surface of the half-space, which allows the base of the beam to move freely but changes the bending angle of the beam depending on the gradient of the surface. This allows the beam to impart a bending moment on the half-space but no horizontal forcing. There is therefore no need to conserve displacements and the horizontal stress must vanish. The reader is reminded that, as discussed in §2.1 on p.2.1, we consider the influence of the resonators under pure bending only. However, in contrast to the case of simple support beams considered in §3, the resonators exert a vertical force on the half-space through rotations in the beams. These rotations can be linked to the gradient of the vertical displacement through the kinematic continuity conditions imposed at the junction points. Whilst this junction condition may initially appear somewhat counter-intuitive, it can be understood in terms of so-called gyroscopic hinges [42].

To construct the junction conditions for the base of the beam, we must first determine the bending moment $M$ in the elastic half-space caused by the flexural deformation of the beam. For a beam with cross-sectional area $A$ the bending moment is

[TABLE]

If the cross-sectional area of the beam is small then the stress can be approximated by the expansion

[TABLE]

which, on evaluation of the integral, leaves only the second term. This therefore represents the leading order contribution in the calculation of the moments. Figure 2 illustrates how this bending moment in the half-space corresponds to the moment at the end of the beam. Furthermore, this figure also shows how, from elementary geometry, the vertical gradient of the half plane is equal to the horizontal gradient of the beam. Together these lead to the junction conditions

[TABLE]

The vertical stress at the junction point is

[TABLE]

where $V(x_{1},x_{2},t)$ is the vertical force exerted on the half-space by the beams and $u_{3,1}$ is the gradient of the vertical displacement along the surface. In particular, solving these in the usual way along with the free end conditions (6) and the beam equation (3) gives the vertical force,

[TABLE]

We will again start with the full unimodal treatment. As we did previously with $H(x_{1},x_{2},t)$, we will use the force distribution (13) and define $\hat{V}$ such that,

[TABLE]

Expressing $u_{3,1}$ in terms of the potentials (1) and making use of (8) the surface stresses on the half-space can be written in the form

[TABLE]

whence, the full unimodal dispersion equation is obtained in the form

[TABLE]

As might be expected, this dispersion relation is very similar to the one obtained for the simply supported case and is likewise difficult to manipulate and interpret. Therefore, in order to provide a clearer understanding of the system, we once again turn to the asymptotic approach.

To apply the asymptotic model to this system, where there is only a vertical stress, we will make use of (19) along the surface of the half-space, whence we obtain

[TABLE]

or equivalently, using (20) and considering time-harmonic waves, in terms of $\psi$

[TABLE]

Combining the surface potential relation (20) with $\hat{V}$ from (43) and the surface junction condition (48), we obtain the asymptotic dispersion relation

[TABLE]

As before in §3, the asymptotic dispersion relation can be obtained by expressing the full unimodal dispersion relation (46) in terms of $R(r)$

[TABLE]

and using the leading order expansion of $R(r)$ for $r\ll 1$, that is, for wave speeds close to that of a pure Rayleigh wave, to obtain

[TABLE]

which is consistent with the asymptotic dispersion relation (49) if $\alpha=\alpha_{R}$ and $\beta=\beta_{R}$, i.e. the decay rates are equal to those of pure Rayleigh waves.

As with the simply supported case, the dispersion relation derived from the asymptotic approach is much simpler than the one produced by the complete analysis. Indeed, the asymptotic dispersion equation is again quadratic and therefore yields an explicit expression for $k$ in terms of $\omega$. A plot of the solutions of the asymptotic dispersion equation along with those of the full unimodal solutions is shown in Figure 4; the parameter values used are again those detailed in Table 1. Figure 4 shows that, once again, there is good agreement between the asymptotic model and the full unimodal solution, particularly near the Rayleigh line. For all of the modes shown in Figure 4, the asymptotic approach predicts the asymptote, the overall shape of the dispersion relation, and the existence of a large band gap for each mode. Indeed, in addition to being a much more elegant approach, it is clear that the asymptotic model is capable of reproducing the main features of the full unimodal solution. Moreover, compared with Figure 3, Figure 4 exhibits significantly larger band gaps.

5 Fully Matched Beams

We now consider the case where both displacements and rotations (surface gradients), in addition to, horizontal and vertical tractions are coupled at the junction points. Although, as discussed in §2.1 on p.2.1 we only consider flexural motion of the resonators and neglect the compressional deformation, the vertical stresses in the half-space couple to the flexural modes in the beams through rotations at the junction point, as in §4; this is fully taken into account here. Predictably, this complex coupling leads to a more involved solution for both the beam motion and the surface wave dispersion equation. For simplicity, we also assume an ideal contact between the beams and half-space, neglecting any local contact strain effects.

In this case, the equations of motion for the beams (3) together with the free end conditions (6) remain unchanged. The junction conditions at the base of the beam $x_{3}=0$ are a combination of those previously considered (23, 40) with each of the quantities illustrated in Figure 2 conserved, giving,

[TABLE]

corresponding to continuity of displacement, rotation, and the balance of moments respectively. The associated stresses at the base of the beam are

[TABLE]

Solving these conditions in the usual way along with (3,6) gives,

[TABLE]

where, for convenience, we have introduced the quantities

[TABLE]

and

[TABLE]

As in the previous sections, we will begin with the full unimodal treatment and then develop the asymptotic model for comparison.

In contrast to the simply supported and rail conditions discussed in §4 & 3, we should not expect the full unimodal solution to intersect with the asymptotic solution at the Rayleigh line. In order for the asymptotic and full unimodal solutions to intersect at the Rayleigh line, the beam must behave as though the coupled end with the half-space is free; this imposes requirements on the displacement and rotation at the base of the beam that are incompatible with the required displacement and gradient for the half-space. We emphasise that this feature does not affect the full unimodal solution but does affect the accuracy of the asymptotic solution since the basic assumption underlying the asymptotic model, that stresses tend to zero at the Rayleigh line, is no longer valid.

To formulate the full unimodal problem we begin by writing the surface stresses in terms of potentials (8), and in terms of the distributed load (13), to obtain

[TABLE]

Writing the displacements in terms of their potentials (1) & (9) and combining (55) with (54), we obtain the dispersion equation

[TABLE]

As one might expect, this dispersion equation is significantly longer and more complicated than those obtained for the previous cases and, correspondingly, more difficult to interpret. Therefore we will, once again, develop an accurate approximation using our asymptotic mode in an effort to obtain a more transparent dispersion equation.

The asymptotic model of Kaplunov et al. [31] consists of two relations which each employ a single surface stress to produce a wave. However, in the present problem, the system has two stresses; one of which is normal to, and the other of which is tangential to, the surface. Therefore, the approach of Kaplunov et al. [31] must be adapted to account for the composition of two stresses into a single wave. In particular, we will assume that there will be two waves, each generated by one of the two stress-states. Since, for the asymptotic model both waves have a speed which is a small perturbation from the Rayleigh speed, it can also be assumed that these two waves travel at the same speed. Each of these two waves will have a corresponding pair of displacement potentials, $\psi_{1}$ and $\phi_{1}$, and $\psi_{2}$ and $\phi_{2}$ such that by linear superposition the overall displacement potentials $\phi$ and $\psi$ are given by a linear combination such that $\phi=\phi_{1}+\phi_{2}$ and $\psi=\psi_{1}+\psi_{2}$. Each pair of potentials then corresponds to one of the surface relations from (17, 48),

[TABLE]

Hence, assuming time-harmonic waves such that the potentials take the form (16), we obtain

[TABLE]

Combining (54) with (59), after writing the displacements in terms of potentials, yields the following dispersion equation

[TABLE]

As with the previous cases, we note that the asymptotic dispersion equation (60) is quadratic in $k$ and can therefore be solved explicitly. Whilst more difficult to interpret than the simply supported asymptotic solution, it is still a substantial improvement on the full unimodal dispersion relation. Figure 5 shows the dispersion curves for corresponding to both the asymptotic (60) and full unimodal (56) solutions. We observe that the asymptotic model predicts the main features of the dispersion curves, including prediction of the existence and general shape of the solution branch for each quasi periodic mode and also giving a close approximation of the intersections with the Rayleigh line. The key feature of the plot, and one present in both the asymptotic and full unimodal results, is the band gaps. This system is notable for the large range of frequencies which do not have a corresponding real wave number solution and will, therefore, not propagate through the system. While it is difficult to locate where the band gaps are from the full unimodal dispersion relation, they can be approximated by the intersections between the asymptotic model and the shear wave line. It is emphasised that, once again, this feature is absent in the previous analyses [22, 23, 37] as the flexural resonances of the beams were neglected.

This system also highlights some of the limitations of this model. In particular, and unlike for the simple supported and rail junctions conditions that involve a single stress, the asymptotics and the full unimodal solution are not guaranteed to intersect at the Rayleigh line. This is due to the way the stresses are treated in the two methods. In the asymptotic model there is no interaction between the two stresses and there will be Rayleigh line intersections only at,

[TABLE]

whereas for the full unimodal dispersion relation there is a cross multiplication of the two stress terms. As discussed earlier, the full unimodal and asymptotic curves will intersect at the Rayleigh line only if the total resultant stress at the junction point is zero. This, perhaps counter-intuitive, point can be seen by assuming that the two stresses have the form,

[TABLE]

Expanding these stresses as before yields,

[TABLE]

If we assume a coincidence at the Rayleigh line then, at this point, we replace $\alpha$ and $\beta$ with $\alpha_{R}$ and $\beta_{R}$ and use the surface relation (20) to obtain,

[TABLE]

Clearly this equation is equivalent to the multi stress junction condition (59) but with an added stress interaction term which the previous asymptotic model does not predict. For this example we can take the surface stresses (54) and put them into the required form which yields,

[TABLE]

Letting $\alpha=\alpha_{R}$ and $\beta=\beta_{R}$ and substituting these into the expansion (64) leads to the dispersion relation,

[TABLE]

which is the same as (60) but with an added term from the stress interaction. Due to cancellation in the stress interaction term this dispersion relation is quadratic on $k$. However, even for simple stresses, the dispersion equation arising from the aforementioned stress interactions will, in general, be a quartic polynomial.

Figure 6 shows the dispersion curves associated with the full unimodal solution (56) and those arising from the asymptotic approach accounting for the stress interaction (60); this should be compared with Figure 5. Figure 6 illustrates that accounting for the stress interaction forces the model to coincide with the full unimodal dispersion relation at the Rayleigh line. That said, comparing Figures 5 and 6 we observe that the dispersion curves arising from the asymptotic approach where the stress interaction is taken into account does not show a significant improvement for the rest of the plot, away from the Rayleigh line, over the standard asymptotic model. This indicates that, for this system, the stress interaction does not have a significant effect on the overall behaviour of the waves away from the Rayleigh line.

Again, comparing Figure 3 with Figures 5 and 6, it is apparent that the fully matched junction condition allows for the creation of significantly larger deep sub-wavelength band gaps compared with the simply supported junction condition. Therefore, the effect of the flexural resonances on the overall behaviour of the system can be increased by an appropriate choice of junction condition.

A subtlety associated with applying the asymptotic model to multi-stress systems is that, the asymptotic approach initial developed in [31] requires that the overall stress is small. It follows that for multi-stress systems, all individual stresses must be small in for the approach to be applicable. If one stress is always large while the other is small then the asymptotic model cannot be valid, even if the result crosses the Rayleigh line.

Finally, away from the Rayleigh line, the dispersion curves approach asymptotes associated with the flexural resonances of the beams; this is a known feature of sub-wavelength resonant arrays (see, for example, [11]). For the previous two cases of simply supported and rail junction conditions, the asymptotic model correctly reproduces the asymptotes. This led to the asymptotic model remaining accurate to a good degree even away from the Rayleigh line. In contrast, the fully matched junction (both uni- and multi-stress) does not accurately reproduce the necessary asymptotes. In particular, the asymptotes of the full unimodel solution are solutions of the transcendental equation

[TABLE]

However, the fully matched asymptotic models produces asymptotes at,

[TABLE]

It is clear that, the fully matched asymptotic solutions and the full unimodal solution will diverge as $k$ becomes large. Therefore, whilst this is a powerful technique for when both stresses are small, the results cannot be relied upon away from this condition.

Figure 6 also shows the mode shapes associated with each branch of the dispersion curves; the colour scheme shows the magnitude of the total displacement. Animated versions of these figures can be found in the supplementary material. The modes were computed, with the Finite Element software COMSOL Multiphysics®, using the exact formulation for an infinite periodic array of resonators as described in §2.1. It can be observed that each mode can be associated with the fundamental modes of the flexural resonators; this is particular apparent from the animated mode shapes that can be found in the supplementary material. Indeed, it is clear that the flexural motion of the beams dominates the behaviour of the system in this regime and emphasises that the flexural resonances can play an important part in such problems where slender beams rest atop an elastic half-space.

Moreover, Figure 7 shows the vertical components of the displacement in the beam ($U_{3}$) and half-space ($u_{3}$). It is clear that, in terms of compressional motion, the beams move as a rigid body in the vertical direction and, therefore, the dynamic effects of the compressional motion of the beams do not materially contribute to the overall behaviour of the system in this regime. This observation is particularly evident in the animated versions of Figure 7, which are available in the supplementary material.

6 Concluding remarks

An array of flexural resonators attached to the surface of an elastic half-space is analysed using an explicit model for the Rayleigh wave. The paper generalises previous considerations using both full unimodal [23] and asymptotic solutions for an elastic half-space in the case of an array of compressional resonators [37].

Three types of junction conditions are considered, and for each the asymptotic solution is verified using full unimodal solutions for the half-space. The first junction condition to be considered treats the resonators as simply supported at the surface of the half-space. Matching of the beam equation yields a horizontal stress only. The resulting dispersion curves exhibit a series of very fine band gaps typical of deep sub-wavelength arrays; the narrow nature of these band gaps means that making use of these stop bands in practical applications such as filtering and mode-conversion is extremely challenging. For example, deep sub-wavelength arrays cannot be used to design broadband filters.

The second considered junction condition treats the resonators as being supported by a freely moving rail at the surface of the half-space. In contrast to the previous system this yields a vertical stress only. The effect of this junction condition is to widen the band gaps in the spectrum. This ability to increase the width of the stop bands through a judicious choice of the junction conditions is an important feature that would allow, for example the implementation of broadband filters in the deep sub-wavelength regime. We also note that this effect will increase the contribution of the dynamic flexural behaviour of the resonators on the overall response of the system. Whilst sub-wavelength resonant arrays have found extensive use in the literature (see, for example, [7]), little attention has thus far been given to the junction condition joining the resonators to the substrate. Indeed, one application of the present paper is in the analysis of the interaction of surfaces waves with dense pine forests and recent experimental data suggests that the flexural resonances as well as the junction conditions may be an important factor [38].

For both the simply supported and rail junction conditions the asymptotic formulation shows good agreement with the full unimodal solution. The asymptotic solution obtained predicts both a coincidence with the full unimodal solution at the Rayleigh line, and also resonances arising from the solution to the beam equation. Moreover, the asymptotic solution has the distinct advantage of being much simpler and easier to manipulate than the full unimodal solution.

The final junction condition considered in this paper fully matches all variables at the junction. Unlike the previous two cases this results in both a vertical stress and a horizontal stress, requiring an alteration to the existing asymptotic formulation to account for both stresses While this asymptotic solution closely matches with the full unimodal solution it does not accurately predict the coincidence with the Rayleigh line or the location of the system resonances. A further addition to the asymptotic model, which accounts for the coupling of stresses at the junction point, produces a solution which does accurately predict coincidence with the full unimodal solution with the Rayleigh line, but still does not accurately predict resonances.

All of the systems considered produce band gaps in the region around the beam resonances, with the largest band gaps produced by the beam on a rail and fully matched junction conditions.

This is an important feature as, usually, sub-wavelength resonant arrays are only capable of producing very narrow band gaps [7] and, as a result, it is challenging to make use of these band gaps in the design of filters and other control devices. However, as we have demonstrated, it is possible to control the width of these deep sub-wavelength band gaps by adjusting the junction condition joining the resonators to the half-space. Indeed, the beam on a rail and full matched junction conditions lend themselves to feasible designs of mechanical structure that exhibit large band gaps and are capable of controlling the propagation of surface waves. In particular, the beam on a rail junction condition is associated with gyroscopic resonators, as discussed in [42], and which can be designed to exhibit the desired properties; the fully match condition is, in some sense, the most natural junction condition to impose as it most closely resembles the conditions of ideal contact.

This demonstrates the potential for such systems to be utilised in controlling the propagation of surface waves.

The explicit asymptotic formulation has the potential to be applied in more sophisticated systems where a full unimodal solution may be difficult or impossible to obtain since it reduces the vector problem in full linear elasticity to a scalar problem along the surface. It has been shown that for the systems considered the asymptotic formulation easily produces simple explicit solutions which match closely with the full unimodal solution.

The formulation discussed involves a surface scalar problem identical to the equation for transverse forces applied to an elastic membrane. Currently, experimental set-ups to model elastic half-spaces can be large and expensive [29] and this formulation gives scope for simpler experiments on membranes which model the behaviour at the surface of an elastic half-space.

This system could also be further developed by considering the near-surface effect of coupled stresses. We also remark that the effect of microstructure, including non-local effects, have previously been considered for near-surface wave motion [43] and, as might be expected, are asymptotically secondary for treating Rayleigh-type waves.

Funding

D. J. C. gratefully acknowledges the financial support of the EPSRC (UK) through programme grant no. EP/L024926/1. P.T.W also acknowledges Keele University for supporting their PhD studies.