Clamped seismic metamaterials: Ultra-low broad frequency stop-bands

The regularity of earthquakes, their destructive power, and the nuisance of ground vibration in urban environments, all motivate designs of defence structures to lessen the impact of seismic and ground vibration waves on buildings. Low frequency waves, in the range $1$ to $10$ Hz for earthquakes and up to a few tens of Hz for vibrations generated by human activities, cause a large amount of damage, or inconvenience, depending on the geological conditions they can travel considerable distances and may match the resonant fundamental frequency of buildings. The ultimate aim of any seismic metamaterial, or any other seismic shield, is to protect over this entire range of frequencies, the long wavelengths involved, and low frequency, have meant this has been unachievable to date. Elastic flexural waves, applicable in the mechanical vibrations of thin elastic plates, can be designed to have a broad zero-frequency stop-band using a periodic array of very small clamped circles. Inspired by this experimental and theoretical observation, all be it in a situation far removed from seismic waves, we demonstrate that it is possible to achieve elastic surface (Rayleigh) and body (pressure P and shear S) wave reflectors at very large wavelengths in structured soils modelled as a fully elastic layer periodically clamped to bedrock. We identify zero frequency stop-bands that only exist in the limit of columns of concrete clamped at their base to the bedrock. In a realistic configuration of a sedimentary basin 15 meters deep we observe a zero frequency stop-band covering a broad frequency range of $0$ to $30$ Hz.


