Acoustic cavities in 2D heterostructures

Two-dimensional (2D) materials offer unique opportunities in engineering the ultrafast spatiotemporal response of composite nanomechanical structures. In this work, we report on high frequency, high quality factor (Q) 2D acoustic cavities operating in the 50–600 GHz frequency (f) range with f × Q up to 1 × 1014. Monolayer steps and material interfaces expand cavity functionality, as demonstrated by building adjacent cavities that are isolated or strongly-coupled, as well as a frequency comb generator in MoS2/h-BN systems. Energy dissipation measurements in 2D cavities are compared with attenuation derived from phonon-phonon scattering rates calculated using a fully microscopic ab initio approach. Phonon lifetime calculations extended to low frequencies (<1 THz) and combined with sound propagation analysis in ultrathin plates provide a framework for designing acoustic cavities that approach their fundamental performance limit. These results provide a pathway for developing platforms employing phonon-based signal processing and for exploring the quantum nature of phonons.

is also one of the modes with the high MoS2 partition index MoS2 (see Supplementary Table 1 below).
The solution of the form ( ) = exp[− ± ] exists for mixed-mode frequencies , and  is anticrossing splitting: where ( ) = sin ( ℎ ) is the eigenmode displacement,  is density, A is area, and h is the plate thickness 11 .
By using the out-of-plane stiffness C33 of MoS2 12 , the effective spring constants were assigned to each of the oscillators , = to match the resonant frequency of the fundamental thickness mode for corresponding plates in the uncoupled state (A, B). The pump-pulse in Supplementary Movie 1 peaks at 10 ps, while the plate is at rest at the start.
Being spatially centered at the monolayer step, the pump provides equal excitation for both stepabutted cavities and initiates an in-phase dilation. Phase differences accumulate due to the difference in thickness (and therefore in eigenfrequency) and are visible at later stages in the movie.

Supplementary Note 3: Frequency comb, h-BN/MoS2 bilayer, time-domain analysis
The bilayer stack is modeled as suspended, while the perimeter treated as low-reflectivity boundary. External strain that imitates the effect of the pump-pulse is applied to the MoS2 layer only. The spatial distribution for the "swelling" is Gaussian in the X-Y plane (rad(1/e) = 0.5 m) and homogeneous across MoS2 thickness (Z direction, normal to film surface). In an effort to reproduce temporal behavior of pump-generated strain, the external strain in the model is given a sharp rise (t  200 fs) at the pulse arrival time (t0 = 10 ps), followed by slow exponential decay ( = 180 ps 16 ). where the effect of boundary layer turned out to be prominent, as it provides significant mechanical isolation between the stacked nanoparticles 10 .
For some MoS2/MoS2 stacked structures, the added inner interface did not have a distinguishable effect on the Q value of the cavity. For others, the Q value decreased by a factor of two. As the difference in performance between the MoS2 stacked cavities and the mono-cavities appears close to the scatter within the mono-cavities themselves, we interpret this as strong evidence for the ability to build intricate 2D material acoustic structures where the interface plays a negligible role.

Supplementary Note 6: Coupled cavities in MoS2/h-BN/MoS2 tri-layer
Coupling, defined as  = (where B is a constant), was used as a fitting parameter to match the calculated frequencies of mixed modes ± with the experimentally observed tri-layer overtones that correspond to the anti-symmetric (2 nd overtone) and the symmetric (3 rd overtone) strain distribution across the stack. of quantized thickness, and is residue free (within our detection limits). We anticipate these features will be critical for low-temperature operation, where a high degree of control and low energy loss is necessary.

