## Abstract

Interfaces impede heat flow in micro/nanostructured systems. Conventional theories for interfacial thermal transport were derived based on bulk phonon properties of the materials making up the interface without explicitly considering the atomistic interfacial details, which are found critical to correctly describing thermal boundary conductance. Recent theoretical studies predicted the existence of localized phonon modes at the interface which can play an important role in understanding interfacial thermal transport. However, experimental validation is still lacking. Through a combination of Raman spectroscopy and high-energy-resolution electron energy-loss spectroscopy in a scanning transmission electron microscope, we report the experimental observation of localized interfacial phonon modes at ~12 THz at a high-quality epitaxial Si-Ge interface. These modes are further confirmed using molecular dynamics simulations with a high-fidelity neural network interatomic potential, which also yield thermal boundary conductance agreeing well with that measured in time-domain thermoreflectance experiments. Simulations find that the interfacial phonon modes have an obvious contribution to the total thermal boundary conductance. Our findings significantly contribute to the understanding of interfacial thermal transport physics and have impact on engineering thermal boundary conductance at interfaces in applications such as electronics thermal management and thermoelectric energy conversion.

## Introduction

Interfaces impede heat flow in micro/nanostructured systems, which highlights the growing importance of thermal management of microelectronics and devices for energy conversion and storage^{1,2}. Thermal boundary conductance (TBC, *G*), the property characterizing thermal transport across an interface, is defined by a temperature drop (Δ*T*) across an interface for a given heat flux (*J*) as *G* = *J*/Δ*T*^{3,4}. For insulator or semiconductor-related interfaces, lattice vibrations (phonons) are the dominant heat carriers for interfacial thermal transport. At the microscale, thermal transport across an interface is usually described by the phonon gas model as phonons impinging the interface, which results in a portion of the energy passing through the interface as characterized by a transmission coefficient. TBC is thus often described as \(G=\left(\mathop{\sum}\limits_{j}\int {v}_{\omega }{c}_{\omega }{t}_{\omega }d{\omega }_{j}\right)/4\), where \({v}_{\omega },{c}_{\omega }{,{t}}_{\omega }\) are, respectively, group velocity, heat capacity, and transmission coefficient of the incident phonon with frequency \(\omega\), and polarization *j*^{4}. A critical task in the interfacial thermal transport research has been calculating the transmission coefficient^{3,5}. The acoustic mismatch model (AMM) and diffuse mismatch model (DMM) are traditional models used for this purpose, and they calculate the transmission coefficient based on the phonon properties of the bulk materials making up the interface^{6,7,8}. However, such treatments are problematic, as they ignore any microscopic details of the interface. It has been proven that factors such as bond strength, interface roughness, and atomic mixing can significantly impact TBC^{1,9,10,11,12}. Atomistic Green’s function (AGF) is able to explicitly include realistic atomistic structures in the TBC calculation, which can yield better predictions compared to DMM or AMM^{13}. However, phonon transmission coefficients calculated using AGF are still expressed as a function of phonons of the bulk materials^{11,14,15}.

In recent years, several molecular dynamics (MD) studies^{16,17,18,19} have shown that interfacial phonon modes can play important roles in interfacial thermal transport. These modes are mainly localized at the interfacial region due to the special interatomic bonding environment not seen in bulk materials. However, these MD simulations either were based on empirical interatomic potentials that were not specifically developed for simulating interfacial thermal transport or employed approximations for interfacial interactions despite the use of first-principles-calculated force constants^{16,17,18,19}. It remains unclear whether such interfacial phonon modes indeed exist at realistic interfaces, which calls for the experimental observation of such modes. The univocal confirmation of the existence of interfacial phonon modes is going to be critical to correctly understanding interfacial thermal transport physics and potentially engineering TBC. However, to date, experiments have been impeded by the technical difficulties related to probing localized modes at relevant size scales and having a properly chosen interface that is free of other factors (e.g., phonon polariton in polar materials, roughness, etc.)^{5,20,21,22}.

In this work, we confirm the existence of interfacial phonon modes by growing high-quality Si-Ge interfaces using molecular beam epitaxy (MBE) and detecting the localized modes with Raman spectroscopy and high-energy-resolution electron energy-loss spectroscopy (EELS) in a scanning transmission electron microscope (STEM). As both Si and Ge are nonpolar materials, there are no strong delocalized phonon polariton modes in the acquired EELS signal, which enables the spatial resolution of the characterization to be atomic scale^{20}. We further confirm that the detected localized modes from these experiments are indeed interfacial phonon modes using MD simulations with a high-fidelity neural network potential (NNP) trained by first-principles calculations specifically for the Si-Ge interface. The calculated TBC values from MD agree well with the experimentally measured values using time-domain thermoreflectance (TDTR). Finally, spectral analysis in MD shows that the interfacial phonon modes can obviously contribute to the overall TBC, despite its limited population and localized nature.

## Results

