Topologically protected elastic waves in phononic metamaterials

Surface waves in topological states of quantum matter exhibit unique protection from backscattering induced by disorders, making them ideal carriers for both classical and quantum information. Topological matters for electrons and photons are largely limited by the range of bulk properties, and the associated performance trade-offs. In contrast, phononic metamaterials provide access to a much wider range of material properties. Here we demonstrate numerically a phononic topological metamaterial in an elastic-wave analogue of the quantum spin Hall effect. A dual-scale phononic crystal slab is used to support two effective spins for phonons over a broad bandwidth, and strong spin–orbit coupling is realized by breaking spatial mirror symmetry. By preserving the spin polarization with an external load or spatial symmetry, phononic edge states are shown to be robust against scattering from discrete defects as well as disorders in the continuum, demonstrating topological protection for phonons in both static and time-dependent regimes.


Supplementary Table 1| Phononic Metamaterials Effective Stiffness Tensor
Stiffness tensor components and mass density of solid and metamaterial aluminum. The second column is extracted from fitting to the dispersion of low-frequency modes in Christoffel model.

Supplementary Note 1 | Effective Medium Theory
The metamaterial with sub-wavelength patterns used in this work was modeled as an effective uniform medium. This approximation holds very accurately as long as the size of the patterns remains much smaller than the wavelength of interest.
To find the effective elasticity parameters of the metamaterial, we studied the dispersion of very-lowfrequency modes of the patterned structure. These modes are described by the Christoffel model of a continuum medium. According to this model, the dispersion of the modes are described by the Christoffel tensor : where is the continuum mass density, ⃗ ⃗ the displacement field, the frequency, and the propagation constant, of the eigenmode. The Christoffel tensor is related to the stiffness tensor and the propagation direction: Note that here we have used the Voigt representation 1 to reduce the number of indices. Supplementary Table 1 lists the nontrivial components of the stiffness tensor for both solid unperforated aluminum and an aluminum metamaterial made by sub-wavelength patterning as shown in Supplementary Figure 1a. Note the strong degree of anisotropy of the perforated structure by comparing the two columns of the table. This strong anisotropy was a key in achieving the degeneracy between the otherwise highly non-degenerate symmetric and anti-symmetric modes of an elastic slab.
Supplementary Figure 1b shows how well the results of the effective medium simulations match the exact solution in the frequency range of interest (<200 kHz). Even though the fitting parameters were obtained by considering the very low frequency modes, as can be seen from the dispersion diagram the match works quite well even for high frequecnies. In very high frequencies (>200 kHz), the perforations start to play signicant role and some resonance features can be seen up in the dispersion diagram. However, as long as we are interested in the wavelengths as large as the unit cell size (here, 9 times the period of the subwavelength holes), no distinction may be drawn between the two approaches. This is why we confidently used the effective medium theory in solving the large-domain simulations, whose computation can be very demanding and frequenctly unfeasible.

Supplementary Note 2 | First-Order Perturbation Theory
First-order perturbation theory can be used to explain the effects of a small perturbation to the eigenmodes of an elastic system. In our case, the system is perturbed by carving out some of the elastic material, resulting in a change of the resonator's free external boundary. The boundary with vacuum is characterized by zero traction ̂⋅ = 0. Below we will explain how we can construct a new set of eigenmodes from the unperturbed ones such that they satisfy this new boundary condition at the shifted boundary.
Both unperturbed and perturbed eigenmodes satisfy the following Cauchy elastic equations: where T is the stress tensor, S the strain tensor, ⃗ the velocity vector, the angular frequency, the density, and = 1/2 ( + ) is the symmetric gradient operator. The only difference between the perturbed and unperturbed cases is in the boundary conditions. To find the projection of the perturbed eigensolutions on the unperturbed ones, we multiply (from the left side) the first line of Supplementary Eq. 3 by ⃗ * and the second line by * . Here, ⃗ and are the velocity and strain fields of the m-th unperturbed eigenmode.
After integrating the result of these dot products over a volume , we find the following projection equations: where is the surface enclosing the volume . Above, we used the vector identity ⃗ ⃗ ⋅ ⋅ = − ⃗ ⃗ ⋅ + ⋅ (⃗ ⃗ ⋅ ). Let us assume that is the volume of the perturbed structure. Therefore, everywhere on the surface , the boundary condition of zero traction holds ̂⋅ = 0. Note that this is only true for the perturbed solution, since the surface is only a virtual boundary for the unperturbed system. Summing the two lines in Supplementary Eq. 4, we obtain the expression for the new eigenvalues: The surface integral can be rewritten into a volume integral proportional to the small change in the volume . Reversing the sign of surface normal so that the volume enclosed by S is Δthe removed volumeand using the above-mentioned vector identity once again, one can rewrite the surface integral as a volume integral: Combing Supplementary Eqs. 5 and 6, we reach at the final expression for the change in the eigenvalue due to the removal of the volume by Δ : Supplementary Eq. 7 is closely related to the Slater theory of perturbed electromagnetic cavities. 2 Assuming that the perturbation is small, the new eigenfunctions can be expanded in terms of the unperturbed ones: ⃗ = ∑ ⃗ , and = ∑ . In this case, Supplementary Eq. 7 can be recast into an eigenvalue problem:

Case of 4 Modes Representing Two Overlaid Dirac Bands
The first-order perturbation theory developed in the previous section can be used to treat degenerate cases as well. In the system considered in this paper, there are 2 overlaid Dirac cones, amounting to 4 degenerate modes at the K point (bands 8 to 11, as shown in Supplementary Figure 3) and quasi-degenerate in its vicinity due to having approximately the same group velocity. Here, in this section, we will develop a lowenergy Hamiltonian describing the hybridization between these modes in the close vicinity of the K point.
We assume that the linear dispersion of the modes is engineered such that they all have the same Dirac In the main text we demonstrated that by adjusting the slab thickness, it is possible to achieve two sets of double degenerate eigenmodes with a bandgap in between such that the slope of the dispersion is zero (i.e., It is clear that the off-diagonal elements of the matrix induce hybridization between the LCP/RCP modes of the A and S class of modes. The hybridized modes that are the new unmixed eigenfunctions of the system can be found from the following unitary transformation: Under this transformation, the new Hamiltonian (obtained by ′ ′ −1 ) exactly matches the Kane-Mele low-energy Hamiltonain, 3 This can be written in a compact form using the Pauli matrices: where ̂, ̂, and ̂ are the band, inter-valley and pseudo-spin subspace Pauli matrices. Although, we didn't explicitly examine the intervalley subspace (K and K′), the Hamiltonian given in Supplementary Eq. 15 can be directly deduced from the time-reversal symmetry of the system. To put it differently, it is easy to show that the Berry phase of each band around the K point is × ( ) and according to the time-reversal symmetry, the overall Berry phase must be zero, hence, at the K′ point, the mass term should reverse sign.
It is important to note that by reversing the perturbation (carving out the other side of the crystal face instead), the mass term changes sign. This can be easily verified from the unperturbed modes profile and the definition of the tensor, given below Supplementary Eq. 8.

Supplementary Note 3 | Topological Character and Spin Chern Number
The dispersion diagram ( ), historically assumed to contain most of the useful information about the system, lack the information about the system's topology. Topological numbers can only be calculated from the eigenfunctions themselves. Examples of such topological measures are Chern number and Berry phase.
Here is the detailed description of how we calculate these numbers for the elastic waves.
Berry connection (or potential) of band n is defined as ⃗ ⃗

Supplementary Note 4 | Effects of Leakage and Attenuation Caused by Interaction with Viscous Fluid
In the main text, we assumed that the external boundaries of the system are free, that is, when the background environment is vacuum. This is justified due to the significant impedance mismatch between our crystal made of Aluminum and air as the background. Below, we will show that loading of the system with air, or even water, has little effect on the dispersion of the topological elastic modes. We found the interaction with a viscous fluid gives rise to rather marginal effect of dissipation of the bulk and edge modes.
By considering realistic values for the fluid viscosity, it was found that the lifetime of the modes (in the unit of the elastic wave period) exceed the value of 10,000 in the spectral range considered (~130 kHz).
Since the viscosity scales as the frequency squared, this results in a greater impedance mismatch for higher frequencies.
One therefore may expect that the modes lifetimes would be even higher at higher frequencies.
Nevertheless, the topological nature of the modes remains intact by these loss channels, which also agrees with recent studies of topological states in open and dissipative condensed matter systems 4 .
The acoustic pressure field in a homogenous medium satisfies the wave equation 5