Supplementary Note 7: Phonon lifetime calculations
Phonon-phonon scattering: Lattice anharmonicity gives rise to phonon-phonon scattering, thus accounting for a variety of vibrational properties of materials including finite thermal conductivities, mechanical dissipation, and temperature-dependent linewidths/lifetimes. Taylor expansion of the crystal potential gives a lowest order anharmonic perturbation corresponding to three-phonon interactions 21,22 (see Supplementary Figure 12).
Each scattering conserves energy and crystal momentum, and its transition probability is determined by Fermi's golden rule: where ⃗ and j label a phonon's wave vector and polarization, ⃗ is the angular frequency, ⃗ are the Bose equilibrium populations, and ± is the anharmonic perturbation constructed from anharmonic interatomic forces constants (IFCs) and creation and annihilation operators 21,22 .
After application of the perturbation and algebraic accounting of combinatorial processes the transition probabilities are given by 23 : where N is the number of unit cells in the system, ⃗ is a reciprocal lattice vector, and the scattering matrix elements are given by: where ⃗ is the eigenvector of phonon ⃗ for the α th Cartesian component of the k th atom in a unit cell, mk is the isotopically averaged mass of the k th atom, ⃗ is the lattice vector locating the l th unit cell, and Φ , , are third-order anharmonic interatomic force constants linking atoms , , and in the origin, , and unit cells, respectively, as depicted Supplementary   Figure 13.
Summing over all possible three-phonon scatterings that conserve energy and momentum then gives the inverse phonon lifetime 21 : The 1/2 factor in the second term accounts for double counting of identical processes. The frequencies and eigenvectors that enter this quantum perturbation formalism are determined by diagonalization of the dynamical matrix: where Phonon-boundary scattering: Attenuation of phonons due to sample boundaries was modeled via the empirical formula 21 : where ⃗ = ⃗ / is the phonon velocity, h is sample thickness, z is the direction perpendicular to the sample surfaces, and is the specularity parameter 21 : Where = 2 / is the wavelength of the incident phonon and η is the root mean square (rms) roughness of the surfaces. To correlate our bulk calculations with the finite sample sizes from experiments we define ℎ = /2.
Other scattering considerations: Phonon lifetimes can also be limited by scattering from higher order anharmonic processes (e.g., four-phonon scattering), isotopic mass variance, and from sample-specific extrinsic defects (e.g., vacancies, grain boundaries). Here we briefly comment on the potential impact of such mechanisms when comparing calculations (which do not consider these) with measurements (which do not quantify these).
Like the anharmonic scattering calculations discussed above, phonon-isotope scattering can be calculated fully using first principles within the context of quantum perturbation theory and density functional theory 24,25 . Phonon-isotope scattering generally scales as ω 4 , thus for the ultralow frequency range considered here this mechanism is unimportant, unlike the situation when considering lifetimes of high frequency optic phonons 26 . We verified that phonon-isotope scattering was negligible in the ω < 1 THz frequency range for each material, and thus only include this scattering mechanism in Figures 5a and 5b of the main text, which considers the entire frequency spectrum of each material.
Similarly, other point defects are expected to have similar low frequency scattering behavior as phonon-isotope scattering, scaling as ω 4 27, 28 . Thus, we expect that vacancies (which are prevalent on Sulfur sites in MoS2 29 ) and other point defects will not play a significant role in the frequency range of our measurements.
Higher order anharmonic interactions have been shown to play an important role in limiting phonon lifetimes, even at room temperature, in strongly anharmonic materials 30,31 and in materials where three-phonon scattering is especially weak 32 . Nonetheless, calculations have demonstrated that at room temperature for ω < 1 THz four-phonon scattering is small compared to three-phonon scattering in a variety of materials 30,31,32,33,34 .