To fabricate interfaces that are clean and sharp for interfacial phonon mode detection, we chose the relatively low lattice-mismatched Si-Ge interface and grew high-quality interfaces using MBE, which enabled a precise control of the sample structure. Two sample architectures were grown, each fabricated for interfacial phonon mode detection and TDTR measurements, respectively. The detailed growth processes can be found in the “Methods” section and Supplementary Information (SI). Sample 1 consists of a 10 nm Si cap on a 300 nm Ge layer grown on a Si substrate and it is used for Raman (Fig. 1a) and EELS measurements (Fig. 2a). Measurement details are included in the “Methods” section and SI. Such a structure of Sample 1 is chosen, as the absorption depth of Si is large (774 nm), tens of times greater than that of Ge, at the Raman laser wavelength (488 nm)^{23}. Thus, the 10 nm Si cap layer can allow the Raman laser to transmit through to reach the Si-Ge interfaces. Sample 2 is a 250 nm Ge layer grown on a Si substrate and is used for TBC measurements, where the Ge layer thickness is optimized to maximize the sensitivity of the TBC of the Si-Ge interface in TDTR measurements (Supplementary Fig. 12). To confirm the quality of the Si-Ge interface, high-angle annular dark-field STEM (HAADF-STEM) imaging technique is used to characterize the interface of Sample 1. In Fig. 1a, the Z-contrast HAADF-STEM image shows that the interface is sharp despite a thin atomic mixing layer (Supplementary Figs. 9 and 10). The atomic structures of the two Si-Ge interfaces in Sample 1 and Sample 2 are similar (Supplementary Fig. 9).

As a comparison to the Raman spectrum of Sample 1 which involves a Si-Ge interface, we have also measured the Raman spectra of a pure Ge wafer and a pure Si wafer. Figure 1b shows that the Raman peak of bulk Ge is at ~9 THz and that of Si is at ~15.6 THz. The Si peak in Sample 1 is redshifted by 0.15 THz compared to pure Si, because the epitaxial Si cap layer is slightly strained on the Ge wafer to accommodate the ~4% lattice mismatch. Other than these two peaks originated from Ge and Si, an additional peak around 11.3–12.2 THz shows up in Sample 1, which may be attributed to the Si-Ge interface.

To further confirm that this peak originates from the interface, we employ high-energy resolution (routinely ~1.7–1.9 THz with an exposure time of 1 s) EELS, in a STEM, with a probe size of 1.5 Å to spatially resolve the phonon signal around the Si-Ge interface (see details in the “Methods” section and SI). Figure 2a shows the schematic diagram of the EELS measurements. The peaks in the EELS signal are the energies of vibrational phonon modes in the sample that inelastically interact with the incident electron probe. Figure 2b shows a line-scan of the vibrational spectra across the interface as illustrated by the arrows in Fig. 2a. Under our experimental setup with a large convergence semi-angle (33 mrad), the as-acquired vibrational spectra contain the momentum-integrated vibrational modes throughout the entire Brillouin zone (BZ) and are comparable with the phonon density of states (PDOS)^{24}.

Three representative EELS vibrational spectra from Si (red), Ge (green), and the interface (blue) taken at the locations marked in Fig. 2a are shown in Fig. 2c. The vibrational spectrum of the interfacial region includes the peaks from Si, Ge, and the interfacial vibration modes. Figure 2d shows the intensities of vibrational spectra at 11.6, 12.0, and 12.4 THz as a function of distance to the interface. A clear peak of the vibrational spectral intensity is observed around the interface at 11.6 THz and there is a weak peak at 12.0 THz. It is notable that there are some similar intensities in the same energy range on the Si side, which arise from the bulk phonon modes in the Si as shown in Fig. 2c. Those are from the optical and acoustic phonon modes near the BZ boundary, according to the phonon dispersion curve. Two results are used to distinguish the interfacial modes from these bulk Si modes. First, the interfacial vibrational spectrum has extra intensity near 12.0 THz (Supplementary Fig. 4), compared to the linear combination of bulk Ge and Si vibrational spectra, indicating the existence of additional phonon modes at the interface. Second, we further performed angle-resolved vibrational EELS to collect the vibrational signals at the BZ center (Γ point, Supplementary Fig. 11a) to exclude the interference of phonon modes near the BZ boundary. We observed a clear peak around 12 THz at the interface while that from the bulk Si at this frequency disappeared. By combining the Raman and high-energy-resolution EELS results above, we can claim that the detected modes around 12 THz are indeed interfacial phonon modes localized at the interface. Leveraging the high spatial resolution (1.5 Å) capability of STEM-EELS, this interfacial vibrational mode is determined to be confined within ~1.2 nm of the interface (Fig. 2d).

To further understand the interfacial modes, we use MD simulations to model the Si-Ge interface with the level of interfacial mixing (~0.7 nm thick) similar to those observed in the HAADF intensity (Fig. 3a and Supplementary Fig. 9). In order to reproduce the interfacial modes, a high-fidelity NNP for the Si-Ge interface is specifically developed by training against first-principles density functional theory calculations (see “Methods” section and SI for details). The NNP-MD scheme has been proved to be capable of simulating semiconductor materials and predicting thermal properties with accuracy comparable to first-principles calculations^{25,26}. Using MD simulations with this NNP, PDOS at different locations throughout the simulation domain along the direction perpendicular to the Si-Ge interface is calculated (see “Methods” section and SI for details).