Introduction
The desire to deflect, absorb or redirect waves is ubiquitous across many fields: electromagnetics, optics, hydrodynamics, acoustics and elasticity. Various techniques have been developed, often in one field but not in the others, and we are going to draw upon advances in optics and electromagnetic wave systems to develop a methodology for reflecting long-wavelength seismic waves in sedimentary basins. The motivations to do so are clear: according to the US Geological Survey there are millions of earthquakes every year worldwide, the vast majority are magnitude 3.9 or lower but more than 1000 measure 5.0 or higher on the Richter scale [1]. Ground vibrations, caused by even minor earthquakes, have an impact upon the structural integrity of buildings and similarly intrusive ground vibrations from urban train systems, subways, machinery such as piledrivers and roads often affect property values or land usage. These vibrations are not simply a nuisance, but small magnitude vibration due to machinery, or nearby railway lines, can cause significant damage to buildings, especially over time [2]. Furthermore, for buildings such as nuclear power plants and oil refineries, even a small level of damage can have disastrous consequences. Seismic waves consist of surface waves (elliptically polarized Rayleigh waves Unfortunately, the main drawback of this type of locally resonant structure is the difficulty in obtaining very large efficient stop bands; there is always a trade-off between the relative bandwidth and the efficiency of the attenuation, which is directly linked to the quality factor of the resonators. The frequency bandwidth of wave protection can be enlarged by considering arrays of resonant cylinders with different eigenfrequency for twodimensional stop-bands [42], or cubic arrays of resonant spheres [43] for 3D stop-bands, but such mechanical metamaterials would be hard to implement at the civil engineering scale. Buried isochronous mechanical oscillators have been also envisaged to filter the S waves of earthquakes [44]. A complementary approach employed in [45,46] is to draw upon the metamaterials literature that utilises subwavelength resonators arranged, in their case, upon the surface of an elastic half-space; this results in surface to bulk wave conversion, and surface wave filters with band-gaps but again at higher frequency than those required for seismic protection.
Identifying that the main objective of a seismic metamaterial must be to achieve a broad low-frequency stopband, or even better a zero-frequency stop-band, we turn to the apparently, distantly connected, field of thin elastic structured plates; these are of interest in terms of flexural waves connected with the vibration of shells. The simplest model, the thin plate model of Kirchhoff and Love is well-known [47,48] and is a fourth-order scalar partial differential equation for the vertical displacement field; this is a dramatic simplification over the full vector elastic system [49] and, for homogeneous plates, this simpler model is, in theory, valid for very thin plates with a thickness less than 1/20th of the typical wavelength [50]. None the less surprising accuracy can be obtained in structured plates, where this thickness restriction can be relaxed, and recent detailed comparison of theory and experiment [51] demonstrate that conclusions from the Kirchhoff-Love (K-L) theory carry across into thin plate modelling even at much higher frequencies where it might naively be thought to be invalid (typically in practice a plate thickness less that 1/5th of the wavelength provides good approximations for homogeneous plates); figure 1 shows dispersion curves from K-L theory alongside those of full-elasticity, both sets calculated using finite element (FE) software [52]. These are shown for a structured plate consisting of small clamped circular regions (1.5 mm radius) arranged on a square lattice of pitch 1 cm and for a 0.5 mm thick plate of duraluminium ( ). An important, indeed critical for our purposes here, observation, is that this system has a zero-frequency stop-band. As noted earlier there is close correspondence of full vector elastic calculations with those from K-L theory far outside the range of frequencies one might usually associate with thin plate K-L theory; the reason, as noted, in [51] is that what actually matters is not a constraint from homogeneous plates, but how the wavelengths in the periodic structured system compare to the plate thickness and they are actually large. Figure 1 takes advantage of the periodicity of the structured elastic plate, as is well-known in solid state physics [53] Bloch's theorem means that, for an infinite array, one need only consider the wavenumbers in the irreducible Brillouin zone (IBZ) which, for a square lattice, are those in the triangle XM shown as the inset to figure 1(a) and for clarity we show the frequency dependence versus wavenumbers going around the edges of the IBZ. We will exclusively use square lattices in this article and comment upon other lattice geometries in section 3 (see also figure 11).
Having successfully identified a scenario, all be it in a different situation of thin elastic plates far removed from the seismic application of thick elastic substrates, giving a zero-frequency band gap we use this to design a seismic metamaterial to have these characteristics [54]. The aim of this article is to investigate structuring soil, as shown in figure 2, to protect a building or portion of the surface. The structuration consists of columns, in a layer of soil, that are clamped to underlying bedrock; the columns are arranged in a periodic fashion. The key difference from all previous studies is that, motivated by thin plate calculations of figure 1 and the resultant zero frequency band gap, we consider the influence of clamping columns to the bedrock.

Results
We note that many seismic applications such as ground reinforcement, concern layers of softer soils overlying more solid bedrock and modern civil engineering processes allow columns to be clamped, that is, rigidly attached, to the bedrock [55]. The layers are no longer thin enough to employ K-L thin plate theory and the effect of depth (equivalently thickness of the plate) is now important. We explore the potential of seismic metamaterial devices shown in figure 2; the side view shows a structure atop a soil layer of finite thickness that overlays the bedrock; columns clamped to the bedrock puncture the soil layer and reach to the surface, or close to the surface, a top view shows the array of columns encircling the building to be protected.
For elastic waves in a finite layer, overlying a deep infinite substrate, two types of wave modes exist [48,49]. Considering shear, SH, polarised waves one finds Love waves and these are characterised by motion concentrated in the layer and out-of-plane displacements as shown in figure 3. Love waves only exist in the layered system and not if one simply has a deep substrate alone.
An elastic substrate, of infinite depth, supports Rayleigh waves that have exponential decay with depth and whose displacements are elliptically polarised in the plane. In a finite layer depth, provided the frequency is high enough, Rayleigh waves exist as their displacement has effectively vanished at the base of the layer. As the frequency decreases the waves interact with the base and one obtains Rayleigh-Lamb modes that arise from considering an elastic waveguide. In the situations we consider we have Rayleigh-like modes that are not dominated by basal interactions and we refer to these as Rayleigh waves. Figure 3 shows the polarizations and displacements of a typical Rayleigh-Lamb mode.
We consider small strain seismicity, [56], with a rate of deformation (γ such that ), in this situation the duration of seismic disturbance, for this type of earthquake, is sufficiently short to accept the validity of the hypothesis of elastic behaviour for the soil. We can therefore take realistic soil, rock and depth parameters from geophysics with the Young's modulus,   rigid material beneath it (i.e. the bedrock interface is modelled with Dirichlet boundary data, which is zero displacement there).
Given the geometric and physical complexity we proceed to numerical simulations. We have three phases: the bedrock, the matrix material and the inclusions. Each of which is modelled as linear, isotropic elastic media represented by the elastic Navier equation [57] with time-harmonic motion considered at fixed frequency i.e. in the absence of any sources we take, in each phase, There is continuity of stress and displacement across interfaces between the media and , l m are Lamé coefficients, ρ the density, ω the angular wave frequency and u u u x x x u x , , , We discretize the weak form of (1) using FE methods and in particular we utilise COMSOL [52]. For band structure calculations, we take advantage of the periodicity of the system in the horizontal x x , 1 2 ( )-plane to consider a single elementary cell [58]. For this we use the Bloch-Floquet theorem i.e. we assume that solutions of (1) are such that = ( ) is the Bloch (or momentum) wavevector in 2D reciprocal space (both a and k lie in the horizontal plane). We then construct the Bloch band spectrum, see for instance [5,9], to obtain dispersion curves of frequency versus Bloch wavenumber. The dispersion surfaces computed are critical in guiding our understanding of the Bloch wave behaviour, we refer to [51] for structured media displaying elliptic, parabolic and hyperbolic shapes of dispersion surfaces leading to extreme elastic wave control at certain frequencies revealed by so-called high-frequency homogenization [59]. Importantly, when k k k 0 , the quasi-static limit of the band spectrum allows one to deduce effective properties of the periodic structure (such as anisotropy) through the slope or curvature of dispersion surfaces in the neighbourhood of Γ point [9,60] depending upon whether they have a conical or parabolic shape. One can then interpret the effective behaviour of the low frequency Bloch waves in terms of their effective group velocity, or effective mass [61]. However, when some Dirichlet data, u 0 = , is set on some domain within the IBZ, corresponding to a clamped (possibly very small) inclusion, one infers from the maximum principle [9,60] that u is null everywhere in the IBZ, and is not an admissible eigenfield. This seemingly minor remark has important practical consequences: it shows that a periodically constrained elastic structure has a zero-frequency stop-band, which can be used to reflect arbitrarily low frequency elastic waves. Making use of the symmetries of the periodic structure, it is common to consider the dispersion properties of Bloch waves by letting k vary only along the edges of the IBZ [4], and so one need only compute dispersion curves, and no longer dispersion surfaces, which greatly reduces computational effort. However, this needs to be done carefully, otherwise edges of stop-bands might be erroneously estimated [62][63][64].
First, simply taking into account just the finite soil layer atop the bedrock, and considering waves in this layer, leads as shown in figure 4(a), to a small zero frequency stop-band with the cut-off frequency at 2.7 Hz. This is to be expected: if one considers this layered system to be a waveguide then fixing one wall to be rigid will automatically shift the modal cut-off frequencies and no modes will propagate for extremely low frequencies.
One point of note is that results, not shown, also demonstrate a weak dependence upon soil depth with the cutoff frequency rising as the depth decreases and vice versa. Although instructive, this case does not lead to the desired wide zero-frequency stop-band. Turning now to another idealised situation, that of pillars (of radius 0.3 m) that are completely rigid, and immovable all along their length, and that are buried in the bedrock, we see in figure 4(b) that an extremely wide zero-frequency stop-band stretching all the way to 30 Hz is achieved. This is unphysical as it is not possible to clamp a column all along its length, but it does demonstrate that there is potential to create a very wide band gap. The other issue is that dispersion curves are reliant upon an infinite array and for realistic systems the array will be finite.
To further motivate the importance of a zero frequency stop-band and illustrate its effect, we simulate normal incidence upon a finite array protecting a region and generate the transmission spectra. The computational domain is shown in figure 5, which consists of the clamped columns (piles) in the soil attached to bedrock with waves incident from the line source. To ensure that there is no spurious effect from the size of the computational domain elastic perfectly matched layers (PMLs) [65] are employed to prevent spurious reflections in the propagation and reflection directions. Since normal incidence is studied we will consider a single row of columns ( figure 5(b)) and then use periodic conditions spanwise. We again proceed numerically in the time harmonic situation and solve the elastic Navier equation using FE methods and in particular utilise COMSOL [52] and use a line source at the surface (invariant across the domain) to initiate either Rayleigh surface waves or SH polarized waves; we use Rayleigh excitation in figure 6.
The effect of the clamped columns is most unequivocally seen by looking at the transmission spectra, that is the transmission computed at the soil/PML computational interface at the opposite side of the columns from the line source normalised by the amplitude of the source. The left hand panel of figure 6 shows very clearly that the fields are highly attenuated beneath the cut-off frequencies of each case, at 2.7 and 30 Hz respectively. The vertical scale is in decibels and so this represents the attenuation that would be seen, and in the clamped column case this is for just four clamped columns and is highly effective at blocking incoming waves. This is further exemplified by examining the physical fields as shown in figures 6(a), (c) for the clamped bedrock and (b), (d) for the idealised clamped, perfectly immovable, columns. The effect of the columns is clearly seen in the fields with almost perfect reflection at the frequencies illustrated.
The results above are not physically feasible as no columns can be completely fixed and rigid all along their length (i.e. zero displacement on all of its boundaries); this provides an ideal upper bound on the possible efficiency of this design, but it is highly encouraging that clamped columnar structures could, at least, achieve part of this behaviour.
We now turn to realistic designs using illustrative columnar materials and geometries. We choose two examples and treat them in detail, first cylindrical columns of steel (density 7850 kg m −3 , Poisson ratio 0.33), 1.2 m in diameter, in a soil layer of 15 m and bedrock of depth 5 m; the columns are buried in the bedrock to a depth of 80 cm but otherwise free to vibrate and interact with the soil and bedrock. For an infinite array placed at the vertices of a square array (of pitch 2 m) the dispersion curves are shown in figure 7(a). It is clear that this has a  zero-frequency stop-band up to 4.5 Hz, which is not dramatic in extent but is none the less an improvement over the unstructured soil. It is, however, far from the optimal case shown in figure 4.
It is well known in the phononic crystal community that very broad band gaps can be designed by the use of cross-shaped inclusions [66], reminiscent of lattice models [67], so we utilise this shape here and further augment this by linking the columns. We consider a cross-shaped inclusion having the same cross-sectional area as the cylinder thereby comparing inclusions that have the same filling fraction. In this second configuration we take a cylindrical column of radius 0.2 m and each column is then linked to its nearest neighbours by struts (long steel plates of 0.2 m in thickness), see figure 7(b), and this also gives the columns a structural link to its neighbours. This additional reinforcement then leads to a dramatic enhancement in the protected frequency range, which is now up to 26 Hz. The dispersion curves of figure 7 are for an infinite array and do not reveal the potential efficiency of the design, as for this we need the transmission spectra.
We consider both shear wave and Rayleigh wave incidence, as they have different character. Shear wave incidence is illustrated in figure 8 that shows the attenuation of the clamped cylinders and the clamped crossshaped inclusions; the latter clearly perform much better and the attenuation is strongly enhanced below 26 Hz. Notably this performance is for normal incidence upon a finite array only four columns deep, and improves further if the number of columns is increased. Movies constructed using the time harmonic solutions to show the time evolution of an incoming wave incident upon the columns are shown in the supplementary material (available online at stacks.iop.org/NJP/19/063022/mmedia). As an aside in figure 8 there is a dip in the transmission spectra for the clamped cylindrical columns at approximately 40 Hz, which is due to the polarisation of the source and its interaction with the phononic crystal; this can be identified in figure 13(A) when taking into account the polarization state of shear waves.  Rayleigh wave incidence is illustrated in figure 9 showing again attenuation of the clamped cylinders and the clamped cross-shaped inclusions and the associated displacement fields. The attenuation is not as dramatic as for shear waves due to the stronger coupling of Rayleigh waves into the structure and the finite structure as a whole can vibrate, i.e. the signal at 15 Hz, none the less strong attenuation is achieved for low frequencies.
We now address the issue of whether there is an optimal filling fraction. In figure 10 we show the upper extent of the zero frequency band-gap for the cross-like inclusions and cylindrical inclusions. In the latter case the upper extent increases monotonically as the filling fraction (also called the soil substitution in the civil engineering literature) increases, as one would expect. For the cross-like inclusion there is a non-monotonic behaviour, all be it very gradual, with a maximum for a filling faction of around 0.6. In practical terms, a soil substitution rate beyond 10%-20% is quite substantial in civil engineering, and in practice the smaller filling fractions are more relevant.
It is natural, of course, given that we have presented in detail a design on a square lattice to consider whether changing the underlying lattice structure to another Bravais lattice creates fundamental changes. It does not, it alters the detail of the extent of the zero frequency band gap but not the physics behind this. As a further illustration we show, in figure 11, the dispersion curves for the steel columns (radius 0.3 m and buried 80 cm in the bedrock) placed on a honeycomb lattice with plate attachments along the honeycomb. Again one notes an extensive zero frequency bandgap extending up to 20 Hz, thus there is unfortunately no specific advantage of using a honeycomb, instead of square, array. Similarly we also considered hexagonal arrays with a band-gap extending to 18 Hz and it is clear that unlike the case of phononic crystals or locally resonant materials, the hexagonal lattice is less efficient compared to the honeycomb or square lattices.
A key point of the design, further emphasised in [68], is the burying of columns in the bedrock and connectivity between columns and to demonstrate how important this is we consider an array of horizontal steel     plates of thickness 0.2 m that splits the soil layer into a laminate structure and a column of radius 0.3 m; the band diagrams for this structure are shown in figure 12 that also shows the polarization p 2 (along x 2 ) defined as Finally, we investigate the sensitivity of the upper-extent of the zero-frequency band gap to attachment in the bedrock and to modifying the square lattice structure of a steel rod attached to its neighbours by steel plates. We do this in figure 13 and show the dispersion curves only in the X G direction, that is for normal propagation to an array, and consider 5 cases: A, is the steel column of radius 0.6 m of figure 7(a) buried 80 cm in bedrock showing a small zero-frequency band-gap, (B), is the reference case of the steel column of radius 0.3 m joined to its neighbours in a square lattice, buried 80 cm of figure 7(b) and sharing the same filling fraction as A. Panels (C) and (D) are not buried in the rock and merely are joined at the soil/bedrock interface to the bedrock and then have the central steel column absent respectively; the upper extent of the band-gap decreases. (E) goes one stage further and decouples the joining plates at the edges of the periodic cell and then the band-gap width decreases dramatically. It is clear that the combination of having the steel column, buried in the bedrock, and the coupling of each column to its neighbours via the steel plates is essential to obtain the broad ultra-low frequency stop-bands.

Discussion
In geotechnics, so-called composite soils made of inclusions inserted in a matrix of soil are mainly developed using a pseudo-static analysis; the dynamic loading is converted into an equivalent static loading and the objective is to obtain a higher value of Youngs modulus or shear modulus G (or μ) in order to reduce the deformation of the soil-deep foundation system [69,70]. Homogenization techniques have been also applied to composite soils [71,72] and recently developed for dynamic loading [73]. The emerging field of seismic metamaterials, based on wave physics, enables us to revisit several longstanding problems of earthquake protection from this fully dynamic point of view.
In this vein, we have considered here how one could use these wave physics-based ideas to protect specific areas from low-frequency vibration. We have demonstrated, conclusively, that it is possible to design realistic seismic metamaterial devices, in the sense of the structures that form the material being subwavelength and in the limit of small strains such that we have linear elastic conditions. The desired performance is achieved by creating zero-frequency stop-bands using clamping of columns to underlying bedrock. Clearly the challenge is now to build such structures in the spirit of [40] and we anticipate that these results will motivate large-scale experiments to verify this since there are potential applications in seismic defence structures [68]. It is notable that zero-frequency stop-bands will arise in other areas of physics whenever periodic arrays of inclusions have zero field, i.e. Dirichlet data, on their boundary see for instance, [74], in the context of electromagnetic waves or [75] for clamped inclusions in full vector elasticity. One can anticipate this, following [74] any Floquet-Bloch eigenvalue problem involving inclusions with zero Dirichlet data on part of their boundary, will display zero frequency stop-bands, since the existence of an eigenfield associated with a zero frequency would imply by the maximum principle that the eigenfield is null everywhere. Therefore, any wave (electromagnetic, acoustic, hydrodynamic, elastodynamic) will be prohibited to propagate within a periodic structure at very low frequencies provided inclusions can be modelled with zero Dirichlet data (e.g. infinite conducting inclusions in the case of electromagnetics and clamped inclusions in the case of elastodynamics).
Although we have concentrated here on the lowest dispersion curve, at high frequencies, the dispersion curves and asymptotics based around generating dynamic effective media [75][76][77] suggest that interesting features will also occur. These include ultra-directivity whereby a periodically pinned plate behaves like an extremely anisotropic effective medium at frequencies associated with flat bands on dispersion diagrams [51].
Finally, we note that we assumed a linear elastic model for soil and columns of concrete, however, the former has visco-elastic features, that could also be implemented e.g. via a Kelvin-Voigt contact condition at the soilcolumn interfaces. Intuitively, we expect that such a model would reduce the local resonances, and the extremely high damping seen in figures 5 and 8 around certain frequencies, would be less pronounced.

Acknowledgments
RC thanks the EPSRC for their support through research grants EP/I018948/1, EP/L024926/1, EP/J009636/1. SG and YA are thankful for an ERC starting grant (ANAMORPHISM) that facilitated the collaboration with Imperial College London.
We also thank A Diatta for his help in the implementation of Cartesian PMLs in elasticity.

Appendix. Material values and parameters
For clarity and convenience we summarize the values used in the text.