Supplementary Note 8: Density functional theory and computational details
Low frequency scattering rates: The harmonic and anharmonic IFCs were computed within the framework of density functional perturbation theory (DFPT) 35 using the generalized gradient approximation as implemented in QUANTUM ESPRESSO 36,37 and the D3Q package 38  Becke-Johnson damping 41 , which has been shown to reproduce measured lattice parameters and interlayer distances 2,42 . A comparison between calculated and measured lattice parameters is given in Supplementary Table 3 One strategy to increase the density of elastic energy within the 2D acoustic cavity is to shrink the cavity footprint. This is accomplished by employing sharply focused laser beams for exciting the vibrations. As the spot size of the diffraction-limited pump beam becomes comparable to the film thickness, the in-plane evolution of elastic waves gains significance. Theoretical analysis 45 shows that spatial variation of the photo-induced elastic strain (defined by the intensity profile of the pump beam) gives rise to symmetric S1 Lamb waves with in-plane wavevectors kL  2/Rbeam that propagate outward, taking out the energy stored in the cavity. The energy leak rate can therefore be governed by the group velocity of the elastic waves (as opposed to internal friction) and is defined by the Lamb waves' dispersion law. In the vicinity of the thickness mode resonance, i.e. in the limit k 0 the S1 wave dispersion exhibits parabolic nature ~ k 2 (see Supplementary Note 10), leading to group velocity proportional to film thickness: A parabolic fit to the FEM-calculated dispersion law 6 in MoS2 plates (Supplementary Figure 18) provides a value of vgr(kL= 2/Rbeam)  320 m/s for a 160L-thick film. The resulting escape time estimated as  = 2Rbeam/vgr  3 ns is consistent with the ring-down time measured for thicker films in Figure 6a of the main text and shows that elastic wave spreading can be a factor affecting cavity performance.
In order to get a more quantitative assessment for the escape time, we use time domain simulations that explicitly model the MoS2 plate response to a pulse-like elastic excitation (Supplementary Figure 19). The validity of this approach for estimating the Lamb waves escape time is limited to structures featuring optical penetration depth that greatly exceeds the film thickness, as the symmetric S1 mode is expected to dominate the vibrational spectrum under that condition.
We note that even though the estimates accounting for the lateral spreading of elastic energy agree well with our experimental data for relaxation times measured in MoS2 cavities in sub-100 GHz frequency range, extra effort is needed to make this comparison fully quantitative. A different experimental approach (e.g., similar to that described in 46 ) might be required in order to circumvent the timing limitations of our pump-probe setup. The accuracy of measuring ringdown times in excess of 1 ns is convoluted by the high repetition rate frep  1 GHz, which as  becomes longer can cause a partial overlap of multiple slowly-decaying traces generated by preceding pump pulses.

Supplementary Note 10: Lamb wave dispersion
Consider the dispersion relation for the symmetric Lamb modes 47 , where h = layer thickness, and k magnitude of the wavenumber parallel to the layer surface. The notation is: where we omit terms of order k 4 and higher. With kT and kL proportional to ω, a solution (albeit a complicated one) can immediately be written for k in terms of a transcendental function of ω. The solutions can be grouped into quasi-longitudinal solutions ωnL and quasi-transverse solutions ωnT.
We seek an analytical approximation to the inverse solutions ωnL(k). Such a solution was given explicitly for the modes ω1L(T)(k) in the case of near degeneracy of the modes ω1L(k) and ω1T(k).
To obtain an analytical solution in general, we follow the same procedure as regards ω: let ω = is the eigenfrequency for the nth longitudinal solution for k = 0. To consistently retain terms of order k 2 we must retain terms of order δω, δω 2 (in addition to those of order k 2 ) because the δω may be proportional to either k or k 2 . Define, where α = cL/cT, the ratio of longitudinal to shear wave speeds. Then we approximate the trigonometric functions appearing in the above dispersion relation as tan(kTh/2) ≈ tanκnT + αδωd sec 2 with m = 2n -1, an odd integer. This result is valid for all symmetric, longitudinal Lamb modes.
For sufficiently small δω (and correspondingly k) one may observe that the dispersion relation is always quadratic in k. The curvature may be positive, negative, or nearly vanish depending on the quantity κnT and α and this in turn is closely related to the degeneracy of S1L and S1T modes.
We note that the considerations above are given for Lamb waves in isotropic materials.
However, our numerical calculations (see example of dispersion curve in Supplementary Figure   18) indicate that the main features -parabolic dispersion at k0 and linear scaling with the film thickness h are preserved for transversely isotropic materials, such as MoS2.