As shown in Fig. 3b, each panel is the PDOS of the atoms at a certain location shown in Fig. 3a. The interfacial region is composed of the mixed region and the first layer away from the mixture on each side. For the regions that are 10 nm away from the interface, the PDOS of both Si and Ge are converged to the ones of bulk materials and are not affected by the interface. The PDOS peaks of the Si and Ge optical phonons are ~15 THz and ~9 THz, which are consistent with the Raman peaks. Near the interface, a peak at ~11.4–12.8 THz emerges for both the Si and Ge atoms, which corresponds well with the interfacial modes detected from the Raman and EELS measurements. This calculated interfacial mode peak agrees better with our experiments compared to other results of MD simulations using empirical potentials. Chalopin et al.^{16} obtained such a peak at ~13.5 THz using a Stinger-Weber potential, while Gordiz et al.^{17} found the peak between 12–13 THz using a Tersoff potential. Our calculated peak is similar to that calculated from first-principles force constants (~12 THz), which, however, used Si force constants to approximate those for Ge and interfacial interactions^{18}. We note that all these previous studies modeled sharp interfaces without considering atomic mixing effects. For comparison, we also performed simulations on an ideally sharp Si-Ge interface and also found a peak around 11.6–12.7 THz (Supplementary Fig. 7).

To understand how these interfacial modes contribute to interfacial thermal transport, we performed non-equilibrium MD (NEMD) simulations to calculate the TBC (see “Methods” section and SI for details). The calculated TBC is ~250 MW m^{−2} K^{−1} at 300 K. We further used spectral analysis to calculate the contributions of phonons with different frequencies to the TBC (see “Methods” section and SI for details). This does not only allow us to directly visualize the role of interfacial modes, but also allows quantum correction to be applied so that the MD-calculated TBC can be directly compared to our experimentally measured TBC. As shown in Fig. 4a, there is an obvious peak of the frequency-dependent TBC around the interfacial phonon mode region (11.4–12.8 THz) from our spectral analysis. As MD are classical simulations, all phonons are equally excited. However, in reality, phonon follows the Bose–Einstein distribution. Such an error can be corrected by applying a quantum correction by weighing the contributions of different frequencies with the ratio of the quantum and classical heat capacities of that frequency (see “Methods” section and SI for details). The quantum correction leads to reduced contributions from high-frequency modes as they are not fully excited, but the contributions from the interfacial modes are still obvious. We integrate the frequency-dependent TBC to obtain the cumulative TBC, and find that the interfacial modes contribute ~5% of the total TBC at 300 K (Fig. 4b). The quantum-corrected TBC (231 MW m^{−2} K^{−1}) agrees very well with our TDTR measurements (244 and 236 MW m^{−2} K^{−1} by different TDTR systems). We note that these calculations are on the mixed interface. For a sharp interface, the calculated TBC is 197 MW m^{−2} K^{−1}, which is lower than the mixed interface and farther away from the measured TBC. This is consistent with previous findings that interfacial mixing can enhance TBC^{11,12}. Figure 4c also show the comparison of TBC values of Si-Ge interfaces calculated from DMM, AGF (both harmonic and anharmonic), and previous MD simulations using different potentials. The NNP-MD results for both the mixed interface and the sharp interface respectively agree well with their counterparts from AGF calculations, highlighting the importance of correctly capturing the microscopic interface conditions for TBC predictions. On the other hand, other models (e.g., DMM) that ignore such microscopic interfacial details show inferior agreement with experimentally measured TBC. The calculation details of non-equilibrium Landauer approach based on DMM can be found in the “Methods” section and SI. The fact that NNP-MD is accurate in predicting TBC also gives confidence on its predicted non-trivial contributions of interfacial modes to the TBC.

We further performed temperature-dependent TDTR measurements and MD simulations from 300 to 500 K. The TDTR setup we used for high-temperature measurements is at University of Illinois at Urbana-Champaign (UIUC). The measured TBC matches the value measured by the TDTR at Georgia Institute of Technology (GT) at 300 K, as shown in Fig. 5a. The measured TBC values also agree well with the MD-calculated TBC values with NNP, which again confirms the accuracy of the present NNP for the Si-Ge interface. Both the measured and calculated TBC values show a weak increase with increasing temperature, indicating that the population of the phonons involved in interfacial thermal transport only increases slightly in this temperature range —a trend commonly seen in the temperature-dependent TBC of many other solid interfaces at high temperatures.

To compare with the interfacial mode contributions to TBC at 300 K, spectral analysis is also conducted for the interface at 460 K with MD and the quantum-corrected results are shown in Fig. 5b, c. The two accumulation curves at different temperatures in Fig. 5c share a similar trend, with an obvious jump around the frequency of the interfacial mode (11.4–12.8 THz). Although the absolute contributions from the interfacial modes increase with increasing temperature, the relative contributions remain around 5% as the overall TBC also increases slightly at 460 K.

It is noted that the frequencies of the interfacial modes are higher than those of the optical phonons of bulk Ge. Elastic processes across interfaces require that bulk phonon modes on both sides have the same frequency. Therefore, for the Si-Ge interface, phonons with frequencies above ~10 THz cannot contribute to TBC via elastic processes while they can contribute to TBC via inelastic processes. Inelastic processes allow phonons with different frequencies to exchange energy and contribute to TBC^{18,28,29,30,31}. Previous simulations have shown that the role of the interfacial modes localized at the interface is to act as a bridge coupling the bulk phonon modes on both sides, which facilitates the inelastic transport across the interfaces for phonons with different frequencies^{17,28,32}. We note that although we still call this an inelastic process, because energy is transferred from one frequency to another across the interface, the underlying mechanism is not originated from a conventional phonon scattering picture. Instead, it is due to the correlation between the localized interfacial phonons and other phonons that extend into the materials making up the interface^{17,32}. It is notable that the MD simulations with NNP in this work intrinsically include both elastic and inelastic processes. However, other theoretical models such as DMM, AMM, and AGF typically only consider the elastic phonon transport even though the inelastic processes can contribute to the total TBC^{3,28,33}. The recent development of AGF enables the consideration of both elastic and inelastic contributions but it is very computationally expensive^{13}.

To explore the physical origin of the interfacial modes, we visualize the eigenmodes of an example interfacial mode at ~12.4 THz near the interface with and without atomic mixing (Fig. 5d). It can be seen that the eigenvectors are mainly localized around interfacial atoms that involve Si-Ge bonds for both the sharp and mixed interfaces. For the mixed interface, the delocalization length of the interfacial modes is around 1.1 nm in the interfacial region, which agrees well with that of the EELS-detected interfacial modes (~1.2 nm around the interface) in Fig. 2d. The origin of the interfacial modes is thus believed to be from the Si-Ge bonds that do not exist in bulk Si or bulk Ge. It has been reported that the Si atoms encapsulated in Ge atomic cages can lead to modes at a frequency that does not exist in bulk Si and Ge, and the observed frequency is similar to that of the interfacial modes^{34}. This suggests that the special modes from the Si-in-Ge cage structure and interfacial modes have the same origin, which is the Si-Ge bonds.

In summary, we experimentally observed the existence of localized interfacial vibrational modes within 1.2 nm at an epitaxial Si-Ge heterointerface by Raman spectra and high-energy-resolution EELS in a STEM. These interfacial modes are found to be reproducible using MD simulations with a high-fidelity NNP, which also yield TBC agreeing favorably to TDTR measurements. The spectral analysis in MD indicates that these interfacial modes, which are mainly localized near the interface, can have non-trivial contributions to TBC. This work paves the way to fundamentally understanding heat transport across realistic interfaces and may stimulate the development of theories of interfacial thermal transport. It will impact applications such as electronics thermal management and thermoelectric energy conversion.

## Methods

### Materials growth

The two samples in this work were grown by MBE. The silicon substrates (boron doped 5 × 10^{14}–1 × 10^{15} atoms/cm^{3}) were given a modified Radio Corporation of America cleaning procedure followed by an HF dip as the last step. Then, the samples were immediately loaded into the vacuum systems and pumped down from atmospheric pressure to <10^{−6} Torr in 15 min. Details can be found in the literature^{35}. For Sample 1 (10 nm Si-300 nm Ge-Si substrate), the temperature of the silicon substrate (boron doped 5 × 10^{14}–1 × 10^{15} atoms/cm^{3}) was first increased to 700 °C in 3 min. A 25 nm Si buffer was grown on the Si substrate at 700 °C, followed by a 300 nm Ge layer in three parts (10 nm Ge at 300 °C, 80 nm Ge at 400 °C, and annealed at 700 °C, then a 210 nm Ge at 400 °C), and finally a 10 nm Si cap was grown on the Ge at 400 °C. The annealing time of the Ge layer at 700 °C is 16.67 min (1000 s). The cooling time after the Ge annealing is ~6 min (~1 °C/s from 700 °C to 400 °C). Similarly, the cooling rate from the Si homo-epitaxial layer to the initial Ge growth layer is ~ 1 °C/s with a 2 min settling time (~8 min in total). Sample 2 (250 nm Ge-Si substrate) was grown in a similar way but without the 700 °C annealing process. The growth temperature of the 250 nm Ge layer was dropped over an 8 min period to ~400 °C, where the 250 nm Ge was grown at a constant temperature. More details can be found in SI.

### Raman spectroscopy

Raman measurements were performed on Sample 1, a Si wafer, and a Ge wafer with a Renishaw InVia Raman system. The Ar^{+} laser wavelength is 488 nm. The objectives are Leica 0.50 NA ×50 and 0.75 NA ×100 with theoretical diffraction-limited spot sizes of 1.19 μm and 0.79 μm, respectively. The ×50 and ×100 objectives were used to measure Sample 1. The ×50 objective was used to measure the reference Si and Ge wafers. The acquisition time of each Raman spectrum is 30 s with five accumulations.

### Thermal characterizations

The TBC of Sample 2 were measured by TDTR. TDTR is an ultrafast-laser-based pump-probe technique which is capable of measuring thermal transport properties of both nanostructured and bulk materials^{36,37}. A modulated pump beam heats the sample surface periodically while a delayed probe beam detects the temperature variation of the sample surface. The signal picked up by a photodetector and a lock-in amplifier was fitted with an analytical heat transfer solution of the sample structure to infer unknown thermal properties. An 80 nm Al layer was deposited on the sample as the TDTR transducer. For the TDTR measurements at GT, a ×10 objective was used for all the measurements with a pump radius of 10.7 μm and a probe radius of 5.8 μm. The modulation frequency is 8.8 MHz. For the TDTR measurements at UIUC, a ×5 objective was used for all the measurements with a pump/probe radius of 10.7 μm. The modulation frequency is 9.3 MHz. More details about thermal measurements can be found in the SI.

### TEM measurements

HAADF-STEM was performed using JEOL JEM-ARM300CF S/STEM, equipped with double aberration correctors with a spatial resolution of ~0.6 Å at 300 kV. HAADF-STEM images were recorded using a convergence semi-angle of 22 mrad and collection semi-angles at 83–165 mrad.

### EELS measurements

The vibrational spectra were performed using a Nion UltraSTEM 200 microscope operating at 60 kV. The convergence semi-angle was 33 mrad and the probe current was about 120 pA, offering a spatial resolution of 1.5 Å. The probe current needs to remain high to warrant adequate vibrational signals on the EELS detector. The EELS collection semi-angle was 25 mrad to allow all phonon branches inside the entire BZ to be recorded by EELS. With the implementation of an alpha-type monochromator and designed magnetic-prism spectrometer, the energy resolution of EELS is routinely 7–8 meV (1.7–1.9 THz) with an exposure time of 1 s^{38}. High spatial-resolution EELS line-scan datasets were automatically performed and assembled by running our custom-developed python script on top of Nion’s Swift software. To optimize the signal-to-noise ratio and mitigate the influence of high voltage instability, the 40-pixel line-scan data in Fig. 2b with a step size of 0.15 nm was acquired by summing 50 frames with 1 s exposure per frame, which were then all aligned by the center of zero-loss peak. Besides, the single point spectra in Fig. 2c were obtained by summing 200 frames with 1 s exposure per frame, which were captured separately from the line-scan. The EELS dispersion was about 0.3 meV/channel. Background subtraction was carried out for each spectrum by fitting an exponential polynomial function.

### Non-equilibrium Landauer approach

Phonons near the interfaces are in strong non-equilibrium conditions because of the difference in modal transmission coefficients and reservoir temperatures. By considering this non-equilibrium effect, a corrected non-equilibrium Landauer approach was developed to calculate TBC^{3,39}. Here, DMM was used to calculate the transmission. A detailed comparison of non-equilibrium Landauer approach and conventional Landauer approach is included in the SI.

### NNP construction

To build a NNP for the Si-Ge interfaces, we first generated a training dataset consisting of snapshots from ab initio MD (AIMD) simulations^{40}. Each snapshot contains the data of potential energy, atomic forces, atomic coordinates, and supercell lattice vectors. An 8 × 3 × 3 Si-Ge interface supercell, with a Si lattice constant (*a* = 5.4307 Å), is simulated in the canonical ensemble (NVT) at temperatures ranging from 50 to 600 K. The system is equilibrated for 1.0 ps with a time step of 1.0 fs and the snapshots are then collected every five steps during 4.0 ps-long production runs. To capture the effect of interfacial mixing, we also generated a dataset for two different 4 × 3 × 3 Si-Ge crystalline alloys with the Si lattice constant adopting the same procedure. AIMD simulations are carried out using the QUICKSTEP algorithm implemented in the CP2K package, which leverages the Gaussian and plane waves approaches^{41,42}. We employ the Goedecker–Teter–Hutter (GTH) double-ζ, single polarization (DZVP-MOLOPT-GTH) basis set, and the GTH-PBE pseudopotentials^{43} to describe the core-valence interactions. The wave-plane cutoff is set to be 300 Ry and the BZ is sampled at the Γ point. Finally, the training dataset contains 9600 snapshots for a perfect Si-Ge interface and 8000 snapshots for the Si-Ge crystalline alloy mixtures.

The NNP is obtained by training a deep neural network potential scheme developed by Zhang et al*.*^{44,45}, which has been proved to be capable of simulating Pt/MoS_{2} interfaces^{46}. This scheme comprises an embedding network and a fitting network. The embedding network first maps the chemical species and atomic coordinates to the embedded features that preserve all the natural symmetries. Then a three-layer feed-forward network learns the potential energies and forces by minimizing the mean squared error loss during training on the ab initio dataset. The interaction cutoff for this potential is set to be 6 Å. More details can be found in the SI.

### NEMD simulations

NEMD simulations were performed using the LAMMPS package^{47} with this NNP for the Si-Ge interfaces. For both Si and Ge sides, the number of conventional cells along the x, y, and z directions are chosen to be equal to 25, 7, and 7, respectively (see SI for details). A 0.7 nm-long mixed region is formed by randomly shuffling the atoms near the interface to mimic our measured Si-Ge interfaces. The interface is a plane perpendicular to the *x*-direction, which is at the middle of the interfacial mixture and perpendicular to the [100] crystallographic direction. Periodic boundary conditions are applied to all three spatial directions, and a time step of 1.0 fs is used for the MD simulations. After relaxing the structure in the NVT ensemble for 0.5 ns, we freeze the two ends of the system and apply a temperature difference using Langevin thermostats in the microcanonical (NVE) ensemble for 2.0 ns. The temperatures applied on the two thermostats are ±25 K of the target mean temperature. The TBC is obtained by

\(\begin{array}{c}{G}_{{{{{\mathrm{Si}}}}}-{{{{{\mathrm{Ge}}}}}}}=Q/\left(A\triangle T\right),\end{array}\) where *Q* is the steady-state heat flux along the *x*-direction, *A* is the cross-sectional area, and ∆*T* is the temperature difference at the interface determined by extrapolating linear fits of the temperatures of the two sides to the interface and calculating the difference. Considering that classical statistics overestimates modal heat capacity at high frequencies^{48}, it is necessary to apply quantum correction to these predictions. To account for the quantum effect, we need to quantify the contributions to the interfacial thermal transport on the mode level. The spectral analysis^{49,50,51} is employed to determine the modal contributions to TBC. With the modal decomposition of TBC, we apply quantum correction by multiplying *G*(*ω*) by the factor \({u}^{2}{e}^{u}/{({e}^{u}-1)}^{2}\), where \(u=\hslash \omega /{k}_{{{{{\mathrm{B}}}}}}T\) represents the ratio of the quantum heat capacity to the classical one. More details can be found in the SI.

### Reporting summary

Further information on research design is available in the Nature Research Reporting Summary linked to this article.

## Data availability

The datasets generated during and/or analyzed during the current study are available from the corresponding authors upon reasonable request.

## Code availability

The code used for calculations, simulations, and data analysis is available from the corresponding authors upon reasonable request.

## References

- 1.
Losego, M. D., Grady, M. E., Sottos, N. R., Cahill, D. G. & Braun, P. V. Effects of chemical bonding on heat transport across interfaces.

*Nat. Mater.***11**, 502–506 (2012). - 2.
Cahill, D. G. et al. Nanoscale thermal transport. II. 2003–2012.

*Appl. Phys. Rev.***1**, 011305 (2014). - 3.
Cheng, Z. et al. Thermal conductance across harmonic-matched epitaxial Al-sapphire heterointerfaces.

*Commun. Phys.***3**, 1–8 (2020). - 4.
Gaskins, J. T. et al. Thermal boundary conductance across heteroepitaxial ZnO/GaN interfaces: assessment of the phonon gas model.

*Nano Lett.***18**, 7469–7477 (2018). - 5.
Giri, A. & Hopkins, P. E. A review of experimental and computational advances in thermal boundary conductance and nanoscale thermal transport across solid interfaces.

*Adv. Funct. Mater*.**30**, 1903857 (2019). - 6.
Swartz, E. T. & Pohl, R. O. Thermal boundary resistance.

*Rev. Mod. Phys.***61**, 605 (1989). - 7.
Swartz, E. & Pohl, R. Thermal resistance at interfaces.

*Appl. Phys. Lett.***51**, 2200–2202 (1987). - 8.
Little, W. The transport of heat between dissimilar solids at low temperatures.

*Can. J. Phys.***37**, 334–349 (1959). - 9.
Hopkins, P. E., Duda, J. C., Petz, C. W. & Floro, J. A. Controlling thermal conductance through quantum dot roughening at interfaces.

*Phys. Rev. B***84**, 035438 (2011). - 10.
Prasher, R. Acoustic mismatch model for thermal contact resistance of van der Waals contacts.

*Appl. Phys. Lett.***94**, 041905 (2009). - 11.
Tian, Z., Esfarjani, K. & Chen, G. Enhancing phonon transmission across a Si/Ge interface by atomic roughness: first-principles study with the Green’s function method.

*Phys. Rev. B***86**, 235304 (2012). - 12.
Fagas, G. et al. Lattice dynamics of a disordered solid-solid interface.

*Phys. Rev. B***60**, 6459 (1999). - 13.
Dai, J. & Tian, Z. Rigorous formalism of anharmonic atomistic Green’s function for three-dimensional interfaces.

*Phys. Rev. B***101**, 041301 (2020). - 14.
Zhang, W., Fisher, T. & Mingo, N. Simulation of interfacial phonon transport in Si–Ge heterostructures using an atomistic Green’s function method.

*J. Heat. Transf.***129**, 483–491 (2007). - 15.
Li, X. & Yang, R. Effect of lattice mismatch on phonon transmission and interface thermal conductance across dissimilar material interfaces.

*Phys. Rev. B***86**, 054305 (2012). - 16.
Chalopin, Y. & Volz, S. A microscopic formulation of the phonon transmission at the nanoscale.

*Appl. Phys. Lett.***103**, 051602 (2013). - 17.
Gordiz, K. & Henry, A. Phonon transport at crystalline Si/Ge interfaces: the role of interfacial modes of vibration.

*Sci. Rep.***6**, 23139 (2016). - 18.
Murakami, T., Hori, T., Shiga, T. & Shiomi, J. Probing and tuning inelastic phonon conductance across finite-thickness interface.

*Appl. Phys. Exp.***7**, 121801 (2014). - 19.
Zhou, Y. & Hu, M. Full quantification of frequency-dependent interfacial thermal conductance contributed by two-and three-phonon scattering processes from nonequilibrium molecular dynamics simulations.

*Phys. Rev. B***95**, 115313 (2017). - 20.
Venkatraman, K., Levin, B. D., March, K., Rez, P. & Crozier, P. A. Vibrational spectroscopy at atomic resolution with electron impact scattering.

*Nat. Phys.***15**, 1237–1241 (2019). - 21.
Krivanek, O. L. et al. Vibrational spectroscopy in the electron microscope.

*Nat***514**, 209 (2014). - 22.
Lagos, M. J., Trügler, A., Hohenester, U. & Batson, P. E. Mapping vibrational surface and bulk modes in a single nanocube.

*Nature***543**, 529 (2017). - 23.
Green, M. A. & Keevers, M. J. Optical properties of intrinsic silicon at 300 K.

*Prog. Photovolt. Res. Appl.***3**, 189–192 (1995). - 24.
Yan, X. et al. Single-defect phonons imaged by electron microscopy.

*Nature***589**, 65–69 (2021). - 25.
Li, R., Lee, E. & Luo, T. A unified deep neural network potential capable of predicting thermal conductivity of silicon in different phases.

*Mater. Today Phys*.**12**, 100181 (2020). - 26.
Li, R. et al. A deep neural network interatomic potential for studying thermal conductivity of β-Ga2O3.

*Appl. Phys. Lett.***117**, 152102 (2020). - 27.
Landry, E. & McGaughey, A. Thermal boundary resistance predictions from molecular dynamics simulations and theoretical calculations.

*Phys. Rev. B***80**, 165304 (2009). - 28.
Feng, T., Zhong, Y., Shi, J. & Ruan, X. Unexpected high inelastic phonon transport across solid-solid interface: Modal nonequilibrium molecular dynamics simulations and Landauer analysis.

*Phys. Rev. B***99**, 045301 (2019). - 29.
Hopkins, P. E., Norris, P. M. & Stevens, R. J. Influence of inelastic scattering at metal-dielectric interfaces.

*J. Heat. Transf.***130**, 022401 (2008). - 30.
Hopkins, P. E. Multiple phonon processes contributing to inelastic scattering during thermal boundary conductance at solid interfaces.

*J. Appl. Phys.***106**, 013528 (2009). - 31.
Kosevich, Y. A. Fluctuation subharmonic and multiharmonic phonon transmission and Kapitza conductance between crystals with very different vibrational spectra.

*Phys. Rev. B***52**, 1017 (1995). - 32.
Rohskopf, A., Li, R., Luo, T. & Henry, A. A computational method for studying vibrational mode dynamics. Preprint at

*arXiv*https://arXiv.org/2108.04795 (2021). - 33.
Lyeo, H.-K. & Cahill, D. G. Thermal conductance of interfaces between highly dissimilar materials.

*Phys. Rev. B***73**, 144301 (2006). - 34.
Kosevich, Y. A. & Strelnikov, I. Extraordinary phonon transmission through hidden lattice-wave nanochannels as resonance quantum phonon tunneling. In

*AIP Conference Proceedings***2241**, 1020023 (AIP Publishing LLC, 2020). - 35.
Jernigan, G. G. & Thompson, P. E. Evolution of Si0.8Ge0.2 quantum dots during Si encapsulation.

*Surf. Sci.***600**, 933–940 (2006). - 36.
Cahill, D. G. Analysis of heat flow in layered structures for time-domain thermoreflectance.

*Rev. Sci. Instrum.***75**, 5119–5122 (2004). - 37.
Cheng, Z. et al. Probing growth-induced anisotropic thermal transport in high-quality CVD diamond membranes by multi-frequency and multi-spot-size time-domain thermoreflectance.

*ACS Appl. Mater. Inter.***10**, 4808–4815 (2018). - 38.
Yan, X. et al. Unexpected strong thermally induced phonon energy shift for mapping local temperature.

*Nano Lett.***19**, 7494–7502 (2019). - 39.
Shi, J., Yang, X., Fisher, T. S. & Ruan, X. Dressed and intrinsic thermal boundary conductance across interfaces from a nonequilibrium Landauer approach. Preprint at

*arXiv*https://arXiv.org/1812.07910 (2018). - 40.
Car, R. & Parrinello, M. Unified approach for molecular dynamics and density-functional theory.

*Phys. Rev. Lett.***55**, 2471 (1985). - 41.
VandeVondele, J. et al. Quickstep: fast and accurate density functional calculations using a mixed Gaussian and plane waves approach.

*Comp. Phys. Commun.***167**, 103–128 (2005). - 42.
Kühne, T. D., Krack, M., Mohamed, F. R. & Parrinello, M. Efficient and accurate Car-Parrinello-like approach to Born-Oppenheimer molecular dynamics.

*Phys. Rev. Lett.***98**, 066401 (2007). - 43.
Goedecker, S., Teter, M. & Hutter, J. Separable dual-space Gaussian pseudopotentials.

*Phys. Rev. B***54**, 1703 (1996). - 44.
Wang, H., Zhang, L., Han, J. & Weinan, E. DeePMD-kit: a deep learning package for many-body potential energy representation and molecular dynamics.

*Comp. Phys. Commun.***228**, 178–184 (2018). - 45.
Zhang, L., Han, J., Wang, H., Car, R. & Weinan, E. Deep potential molecular dynamics: a scalable model with the accuracy of quantum mechanics.

*Phys. Rev. Lett.***120**, 143001 (2018). - 46.
Zhang, L. et al. End-to-end symetry preserving inter-atomic potential energy model for finite and extended systems.

*Adv. Neural Inform. Process. Sys*. 4436–4446 (2018). - 47.
Plimpton, S. Fast parallel algorithms for short-range molecular dynamics.

*J. Comput. Phys.***117**, 1–19 (1995). - 48.
Turney, J., McGaughey, A. & Amon, C. Assessing the applicability of quantum corrections to classical thermal conductivity predictions.

*Phys. Rev. B***79**, 224305 (2009). - 49.
Sääskilahti, K., Oksanen, J., Volz, S. & Tulkki, J. Frequency-dependent phonon mean free path in carbon nanotubes from nonequilibrium molecular dynamics.

*Phys. Rev. B***91**, 115426 (2015). - 50.
Fan, Z. et al. Thermal conductivity decomposition in two-dimensional materials: application to graphene.

*Phys. Rev. B***95**, 144309 (2017). - 51.
Li, R. et al. Effect of light atoms on thermal transport across solid–solid interfaces.

*Phys. Chem. Chem. Phys.***21**, 17029–17035 (2019).

## Acknowledgements

We acknowledge the financial support from Office of Naval Research MURI Grant Number N00014-18-1-2429. We thank helpful discussions with Asegun Henry, Andrew Rohskopf, and Tianli Feng. The TEM work was supported by the Department of Energy (DOE), Office of Basic Energy Sciences, Division of Materials Sciences and Engineering under Grant DE-SC0014430, and partially supported by the Center for Nanophase Materials Sciences, which is a DOE Office of Science User Facility (JCI). We acknowledge the use of facilities and instrumentation at the UC Irvine Materials Research Institute (IMRI), which is supported in part by the National Science Foundation through the UC Irvine Materials Research Science and Engineering Center (DMR-2011967).

## Author information

### Affiliations

### Contributions

Z.C., S.G. and T.L. initialized the project. Z.C. and N.J.H. did the Raman measurements. R.L., E.L. and T.L. finished the NNP development and NEMD simulations. X.Y., C.G., X. P. and J.C.I. performed the TEM and EELS measurements. Z.C. performed the TDTR measurements. G.J. and K.D.H. grew the samples. J.S. did the non-equilibrium Landauer calculations. M.E.L. and M.S.G. performed the XRD measurements. Z.C. coordinated the project and wrote the manuscript with inputs from all authors. X.Y. and R.L. assisted with manuscript preparation, prepared figures, and reviewed the manuscript. T.L. provided insights in understanding the data and edited the manuscript. S.G. initially motivated the investigation, provided overall guidance to the project, and reviewed the manuscript.

### Corresponding authors

## Ethics declarations

### Competing interests

The authors declare no competing interest.

## Additional information

**Peer review information** *Nature Communications thanks* Guang-Kun Ren and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available.

**Publisher’s note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

## Supplementary information

## Rights and permissions

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

## About this article

### Cite this article

Cheng, Z., Li, R., Yan, X. *et al.* Experimental observation of localized interfacial phonon modes.
*Nat Commun* **12, **6901 (2021). https://doi.org/10.1038/s41467-021-27250-3

Received:

Accepted:

Published:

DOI: https://doi.org/10.1038/s41467-021-27250-3

## Comments

By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.