Improving interfacial thermal conductivity by constructing covalent bond between Ga₂O₃ and SiC

Domain: physics
Published: 2025-11-28
Editorial outcome: major_revision

Abstract
Gallium oxide (Ga₂O₃) is emerging as a promising semiconductor for next-generation power and radio-frequency electronics due to its ultra-wide bandgap and high breakdown field. Yet, its intrinsic thermal conductivity is extremely low, which causes severe self-heating and limits reliable device operation. A common approach to overcome this challenge is to integrate Ga₂O₃ with substrates of high thermal conductivity such as silicon carbide (SiC). However, weak bonding across the heterojunction interface creates large thermal resistance, preventing efficient heat removal. Here we show that constructing strong covalent bonds between Ga₂O₃ and SiC through an engineered interlayer enables both structural compatibility and efficient phonon transport. This bonding strategy increases the interface thermal conductivity to 162 MW/m²·K, the highest value reported for Ga₂O₃ heterostructures. Infrared thermography confirms that the bonded devices exhibit a temperature reduction of up to 29 °C under high power densities, demonstrating significant mitigation of self-heating. These findings establish a practical route for enhancing thermal management in Ga₂O₃ electronics and highlight the importance of interfacial bonding design. Beyond Ga₂O₃, this approach may be extended to other wide-bandgap semiconductors where thermal bottlenecks constrain device performance.

Introduction
Since modern electronic devices have evolved toward the higher power, frequency and temperature capabilities, wide-bandgap semiconductor materials have garnered significant research interest due to their exceptional physical properties1–3. Silicon carbide (SiC) is widely utilized in high-power applications because of its wide bandgap of 3.2 eV and low conduction loss. However, the low saturated electron mobility still limit the effectiveness of SiC in high-frequency domains. In contrast, gallium nitride (GaN) offers notable advantages under high-frequency conditions but lacks the ability to produce self-supporting substrates, which restricts its use in high-power scenarios. Gallium oxide (Ga₂O₃), a semiconductor with an ultra-wide bandgap (approximately 4.8 eV), has attracted considerable attention for its potential in high-voltage power and high-frequency electronic devices4. Ga₂O₃ exhibits exceptionally high critical breakdown field (projected to be up to 8 MV/cm) and low on-resistance, which can significantly reduce losses in power devices. Furthermore, it has a high saturation electron velocity (∼2×107 cm/s) and an impressive Johnson’s figure of merit (JFOM), thus being an ideal candidate for high-frequency applications5–7. Beyond its outstanding electronic properties, Ga₂O₃ possesses a combination of unique material advantages that make it a highly promising candidate for next-generation semiconductor applications8. It demonstrates excellent thermal and chemical stability, maintaining performance under extreme conditions such as ultra-low temperatures and strong radiation environments9. Additionally, Ga₂O₃ exhibits strong absorption in the deep-ultraviolet (DUV) spectral region (~250–280 nm), making it highly suitable for solar-blind UV photodetectors and space-based optoelectronic applications10. From a material perspective, Ga₂O₃ is advantageous in both bulk and epitaxial growth. It supports melt-based bulk crystal growth techniques such as edge-defined film-fed growth (EFG) and Czochralski (CZ), enabling the production of large-area, low-cost single crystals—unlike SiC, which requires expensive high-temperature vapor-phase growth11. At the same time, Ga₂O₃ shows excellent compatibility with various heteroepitaxial growth methods, including MOCVD (Metal-Organic Chemical Vapor Deposition), MBE (Molecular Beam Epitaxy), and PLD (Pulsed Laser Deposition), on substrates such as sapphire, Si, and SiC12. In addition, Ga₂O₃ has relatively low hardness compared to SiC, reducing wafer slicing and processing costs. Given the breakthrough in manufacturing, including the production of 6-inch substrates, Ga₂O₃ shows great promise for the development of next-generation high-power and high-frequency devices13,14.

However, one significant disadvantage of Ga₂O₃ is the low thermal conductivity (approximately 10–27 W/m·K)15,16, which is inferior to that of other semiconductors, such as Si (150 W/m·K), AlN (180 W/m·K), GaN (230 W/m·K), SiC (490 W/m·K), and diamond (2000 W/m·K)5. Such a low thermal conductivity value can lead to serious thermal management issues when operating at high temperatures and powers. Energy losses during conduction and switching generate heat, raising the temperature of the device’s active junction and package. In turn, excessive junction temperatures (Tj) negatively affect the performance and reliability of the device17. To ensure the long-term reliability, Tj in commercial power devices is typically kept18. The “10 K rule” indicates that the failure rate of electronic components can double for every 10 K increase19. Chatterjee et al. have found that Tj in Ga₂O₃ devices can reach up to 1500 °C at a power density of 10 W/mm20. Therefore, this thermal management is a major bottleneck in the use of Ga₂O₃ in high-power and high-frequency applications.

In recent years, various strategies have been proposed to enhance the heat dissipation performance of Ga₂O₃-based devices9. As illustrated in Supplementary Fig. 1, heat generated during operation is typically transferred through a thermal interface material (TIM) to the substrate and then dissipated via an external heatsink21,22. While adding heatsink modules can improve thermal efficiency, this strategy often leads to bulkier device structures and limited heat removal from localized hotspots, which hampers device miniaturization and reliability. A possible strategy to circumvent the low thermal conductivity of β-Ga₂O₃ is to thin the Ga₂O₃ substrate23. However, while improving heat dissipation, this route may compromise the mechanical strength under harsh operating conditions and reduce reliability of the appliance. A more effective approach is near-junction thermal management (NJTM), using a highly conductive substrate in direct contact with the heat source so as to ensure the direct heat dissipation, e.g., Si24, SiC25–27, Ge28 and diamond29–31 substrates. NJTM is a highly effective heat dissipation method, enabling to position the heat-conducting materials or heat-sinking structures directly near the heat source of a chip to rapidly transfer heat outward, thereby enhancing the overall thermal efficiency32. The effective NJTM is essential to ensure the optimal heat dissipation, which has been demonstrated in GaN high-electron-mobility transistors through the use of substrates with high thermal conductivity33,34. Among high-thermal conductivity substrates, SiC is considered one of the most suitable substrates for heat dissipation. Yuan et al. demonstrated that although SiC, with a thermal conductivity of 490 W/m·K, has a lower value than diamond, a substrate thermal conductivity exceeding 400 W/m·K is sufficient to achieve comparable heat dissipation performance to that of diamond35. Besides that, SiC possesses excellent thermal stability and lattice matching with β-Ga₂O₃. At the same time, compared with diamond, SiC is more suitable for the preparation of large-size substrates at low cost, which is conducive to the large-scale industry of relevant devices. Numerous studies have been reported using composite substrates bonded with SiC and Ga₂O₃ so as to prepare high-frequency devices36 and high-power devices37–39.

Although SiC possesses noticeable thermal management advantages, there still exists the problem of low interface thermal conductivity (ITC), resulting in the ineffective heat conduction to the high-thermal conductivity substrate, limiting its applicability in high-power systems40. Song et al. introduced SiNx into the interface to enhance the ITC of Ga₂O₃ and SiC bonding, but the value was still only 22 MW/m²·K25. Chen et al. bonded a Ga₂O₃ single crystal to a SiC substrate via the ion-cutting technique27. As a result, the ITC was increased from 72 to 100 MW/m²·K once the thickness of the Al₂O₃ interlayer was decreased from 30 to 9 nm. Xu et al. found that the introduction of Al₂O₃ between SiC and Ga₂O₃ lowered the ITC to 70 MW/m²·K24. Moreover, the interfacial microstructure has a pronounced impact on the thermal conduction efficiency41. Good lattice matching and low interfacial roughness enable to reduce phonon scattering and thus increase the ITC42. On the contrary, defects and impurities within the interface increase the thermal resistance and decrease the heat transfer efficiency. Research indicates that for diamond-based Ga₂O₃ devices, the reduction in junction temperature can be optimized by increasing the ITC to 200 MW/m²·K35. The researchers found that covalent binding can increase ITC by 1–2 orders of magnitude compared to Van der Waals force (VdW) binding43,44. In addition, the local molecular environment between heterogeneous materials can regulate the ITC45. In particular, the enhancement of phonon transport occurs when the chemical bond changes from a weak VdW interaction to a strong covalent bond. M-D. Losego et al. also found through research that as the number of covalent bonding sites increases, the coupling of phonons at the interface will be stronger, resulting in a significant increase in interfacial thermal conductivity46. Therefore, properly adjusting the interface microstructure and constructing covalent bonds are among the key strategies to furtherly improve the ITC of SiC-based Ga₂O₃.

In view of the above, this work achieved high ITC by constructing a covalent bond between Ga₂O₃ and SiC. Through precise interfacial engineering, a stable and flat SiO₂ interlayer was successfully formed, facilitating the formation of strong covalent bonds between the SiO₂ layer and both adjacent layers—Ga₂O₃ and SiC. As a result, the ITC increased ninefold, reaching 162 MW/m²·K, which is the highest value ever reported so far. The NJTM approach effectively enhanced the dissipation of heat generated under light source irradiation, thereby reducing the operating temperature of 29 °C. Therefore, these results highlight the potential of interfacial covalent bonding as a cutting-edge method for improving thermal management in wide-bandgap semiconductor materials, with substantial academic and practical significance.

Results
Sample preparation and interfacial characterization
In this work, an approximately 70-nm-thick Ga₂O₃ film was deposited on a SiC substrate using atomic layer deposition (ALD), as depicted in Fig. 1a. The microstructural evolution of the interface was monitored using transmission electron microscopy (TEM) in coupling with energy-dispersive X-ray spectroscopy (EDS) (see Fig. 1b–g). Figure 1b displayed a representative region of the Ga₂O₃/SiC interface for the deposited sample, which appeared very flat and exhibited the Ga₂O₃ thickness of about 71.4 nm. The corresponding EDS elemental maps and vertical line-scan data indicated that the Ga₂O₃ film was homogeneous and had a single layer structure. The high-resolution TEM (HRTEM) image and selected area electron diffraction (SAED) results, as shown in Fig. 1d, indicate that the as-deposited Ga₂O₃ layer is amorphous, while the underlying SiC substrate exhibits high-quality single-crystal characteristics. A cross-sectional TEM image of the sample after 120 min of annealing at 1173 K is depicted in Fig. 1e. Compared to the untreated samples, the annealed ones clearly demonstrated the formation of a reaction layer between Ga₂O₃ and SiC. Simultaneously, the Ga₂O₃ film became thinner, and the interface became rougher. According to the EDS elemental maps, the intermediate reaction layer was rich in Si and O, while being poor in Ga and C. Figure 1f depicts the elemental content variation in the vertical direction, confirming that the Ga₂O₃ film and SiC substrate maintained their typical elemental percentages. A zoomed-in view of the interface is shown in the inset image of Fig. 1f. The results suggested that the intermediate reaction layer had an oxygen-to-silicon ratio of approximately 2:1, which corresponded to a SiO₂ structure. A HRTEM image of the Ga₂O₃/SiC interface is displayed in Fig. 1g, in which the crystalline regions were ascribed to Ga₂O₃ and SiC, while the SiO₂ layer had a structurally amorphous state. The SAED patterns in Fig. 1g, taken from Ga₂O₃ and SiC regions, confirmed their β-Ga₂O₃ and 4H-SiC structures, corresponding to the [132] and [110] zone axes, respectively. The evolution of microstructure under different annealing conditions was further investigated by TEM, as shown in Supplementary Figs. 2–4. The data show that the thickness of the reactive layer and the interface roughness gradually increase with increasing annealing time and annealing temperature. When annealed at the lower temperature of 1073 K, the thickness of the SiO₂ was slightly thinner, but the overall interface was flatter. When annealed at 1273 K, the growth rate of SiO₂ was significantly higher than that at 1073 K and 1173 K. However, deterioration of the interface and the Ga₂O₃ film, as well as an increase in roughness, was also clearly observed.

To gain deeper insight into the structural composition of each film layer, X-ray diffraction (XRD) and X-ray reflectivity (XRR) measurements were performed on the samples annealed at 1173 K. In addition to the prominent SiC-related peak around 35.6° (Fig. 2a), the annealed Ga₂O₃ exhibited a distinct crystalline peak, while a broad hump appeared near 22°. This was due to the diffuse scattering from the amorphous material, confirming the formation of a SiO₂ layer in the middle of the structure. Figure 2b shows measured and fitted XRR curves of the as-deposited and annealed samples and corresponding information of different layers was plotted in Fig. 2c–e. For the untreated samples, the 5-degree oscillatory cycles were observed, indicating that the interface and film surface were flat and possessed the low roughness. Once the annealing time increased, there was the decay in the oscillations, accompanied by another periodic oscillation in the curve, which corresponded to a thinner SiO₂ film layer. As soon as Ga₂O₃ reacted with SiC during the annealing process, the thickness of the Ga₂O₃ film decreased, while there was little change in the roughness of the SiC substrate. With the further annealing, the crystallinity of the sample increased, enhancing the roughness of the Ga₂O₃ surface. More details regarding the atomic force microscopy (AFM) analysis of the Ga₂O₃ surface morphology are described in Supplementary Fig. 5 of the Supporting Information. However, after 240 min of annealing, a sharp increase in the roughness of both the substrate and the SiO₂ interface was observed, indicating the deterioration of the interfacial structure. A representative region of the Ga₂O₃/SiO₂/SiC interface after annealing for 240 min is displayed in Supplementary Fig. 6. The TEM maps revealed humps of approximately 3–8 nm in height, which were distributed across the entire interface, and the EDS maps suggest that these were SiO₂. Moreover, the TEM characterization confirmed that the prolonged heat treatment led to the formation of an undulating and uneven interface, which coincided with the larger roughness observed in the XRR results.

To further elucidate the bonding conditions at the interface and the changes in chemical composition, the X-ray photoelectron spectroscopy (XPS) experiments were conducted across the Ga₂O₃/SiC heterojunction at various etching depths. The depth-profiled XPS spectral maps for the Ga 3s, Si 2s, and O 1s core levels of the Ga₂O₃/SiC heterojunction are displayed in Fig. 3 along with the fitted XPS spectra. Figure 3a depicts a pronounced interface between the Ga₂O₃ film and SiC substrate. As seen from Stage I of Fig. 3a, the representative O 1s XPS spectrum within the Ga₂O₃ film consisted solely of the O-Ga peak47. At the interface, the Ga 3s core-level spectrum exhibited a modest downward band bending, whereas the Si 2s component shifted toward the higher binding energies. Additionally, the Si 2s core-level spectrum revealed the presence of a Si-O bond at 153.93 eV alongside the Si-C bond at 151.92 eV, which were interpreted based on reports48,49, indicating the weak Si-O bonding at the interface. In Stage VI, a low-intensity O-Ga peak remained still detectable despite etching up to the SiC substrate. This was likely due to a small amount of the product sputtering onto the fresh surface during the etching process. From the spectral maps of annealed samples, it was concluded that the surface layer had a pure Ga₂O₃ structure similar to the untreated specimen (Fig. 3c). As soon as the etching depth increased to 65 nm, the emergence of signals related to the Si-O and O-Si bonds was observed in Si 2s and O 1s spectra, respectively. This indicated that the reaction between Ga₂O₃ and SiC at high temperatures facilitated the formation of SiO₂. To analyze the chemical bonding in detail, a split-peak fit of the core-energy levels was performed at specific locations, as shown in Fig. 3d. The signal associated with the Si-O bond was located at a binding energy of approximately 154.55 eV, while that ascribed to the Si-C bond was centered at 151.05 eV. Taking the binding energy of the Si-C bond as a reference, the difference was about 3.50 eV. For the untreated sample, the binding energy of the Si-O bond was only 2.01 eV higher than that of the Si-C bond, providing evidence that the heat-treated samples possessed a higher binding energy at the interface. This suggests that strong covalent bonding between Ga₂O₃ and SiC was enabled by the SiO₂ layer, as opposed to the VdW force connections formed directly in the untreated samples. Similarly, from the O 1s core-level mapping, it can be found that the covalent bond between Ga₂O₃ and SiO₂ is constructed by O atoms. With a further increase in etching depth, the Si-O bonds were transformed to the Si-C ones.

To describe the reaction between Ga₂O₃ and SiC, the interaction between both components was simulated by combining density functional theory (DFT) and molecular dynamics (MD) methods. Calculating the reaction energy enables to predict whether the reaction will proceed spontaneously or not. The heat released during the reaction between Ga₂O₃ and SiC was calculated to be −460 kJ/mol at 0 K, indicating its exothermic nature under standard conditions. The corresponding reaction can be expressed as follows:

Ga₂O₃ + SiC → SiO₂ + CO(g) + 2Ga (1)

where Ga₂O₃ reacted with SiC to produce Ga, SiO₂, and CO gas. At 0 K, the negative heat value indicated that the reaction occurred spontaneously and was accompanied by the heat release50. At 1173 K a continuous heat input from the environment promoted the ongoing reactions. The thermodynamic analysis confirmed that Ga₂O₃ reacted with SiC to form SiO₂. To gain a deeper understanding of the reaction mechanism and kinetic behavior at the molecular level, an MD simulation was conducted at 2000 K. The results are depicted in Fig. 3e–h. First, the β-Ga₂O₃ is annealed at a high temperature of 2600 K to obtain a disordered amorphous Ga₂O₃ structure, as shown in Fig. 3e. From Fig. 3f, the radial distribution function has multiple peaks within 4 Å, and gradually approaches 1 after it is greater than 4 Å, indicating that the system exhibits a short-range ordered and long-range disordered structure, corresponding to the amorphous state. The obtained amorphous Ga₂O₃ and 4H-SiC are used to construct a heterostructure as shown in Fig. 3g. After a thermal reaction at 2000 K for 1 ns, it can be found from the interface that some of these Si atoms will completely detach from SiC and form the structure of SiO₂ tetrahedra inside Ga₂O₃ after a reaction of 1 ns. In the local magnified image, SiO₂ crystal structure can be clearly found, proving that Ga₂O₃ and 4H-SiC can form SiO₂ at the interface at high temperature. At the same time, the change curve of the number of Si-O bonds over time is also statistically analyzed as shown in Fig. 3h. It can be observed that in the early stage, with the contact between amorphous Ga₂O₃ and 4H-SiC, the number of Si-O bonds increases rapidly. Afterwards, the bonding rate slows down, indicating that the reaction of Ga₂O₃ and SiC to form SiO₂ is not a slow process. At the same time, the number of Si-O bonds increases steadily over time with no saturation trend, indicating that SiC will continue to be oxidized.

To validate the accuracy of ReaxFF in describing the energy variations during MD simulations, the formation energy of SiO₂ obtained from DFT was compared with that from ReaxFF. As shown in Supplementary Fig. 10a, the two values are in close agreement, with only a ~1.3% difference, confirming that the ReaxFF potential can reliably capture the energy change during the SiO₂ formation process. In addition, to further verify the reliability of the employed potential, we calculated the binding energies of the O–SiC, O–Ga₂O₃, C–SiC, and Si–Ga₂O₃ systems using DFT, and ReaxFF-based MD, as shown in Supplementary Fig. 10b. Both methods exhibit consistent trends. The Si–Ga₂O₃ binding energy is positive, suggesting that Si atoms released from SiC do not tend to bond directly with Ga₂O₃ but are more likely to undergo alternative reactions such as SiO₂ formation. This observation is consistent with the reaction pathway identified in our MD simulations of Ga₂O₃ and SiC.

Characterization of interfacial heat transfer properties
To investigate the covalent bond formation within the Ga₂O₃/SiC interface based on the heat transfer characteristics, the ITC values were measured via the time-domain thermoreflectance (TDTR) experiment, as depicted in Fig. 4a. The experimental details of the TDTR measurements are provided in Supplementary Section 2. The thickness of SiO₂ layer as a function of time at different temperatures is shown in Fig. 4b and corresponding ITC values as a function of thickness of SiO₂ layer are presented in Fig. 4c. The amplitudes of TDTR signals versus the delay time and their best fits are shown in Supplementary Figs. 7–9. The corresponding results are summarized in Supplementary Table 4. The ITC value of the untreated sample was found to be 17.5 MW/m²·K. After covalently bonded treatment in 1173 K, the emergence of covalent bonds led to an exponential increase in ITC to 100.9 MW/m²·K. Once the thickness of SiO₂ layer increased to 5.8 nm, the ITC achieved a value of 162.1 MW/m²·K, which was nine times higher than that of the sample undergoing the VdW interactions. At 1073 K, SiO₂ can be formed after 60 minutes, but at a much slower rate compared to the 1173 K conditions. For similar thickness, its ITC is also slightly lower than that of the 1173 K annealing conditions. When the temperature is raised to 1273 K, the rate of SiO₂ formation increases, but the quality significantly decreases, resulting in a final ITC much lower than that of the 1173 K annealing conditions. This result demonstrated that the low interface thermal resistance (ITR) can be achieved through the formation of covalent bonds, enhancing heat transport compared to the VdW interactions, which was consistent with the previous research51,52. SiO₂ itself has a low thermal conductivity, and an excessively thick layer of SiO₂ significantly increases the interfacial thermal resistance. However, since the phonon mean free path of SiO₂ is about 6 nm, the thermal conductivity decreases further when its thickness is less than 5 nm. In addition, the interfacial mass also has an important effect on the thermal conductivity. Therefore, under the premise of maintaining covalent bonding, it is experimentally measured that when the thickness of the SiO₂ layer is 5.8 nm, an optimal balance between thickness and interfacial mass can be achieved, resulting in the highest interfacial thermal conductivity. Meanwhile, the thermal conductivity of the Ga₂O₃ film before treatment was 5.7 W/m·K. The heat treatment eliminated defects within the Ga₂O₃ film, causing thermal conductivity to achieve to 6.2 W/m·K. Nevertheless, the thermal conductivity value of the defect-optimized film remained lower than the levels of 10–27 W/m·K reported for bulk materials15,16. This discrepancy is due to the fact that the phonon mean free path of Ga₂O₃ is on the order of hundred nanometers (up to 700 nm) in the [010] direction53. Thus, in the layers thinner than the phonon mean free path, the phonon-boundary scattering significantly affects the thermal conductivity values, especially for films thinner than 100 nm27. The room-temperature ITC values of Ga₂O₃ with different high-thermal conductivity substrates (GaN, Sapphire, SiC, and Diamond) are summarized in Table 1, indicating that the ITC values obtained in this work exceeded those from other reports. However, when the thickness of SiO₂ layer increased to 6 nm, the ITC exhibited a slight decrease to 148.4 MW/m²·K, primarily due to the deterioration of the interfacial structure caused by the prolonged heat treatment. This led to a significant increase in interface defects, such as atomic mixing, disorder, dislocations and roughness, which altered the phonon transmission path and reduced thermal transfer efficiency54,55. For instance, multiple reflections at the interface impaired phonon transmission, while increasing the interface roughness and enhancing diffuse reflection of phonons.

An in-depth analysis was conducted to elucidate the interfacial heat transfer mechanisms of covalent bonding and VdW force contacts using molecular dynamics (MD) simulations, with the computational details provided in Supplementary Section 1. Three distinct interface structure models for the Ga₂O₃/SiC and Ga₂O₃/SiO₂/SiC trilayer composite structure were constructed, as depicted in Fig. 4d, g, j. Additionally, to validate the accuracy of the constructed model, binding energies were calculated in non-reactive molecular dynamics simulations. The results are shown in Supplementary Fig. 10b. It can be observed that these results are consistent with those obtained from ab initio and MD calculations with the ReaxFF potential. The relevant heat flux results are displayed in Fig. 4e, h, k. The heat flow value of 400 eV/ns associated with covalent bonding was found to be higher by an order of magnitude than that relying on VdW forces. In the covalently bonded interface, the thermal resistance was drastically reduced from the hot source to the cold end. The corresponding temperature profiles are illustrated in Fig. 4f, i, l. Compared to the direct Vander Waals contact between Ga₂O₃ and SiC, the insertion of a SiO₂ layer introduces an additional van der Waals-bonded interface, resulting in a higher ITR due to the presence of this extra interface and the low thermal conductivity of SiO₂. Figure 4i also indicated that the ITR between SiC and SiO₂ was smaller than the respective value between Ga₂O₃ and SiC. This discrepancy is attributed to stronger VdW forces between SiC and SiO₂, resulting in a tighter bond at the interface, which promoted heat transfer and consequently decreased the ITR. It was evident that the temperature difference between the Ga₂O₃ and SiC interfaces decreased dramatically, following the formation of covalent bonds. Furthermore, the temperature distribution suggested that the ITC of SiO₂ with Ga₂O₃ was greater than that with SiC, which was due to the better thermal property matching between Ga₂O₃ and SiO₂. Some researchers have also found through simulations of different intermediate layer materials that SiO₂ can act as a phonon bridge, transmitting phonon vibrations on both sides of Ga₂O₃ and SiC to improve the ITC56. This work further constructed a covalent bond connection through SiO₂, which greatly improved the thermal conductivity of the interface. Using the heat flow and temperature difference at the interface, the ITC can be calculated as:

G = J / (AΔT) (2)

where J represents the heat flux continuity at the interface, A is the cross-sectional area for energy transfer (900 Å² in this model), and ΔT is the temperature jump at the interface. For the Ga₂O₃, SiO₂, and SiC layers connected through VdW interactions, the ITC values were 20.6 MW/m²·K between Ga₂O₃ and SiO₂ and 27.7 MW/m²·K between SiO₂ and SiC. Once the covalent bonding formed within the Ga₂O₃/SiO₂/SiC structure, there was an order of magnitude enhancement in thermal conductivity. The ITC value between Ga₂O₃ and SiO₂ reached 919.9 MW/m²·K, while that between SiO₂ and SiC was measured to be 773.9 MW/m²·K.

Characterization of thermal management capabilities
The effects of covalent bonding and interfacial thermal conductivity enhancement were demonstrated at the atomic scale using MD. In order to further elaborate the intrinsic mechanism, finite element simulation was carried out on the micrometer scale, and the corresponding simulation structure was shown in Fig. 5a, and more detailed simulation details are shown in Supplementary Fig. 11 and Supplementary Tables 5–6 of the Supporting Information. In the computational results of Fig. 5b–d, it can be seen that when 5.8 nm of SiO₂ forms covalent contacts with Ga₂O₃ and SiC, respectively, the corresponding interfacial temperature changes are consistent with that of Fig. 4, and the interfacial temperature difference (temperature difference = 72.5 K) is significantly lower than that in the case of VdW contact (temperature difference = 90.4 K). In addition, we calculated the heat flux of the covalent contact based on this steady-state temperature distribution, and it is 2.23×109 pW/μm2, which is 187% higher than that of the VdW contact (7.77×109 pW/μm2), and the excellent heat transfer performance of the covalent contact is directly proved at the micrometer scale.

To evaluate the heat dissipation ability of the Ga₂O₃/SiC heterojunction, the VdW and covalently bonded samples were further compared conforming to the experimental design in Fig. 5e. First, samples of the same size were placed on a cold stage and maintained at room temperature. The temperatures of both samples were then increased rapidly upon illumination at 390 mW/cm², reaching thermal equilibrium. Eighty-five seconds after the light was turned off, the temperature dropped rapidly and gradually returned to room temperature, as shown in Fig. 5f. Corresponding temperature distributions of the samples surface captured in situ at different time were shown in Fig. 5g. Moreover, the surface temperature of the as-deposited specimen increased more rapidly than that of the covalent bond sample. Notably, the peak temperature of the as-deposited sample reached a value of approximately 76 °C, which was about 29 °C higher than that of the thermally treated sample. This confirmed that the covalent bonded sample possessed a superior ITC. Supplementary Fig. 12 of Supporting Information presents a histogram illustrating the measured heating and cooling rates for the samples. During the heating process, sample bonded via VdW forces exhibits a rapid rise in surface temperature, primarily due to the inefficient thermal conduction through the SiC substrate. In contrast, covalently bonded sample demonstrates a significantly lower surface heating rate, attributed to its higher ITC. Moreover, the surface temperature of the covalently bonded sample reaches equilibrium more quickly compared to their VdW bonded counterparts. The results highlight significantly improved thermal transport efficiency during both the heating and cooling processes in the covalently bonded sample.

In order to quantitatively evaluate the impact of heat dissipation on the electrical properties of the device, the on-state characteristics of Ga₂O₃ metal-oxide-semiconductor field-effect transistors (MOSFETs) with homogeneous and heterogeneous SiC substrates were further investigated via a technology computer-aided design (TCAD) approach20, as shown in Fig. 5h. In the heterogeneous MOSFET structure, different ITC values between Ga₂O₃ and SiC were set. To observe the effect of ITC on heat dissipation, the lattice temperature maps and zoomed-in views of localized hotspots were acquired, as shown in Fig. 5i. Because of the low thermal conductivity value of Ga₂O₃, the homogeneous substrate exhibited a large high-temperature area within the MOSFET. At the low ITC (17.5 MW/m²·K), heat remained trapped within the Ga₂O₃ layer and was not efficiently conducted to the SiC substrate. In contrast, once the ITC increased to 162.1 MW/m²·K, the SiC substrate effectively conducted heat away from the device, reducing the internal temperatures. The output characteristics of various MOSFET structures are displayed in Fig. 5j at a gate voltage (VGS) of 10 V. Then, the maximum ID as a function of ITC values was summarized in Fig. 5k. The SiC substrate enhanced the peak current of the Ga₂O₃ MOSFET, and the output current degradation of the Ga₂O₃ device on the SiC substrate with 17.5 MW/m²·K of ITC was less pronounced than that of the homogeneous Ga₂O₃ device with the increase in drain voltage. Once the ITC rose from 17.5 to 162.1 MW/m²·K, the output current of the MOSFET increased from 321 to 425 mA/mm, representing a roughly 104 mA/mm enhancement. The maximum output current rises significantly, eventually approaching its limit at around 200 MW/m²·K. Additionally, the ID of degradation as a function of ITC values was summarized in Fig. 5l. When the drain voltage reaches 100 V, the current degradation decreases from 147 to 87 mA/mm as the ITC increases from 17.5 to 162.1 MW/m²·K, indicating that enhancing the ITC significantly alleviates the self-heating effect. When the ITC exceeds 300 MW/m²·K, the performance begins to stabilize near its maximum value. Furthermore, the device exhibits low heat generation during low-voltage operation at 50 V. For devices with an ITC of 162.1 MW/m²·K, the current only decreased by 9 mA/mm from the maximum value, while samples with an ITC of 17.5 MW/m²·K decreased by 67 mA/mm. It is well known that the carrier mobility decreases exponentially with the increasing temperature35,57. Therefore, increasing the ITC of the SiC/Ga₂O₃ heterojunction can effectively improve performance, stability, and reliability of the device. Additionally, the large bandgap of SiO₂ (approximately 8 eV) effectively prevents the leakage current conduction to the substrate, enhancing energy efficiency and thus reducing power consumption, extending device lifetime, and improving electrical performance, particularly in high-power applications58,59.

In summary, the heat dissipation of a Ga₂O₃ film on the SiC substrate was successfully enhanced through a simple route consisting of the formation of covalent bonds, which is compatible with existing production lines. Furthermore, we conducted TCAD simulations by treating the ITC derived from the experimental tests on both covalent and VdW contacts as a variable. The results showed that the increase of ITC from 17.5 to 162.1 MW/m²·K significantly reduced the impact of the self-heating effect, leading to an 83 mA/mm improvement in the maximum output current of the MOSFET along with a notable drop in the local junction temperature. By optimizing the heat treatment process, the microstructure at the interface was adjusted so as to facilitate the formation of robust covalent bonds and the creation of effective phonon transport channels. Time-domain thermoreflectance (TDTR) measurements revealed that the ITC between Ga₂O₃ and SiC after the covalent bond formation achieved a remarkable value of 162 MW/m²·K, representing a nine-fold increase compared to the van der Waals (VdW) force connection (17.5 MW/m²·K). Furthermore, the NJTM approach was applied to enhance the ITC whereby the bonding temperature of samples exposed to intense illumination was reduced by approximately 29 °C relative to those with VdW connections. This significantly enhanced performance and reliability of the appliance. Moreover, in-situ XPS studies during the deep etching confirmed the formation of SiO₂ upon heat treatment. The Si-O and Si-C covalent bonds at the interface ensured the effective phonon transport, in which the SiO₂ layer acted as a phonon bridge between Ga₂O₃ and SiC. MD simulations further confirmed that the formation of covalent bonds considerably increased the heat flow at the interface, thereby reducing its thermal resistance. The ITC value between Ga₂O₃ and SiO₂ was enhanced by approximately 40 times compared to that in the case of VdW forces, reaching an impressive value of 919.9 MW/m²·K. Through the comprehensive experimental and theoretical analyses, it was concluded that introducing a SiO₂ layer between Ga₂O₃ and SiC was conducive to the formation of covalent bonding, significantly enhancing heat diffusion and providing a promising strategy to address the low thermal conductivity of Ga₂O₃.

Methods
Sample preparation
The Ga₂O₃ films were deposited onto 4H-SiC substrates by means of a BENEQ TFS200 ALD system (BENEQ Inc., Espoo, Finland). The growth temperature was maintained at 473 K using Trimethylgallium (TMGa) and high-purity O₂ (99.995%) as Ga and O precursors and N₂ as the carrier gas. The precursor for the oxygen was O₂ plasma, which was activated at a power of 200 W. For every cycle, TMGa and O₂ plasma were applied for 50 ms and 8 s, respectively. After that, the chamber was purged for 2 s with Ar gas at the end of each step. The thickness of Ga₂O₃ films was kept at around 75 nm. The fabricated Ga₂O₃ films were annealed in a tube furnace under an O₂ atmosphere. To investigate the effect of annealing conditions, post-annealing was conducted at 1073 K for 60, 120, and 240 minutes; at 1173 K for 20, 60, 120, and 240 minutes; and at 1273 K for 20, 60, and 120 minutes.

Characterization
The XRD measurements (Bruker D8 ADVANCE) were carried out using Cu Kα radiation (λ = 1.54 Å) at a scanning speed of 2.5°/min to obtain the information about the crystal structure, surface morphology, and microstructure of Ga₂O₃ thin films. The X-ray reflection (XRR) was tested by D8 Discover from Bruker. The cross-sectional microstructures and elements distributions were observed via the transmission electron microscopy (TEM, Talos F200x) coupled with energy-dispersive X-ray spectroscopy (EDS). The cross-sectional TEM lamellar samples were fabricated by means of a dual-beam focused ion beam (FIB) system (Thermo Scientific, Helios 5 CX). X-ray photoelectron spectroscopy (XPS, Shimadzu AXIS SUPRA+) was employed to characterize the chemical state and the interior of samples after exposure to Ar+ ion etching. The XPS spectra were calibrated with reference to the C 1s peak (284.8 eV) from adventitious carbon. This work employs the Austin Pioneer ONE TDTR system for measurement for the thermal characteristics, as described in detail in Supplementary Section 2. All characterizations were conducted at room temperature. Infrared thermography (FLIR A325sc) was utilized to in situ monitor the temperature evolution of the samples.

Simulation and calculation
Density functional theory (DFT) calculations were carried out to evaluate the reaction energies between Ga₂O₃ and SiC. The DFT simulations were done using the projector augmented plane-wave basis embedded in Vienna ab-initio simulation package60. During simulations, the plane waves were cut-off at 400 eV. The exchange-correlation functionals of electrons were described in terms of generalized gradient approximations proposed by Perdew, Burke, and Ernzerhof61. A DFT-D3 approach was utilized to capture the long-range VdW interactions62. The Brillouin zone was sampled by employing resolutions better than 0.05 Å−1, using the Monkhorst-Pack scheme63. For the heterostructure model, a vacuum of 15 Å was applied along the thickness direction so as to eliminate the interactions between periodic images. The energy converge criterion for solving the self-consistent Kohn-Sham equations was kept at 10−5 eV. All the structures were fully optimized until the residual Hellman-Feynman forces and stresses were below 0.05 eV/Å. The reaction energy was calculated as follows:

E_reax = E[CO(g)] + E[SiO₂] + 2E[Ga] - E[Ga₂O₃] - E[SiC]

where E[CO(g)] is the energy of CO gas, and E[Ga] is the energy of Cmce-Ga. The Langevin dynamical algorithm was utilized to simulate the NpT ensemble64,65 at the time interval of 2.0 fs.

Molecular dynamics (MD) simulations with LAMMPS package66 were performed to examine the reaction process between Ga₂O₃ and SiC that leads to the formation of SiO₂. The reactive force fields (ReaxFF) are utilized to describe the interaction potentials between atoms67,68:

E_system = E_bond + E_over + E_under + E_lp + E_val + E_tor + E_vdW + E_Coulomb + E_qeq

in which E_bond denotes the bonding energy, E_over and E_under are energy penalty to prevent over- and under coordination, E_lp is the energy of long-pairs of atoms, E_val and E_tor are three-body energy from valence angle and four-body torsions, E_vdW is van der Waals dispersion energy, E_Coulomb is electrostatic energy and E_qeq is the charge equilibrium energy67,69,70. The force field parameters are in line with former study, which is suitable to describe the interations between SiC, Ga and oxygen71. The timestep for integration of equation of motion is 0.1 fs. For the simulation of annealing to obtain amorphous Ga₂O₃ film structure, the NpT ensemble with temperature of 2600 K is utilized and 50 ps is to conducted to drive the system into liquid phase. Then, another 100 ps simulations is utilized to quench the temperature to 300 K to obtain amorphous. For the simulation of SiC/amorphous Ga₂O₃, 1 ns simulation is conducted in total, under NpT ensemble. To accelerate the chemical processes, temperature is set to 2000 K.

The MD simulation also was analyzed the heat transfer behavior in Ga₂O₃/SiC and Ga₂O₃/SiO₂/SiC heterostructures. The whole process was done using LAMMPS molecular dynamics software. The cvff parameter was chosen to describe the force field, and periodic boundary conditions were applied. The equilibrium molecular dynamics simulations were carried out within a certain system, and the NVT equilibrium system was employed. The kinetic time was set to 10 ns and 1 fs as the time step. The NPT simulations were then performed under a standard atmospheric pressure of 1 kPa and a temperature T of 300 K. The running kinetic time was 10 ns. The thermal conductivity along the x-direction was assessed by means of the non-equilibrium molecular dynamics. The right-end region of the model was set as a cold source (200 K), while the left-end of the module was set as a hot source (370 K). After 1 ns of structural relaxation, the temperature distribution tended to achieve dynamic equilibrium over the space along the x-direction. The cross-sectional area of the model was 30×30 Å. The formulas and parameters for non-bonded interactions and covalent bond interactions were shown in Supplementary Section 1.

The formation energy of SiO₂ and the binding energies with different systems were calculated. The calculations of density functional theory are done with the projector augmented plane-wave basis, which is implemented in Vienna ab-initio simulation package60. And the plane-waves are cut-off at 400 eV. The exchange-correlations of electrons are described by the generalized gradient approximations with the form proposed by Perdew, Burke, and Ernzerhof61. DFT-D3 approach is utilized to capture the long-range van der Waals interactions62. The Brillouin zone is sampled with resolutions better than 0.05 Å−1, using the scheme of Monkhorst-Pack63. A vacuum of 15 Å is placed at thickness direction to eliminate the interactions between periodic images. The energy converge criterion for solving self-consistent Kohn-Sham equations is 10−5 eV. All the structures all fully optimized until the residual Hellman-Feynman Forces and stresses below 0.05 eV/Å. The force field calculations are done with LAMMPS package66. The reactive force fields (ReaxFF) are utilized to describe the interaction potentials between atoms67,68:

E_system = E_bond + E_over + E_under + E_lp + E_val + E_tor + E_vdW + E_Coulomb + E_qeq

in which E_bond denotes the bonding energy, E_over and E_under are energy penalty to prevent over- and under coordination, E_lp is the energy of long-pairs of atoms, E_val and E_tor are three-body energy from valence angle and four-body torsions, E_vdW is van der Waals dispersion energy, E_Coulomb is electrostatic energy and E_qeq is the charge equilibrium energy67,69,70. The force field parameters are in line with former study, which is suitable to describe the interations between SiC, Ga and oxygen71. The binding energy is calculated via

E_b = E[AB] - E[A] - E[B]

where E[AB] is the energy of linked system, and E[A=B] is the energy of individual system. The formation energy of SiO₂ is defined by:

E_form = E[SiO₂] - E[Si] - E[O₂(g)]

where E[SiO₂] is the energy of α-SiO₂, E[Si] is the energy of Si, E[O₂(g)] is the energy of oxygen molecule.

In addition, TCAD simulations were used to model two types of Ga₂O₃-based MOSFETs, namely Ga₂O₃ and Ga₂O₃/SiC devices, in order to investigate their electrical characteristics. Both kinds of devices had the same structure which consisted of a 20-nm-thick Al₂O₃ dielectric, a ~ 300-nm-thick Ga₂O₃ drift layer with the low dopant content of 5×1017 cm−3, and a semi-insulating Ga₂O₃ or SiC substrate with the lower dopant of 1×1013 cm−3. The total device length, source-to-gate spacing (LGS), gate length (LG), and gate-to-drain spacing (LGD) were set to 40, 11, 2, and 11 µm, respectively. The source and drain had highly dopant content of 1×1019 cm−3.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "The article presents an interesting study on improving the thermal performance of Ga₂O₃-based devices. The authors provide a clear motivation for their work and conduct experimental measurements and simulations. However, in its current form, I believe the article should not be accepted at this stage; it requires significant revision. Below is a brief list of key issues that need to be addressed:\n- The paper is challenging to read. It would be beneficial to separate the methods and results sections.\n- While the motivation for thermal transport management in the studied structures is well described, the significance of Ga₂O₃ itself remains unclear. The authors should elaborate on why this material is crucial for the study.\n- Figure 1 should be entirely revised. Some data is more suitable for supporting information than the main text; it is unclear if it is simulations, experimental results, etc. This can probably partially be correlated with the paper's final results.\n- The authors should perform a thorough thermal characterization of the fabricated structures by measuring ITC. Without this, the study lacks logical coherence.\n- The details of the atomistic simulations are insufficiently presented. The paper briefly mentions the methods used (ab initio, ab initio MD, and classical MD), but the specific simulation parameters and setup are missing. Without this information, it is impossible to assess the validity of the results. After providing this data, the paper should be reviewed again\n- Despite the main context of the paper being thermal management of the Ga₂O₃-based devices, the final thermal characterization is relatively pure. The authors should analyze the fabricated systems more deeply under different conditions. On the contrary, the sections connected with fabrication and characterization are somewhat oversaturated.\n- A more direct connection between experimental and atomistic results should be added, for instance, by the use of FEM simulations.\nIn summary, while the study is relevant and well-motivated, the manuscript requires substantial improvements before being considered for publication. Specifically, the abovementioned issues make the scientific evaluation of the paper quite challenging.", "reviewer_2": "This manuscript focuses on the thermal management of Ga2O3—a wide-bandgap semiconductor anticipated for next-generation high-power, high-frequency devices—and reports a significant improvement in thermal transport properties by controlling its interaction with the substrate. The combination with SiC, which has traditionally exhibited limited heat dissipation due to high interfacial thermal resistance, is demonstrated to improve its device performance. The developed sample synthesis method succeeded to form a SiO2 layer at the interface, altering the bonding state and dramatically improving the interfacial thermal conductance (ITC); indeed, the reported ITC is the highest among those documented for Ga2O3/high thermal conductivity substrate systems.\nWhile the reported heat dissipation performance is impressive, the concept of enhancing ITC by transforming the interface interaction from a weak vdW type to a stronger covalent bond is not novel. In its current form, the scientific impact does not appear sufficient for publication in a high-impact journal such as Nature Communications. The following points require further discussion:\n1. Ref. 59 discusses, based on theoretical calculations, that introducing a SiO2 intermediate layer between Ga2O3 and SiC can reduce interfacial thermal resistance. The authors should provide a more detailed comparison between the results of Ref. 59 and their own findings.\n2. The impact of chemical bonding at the interface on thermal transport properties has been widely discussed, particularly regarding the enhancement of phonon transport when the bonding shifts from weak vdW interactions to strong covalent bonds. For example, this has been addressed in studies such as M-D. Losego et al., Nat. Mater. 11, 502 (2012) and S. Majumdar et al., Nano Lett. 17, 220 (2017). As these works clearly predate the authors’ approach, the manuscript should clearly position this study within the existing literature while citing these references.\n3. The authors frequently use expressions implying that the formation of a SiO2 layer leads to the formation of covalent bonds between Ga2O3 and SiC; however, this lacks precision. It should be clearly stated that the covalent bonds occur between the intermediate SiO2 layer and the adjacent layers on both sides.\n4. In the AIMD simulation presented in Fig. 3(e), the initial structure appears to be a stacking of crystalline Ga2O3 and SiC. However, according to Fig. 2(d), the untreated sample should exhibit an amorphous Ga2O3 layer. Therefore, it should be explicitly stated that this simulation does not replicate the experimental conditions. Moreover, the simulation does not seem to demonstrate the formation of a SiO2 intermediate layer.\n5. The manuscript should include a more detailed description of the TDTR experimental conditions, such as the thickness of the transducer layer and the methods used to calculate the heat capacities of each layer. Additionally, the TDTR signals for all samples should be plotted together in a single figure to qualitatively assess whether there are significant differences in the cooling curves.\n6. If the calculations indicating that the reaction between Ga2O3 and SiC proceeds spontaneously are valid, it suggests that SiO2 could form even at lower annealing temperatures. Then, wouldn't it be possible to optimize the annealing conditions by adjusting both temperature and time?\n7. In the MD simulations shown in Fig. 4(d)-(i), the authors state that they have simulated two types of structures—one with vdW bonding and one with covalent bonding. The authors should clarify how these bonding types are differentiated in their calculations, including a description of the structural parameters used. Moreover, since the authors consider the Ga2O3/SiC interface to be governed by vdW interactions, shouldn't the comparison be between a Ga2O3/SiC structure and a Ga2O3/SiO2/SiC structure?\n8. In the Conclusion, immediately after discussing covalent bonding, thee TCAD simulation is mentioned. However, it should be clarified that this simulation is not related to the context of covalent bonding and is merely a simulation where the ITC is treated as a variable.", "reviewer_3": "Shen and coworkers present a study of the thermal boundary resistance in Ga2O3/SiC interfaces, proposing this system for its efficient thermal management characteristic for Ga2O3 electronics.\nIn a nutshell, the authors show that a proper post-growth annealing cycle improve the quality of the Ga2O3/SiC interface, which is not by itself a very surprising result. Nevertheless, the meticulous analysis performed of how the progressive formation of a quasi-epitaxial interface (i.e., the formation of covalent bonds) determines the thermal boundary resistance is interesting and suggest that the manuscript perhaps may be considered for publication in Nature Communications.\nThere are, however, some important issues that must be convincingly addressed.\n1) The authors show that the quality of the interface improves and the thermal boundary resistance decreases with a suitable annealing cycle. This thermal treatment promotes the formation of covalent bonds and the creation of a thin SiO2 buffer layer, where O from the Ga2O3 and Si from the SiC react.\nHowever, SiO2 thin films, even in crystalline form, should have a rather low thermal conductivity (see e.g. https://www.nature.com/articles/s41598-018-28925-6). Therefore, there is a trade off in the proposed strategy that is neither discussed, nor mentioned or considered. Promoting the formation of a SiO2 buffer layer is beneficial for thermal management (it provides a dissipation channel trough SiC with a low thermal boundary resistance), but if this layer becomes too thick the whole system becomes too resistive (the contribution of the intrinsic thermal resistance of SiO2 to the overall thermal resistance cannot be neglected).\nWhat is the optimal thickness and how can it be controlled? The tradeoff between thin and thick SiO2 (low thermal boundary resistance, but high SiO2 film resistance) layer must be addressed.\n2) The authors state that: “Meanwhile, the heat treatment eliminated defects within the Ga2O3 film, causing thermal conductivity to achieve to 6.2 W/m·K”.\nHowever, it is not specified what was this value before the treatment. Is 6.2 W/mK large? Is it small?\n3) The AIMD results do not seem to be very significative. At least judging from the snapshots of Figure 3(d) and 3(g), the AIMD simply shows that a thin SiO2 layer sandwiched between SiC and Ga2O3 forms covalent bonds, starting form the initial vdW configuration. This appears to be a rather trivial and largely expected result. The challenge would have been demosntrating the formation of the SiO2 layer bewteen SiC and Ga2O3. While I do understand that this probably too challenging, I don0t think that “the simulation data confirmed that Ga₂O₃ and SiC can form a Si-260 O covalent bond at the interface under high temperatures”.\nOther remarks:\na) The authors say that their results “confirm that the Ga₂O₃ layer was amorphous”, implying that this was expected. Should not it be already crystalline, after the ALD?\nb) Beware of “Colding” in place of “cooling” in Figure 5(b)."}, "author_response": "Dear Reviewers and Editor,\nWe are very grateful to your critical comments and thoughtful suggestions. Based\non these comments and suggestions, we have made careful revision on the original\nmanuscript entitled “Advanced covalent bond engineering for improving interfacial\nthermal conductivity between Ga₂O₃ and SiC” (Manuscript ID: NCOMMS-25-\n15608A). We responded point by point to the reviewer comments as listed below, along\nwith an indication of the location of the revision (all marked with red color in the main\ntext of the revised manuscript).\nIn specific response to the points raised by reviewer 1:\nGeneral comment: The article presents an interesting study on improving the\nthermal performance of Ga O -based devices. The authors provide a clear motivation\n2 3\nfor their work and conduct experimental measurements and simulations. However, in\nits current form, I believe the article should not be accepted at this stage; it requires\nsignificant revision.\nResponse: Thanks for your careful reading of our manuscript and the recognition\nof the experimental contribution of our work. We have responded point by point to the\ncomments listed below.\nComment 1: The paper is challenging to read. It would be beneficial to separate\nthe methods and results sections.\nResponse: We would like to thank the reviewer for carefully reading our\nmanuscript and providing constructive suggestions. Based on the reviewer's comment\nand the formatting requirements of Nature Communications, we have placed the\nMethods section at the end of the manuscript. Additionally, to improve the readability\nof the manuscript, we have divided the Results and Discussion section into three parts\nand added subheadings. These parts are:\n(2.1) Sample Preparation and Interface Characterization;\n(2.2) Analysis of Interface Thermal Conductivity (ITC);\n(2.3) Interface Heat Transfer Performance Verification.\nComment 2: While the motivation for thermal transport management in the\nstudied structures is well described, the significance of Ga O itself remains unclear.\n2 3\nThe authors should elaborate on why this material is crucial for the study.\nResponse: Thank you for your suggestion, and I greatly admire your deep insight\ninto this specific subfield. Based on your comment, we have added a section in the\nintroduction emphasizing the importance of Ga O as a semiconductor material to\n2 3\nsupport the significance of our work. The added section is highlighted in red in the\nmanuscript and can also be found in the Change 1 section below.\n(Change 1) Line 13-25, Page 3 in the Manuscript:\nBeyond its outstanding electronic properties, Ga O possesses a combination of\n2 3\nunique material advantages that make it a highly promising candidate for next-\ngeneration semiconductor applications1. It demonstrates excellent thermal and\nchemical stability, maintaining performance under extreme conditions such as ultra-low\ntemperatures and strong radiation environments2. Additionally, Ga O exhibits strong\n2 3\nabsorption in the deep-ultraviolet (DUV) spectral region (~250–280 nm), making it\nhighly suitable for solar-blind UV photodetectors and space-based optoelectronic\napplications3. From a materials perspective, Ga₂O₃ is advantageous in both bulk and\nepitaxial growth. It supports melt-based bulk crystal growth techniques such as edge-\ndefined film-fed growth (EFG) and Czochralski (CZ), enabling the production of large-\narea, low-cost single crystals—unlike SiC, which requires expensive high-temperature\nvapor-phase growth4. At the same time, Ga₂O₃ shows excellent compatibility with\nvarious heteroepitaxial growth methods, including MOCVD, MBE, and PLD, on\nsubstrates such as sapphire, Si, and SiC5. In addition, Ga O has relatively low hardness\n2 3\ncompared to SiC, reducing wafer slicing and processing costs.\nComment 3: Figure 1 should be entirely revised. Some data is more suitable for\nsupporting information than the main text; it is unclear if it is simulations, experimental\nresults, etc. This can probably partially be correlated with the paper's final results.\nResponse: Thank you for the reviewer’s suggestion. We have divided original\nFigure 1 into two parts based on the content it expresses. Figure 1a has been moved to\nthe Supporting Information as Figure S1 to show the structure of the Ga O device.\n2 3\nThe original Figure 1b–f presents the simulation results. We performed TCAD\nsimulations on Ga O and SiC heterostructures under different interface thermal\n2 3\nresistances, proving the importance of improving the interface thermal conductivity\nfrom the perspective of device electrical characteristics. When the interface thermal\nconductivity is increased, the maximum junction temperature of the device can be\ngreatly reduced, and the maximum output current of the device can be increased.\nTherefore, we have presented this content as Figure 5h–l in the revised manuscript to\nsupport our experimental results. For detailed revisions, please refer to the Change\nsection below.\n(Change 1) Figure S1 in the Supporting Information (Page 9).\nFigure S1: The schematic diagram of the heat dissipation chip structure.\n(Change 2) Figure 5 in the Manuscript (Page 18).\nFigure 5. (a) Schematic diagram of Finite Element Method simulations. (b-c) The\ncomparison of heat dissipation capability based on simulated profiles. (d) The\ntemperature profiles of the different interface. (e) Schematic diagram of the infrared\nthermography measurements. (f) The surface temperature evolution of VdW force and\ncovalent bond samples as a function of time, captured with an infrared camera. (g)\nCorresponding temperature distributions of the samples surface captured in situ at\ndifferent time. (h) Schematic illustrations of the MOSFET device. (i) Simulated 2D heat\nmaps of the cross-section of the MOSFET device. (j) I as a function of V at different\nD D\nITC values. (k) Maximum I as a function of ITC values. (l) I of degradation as a\nD D\nfunction of ITC values.\nComment 4: The authors should perform a thorough thermal characterization of\nthe fabricated structures by measuring ITC. Without this, the study lacks logical\ncoherence.\nResponse: We thank the reviewers for their interest in this issue. We came in a\ncomplete thermal characterization of the prepared structures and measured the overall\nthermal resistance of the samples using time-domain thermoreflectance (TDTR), a\nfemtosecond laser-based non-contact testing technique, in which a metallic film on the\nsurface of the sample (in this paper, we used an Al film) is heated by pumping light,\nand the probe light monitors the change of its reflectivity over time, so as to obtain the\ntemperature decay curves and, by fitting a heat conduction model to extract key\nparameters such as thermal conductivity (κ) or interfacial thermal conductivity (ITC).\nAlthough TDTR measures the overall thermal transport properties of a multilayer\nstructure, we excluded the effects of other variables as much as possible in our\nexperimental design. We chose the same deposition conditions for the Al metal layer to\nensure consistent thermal resistance on its own; we also measured the thicknesses of\nthe Ga₂O₃ film and the Al layer by means of TEM, etc., to minimize the accuracy of the\nfit. Therefore, the difference in thermal resistance between the samples mainly\noriginates from the ITC between Ga₂O₃ and the substrate, which is the core of our study.\nThe thermal properties of each part are listed in Table R1.\nIn high-voltage, high-power devices, Ga₂O₃ faces a serious heat dissipation bottleneck due to its own low thermal\nconductivity, and excessive interfacial thermal resistance will further inhibit its thermal\ndiffusion capability. Therefore, enhancing the ITC is of key significance to enhance the\nheat dissipation capability of thin films and thus improve the reliability and\nperformance of devices. In this study, the positive effects of introducing SiO₂ interlayer\nand constructing covalent bonds to enhance the interfacial thermal conductivity are\neffectively verified through experimental measurements and combined with structural\nmodulation, which has a clear physical basis and application prospect.\nTable R1-1 Thermal characteristics of each part of the fabricated structures\nITC of SiC and Ga₂O₃ κ of Ga₂O₃ ITC of Al and Ga₂O₃\nSample\n(MW/m²·K) (W/m·K) (MW/m²·K)\nAs-deposited 17.5 5.7 159.2\n60 min 90 5.9 151.3\n120 min 97 6.1 158.3\n240 min 105 6.2 155.4\n20 min 100.9 5.9 155.3\n60 min 115 6.2 153.2\n120 min 162.1 6.2 157.4\n240 min 148.4 6.2 152.7\n20 min 73 6.1 151.4\n60 min 75 6.1 147.8\n120 min 98 5.9 143.7\nBased on your suggestions, we have also obtained the ITC of the Al layer and the\nGa O layer. The κ of Ga O is supplemented in Table S5 of the supporting information.\n2 3 2 3\n(Change 1) Table S5 in the Supporting Information (Page 14-15).\nTable S5 Thermal characteristics of each part of the fabricated structures.\nComment 5: The details of the atomistic simulations are insufficiently presented.\nThe paper briefly mentions the methods used (ab initio, ab initio MD, and classical\nMD), but the specific simulation parameters and setup are missing. Without this\ninformation, it is impossible to assess the validity of the results. After providing this\ndata, the paper should be reviewed again.\nResponse: Thank you for carefully reading our work and helping to improve our\nmanuscript. Additionally, we apologize for not including the detailed molecular\ndynamics simulation process in the original manuscript. In the original manuscript,\nthere are two sections involving simulation analysis: one is in Figure 3e–h, where\nmolecular dynamics simulations are used to simulate the reactions between Ga O and\n2 3\nSiC, and the other is in Figure 4, where molecular dynamics is used to analyze the\ninterfacial thermal transport properties of Ga O /SiO /SiC heterostructures. Based on\n2 3 2\nyour suggestion, we have now provided the detailed process of these two simulations\nhere.\n(1) Molecular dynamics simulation of the reaction between gallium oxide and\nsilicon carbide\nThe molecular dynamic simulations are done with LAMMPS package6. The\nreactive force fields (ReaxFF) are utilized to describe the interaction potentials\nbetween atoms7, 8:\n𝐸𝑠𝑦𝑠𝑡𝑒𝑚 = 𝐸𝑏𝑜𝑛𝑑 +𝐸𝑜𝑣𝑒𝑟 +𝐸𝑢𝑛𝑑𝑒𝑟 +𝐸𝑙𝑝 +𝐸𝑣𝑎𝑙 +𝐸𝑡𝑜𝑟 +𝐸𝑣𝑑𝑊 + ·𝐸·· 𝐶 · 𝑜 · 𝑢 ·· 𝑙𝑜 ·· 𝑚 ·· 𝑏 ··(1)\nIn Equation (1), denotes the bonding energy;\n+𝐸𝑞𝑒 𝑞\nand are\nenergy penalty to prev𝐸e 𝑏 n 𝑜 t 𝑛 o 𝑑 ver- and undercoordination; 𝐸i 𝑜 s 𝑣 𝑒 th 𝑟 e energ𝐸y 𝑢 𝑛 o 𝑑 f 𝑒 𝑟 long-\npairs of atoms; and are three-body energy fro𝐸m 𝑙𝑝 valence angle and four-\nbody torsions; 𝐸𝑣𝑎𝑙 is va𝐸n 𝑡 𝑜 d 𝑟 er Waals dispersion energy; is electrostatic\n𝐸𝑣𝑑𝑊 𝐸𝐶𝑜𝑢𝑙𝑜𝑚𝑏\nenergy and is the charge equilibrium energy.7, 9, 10 The force field parameters\nare in line wi𝐸th 𝑞 𝑒 f 𝑞 ormer study, which is suitable to describe the interations between SiC,\nGa and oxygen.11 The time step for integration of equation of motion is 0.1 fs. For the\nsimulation of annealing to obtain amorphous Ga O film structure, the NpT ensemble\n2 3\nwith temperature of 2600 K is utilized and 50 ps is to conducted to drive the system\ninto liquid phase. Then, another 100 ps simulations is utilized to quenching the\ntemperature to 300 K to obtain amorphous. For the simulation of SiC/amorphous\nGa O , 1 ns simulation is conducted in total, under NpT ensemble. To accelerate the\n2 3\nchemical processes, temperature is set to 2000 K.\n(2) The molecular dynamics (MD) simulation of heat transport was performed\nfor the Ga O /SiO /SiC heterostructures\n2 3 2\nThe following are the formulas and parameters for non-bonded interactions and\ncovalent bond interactions. In van der Waals connection, it actually means that there\nare only non-bonded interactions between interfaces; while in covalent bond connection,\nthere are both covalent bonds and non-bonded interactions at the interface.\nIn LAMMPS, pair potentials are defined between pairs of atoms that are within a\ncutoff distance. The standard 12/6 Lennard-Jones potential was adopted, and given by\nEquation (2), where is the cutoff radius.\n𝑟c\n(2)\n𝜎 12 𝜎 6\n𝐸 = 4𝜖[(𝑟) − (𝑟) ], 𝑟 < 𝑟𝑐\nTable R1-2 The σ and ϵ parameters for different atoms.\n12/6 LJ potential (angstrom) (kcal/mol)\nGa 3.87 0.217\n𝝈 𝝐\nO 2.859 0.228\nSi 4.053 0.040\nC 3.474 0.159\nThe harmonic bond style was adopted as the covalent bonding, given by Equation\n(3), where is the equilibrium bond distance, and is the bonding factor.\n𝑟0 𝐾\n(3)\n2\n𝐸 = 𝐾(𝑟−𝑟0)\nTable R1-3 The and parameters for different covalent bonding.\nCovalent bon𝑟0 ding 𝐾 (angstrom) (kcal/mol)\nGa-O 𝑟0 1.932 𝐾 334.5\nSi-O 1.615 392.8\nSi-C 1.809 238.0\nGa-Si 2.106 207.3\nThe simulation details of these two sections have been included as Supplementary\nTexts (Section S2 and Section S3) in the revised supporting information. (Line 9 – 15,\nPage 29).\n(Change 1) Section S1 in the Supporting Information (Line 2 – 19, Page 4).\nSection S2: The details on the molecular dynamics simulation of the reaction\nbetween Ga O and SiC.\n2 3\n(Change 2) Section S2 in the Supporting Information (Line 21 – 27, Page 4\nand Line 1 – 13, Page 5).\nSection S3: The details on the molecular dynamics simulation of interfacial heat\ntransfer.\nComment 6: Despite the main context of the paper being thermal management of\nthe Ga₂O₃-based devices, the final thermal characterization is relatively pure. The\nauthors should analyze the fabricated systems more deeply under different conditions.\nOn the contrary, the sections connected with fabrication and characterization are\nsomewhat oversaturated.\nResponse: We would like to thank the reviewers for their constructive comments.\nBased on your comments, we have added different conditions to the article include 1073\nK (60 min, 120 min, and 240 min) and 1273 K (20 min, 60 min, and 120 min). Figures\nR1-1a and R1-1b show the thickness of the SiO layers and the final ITC under these\n2\nconditions. The corresponding morphological analyses is shown in Figures R1-2 and\nR1-4, and the corresponding TDTR test results are shown in Figures R1-3 and R1-5. At\n2\nthe 1173 K conditions. For similar thickness, its ITC is also slightly lower than that of\nthe 1173 K annealing conditions. When the temperature is raised to 1273 K, the rate of\nSiO formation increases, but the quality significantly decreases, resulting in a final\nannealing condition we selected strikes a balance between the rate and quality of SiO\n2\nformation.\nBased on your suggestion, we have added the results from these different\ntemperatures and times into the revised version. Figures R1-1a and R1-1b have been\nincorporated as Figures 4b and 4c in the revised manuscript, with the corresponding\ndiscussion added. Figures R2 to R5 have been included in the supporting information.\nFigure R1-1. (a) The thickness of SiO layer as a function of time at different\n2\ntemperature. (b) The ITC values established via the TDTR measurements at different\nthickness of SiO layer.\nimages of the annealed sample with (a) 60 min, (c) 120 min, and (e) 240 min at 800 ℃.\nThe EDS line scanning profiles of cross-section under annealing (b) 60 min, (d) 120\nmin, and (f) 240 min at 800 ℃.\nFigure R1-3. (a) TDTR amplitude signals of different samples annealed at 1073 K as\na function of delay time. The best fits as a function of delay time: (b) 60 min, (c) 120\nmin, and (d) 240 min.\nFigure R1-4 The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 20 min, (c) 60 min, and (e) 120 min at 1000 ℃.\nThe EDS line scanning profiles of cross-section under annealing (b) 20 min, (d) 60 min,\nand (f) 120 min at 1000 ℃.\nFigure R1-5. (a) TDTR amplitude signals of different samples annealed at 1273 K as\na function of delay time. The best fits as a function of delay time: (b) 20 min, (c) 60\nmin, and (d) 120 min.\nDetailed changes can be found in the following Change section.\n(Change 1) The corresponding discussion in the Manuscript (Line 25 – 28,\nPage 7 and Line 1 – 2, Page 8).\nThe data show that the thickness of the reactive layer and the interface roughness\ngradually increase with increasing annealing time and annealing temperature. When\nannealed at the lower temperature of 1073 K, the thickness of the SiO₂ thickness was\nslightly thinner, but the overall interface was flatter. When annealed at 1273 K, the\ngrowth rate of SiO₂ was significantly higher than at 1073 K and 1173 K. However,\ndeterioration of the interface and the Ga₂O₃ film, as well as an increase in roughness,\nwas also clearly observed.\n(Change 2) Figure 4b – c in the Manuscript (Line 7 – 9, Page 12).\nFigure 4 (a) Schematic diagram of the sample structure and the TDTR experiment.\n(b) The thickness of SiO layer as a function of time at different temperature. (c) The\n2\n(d, g, f) Structural modeling for molecular dynamics calculations at the Ga₂O₃/SiC\ninterface bonded by (d) VdW and Ga₂O₃/SiO /SiC interface bonded by (g) VdW and (f)\n2\ncovalent interactions. (e, h, k) The accumulative energy of the heat source and heat sink\nof the Ga₂O₃/SiC interface bonded by (e) VdW and Ga₂O₃/SiO /SiC interface bonded\n2\nby (h) VdW and (k) covalent interactions. (f, i, l) Typical temperature profiles of the\nGa₂O₃/SiC interface bonded by (f) VdW and Ga₂O₃/SiO /SiC interface bonded by (i)\n(Change 3) Line 15 - 19, Page 13 in the Manuscript.\nAt 1073 K, SiO can be formed after 60 minutes, but at a much slower rate\n2\ncompared to the 1173 K conditions. For similar thickness, its ITC is also slightly lower\nthan that of the 1173 K annealing conditions. When the temperature is raised to 1273\nK, the rate of SiO formation increases, but the quality significantly decreases, resulting\n2\nin a final ITC much lower than that of the 1173 K annealing conditions.\n(Change 4) Figure S2 and Figure S4 in the Supporting Information (Page 9\nand Page 11).\nFigure S2. The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 60 min, (c) 120 min, and (e) 240 min at 1073\nK. The EDS line scanning profiles of cross-section under annealing (b) 60 min, (d) 120\nmin, and (f) 240 min at 1073 K.\nFigure S4. The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 20 min, (c) 60 min, and (e) 240 min at 1273 K.\nThe EDS line scanning profiles of cross-section under annealing (b) 20 min, (d) 60 min,\nand (f) 240 min at 1273 K.\n(Change 5) Figure S8 and Figure S9 in the Supporting Information (Page 13-\n14);\nFigure S8. (a) TDTR amplitude signals of different samples annealed at 1073 K as a\nfunction of delay time. The best fits as a function of delay time: (b) 60 min, (c) 120 min,\nand (d) 240 min.\nFigure S9. (a) TDTR amplitude signals of different samples annealed at 1273 K as a\nfunction of delay time. The best fits as a function of delay time: (b) 20 min, (c) 60 min,\nand (d) 120 min.\nComment 7: A more direct connection between experimental and atomistic results\nshould be added, for instance, by the use of FEM simulations.\nResponse: Thank you for your suggestion. In the original manuscript, we used\natomic-scale molecular dynamics to analyze the interfacial thermal conduction of\nGa₂O₃/SiO₂/SiC heterostructures in Figure 4. In Figure 5, we employed actual prepared\nsamples and compared the interfacial thermal performance of SiO₂ after forming\ncovalent bonds with Ga₂O₃ and SiC, as well as the differences with the van der Waals\ncontact interface through heating tests. Based on your suggestion, we have added a\nmicro-scale finite element analysis between the two to further support our argument.\nFigure R1-6. (a) Ga O -SiO -SiC three-layer structure thermal simulation model,\n2 3 2\nwith the corresponding interface schematic diagram shown in (b); (c-e) Mesh\ndivision: (c) full mesh, (d) local enlargement, (e) interface mesh.\nThe simulation model used is shown in Figure R1-6a. This simulation model\nconsists of a three-layer structure of Ga₂O₃-SiO₂-SiC, with the corresponding material\nand dimensional parameters listed in Table R1-4. The thicknesses of Ga₂O₃ and SiO₂\nare 68.4 nm and 5.8 nm, respectively, which are derived from experimental testing data;\nthe thickness of SiC is set to be the same as Ga₂O₃, 68.4 nm. The lateral size of the\nmodel is set to 1 × 1 μm, with the plane dimensions much larger than the thickness of\nthe Ga O and SiO films, effectively reducing the boundary effects. Additionally, based\n2 3 2\non the ITC data obtained from atomic-scale molecular dynamics simulation results in\nFigure 4, we assigned the Ga O /SiO interface and SiO /SiC interface in the model.\n2 3 2 2\nThe parameters for the covalent bonded interface and van der Waals (VdW) interface\nare shown in Table S4. Furthermore, referring to the simulation boundary conditions in\nFigure 4, in this thermal analysis, we applied isothermal loads of 373 K and 273 K to\nthe top surface of the Ga O layer and the bottom surface of the SiC layer, respectively,\n2 3\nand calculated the temperature distribution and heat flux through steady-state thermal\nanalysis.\nTable R1-4 The detailed parameters of the components in the simulated system.\nLateral size (μm) Thickness (nm) κ (W/m·K)\nGa O 1 × 1 68.4 6.2\n2 3\nSiO 1 × 1 5.8 1.0\nTable R1-5 The detailed parameters of the interface in the simulated system.\nITC of Ga O /SiO ITC of SiO /SiC\n2 3 2 2\ninterface interface\n(MW/m²·K) (MW/m²·K)\nCovalent bonded interface 919.9 773.9\nVdW interface 20.6 27.7\nIn the computational results of Figure 1-7 (b-d), it can be seen that when 5.8 nm of\nSiO₂ forms covalent contacts with Ga₂O₃ and SiC, respectively, the corresponding\ninterfacial temperature changes are consistent with that of Figure 4, and the interfacial\ntemperature difference (temperature difference = 72.5 K) is significantly lower than\nthat in the case of VdW contact (temperature difference = 90.4 K). In addition, we\ncalculated the heat flux of the covalent contact based on this steady-state temperature\ndistribution, and it is 2.23 × 109 pW/μm2, which is 187% higher than that of the VdW\ncontact (7.77 × 109 pW/μm2), and the excellent heat transfer performance of the\ncovalent contact is directly proved at the micrometer scale.\nFigure R1-7. (a) Schematic diagram of Finite Element Method simulations. (b–c) The\ncomparison of heat dissipation capability based on simulated profiles. (d) The\ntemperature profiles of the different interface.\nIn the above content, we have added Figure R1 as Figure S10 to the revised\nSupporting Information and provided a detailed description of the finite element\nsimulation process. Additionally, Figure R1-7 has been added as Figure 5a-d in the\nrevised Manuscript, with the corresponding discussion included. Detailed changes can\nbe found in the following Change section.\n(Change 1) Figure S10 in the Supporting Information (Page 17);\nFigure S10. (a) Ga O -SiO -SiC three-layer structure thermal simulation model,\n2 3 2\nwith the corresponding interface schematic diagram shown in (b); (c-e) Mesh division:\n(c) full mesh, (d) local enlargement, (e) interface mesh.\n(Change 2) Figure 5a – d in the Manuscript (Page 16).\nFigure 5. (a) Schematic diagram of Finite Element Method simulations. (b-c) The\ncomparison of heat dissipation capability based on simulated profiles. (d) The\ntemperature profiles of the different interface. (e) Schematic diagram of the infrared\nthermography measurements. (f) The surface temperature evolution of VdW force and\ncovalent bond samples as a function of time, captured with an infrared camera. (g)\nCorresponding temperature distributions of the samples surface captured in situ at\ndifferent time. (h) Schematic illustrations of the MOSFET device. (i) Simulated 2D heat\nmaps of the cross-section of the MOSFET device. (j) I as a function of V at different\nD D\nITC values. (k) Maximum I as a function of ITC values. (l) I of degradation as a\nD D\nfunction of ITC values.\n(Change 3) Line 9-19, Page 17 in the Manuscript\nThe effects of covalent bonding and interfacial thermal conductivity enhancement\nwere demonstrated at the atomic scale using MD. In order to further elaborate the\nintrinsic mechanism, finite element simulation was carried out on the micrometer scale,\nand the corresponding simulation structure was shown in Figure 4(a), and more detailed\nsimulation details are shown in Figure S10 and Table S5-6 of the Supporting\nInformation. In the computational results of Figure 4 (b-d), it can be seen that when 5.8\nnm of SiO₂ forms covalent contacts with Ga₂O₃ and SiC, respectively, the\ncorresponding interfacial temperature changes are consistent with that of Figure 3, and\nthe interfacial temperature difference (temperature difference = 72.5 K) is significantly\nlower than that in the case of VdW contact (temperature difference = 90.4 K). In\naddition, we calculated the heat flux of the covalent contact based on this steady-state\ntemperature distribution, and it is 2.23 × 109 pW/μm2, which is 187% higher than that\nof the VdW contact (7.77 × 109 pW/μm2), and the excellent heat transfer performance\nof the covalent contact is directly proved at the micrometer scale.\nIn specific response to the points raised by reviewer 2:\nGeneral comment: This manuscript focuses on the thermal management of\nGa O —a wide-bandgap semiconductor anticipated for next-generation high-power,\n2 3\nhigh-frequency devices—and reports a significant improvement in thermal transport\nproperties by controlling its interaction with the substrate. The combination with SiC,\nwhich has traditionally exhibited limited heat dissipation due to high interfacial thermal\nresistance, is demonstrated to improve its device performance. The developed sample\nsynthesis method succeeded to form a SiO layer at the interface, altering the bonding\n2\nstate and dramatically improving the interfacial thermal conductance (ITC); indeed, the\nreported ITC is the highest among those documented for Ga O /high thermal\n2 3\nconductivity substrate systems. While the reported heat dissipation performance is\nimpressive, the concept of enhancing ITC by transforming the interface interaction\nfrom a weak vdW type to a stronger covalent bond is not novel. In its current form, the\nscientific impact does not appear sufficient for publication in a high-impact journal such\nas\nNature Communications.\nResponse: We sincerely appreciate the reviewer’s thoughtful comments and\nrecognition of the impressive heat dissipation performance demonstrated in our work.\nWe fully understand the concern regarding the novelty of the concept, and in response,\nwe would like to elaborate further on the innovation, technical realization, and scientific\nsignificance of our study. While the general concept of enhancing interfacial thermal\nconductance (ITC) through covalent bonding has been explored in previous theoretical\nstudies and demonstrated in a limited number of material systems, our work presents\nthe first experimental realization of this approach at the Ga₂O₃/SiC heterointerface\nthrough the introduction of a SiO₂ interlayer. By employing a self-limiting interfacial\nreaction strategy during high-temperature growth/annealing, we successfully achieved\nthe controlled formation of an ultrathin (~5.8 nm) amorphous SiO₂ layer at the interface.\nThis SiO₂ layer serves not only as a structural buffer but, more importantly, facilitates\nthe formation of strong covalent bonds between Ga₂O₃ and SiC via C–Si–O and Ga–\nO–Si bonds, effectively transforming the original van der Waals-dominated interface\ninto a covalently bonded one.\nAs a result, we achieved an ITC of ~162.1 MW/m²·K, the highest reported value\namong Ga₂O₃/high thermal conductivity substrate systems to date. Comprehensive\nstructural (HRTEM), chemical (XPS), and thermal (TDTR) characterizations clearly\ncorrelate the ITC enhancement with the modified interfacial bonding state, offering a\nmechanistic understanding that informs future thermal interface designs for ultra-wide\nbandgap semiconductors. Furthermore, our device-level TCAD simulations\ndemonstrate that the enhanced ITC leads to a substantial reduction in self-heating\neffects. This study provides a scalable and generalizable interface engineering strategy,\npushing the boundaries of thermal management in Ga₂O₃-based power devices and\noffering valuable guidance for next-generation high-power electronic systems.\nIn addition, taking your advices, we have considered your proposed comments and\nmade corresponding revisions of our manuscript. We have answered your comments\npoint by point below.\nComment 1: Ref. 59 discusses, based on theoretical calculations, that introducing\na SiO intermediate layer between Ga O and SiC can reduce interfacial thermal\n2 2 3\nresistance. The authors should provide a more detailed comparison between the results\nof Ref. 59 and their own findings.\nResponse: Thank you for your thorough and careful review of our manuscript. We\nnote that Ref. 59 you mentioned represents a preliminary simulation-based study\nconducted by our collaborating team prior to this work. Indeed, in current research on\nenhancing the thermal performance of Ga O devices, the selection of substrate\n2 3\nmaterials and compatible interlayer materials has attracted significant attention.\nTherefore, in earlier work, our collaborators conducted extensive simulations on\ndifferent substrate and interlayer materials. Based on theoretical modeling, they\npredicted that introducing a silicon oxide SiO interlayer between Ga O and SiC could\n2 2 3\neffectively improve interfacial thermal conductance.\nDuring subsequent experimental validation, our initial approach involved directly\ndepositing SiO onto the SiC substrate followed by Ga O deposition. However, the\n2 2 3\nmeasured interface thermal conductivity (ITC) fell short of expectations; it was even\nlower than that of the simple bilayer structure with Ga O deposited directly on SiC.\n2 3\nIn subsequent work, we focused on exploring what type of SiO interlayer would\n2\neffectively enhance the ITC between Ga O and SiC. Through extensive experimental\n2 3\nvalidation, we identified the strategy presented in this work: generating Ga O in-\n2 3\nsitu between Ga O and SiC via direct annealing. The SiO interlayer formed this way\n2 3 2\nestablishes covalent bonds with both Ga O and SiC, increasing theITC from 17.5\n2 3\nMW/m2·K to 162.1 MW/m2·K—an 8.25-fold improvement. This significantly exceeds\nthe 3.7-fold enhancement predicted in our collaborators' earlier simulation (Ref. 59).\nBased on this result, we have also revised the original simulation model for the Ga O -\n2 3\nSiO -SiC trilayer structure discussed in this manuscript.\ninsights screened from prior simulations. Our experimental results further demonstrate\nthat when SiO serves as the interlayer between Ga O and SiC, effective increase of\n2 2 3\nITC requires the formation of covalent bonds at both interfaces—not merely van der\nWaals contact. Compared to the theoretical predictions in Ref. 59, this study represents\nsubstantial progress in experimental design, material fabrication, and actual thermal\nperformance measurement. Furthermore, the experimental conclusions provide\nvaluable supporting information for subsequent simulation analyses. The added section\nis highlighted in red in the manuscript and can also be found in the Change 1 section\nbelow.\n(Change 1) Line 17 – 20, Page 15 in the Manuscript:\nSome researchers have also found through simulations of different intermediate\nlayer materials that SiO₂ can act as a phonon bridge, transmitting phonon vibrations on\nboth sides of Ga₂O₃ and SiC to improve the ITC12. This work further constructed a\ncovalent bond connection through SiO₂, which greatly improved the thermal\nconductivity of the interface.”\nComment 2: The impact of chemical bonding at the interface on thermal transport\nproperties has been widely discussed, particularly regarding the enhancement of\nphonon transport when the bonding shifts from weak vdW interactions to strong\ncovalent bonds. For example, this has been addressed in studies such as M-D. Losego\net al., Nat. Mater. 11, 502 (2012) and S. Majumdar et al., Nano Lett. 17, 220 (2017). As\nthese works clearly predate the authors’ approach, the manuscript should clearly\nposition this study within the existing literature while citing these references.\nResponse: We thank the reviewer for the helpful comment. The suggestion has\ngreatly improved the clarity and depth of the background analysis in our manuscript.\nThe added section is highlighted in red in the manuscript and can also be found in\nthe Change 1 section below.\n(Change 1) Line 21-27, Page 5 in the Manuscript:\nIn addition, the local molecular environment between heterogeneous materials can\nregulate the interfacial thermal conductivity13. In particular, the enhancement of phonon\ntransport occurs when the chemical bond changes from a weak vdW interaction to a\nstrong covalent bond. M-D. Losego et al. also found through research that as the number\nof covalent bonding sites increases, the coupling of phonons at the interface will be\nstronger, resulting in a significant increase in interfacial thermal conductivity14.\nTherefore, properly adjusting the interface microstructure and constructing covalent\nbonds are among the key strategies to furtherly improve the ITC of SiC-based Ga₂O₃.\nComment 3: The authors frequently use expressions implying that the formation\nof a SiO layer leads to the formation of covalent bonds between Ga O and SiC;\n2 2 3\nhowever, this lacks precision. It should be clearly stated that the covalent bonds occur\nbetween the intermediate SiO layer and the adjacent layers on both sides.\nof precision in our expression. Based on your suggestions, we have made revisions to\nthe relevant sections throughout the text. Please refer to the Change section below for\ndetails.\n(Change 1) Line 12-14, Page 2 in the Manuscript\nThese results indicate that covalent bonding between the SiO₂ layer and both Ga₂O₃\nand SiC significantly improves interfacial heat dissipation.\n(Change 2) Line 2-4, Page 6 in the Manuscript\nThrough precise interfacial engineering, a stable and flat SiO₂ interlayer was\nsuccessfully formed, facilitating the formation of strong covalent bonds between the\nSiO₂ layer and both adjacent layers—Ga₂O₃ and SiC.\n(Change 3) Line 2-3, Page 11 in the Manuscript\nSimilarly, from the O 1s core-level mapping, it can be found that the covalent bond\nbetween Ga₂O₃ and SiO is constructed by O atoms.\n2\n(Change 4) Line 1-2, Page 20 in the Manuscript\nThe Si-O and Si-C covalent bonds at the interface ensured the effective phonon\ntransport, in which the SiO₂ layer acted as a phonon bridge between Ga₂O₃ and SiC.\nComment 4: In the AIMD simulation presented in Fig. 3(e), the initial structure\nappears to be a stacking of crystalline Ga O and SiC. However, according to Fig. 2(d),\n2 3\nthe untreated sample should exhibit an amorphous Ga O layer. Therefore, it should be\n2 3\nexplicitly stated that this simulation does not replicate the experimental conditions.\nMoreover, the simulation does not seem to demonstrate the formation of a SiO\n2\nintermediate layer.\nResponse: We sincerely appreciate the reviewer’s valuable comments, which have\nhelped us further improve and enrich this work. Indeed, as you mentioned, simulating\ncrystalline Ga O does not align well with our experimental conditions. Based on your\n2 3\nsuggestion, we have optimized the computational model, firstly by replacing the\nprevious Ab initio molecular dynamics (AIMD) with molecular dynamics (MD), which\nis more suitable for large systems and prolonged simulations. Secondly, therefore to\nobtain the formation of silicon oxide, we extended the reaction time until the formation\nof SiO at the interface. At the same time, the previous crystalline state of Ga O was\n2 2 3\nreplaced by amorphous Ga O in order to be more compatible with the experimental\n2 3\nconditions. The corresponding results are shown in Figure R2-1.\nFigure R2-1 (a) Schematic crystal structures of β-Ga O and amorphous Ga O after\n2 3 2 3\nhigh-temperature thermalization. (b) Radial distribution function of amorphous Ga O .\n2 3\nHeterogeneous interface of Ga O and 4H-SiC (c) before and (d) after the reaction. (e)\n2 3\nSi-O bond ratios varying with reaction time.\nFirst, to obtain the amorphous Ga O model, crystalline β-Ga O was annealed at\n2 3 2 3\na high temperature of 2600 K, as shown in Figure R2-1 (a). From Figure R2-1 (b), the\nradial distribution function has multiple peaks within 4 Å, and gradually approaches 1\nafter it is greater than 4 Å, indicating that the system exhibits a short-range ordered and\nlong-range disordered structure, corresponding to the amorphous state. The obtained\namorphous Ga O and 4H-SiC are used to construct a heterostructure as shown in\n2 3\nFigure R2-1 (c). After a thermal reaction at 2000 K for 1ns, the structure shown in\nFigure R2-1 (d) is obtained. It can be found from the interface that some of these Si\natoms will completely detach from SiC and form the structure of SiO tetrahedra inside\n2 3\nbe clearly found, proving that Ga O and 4H-SiC can form SiO₂ at the interface at high\n2 3\ntemperature. At the same time, the change curve of the number of Si-O bonds over time\nis also statistically analyzed. It can be observed that in the early stage, with the contact\nbetween amorphous Ga O and 4H-SiC, the number of Si-O bonds increases rapidly.\n2 3\nAfterwards, the bonding rate slows down, indicating that the reaction of Ga O and SiC\n2 3\nto form SiO₂ is not a slow process. At the same time, the number of Si-O bonds increases\nsteadily over time with no saturation trend, indicating that SiC will continue to be\noxidized.\nBased on your suggestion and the results mentioned above, we have revised Figure\n3 in the revised manuscript and added the corresponding discussion, as detailed in the\nChanges section below.\n(Change 1) Line 15-27, Page 11 and Line 1-4, Page 12 in the Manuscript\nTo gain a deeper understanding of the reaction mechanism and kinetic behavior at\nthe molecular level, an MD simulation was conducted at 2000 K. The results are\ndepicted in Figure 3(e-h), and the computational details are provided in Section 2 of the\nSupporting Information. First, the β-Ga O is annealed at a high temperature of 2600 K\n2 3\nto obtain a disordered amorphous Ga O structure, as shown in Figure 3 (e). From\n2 3\nFigure 3 (f), the radial distribution function has multiple peaks within 4 Å, and\ngradually approaches 1 after it is greater than 4 Å, indicating that the system exhibits a\nshort-range ordered and long-range disordered structure, corresponding to the\namorphous state. The obtained amorphous Ga O and 4H-SiC are used to construct a\n2 3\nheterostructure as shown in Figure 3 (g). After a thermal reaction at 2000 K for 1ns, it\ncan be found from the interface that some of these Si atoms will completely detach from\nSiC and form the structure of SiO tetrahedra inside Ga O after a reaction of 1ns. In\n2 2 3\nthe local magnified image, SiO₂ crystal structure can be clearly found, proving that\nGa O and 4H-SiC can form SiO₂ at the interface at high temperature. At the same time,\n2 3\nthe change curve of the number of Si-O bonds over time is also statistically analyzed.\nIt can be observed that in the early stage, with the contact between amorphous Ga O\n2 3\nand 4H-SiC, the number of Si-O bonds increases rapidly. Afterwards, the bonding rate\nslows down, indicating that the reaction of Ga O and SiC to form SiO₂ is not a slow\n2 3\nprocess. At the same time, the number of Si-O bonds increases steadily over time with\nno saturation trend, indicating that SiC will continue to be oxidized.\nComment 5: The manuscript should include a more detailed description of the\nTDTR experimental conditions, such as the thickness of the transducer layer and the\nmethods used to calculate the heat capacities of each layer. Additionally, the TDTR\nsignals for all samples should be plotted together in a single figure to qualitatively\nassess whether there are significant differences in the cooling curves.\nResponse: Thank you for your comment. Based on your two suggestions, we have\nmade two revisions regarding the TDTR section of the manuscript:\n(1) In the supporting information, we have provided a detailed description of the\nTDTR measurement method and presented the parameters used for each test in a list\nformat.\n(2) The TDTR signals for all samples have been plotted together in a single figure,\nwith separate plots showing the fitting for each signal.\nThe above changes are detailed in the Change section below.\n(Change 1) Section S4 in the Supporting Information (Page 7-8)\nThe time-domain thermoreflectance (TDTR) method is characterized by non-\ncontact measurement and high - time resolution, which can avoid the error caused by\ncontact thermal resistance and accurately capture the rapid thermal processes of\nmaterials. Its advantage lies in the ability to perform in-situ measurement on micro-\nsamples, thin films and complex-structured samples with minimal damage to the\nsamples, and it is widely used in the study of material thermal properties. The basic\nprinciple of the TDTR method is based on the thermoreflectance effect. When the\nsurface temperature of a material changes, its reflectivity will also change accordingly.\nIn the TDTR experiment, an ultrashort-pulse laser is first used as the pump light to\nirradiate the sample surface, causing a transient temperature rise on the surface. Then,\nanother probe laser irradiates the same area with a certain time delay to measure the\nchange in reflectance. By precisely controlling the time delay between the pump and\nprobe light, the temperature information of the sample at different time points can be\nobtained, and then the thermal response curve of the sample can be constructed. By\nfitting and analyzing the thermal response curve in combination with the geometric size\nof the sample and the material properties, the key thermal physical parameters of the\nsample, such as thermal conductivity and interfacial thermal conductance, can be\naccurately calculated. This work employs the Autinst Pioneer ONE TDTR system for\nmeasurement. The system utilizes a femtosecond laser with an output power of 2 W as\nthe light source, featuring a pulse width of 100 fs and a repetition rate of 80 MHz. The\npump laser is intensity - modulated by an electro-optic modulator (EOM) at a frequency\nof 1 MHz. The radii of the pump and probe spots are 12 and 4 micrometers, respectively.\nThe thickness of the transducer layer of TDTR samples is about 70 nm, which\ndetermined jointly by a profilometer and the TDTR acoustic echo signal. The thickness\nof Ga O film is measured by TEM. The volumetric heat capacity data of the materials\n2 3\ncontained in the tested samples are listed in Table S3.\n(Change 3) Table S3 in the Supporting Information (Page 8)\nTable S3 The parameters used for TDTR experiment.\nMaterial Thickness (nm) K (W/m·K) K (W/m·K) C (J/m3·K)\nz r p\nAl 85 237 237 2420000\n∽\nSiC / 360 360 2160000\nGa O To be tested 284000\n2 3\nBy TEM\n(Change 3) Figure S7 – S9 in the Supporting Information (Page 13-14)\nFigure S7. (a) TDTR amplitude signals of different samples annealed at 1173 K as a\nfunction of delay time. The best fits as a function of delay time: (a) 0 min, (b) 20 min,\n(c) 60 min, (d) 120 min, and (e) 240 min.\nFigure S8. (a) TDTR amplitude signals of different samples annealed at 1073 K as a\nfunction of delay time. The best fits as a function of delay time: (b) 60 min, (c) 120 min,\nand (d) 240 min.\nFigure S9. (a) TDTR amplitude signals of different samples annealed at 1273 K as a\nfunction of delay time. The best fits as a function of delay time: (b) 20 min, (c) 60 min,\nand (d) 120 min.\nComment 6: If the calculations indicating that the reaction between Ga O and\n2 3\nSiC proceeds spontaneously are valid, it suggests that SiO could form even at lower\n2\nannealing temperatures. Then, wouldn't it be possible to optimize the annealing\nconditions by adjusting both temperature and time?\nResponse: We thank the reviewer for their attention to this issue. Indeed, as you\nmentioned, SiO could form even at lower annealing temperatures in our system. In our\n2\nprevious experiments, we also studied the effect of different annealing temperatures\nand times on the formed SiO layer, and measured the thickness and ITC of the SiO\n2 2\nlayers formed under these conditions. Based on these experiments, we found that SiO\n2\nformed at 900°C for 2 hours exhibited the highest ITC. Therefore, we used this\noptimized condition as the basis for a more systematic analysis, which is presented in\nthis manuscript. In response to your interest in this issue, we have further refined the\nmorphology and composition analysis under other experimental conditions and present\nthese results here.\nFigure R2-2. (a) The thickness of SiO layer as a function of time at different\n2\ntemperature. (b) The ITC values established via the TDTR measurements at different\nthickness of SiO layer.\nmin), 900 °C (20 min, 60 min, 120 min, and 240 min), and 1000 °C (20 min, 60 min,\nand 120 min). Figures R2-2a and R2-2b show the thickness of the SiO layers and the\n2\nfinal ITC under these conditions, with corresponding morphological analyses shown in\nFigures R2-3 to R2-4. At 800°C, SiO can be formed after 60 minutes, but at a much\n2\nslower rate compared to the 900 °C conditions. For similar thickness, its ITC is also\nslightly lower than that of the 900°C annealing conditions. When the temperature is\nraised to 1000 °C, the rate of SiO formation increases, but the quality significantly\n2\ndecreases, resulting in a final ITC much lower than that of the 900°C annealing\nconditions. In summary, the 900°C annealing condition we selected strikes a balance\nbetween the rate and quality of SiO formation.\ntemperatures and times into the revised version. Figures R2-2 have been incorporated\nas Figures 4b and 4c in the revised manuscript, with the corresponding discussion added.\nFigures R2-3 and R2-4 have been included in the supporting information.\nFigure R2-3 The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 60 min, (c) 120 min, and (e) 240 min at 800 ℃.\nThe EDS line scanning profiles of cross-section under annealing (b) 60 min, (d) 120\nmin, and (f) 240 min at 800 ℃.\nFigure R2-4 The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 20 min, (c) 60 min, and (e) 120 min at 1000 ℃.\nThe EDS line scanning profiles of cross-section under annealing (b) 20 min, (d) 60 min,\nand (f) 120 min at 1000 ℃.\nDetailed changes can be found in the following Change section.\n(Change 1) Figure 4b – c in the Manuscript (Line 7 – 9, Page 12).\nFigure 4 (a) Schematic diagram of the sample structure and the TDTR experiment.\n(b) The thickness of SiO layer as a function of time at different temperature. (c) The\n2\n(d, g, f) Structural modeling for molecular dynamics calculations at the Ga₂O₃/SiC\ninterface bonded by (d) VdW and Ga₂O₃/SiO /SiC interface bonded by (g) VdW and (f)\n2\ncovalent interactions. (e, h, k) The accumulative energy of the heat source and heat sink\nof the Ga₂O₃/SiC interface bonded by (e) VdW and Ga₂O₃/SiO /SiC interface bonded\n2\nby (h) VdW and (k) covalent interactions. (f, i, l) Typical temperature profiles of the\nGa₂O₃/SiC interface bonded by (f) VdW and Ga₂O₃/SiO /SiC interface bonded by (i)\n(Change 2) The corresponding discussion in the Manuscript (Line 25 - 28,\nPage 7 and Line 1-2, page 8).\nThe data show that the thickness of the reactive layer and the interface roughness\ngradually increase with increasing annealing time and annealing temperature. When\nannealed at the lower temperature of 1073 K, the thickness of the SiO₂ thickness was\nslightly thinner, but the overall interface was flatter. When annealed at 1273 K, the\ngrowth rate of SiO₂ was significantly higher than at 1073 K and 1173 K. However,\ndeterioration of the interface and the Ga₂O₃ film, as well as an increase in roughness,\nwas also clearly observed.\n(Change 3) Line 15-19, Page 13 in the Manuscript.\nAt 1073 K, SiO can be formed after 60 minutes, but at a much slower rate\n2\ncompared to the 1173 K conditions. For similar thickness, its ITC is also slightly lower\nthan that of the 1173 K annealing conditions. When the temperature is raised to 1273\nK, the rate of SiO formation increases, but the quality significantly decreases, resulting\n2\nin a final ITC much lower than that of the 1173 K annealing conditions.\n(Change 4) Figure S2 and Figure S4 in the Supporting Information (Page 9\nand 11);\nFigure S2. The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 60 min, (c) 120 min, and (e) 240 min at 1073\nK. The EDS line scanning profiles of cross-section under annealing (b) 60 min, (d) 120\nmin, and (f) 240 min at 1073 K.\nFigure S4. The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 20 min, (c) 60 min, and (e) 240 min at 1273 K.\nThe EDS line scanning profiles of cross-section under annealing (b) 20 min, (d) 60 min,\nand (f) 240 min at 1273 K.\nComment 7: In the MD simulations shown in Fig. 4(d)-(i), the authors state that\nthey have simulated two types of structures—one with vdW bonding and one with\ncovalent bonding. The authors should clarify how these bonding types are differentiated\nin their calculations, including a description of the structural parameters used. Moreover,\nsince the authors consider the Ga O /SiC interface to be governed by vdW interactions,\n2 3\nshouldn't the comparison be between a Ga O /SiC structure and a Ga O /SiO /SiC\n2 3 2 3 2\nstructure?\nResponse: We thank the reviewer for their constructive suggestions, and we refine\nthe formulas and parameters for non-bonding interactions and covalent bonding\ninteractions in the supplementary material. In van der Waals connection, it actually\nmeans that there are only non-bonded interactions between interfaces; while in covalent\nbond connection, there are both covalent bonds and non-bonded interactions at the\ninterface.\nIn LAMMPS, pair potentials are defined between pairs of atoms that are within a\ncutoff distance. The standard 12/6 Lennard-Jones potential was adopted, and given by\nEquation (2), where is the cutoff radius.\n𝑟c\n(2)\n𝜎 12 𝜎 6\n𝐸 = 4𝜖[(𝑟) − (𝑟) ], 𝑟 < 𝑟𝑐\nTable R2-1. The σ and ϵ parameters for different atoms.\n12/6 LJ potential (angstrom) (kcal/mol)\n𝝈 𝝐\nGa 3.87 0.217\nO 2.859 0.228\nSi 4.053 0.040\nC 3.474 0.159\nThe harmonic bond style was adopted as the covalent bonding, given by Equation\n(3), where is the equilibrium bond distance, and is the bonding factor.\n𝑟0 𝐾\n(3)\n2\n𝐸 = 𝐾(𝑟−𝑟0)\nTable R2-2. The and parameters for different covalent bonding.\n𝑟0 𝐾\nCovalent bonding (angstrom) (kcal/mol)\n𝒓𝟎 𝑲\nGa-O 1.932 334.5\nSi-O 1.615 392.8\nSi-C 1.809 238.0\nGa-Si 2.106 207.3\nFigure R2-5 (a, d, g) Structural modeling for molecular dynamics calculations at the\nGa₂O₃/SiC interface bonded by (a) VdW and Ga₂O₃/SiO /SiC interface bonded by (d)\nand heat sink of the Ga₂O₃/SiC interface bonded by (b) VdW and Ga₂O₃/SiO /SiC\n2\ninterface bonded by (e) VdW and (h) covalent interactions. (c, f, i) Typical temperature\nprofiles of the Ga₂O₃/SiC interface bonded by (c) VdW and Ga₂O₃/SiO /SiC interface\n2\nbonded by (f) VdW and (i) covalent interactions.\nMoreover, based on your further suggestions, we have added a simulation analysis\nof a Ga O /SiC bilayer structure by vdW interactions, as the results shown in Figure\n2 3\nR2-5. When a vdW contact is formed between Ga O and SiC, the calculated ITC is\n2 3\n13.6 MW/m²·K. If a vdW contact SiO₂ interlayer is directly introduced between the two,\nit leads to a decrease in ITC, primarily due to the low thermal conductivity of SiO₂ itself\nand the creation of a poorly conductive interface. In contrast, only when SiO₂ forms\ncovalent contacts with both Ga O and SiC can the thermal performance between them\n2 3\nbe effectively enhanced.\nBased on your suggestion, the above simulation details have been included as\nSupplementary Texts (Section S3) in the revised supporting information. (Line 15 – 29,\nPage 5 and Line 1 – 15, Page 6). At the same time, we have revised Figure 4 in the text\nbased on Figure R2-5 and added the corresponding discussion, as detailed in the\nChange section below.\n(Change 1) Section S3 in the Supporting Information (Line 15 – 29, Page 5\nand Line 1 – 15, Page 6).\nSection S2: The details on the molecular dynamics simulation of interfacial heat\ntransfer.\n(Change 2) Figure 4 in the Manuscript (Page 12).\nFigure 4 (a) Schematic diagram of the sample structure and the TDTR experiment.\n(b) The thickness of SiO layer as a function of time at different temperature. (c) The\n2\n(d, g, f) Structural modeling for molecular dynamics calculations at the Ga₂O₃/SiC\ninterface bonded by (d) VdW and Ga₂O₃/SiO /SiC interface bonded by (g) VdW and (f)\n2\ncovalent interactions. (e, h, k) The accumulative energy of the heat source and heat sink\nof the Ga₂O₃/SiC interface bonded by (e) VdW and Ga₂O₃/SiO /SiC interface bonded\n2\nby (h) VdW and (k) covalent interactions. (f, i, l) Typical temperature profiles of the\nGa₂O₃/SiC interface bonded by (f) VdW and Ga₂O₃/SiO /SiC interface bonded by (i)\n(Change 3) Line 8 - 11, Page 15 in the Manuscript.\nThe corresponding temperature profiles are illustrated in Figures 3(f, i, and l).\nCompared to the direct Van der Waals contact between Ga₂O₃ and SiC, the insertion of\na SiO₂ layer introduces an additional van der Waals-bonded interface, resulting in a\nhigher ITR due to the presence of this extra interface and the low thermal conductivity\nof SiO₂.\nComment 8: In the Conclusion, immediately after discussing covalent bonding,\nthe TCAD simulation is mentioned. However, it should be clarified that this simulation\nis not related to the context of covalent bonding and is merely a simulation where the\nITC is treated as a variable.\nResponse: Thank you very much for your comments, and we greatly admire your\nin-depth insight into this field. Indeed, as you mentioned, in our TCAD simulation,\nwe treated ITC as an adjustable variable. This is mainly because ITC testing is at the\nmaterial level, while TCAD simulation is at the device level. We cannot distinguish\nbetween covalent and van der Waals bonding in the TCAD simulation model, so we\nused the ITC obtained from experimental testing as a variable in this indirect manner.\nBased on your suggestion, we have revised the description in the conclusion for greater\nprecision, as detailed in the Change section below.\n(Change 1) Line 17-21, Page 19 in the Manuscrtpt.\nIn summary, the heat dissipation of a Ga₂O₃ film on the SiC substrate was\nsuccessfully enhanced through a simple route consisting in the formation of covalent\nbonds, which is compatible with existing production lines. Furthermore, we conducted\nTCAD simulations by treating the ITC derived from the experimental tests on both\ncovalent and vdW contacts as a variable. The results showed that the increase of ITC\nfrom 17.5 to 162.1 MW/m²·K significantly reduced the impact of the self-heating effect,\nleading to an 83 mA/mm improvement in the maximum output current of the MOSFET\nalong with a notable drop in the local junction temperature.\nIn specific response to the points raised by reviewer 3:\nGeneral comment: Shen and coworkers present a study of the thermal boundary\nresistance in Ga O /SiC interfaces, proposing this system for its efficient thermal\n2 3\nmanagement characteristic for Ga O electronics. In a nutshell, the authors show that a\n2 3\nproper post-growth annealing cycle improve the quality of the Ga O /SiC interface,\n2 3\nwhich is not by itself a very surprising result. Nevertheless, the meticulous analysis\nperformed of how the progressive formation of a quasi-epitaxial interface (i.e., the\nformation of covalent bonds) determines the thermal boundary resistance is interesting\nand suggest that the manuscript perhaps may be considered for publication in Nature\nCommunications.\nResponse: We sincerely appreciate your professional comments and valuable\nsuggestions, which have greatly helped us improve our manuscript. In response to your\nadvice, we have carefully considered each of your comments and made the\ncorresponding revisions. Below, we provide detailed point-by-point responses for your\nreference.\nComment 1: The authors show that the quality of the interface improves and the\nthermal boundary resistance decreases with a suitable annealing cycle. This thermal\ntreatment promotes the formation of covalent bonds and the creation of a thin SiO\n2\nbuffer layer, where O from the Ga O and Si from the SiC react.\n2 3\nHowever, SiO thin\n2\nfilms, even in crystalline form, should have a rather low thermal conductivity (see e.g.\nhttps://www.nature.com/articles/s41598-018-28925-6). Therefore, there is a trade off in\nthe proposed strategy that is neither discussed, nor mentioned or considered. Promoting\nthe formation of a SiO buffer layer is beneficial for thermal management (it provides\n2\na dissipation channel trough SiC with a low thermal boundary resistance), but if this\nlayer becomes too thick the whole system becomes too resistive (the contribution of the\nintrinsic thermal resistance of SiO to the overall thermal resistance cannot be\n2\nneglected). What is the optimal thickness and how can it be controlled? The tradeoff\nbetween thin and thick SiO (low thermal boundary resistance, but high SiO film\n2 2\nresistance) layer must be addressed.\nResponse: Thank you to the reviewers for their constructive suggestions and we\nwill add a discussion of this in the text. As you said regardless of the form of SiO ,\n2\nit has a low thermal conductivity and too thick a layer of SiO will increase the interfacial\n2\nthermal resistance and reduce the thermal efficiency 15. In this respect the thinner the\nthickness of the SiO layer the better the heat dissipation based on the formation of\n2\ncovalent bonds. On the other hand, it has been shown that the phonon free path of SiO\n2\nis around 6 nm, and when the film thickness is lower than 5 nm, the thermal\nconductivity of the SiO layer will be significantly reduced 16. In addition, the quality\n2\nof the interface needs to be taken into account, as a rough and deteriorated interface\nwill increase the phonon heat dissipation at the interface and also significantly reduce\nthe interfacial thermal conductivity. Therefore, as stated by the reviewer there is a trade-\noff in the proposed strategy. In this work, we have tuned the thickness of the silicon\noxide layer and the interface by varying the temperature and time of annealing. The\noptimal interfacial thermal conductivity was obtained by combining factors such as\ninterfacial roughness and thickness, which reached 162.1 MW/m²·K. Figures R3-1a and\nR3-1b show the thickness of the SiO layers and the final ITC under these conditions.\nresults of TEM and TDTR at 800°C and 1000°C in the Supporting Information. Figures\nR3-1a and R3-1b have been incorporated as Figures 4b and 4c in the revised manuscript,\nwith the corresponding discussion added.\nFigure R3-1. (a) The thickness of SiO layer as a function of time at different\n2\ntemperature. (b) The ITC values established via the TDTR measurements at different\nthickness of SiO layer.\nthe Change 1 section below.\n(Change 1) The corresponding discussion in the Manuscript (Line 25 - 28,\nPage 7 and Line 1 - 2, Page 8).\nThe data show that the thickness of the reactive layer and the interface roughness\ngradually increase with increasing annealing time and annealing temperature. When\nannealed at the lower temperature of 1073 K, the thickness of the SiO₂ thickness was\nslightly thinner, but the overall interface was flatter. When annealed at 1273 K, the\ngrowth rate of SiO₂ was significantly higher than at 1073 K and 1173 K. However,\ndeterioration of the interface and the Ga₂O₃ film, as well as an increase in roughness,\nwas also clearly observed.\n(Change 2) Line 15-19, Page 13 in the Manuscript.\nAt 1073 K, SiO can be formed after 60 minutes, but at a much slower rate\n2\ncompared to the 1173 K conditions. For similar thickness, its ITC is also slightly lower\nthan that of the 1173 K annealing conditions. When the temperature is raised to 1273\nK, the rate of SiO formation increases, but the quality significantly decreases, resulting\n2\nin a final ITC much lower than that of the 1173 K annealing conditions.\n(Change 3) Line 21 - 27, Page 13 and Line 1 - 2, Page 14 in the Manuscript.\nSiO itself has a low thermal conductivity, and an excessively thick layer of SiO\n2 2\nsignificantly increases the interfacial thermal resistance. However, since the phonon\nmean free path of SiO is about 6 nm, the thermal conductivity decreases further when\n2\nits thickness is less than 5 nm. In addition, the interfacial mass also has an important\neffect on the thermal conductivity. Therefore, under the premise of maintaining\ncovalent bonding, it is experimentally measured that when the thickness of the SiO\n2\nlayer is 5.8 nm, an optimal balance between thickness and interfacial mass can be\nachieved, resulting in the highest interfacial thermal conductivity.\n(Change 4) Figure S2 and S4 in the Supporting Information (Page 9 and 11)\nFigure S2. The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 60 min, (c) 120 min, and (e) 240 min at 1073\nK. The EDS line scanning profiles of cross-section under annealing (b) 60 min, (d) 120\nmin, and (f) 240 min at 1073 K.\nFigure S4. The cross-sectional TEM image and the corresponding EDS mapping\nimages of the annealed sample with (a) 20 min, (c) 60 min, and (e) 240 min at 1273 K.\nThe EDS line scanning profiles of cross-section under annealing (b) 20 min, (d) 60 min,\nand (f) 240 min at 1273 K.\n(Change 4) Figure S8 and S9 in the Supporting Information (Page 13-14)\nFigure S8. (a) TDTR amplitude signals of different samples annealed at 1073 K as a\nfunction of delay time. The best fits as a function of delay time: (b) 60 min, (c) 120 min,\nand (d) 240 min.\nFigure S9. (a) TDTR amplitude signals of different samples annealed at 1273 K as a\nfunction of delay time. The best fits as a function of delay time: (b) 20 min, (c) 60 min,\nand (d) 120 min.\nComment 2: The authors state that: “Meanwhile, the heat treatment eliminated\ndefects within the Ga O film, causing thermal conductivity to achieve to 6.2 W/m·K”.\n2 3\nHowever, it is not specified what was this value before the treatment. Is 6.2 W/mK large?\nIs it small?\nResponse: We are grateful for the reviewer’s careful comment on this issue. The\nthermal conductivity of the Ga O film before treatment was 5.7 W/m·K, which was\n2 3\nslightly lower than that of the treated sample. Based on your comments, we have revised\nthe relevant descriptions in the manuscript. Please refer to the Change section below\nfor details.\n(Change 1) Line 1 – 3, Page 14 in the Manuscrtpt.\nMeanwhile, the thermal conductivity of the Ga₂O₃ film before treatment was 5.7\nW/m·K, which is consistent with reports in other literature 17 18, 19. The heat treatment\neliminated defects within the Ga₂O₃ film, causing thermal conductivity to achieve to\n6.2 W/m·K.\nComment 3: The AIMD results do not seem to be very significative. At least\njudging from the snapshots of Figure 3(d) and 3(g), the AIMD simply shows that a thin\nSiO layer sandwiched between SiC and Ga O forms covalent bonds, starting from the\n2 2 3\ninitial vdW configuration. This appears to be a rather trivial and largely expected result.\nThe challenge would have been demonstrating the formation of the SiO2 layer between\nSiC and Ga O . While I do understand that this probably too challenging, I don0t think\n2 3\nthat “the simulation data confirmed that Ga O and SiC can form a Si-260 O covalent\n2 3\nbond at the interface under high temperatures”\nResponse: Thank you for your comment, which has been extremely helpful in\nimproving the quality of our work. Based on this, we re-conducted the molecular\ndynamics simulation following your approach. Additionally, in response to suggestions\nfrom other reviewers and our own microscopic morphology analysis, we replaced the\ncrystalline Ga O in the original model with amorphous Ga O . On this basis, we\n2 3 2 3\nextended the reaction time to 1 ns to observe the formation process of SiO .\n2 3 2 3\nhigh-temperature thermalization. (b) Radial distribution function of amorphous Ga O .\n2 3\nHeterogeneous interface of Ga O and 4H-SiC (c) before and (d) after the reaction. (e)\n2 3\nSi-O bond ratios varying with reaction time.\nThe simulation results are shown in Figure R3-2. First, to obtain the amorphous\nGa O model, crystalline β-Ga O was annealed at a high temperature of 2600 K,\n2 3 2 3\nas shown in Figure R3-2 (a). From Fig. R3-2 (b), it can be seen that the radial distribution\nfunction has multiple peaks in the range of 4 Å, and gradually converges to 1 after it is\nlarger than 4 Å, indicating that the system exhibits proximally ordered and remotely\ndisordered structures, which corresponds to an amorphous state. The heterostructure\nwas constructed with the obtained amorphous Ga O and 4H-SiC as shown in Fig. R3-\n2 3\n2 (c). After thermal reaction at 2000 K for 1 ns, the structure shown in Fig. R3-2 (d)\nwas obtained. From the interface, it can be found that after the reaction of 1 ns, some\nSi atoms will be completely separated from SiC, and the SiO tetrahedral structure will\n2\nbe formed inside Ga O .\n2 3\nIn the local magnified image, the SiO₂ crystal structure can be clearly found, which\nproves that Ga O and 4H-SiC can form SiO at the interface at high temperature.\n2 3 2\nMeanwhile, the change curve of Si-O bond number with time was statistically analyzed.\nIt can be observed that at the initial stage, the number of Si-O bonds increases rapidly\nwith the contact between amorphous Ga O and 4H-SiC. After that, the bonding rate\n2 3\nbecomes slower, indicating that the reaction between Ga O and SiC to form SiO is\n2 3 2\nnot a slow process. Meanwhile, the number of Si-O bonds increases steadily with time\nwithout a saturation trend, indicating that SiC will continue to be oxidized for further\nformation of more SiO .\n3 in the revised manuscript and added the corresponding discussion, as detailed in the\nChanges section below.\n(Change 1) Line 15 – 27, Page 11 and Line 1 – 4, Page 12 in the Manuscript\nTo gain a deeper understanding of the reaction mechanism and kinetic behavior at the\nmolecular level, an MD simulation was conducted at 2000 K, and the results are\ndepicted in Figure 3(e-h). First, the β-Ga O is annealed at a high temperature of 2600\n2 3\nK to obtain a disordered amorphous Ga O structure, as shown in Figure 3 (e). From\n2 3\nFigure 3 (f), the radial distribution function has multiple peaks within 4 Å, and\ngradually approaches 1 after it is greater than 4 Å, indicating that the system exhibits a\nshort-range ordered and long-range disordered structure, corresponding to the\namorphous state. The obtained amorphous Ga O and 4H-SiC are used to construct a\n2 3\nheterostructure as shown in Figure 3 (g). After a thermal reaction at 2000 K for 1ns, it\ncan be found from the interface that some of these Si atoms will completely detach from\nSiC and form the structure of SiO tetrahedra inside Ga O after a reaction of 1ns. In\n2 2 3\nthe local magnified image, SiO₂ crystal structure can be clearly found, proving that\nGa O and 4H-SiC can form SiO₂ at the interface at high temperature. At the same time,\n2 3\nthe change curve of the number of Si-O bonds over time is also statistically analyzed.\nIt can be observed that in the early stage, with the contact between amorphous Ga O\n2 3\nand 4H-SiC, the number of Si-O bonds increases rapidly. Afterwards, the bonding rate\nslows down, indicating that the reaction of Ga O and SiC to form SiO₂ is not a slow\n2 3\nprocess. At the same time, the number of Si-O bonds increases steadily over time with\nno saturation trend, indicating that SiC will continue to be oxidized.\nComment 4. The authors say that their results “confirm that the Ga O layer was\n2 3\namorphous”, implying that this was expected. Should not it be already crystalline, after\nthe ALD?\nResponse: Thank you for your insightful comment. We agree that the original\nphrasing may have caused some misunderstanding. In fact, the amorphous nature of the\nGa₂O₃ layer was not assumed beforehand, but rather determined based on the\nexperimental results from HRTEM and SAED analyses. Atomic layer deposition (ALD)\nof Ga₂O₃ can result in either amorphous or crystalline films depending on specific\ngrowth conditions, such as temperature, precursor chemistry, and post-deposition\ntreatment. In our case, under the deposition conditions used, the as-deposited Ga₂O₃\nwas found to be amorphous, as evidenced by the absence of lattice fringes in HRTEM\nand the diffuse halo rings in the SAED pattern. Based on your comments, we have\nrevised the sentence in the manuscript.\n(Change 1) Line 10 – 12, Page 7 in the Manuscript\nThe high-resolution TEM (HRTEM) image and selected area electron diffraction\n(SAED) results, as shown in Figure 1(d), indicate that the as-deposited Ga₂O₃ layer is\namorphous, while the underlying SiC substrate exhibits high-quality single-crystal\ncharacteristics.\nComment 5: Beware of “Colding” in place of “cooling” in Figure 5(b).\nResponse: We thank the reviewer for carefully reading our manuscript and pointing\nout our typographical error. Based on your comments, we have replaced “Colding” with\n“Cooling” in this Figure.\n(Change 1) Figure 5 in the Manuscript (Page 17).\nWe appreciate your earnest work, and hope that the correction will meet with\napproval. The manuscript has been overall checked, and the changes have been marked\nin red. We hope that this revision is sufficient to make our manuscript acceptable for\npublication in Nature Communications. If you think any additional clarification needs\nto be addressed, I will be happy to include them. Thank you very much for your\ncomments and suggestions.\nYours sincerely,\nHong-Ping Ma\nFudan University\nE-mail: hpma@fudan.edu.cn\nReferences\n1. Higashiwaki M, Sasaki K, Kuramata A, Masui T, Yamakoshi S. Development of gallium oxide\npower devices. physica status solidi (a) 211, 21-26 (2014).\n2. Qin Y, Wang Z, Sasaki K, Ye J, Zhang Y. Recent progress of Ga O power technology: large-\n2 3\narea devices, packaging and applications. Jpn J Appl Phys 62, Sf0801 (2023).\n3. Xu J, Zheng W, Huang F. Gallium oxide solar-blind ultraviolet photodetectors: a review. J Mater\nChem C 7, 8753-8770 (2019).\n4. Pearton SJ, et al. A review of Ga O materials, processing, and devices. Appl Phys Rev 5, 011301\n2 3\n(2018).\n5. Higashiwaki M, Wong MH. Beta-Gallium Oxide Material and Device Technologies. Annual\nReview of Materials Research 54, 175-198 (2024).\n6. Thompson AP, et al. LAMMPS - a flexible simulation tool for particle-based materials\nmodeling at the atomic, meso, and continuum scales. Computer Physics Communications 271,\n108171 (2022).\n7. Aktulga HM, Fogarty JC, Pandit SA, Grama AY. Parallel reactive molecular dynamics:\nNumerical methods and algorithmic techniques. parallel computing 38, 245-259 (2012).\n8. Chenoweth K, van Duin ACT, Goddard WA. ReaxFF Reactive Force Field for Molecular\nDynamics Simulations of Hydrocarbon Oxidation. The Journal of Physical Chemistry A 112,\n1040-1053 (2008).\n9. Rappe AK, Goddard WA, III. Charge equilibration for molecular dynamics simulations. The\nJournal of Physical Chemistry 95, 3358-3363 (1991).\n10. Nakano A. Parallel multilevel preconditioned conjugate-gradient approach to variable-charge\nmolecular dynamics. Computer Physics Communications 104, 59-69 (1997).\n11. Niefind F, et al. Watching (De)Intercalation of 2D Metals in Epitaxial Graphene: Insight into\nthe Role of Defects. Small 20, 2306554 (2024).\n12. Li Y, Sun F, Feng Y. Thermal boundary conductance in heterogeneous integration between β-\nGa O and semiconductors. Ceram Int 50, 18787-18796 (2024).\n2 3\n13. Majumdar S, Malen JA, McGaughey AJH. Cooperative Molecular Behavior Enhances the\nThermal Conductance of Binary Self-Assembled Monolayer Junctions. Nano Lett 17, 220-227\n(2016).\n14. Losego MD, Grady ME, Sottos NR, Cahill DG, Braun PV. Effects of chemical bonding on heat\ntransport across interfaces. Nature Materials 11, 502-506 (2012).\n15. Zhu W, Zheng G, Cao S, He H. Thermal conductivity of amorphous SiO(2) thin film: A\nmolecular dynamics study. Sci Rep 8, 10537 (2018).\n16. Martin E, Ori G, Duong T-Q, Boero M, Massobrio C. Thermal conductivity of amorphous SiO2\nby first-principles molecular dynamics. J Non·Cryst Solids 581, (2022).\n17. Cheng Z, et al. Thermal Transport across Ion-Cut Monocrystalline β-Ga O Thin Films and\n2 3\nBonded β-Ga O -SiC Interfaces. ACS Appl Mater Interfaces 12, 44943-44951 (2020).\n2 3\n18. Song Y, et al. Characterization of the Thermal Boundary Resistance of a Ga O /4H-SiC\n2 3\nComposite Wafer. In: 2020 19th IEEE Intersociety Conference on Thermal and\nThermomechanical Phenomena in Electronic Systems (ITherm)). IEEE (2020).\n19. Xu W, et al. Thermal Transport Properties of β-Ga O Thin Films on Si and SiC Substrates\n2 3\nFabricated by an Ion-Cutting Process. ACS Appl Electron Mater 6, 1710-1717 (2024)."}, {"round": 2, "reviewer_comments": {"reviewer_1": "From my point of view, the authors addressed a lot of my comments, and the quality of the paper improved. Nevertheless, the paper still requires some revision:\n- The methodology section should contain all methods that the authors used to achieve the objectives and to formulate the conclusions. Currently, it has only material fabrications and structural characterization subsections.\n- It would be nice to correlate the results of ab-initio and MD calculations with the ReaxFF potential to validate the potential for the considered case. Also, it would be nice to correlate ReaxFF and non-reactive MD results. This can be done, for example, by calculating the vibration frequencies or binding energies.\n- The TDTR requires some additional details: acquisition time, properties of the top layer (material, thickness, thermal properties) in the text. The quality of the paper will also be improved if the model for the treatment of TDTR is added.\n- It is not clear what the distance is in Fig. 5d, the corresponding sketch will be helpful.\n- It would be nice to reproduce Fig. 5f with the FEM method and the ITC evaluated with MD for cross-validation of the experiment and multiscaling simulations.", "reviewer_2": "I think that the authors have responded to review comments and revised in general satisfactorily, and the manuscript looks like it can be published.", "reviewer_3": "I thank the authors for the careful consideration given to my criticism.\nIn its revised version, I can recommend the paper for publication in Nature Communication."}, "author_response": "In specific response to the points raised by reviewer 1:\nGeneral comment: From my point of view, the authors addressed a lot of my\ncomments, and the quality of the paper improved. Nevertheless, the paper still requires\nsome revision:\nResponse: We sincerely thank you for carefully reviewing our revised manuscript and\nresponses, and for providing valuable comments and suggestions for further\nimprovement. We have carefully considered each of your points and made\ncorresponding revisions to the manuscript. A detailed point-by-point response is\nprovided below for your kind review, together with the corresponding revisions in the\nmanuscript, which are highlighted in red in the revised text.\nComment 1: The methodology section should contain all methods that the authors used\nto achieve the objectives and to formulate the conclusions. Currently, it has only\nmaterial fabrications and structural characterization subsections.\nResponse: We appreciate the reviewer’ suggestion, which have helped us improve the\nreadability and completeness of the paper. As your recommendation, we have\nincorporated a description of the simulation calculations into the methodology section.\nThe above changes are detailed in the Change section below.\n(Change 1) The corresponding description in the Manuscript (Line 15 - 23, Page\n21).\nSimulation and calculation: Density functional theory (DFT) calculations were\ncarried out to evaluate the reaction energies between Ga₂O₃ and SiC (Supporting\nInformation, Section S1). Molecular dynamics (MD) simulations were performed to\nexamine the reaction process between Ga₂O₃ and SiC that leads to the formation of SiO₂\n(Supporting Information, Section S2), as well as to analyze the heat transfer behavior\nin Ga₂O₃/SiC and Ga₂O₃/SiO₂/SiC heterostructures (Supporting Information, Section\nS3). The formation energy of SiO and the binding energies with different systems were\n2\ncalculated (Supporting Information, Section S4). In addition, TCAD simulations were\nused to model two types of Ga₂O₃-based MOSFETs, namely Ga₂O₃ and Ga₂O₃/SiC\ndevices, in order to investigate their electrical characteristics (Supporting Information,\nSection S6).\nComment 2: It would be nice to correlate the results of ab-initio and MD calculations\nwith the ReaxFF potential to validate the potential for the considered case. Also, it\nwould be nice to correlate ReaxFF and non-reactive MD results. This can be done,\nfor example, by calculating the vibration frequencies or binding energies.\nResponse: We sincerely appreciate the reviewer’s constructive suggestion.\nFollowing\nyour advice, we have performed additional calculations to further validate the accuracy\nof the potential used in our simulations. First, we compared the SiO₂ formation energies\nobtained from ab initio (DFT) calculations and molecular dynamics (MD) simulations\nwith the ReaxFF potential. As shown in Figure R1 (a), the two results are in excellent\nagreement, with only a ~1.3% difference, confirming that the ReaxFF potential can\nreliably capture the energy variation during the SiO₂ formation process.\nIn addition, to further verify the reliability of the employed potential, we calculated\nthe binding energies of the O–SiC, O–Ga₂O₃, C–SiC, and Si–Ga₂O₃ systems using DFT,\nReaxFF-based MD, and non-reactive MD, as shown in Figure R1 (b). All three methods\nexhibit consistent trends. Specifically, the DFT results show that the binding energies\nof O–SiC, O–Ga₂O₃, and C–SiC systems are all negative, indicating that these bonds\nare thermodynamically stable, with their magnitudes decreasing in sequence. By\ncontrast, the Si–Ga₂O₃ binding energy is positive, suggesting that Si atoms released\nfrom SiC do not tend to bond directly with Ga₂O₃ but are more likely to undergo\nalternative reactions such as SiO₂ formation. This observation is consistent with the\nreaction pathway identified in our MD simulations of Ga₂O₃ and SiC.\nSince the ReaxFF results reproduce the exact same trend as the DFT calculations,\nwe conclude that the ReaxFF potential can reliably describe the interfacial reaction\nbetween Ga₂O₃ and SiC leading to SiO₂ formation. Moreover, the non-reactive MD\nresults also follow the same trend, further supporting the validity of both our simulation\nstructures and methodologies.\nFigure R1. (a) The formation energy of SiO by calculations of DFT and ReaxFF.\n2\n(b)\nThe binding energy of the O-SiC, O-Ga\n3\n, Si-SiC, and Si-Ga\n3\nsystems.\nThe added section is highlighted in red in the manuscript and can also be found in\nthe Change 1 section below.\n(Change 1) Figure S10 in the Supporting Information (Page 17).\nFigure S10. (a) The formation energy of SiO by calculations of DFT and ReaxFF.\n2\n(b) The binding energy of the O-SiC, O-Ga\n3\n, Si-SiC, and Si-Ga\n3\nsystems.\n(Change 2) Section 4 in the Supporting Information (Page 8-9).\nSection S4. The details on of the formation energy of SiO and the binding energy.\nplane-wave basis, which is implemented in Vienna ab-initio simulation package1.\nAnd the plane-waves are cut-off at 400 eV. The exchange-correlations of electrons\nare described by the generalized gradient approximations with the form proposed by\nPerdew, Burke, and Ernzerhof2. DFT-D3 approach is utilized to capture the long-\nrange van der Waals interactions3. The Brillouin zone is sampled with resolutions\nbetter than 0.05 Å-1, using the scheme of Monkhorst-Pack4. A vacuum of 15 Å is\nplaced at thickness direction to eliminate the interactions between periodic images.\nThe energy converge criterion for solving self-consistent Kohn-Sham equations is\n10-5 eV. All the structures all fully optimized until the residual Hellman-Feynman\nForces and stresses below 0.05 eV/Å. The forcefield calculations are done with\nLAMMPS package.5 The reactive force fields (ReaxFF) are utilized to describe the\ninteraction potentials between atoms:6, 7\n𝐸𝑠𝑦𝑠𝑡𝑒𝑚 = 𝐸𝑏𝑜𝑛𝑑 +𝐸𝑜𝑣𝑒𝑟 +𝐸𝑢𝑛𝑑𝑒𝑟 +𝐸𝑙𝑝 +𝐸𝑣𝑎𝑙 +𝐸𝑡𝑜𝑟 +𝐸𝑣𝑑𝑊 +𝐸𝐶𝑜𝑢𝑙𝑜𝑚𝑏 +𝐸𝑞𝑒𝑞\nin which denotes the bonding energy, and are energy penalty\n𝐸𝑏𝑜𝑛𝑑 𝐸𝑜𝑣𝑒𝑟 𝐸𝑢𝑛𝑑𝑒𝑟\nto prevent over- and undercoordination, is the energy of long-pairs of atoms,\n𝐸𝑙𝑝\nand are three-body energy from valence angle and four-body torsions,\n𝐸𝑣𝑎𝑙 𝐸𝑡𝑜𝑟\nis van der Waals dispersion energy, is electrostatic energy and\n𝐸𝑣𝑑𝑊 is the charge equilibrium energy.6, 8, 9 The 𝐸 f 𝐶 o 𝑜 r 𝑢 c 𝑙 e 𝑜 𝑚 fi 𝑏 eld parameters are in line with\n𝐸𝑞𝑒𝑞\nformer study, which is suitable to describe the interations between SiC, Ga and\noxygen.10 The binding energy is calculated via\n𝐸𝑏 = 𝐸[𝐴𝐵]−𝐸[𝐴]−𝐸[𝐵]\nwhere is the energy of linked system, and is the energy of\n𝐸[𝐴𝐵] 𝐸[𝐴/𝐵]\nindividual system. The formation energy of SiO2 is defined by:\n𝐸𝑓𝑜𝑟𝑚 = 𝐸[𝑆𝑖𝑂2]−𝐸[𝑆𝑖]−𝐸[𝑂2(𝑔)]\nwhere is the energy of α-SiO , is the energy of Si, is the\n2\n𝐸[𝑆𝑖𝑂2] 𝐸[𝑆𝑖] 𝐸[𝑂2(𝑔)]\nenergy of oxygen molecule.\n(Change 3) The corresponding description in the Manuscript (Line 27 - 29,\nPage 11 and 1-7, Page 12).\nTo validate the accuracy of ReaxFF in describing the energy variations during MD\nsimulations, the formation energy of SiO₂ obtained from DFT was compared with that\nfrom ReaxFF. As shown in Figure S10(a), the two values are in close agreement, with\nonly a ~1.3% difference, confirming that the ReaxFF potential can reliably capture the\nenergy change during the SiO₂ formation process. In addition, to further verify the\nreliability of the employed potential, we calculated the binding energies of the O–SiC,\nO–Ga₂O₃, C–SiC, and Si–Ga₂O₃ systems using DFT, and ReaxFF-based MD, as shown\nin Figure S10 (b). All two method exhibit consistent trends. The Si–Ga₂O₃ binding\nenergy is positive, suggesting that Si atoms released from SiC do not tend to bond\ndirectly with Ga₂O₃ but are more likely to undergo alternative reactions such as SiO₂\nformation. This observation is consistent with the reaction pathway identified in our\nMD simulations of Ga₂O₃ and SiC.\n(Change 4) The corresponding description in the Manuscript (Line 5 - 8, Page\n15).\nAdditionally, to validate the accuracy of the constructed model, binding energies\nwere calculated in non-reactive molecular dynamics simulations. The results are shown\nin Figure S10(b). It can be observed that these results are consistent with those obtained\nfrom ab initio and MD calculations with the ReaxFF potential.\nComment 3: The TDTR requires some additional details: acquisition time, properties\nof the top layer (material, thickness, thermal properties) in the text. The quality of the\npaper will also be improved if the model for the treatment of TDTR is added.\nResponse: Thank you for your comment. Based on your suggestions, we have made\nrevisions regarding the TDTR section of the manuscript. In the supporting information,\nwe have provided a detailed description of the TDTR measurement method and\npresented the parameters used for each test in a list format. Our study has been\nperformed with reference to the well-established TDTR modeling framework 11, which\nprovides a solid theoretical basis for the thermal conductivity and interfacial thermal\nconductance results reported in the manuscript. The model for the treatment of TDTR\nis also added in the supporting information.\n(Change 1) Section S5 in the Supporting Information (Page 9-10)\nThe time-domain thermoreflectance (TDTR) method is characterized by non-\ncontact measurement and high - time resolution, which can avoid the error caused by\ncontact thermal resistance and accurately capture the rapid thermal processes of\nmaterials. Its advantage lies in the ability to perform in-situ measurement on micro-\nsamples, thin films and complex-structured samples with minimal damage to the\nsamples, and it is widely used in the study of material thermal properties. The basic\nprinciple of the TDTR method is based on the thermoreflectance effect. When the\nsurface temperature of a material changes, its reflectivity will also change accordingly.\nIn the TDTR experiment, an ultrashort-pulse laser is first used as the pump light to\nirradiate the sample surface, causing a transient temperature rise on the surface. Then,\nanother probe laser irradiates the same area with a certain time delay to measure the\nchange in reflectance. By precisely controlling the time delay between the pump and\nprobe light, the temperature information of the sample at different time points can be\nobtained, and then the thermal response curve of the sample can be constructed. To\nquantitatively extract thermal properties, a three-dimensional heat diffusion model for\nmultilayer systems is solved, and the Green’s function is constructed using the\nquadrupole method 11. On this basis, a complete signal modeling framework is\nestablished. In this framework, the pump beam is treated as a train of modulated pulses,\nand the surface temperature response is obtained through convolution with the Green’s\nfunction. By further incorporating probe beam delay modeling and applying lock-in\ndetection, the system response signal at the modulation frequency is derived. This\nsystem response is then fitted to extract the thermal conductivity and interfacial thermal\nconductance of the sample with high accuracy.\nThis work employs the Autinst Pioneer ONE TDTR system for measurement. The\nsystem utilizes a femtosecond laser with an output power of 2 W as the light source,\nfeaturing a pulse width of 100 fs and a repetition rate of 80 MHz. The pump laser is\nintensity - modulated by an electro-optic modulator (EOM) at a frequency of 1 MHz.\nThe radii of the pump and probe spots are 12 and 4 micrometers, respectively. The\nthickness of the transducer layer of TDTR samples is about 70 nm, which determined\njointly by a profilometer and the TDTR acoustic echo signal. The thickness of Ga O\n2 3\nfilm is measured by TEM. The volumetric heat capacity data of the materials contained\nin the tested samples are listed in Table S3.\n(Change 2) Table S3 in the Supporting Information (Page 10)\nTable S3 The parameters used for TDTR experiment.\nMaterial Thickness (nm) K (W/m·K) K (W/m·K) C (J/m3·K)\nz r p\nAl 85 237 237 2420000\n∽\nSiC / 360 360 2160000\nGa O To be tested 284000\n2 3\nBy TEM\nComment 4: It is not clear what the distance is in Fig. 5d, the corresponding sketch\nwill be helpful.\nResponse: We thank the reviewer for carefully reading our manuscript. This distance\nrepresents the relative distance within the simulation model. Based on your comments,\ncorresponding color blocks have been added to the diagram for better reader\ncomprehension.\nFigure R2. (a) The temperature profiles of the different interface. (b) The\ncorresponding magnified view of the interfaces.\n(Change 1) Figure 5 in the Manuscript (Page 17).\nComment 5: It would be nice to reproduce Fig. 5f with the FEM method and the ITC\nevaluated with MD for cross-validation of the experiment and multiscaling\nsimulations.\nResponse: Thank you for your professional suggestions. Regarding the multiscale\nsimulation method combining FEM and MD that you mentioned, we currently lack the\ncapability to accomplish it. We have thoroughly reviewed extensive literature and\nsought advice from relevant experts, but consistently, it has been pointed out that it is\nnot feasible to establish a model that simultaneously supports both FEM and MD\nsimulations. We sincerely apologize for this. However, in response to your concerns,\nwe have provided a comprehensive explanation in the previously submitted manuscript,\nutilizing a stepwise cross-validation approach involving experiments and simulations\nacross different scales.\nFirst, in Figures 4a–c, we experimentally demonstrated at the nanoscale using the\nTDTR testing method that introducing a SiO₂ transition layer between Ga₂O₃ and SiC\nsignificantly enhances the ITC compared to direct contact between Ga₂O₃ and SiC.\nBased on these experimental results at the micro-nano scale, we proposed that the\nunderlying reason for the improved ITC by introducing a low-thermal-conductivity\nSiO₂ transition layer is that the in-situ formed SiO₂ establishes covalent contact\ninterfaces with both Ga₂O₃ and SiC.\nBuilding on this hypothesis, we compared the heat transfer capabilities of two\nmodels: a covalently bonded Ga₂O₃/SiO₂/SiC model and a VdW interface-formed\nGa₂O₃/SiO₂/SiC model, using nanoscale MD simulations. The corresponding results are\nshown in Figures 4g–i. In the covalently bonded Ga₂O₃/SiO₂/SiC model, the ITC values\nfor the Ga₂O₃/SiO₂ and SiO₂/SiC interfaces are 919.9 MW/m²·K and 773.9 MW/m²·K,\nrespectively. In contrast, for the VdW interface-formed Ga₂O₃/SiO₂/SiC model, the ITC\nvalues for the Ga₂O₃/SiO₂ and SiO₂/SiC interfaces are 20.6 MW/m²·K and 27.7\nMW/m²·K, respectively. At the nanoscale, we cross-validated our experimental results\nusing both TDTR testing and MD simulations.\nNext, at the microscale, we performed validation through FEM simulations, with\nthe results shown in Figures 5a–d and Figure S10. For this scale, we established a heat\ntransfer model for Ga₂O₃/SiO₂/SiC with planar dimensions of 1 × 1 μm, where the\nthicknesses of Ga₂O₃ and SiO₂ were based on experimentally measured values of 68.4\nnm and 5.8 nm, respectively. In this model, the ITC values for the Ga₂O₃/SiO₂ and\nSiO₂/SiC interfaces were directly adopted from the MD calculation results (covalent\ncase: 919.9 MW/m²·K and 773.9 MW/m²·K for the Ga₂O₃/SiO₂ and SiO₂/SiC interfaces,\nrespectively; VDW case: 20.6 MW/m²·K and 27.7 MW/m²·K for the Ga₂O₃/SiO₂ and\nSiO₂/SiC interfaces, respectively). By incorporating the parameters obtained from MD\ncalculations into the FEM model, we calculated that the heat flux of the covalently\ncontacted Ga₂O₃/SiO₂/SiC (2.23 × 10⁹ pW/μm²) is 187% higher than that of the van der\nWaals contact case (7.77 × 10⁸ pW/μm²). At this stage, we demonstrated the superior\nheat transfer performance of the covalent interface at the microscale through\nsimulations.\nFinally, at the millimeter scale, we further validated the findings experimentally,\nwith the corresponding results shown in Figures 5e–g. At this scale, we prepared both\ncovalently contacted and van der Waals contacted Ga₂O₃/SiO₂/SiC samples with planar\ndimensions of 10 × 10 mm. Both samples were subjected to heating and cooling under\nthe same conditions. The experimental results showed that the maximum temperature\nof the covalently contacted sample was 29°C lower than that of the van der Waals\ncontacted sample, directly demonstrating its superior heat transfer performance.\nIn summary, throughout the validation process, we conducted analyses using both\nexperiments and simulations at the nanoscale and derived key thermal resistance\nparameters. These parameters were then incorporated into microscale finite element\nanalysis for further simulation validation. Finally, we performed experimental\nvalidation at the millimeter scale. We hope that this stepwise cross-validation approach,\ninvolving experiments and simulations across different scales, adequately addresses\nyour concerns. In our future work, we will continue to explore and advance the direct\ncoupling of FEM and MD, with the ultimate goal of realizing this methodology.\nReference\n1. Kresse G, Joubert D. From ultrasoft pseudopotentials to the projector augmented-wave\nmethod. Physical Review B 59, (1999).\n2. Perdew JP, Burke K, Ernzerhof M. Generalized Gradient Approximation Made Simple.\nPhysical Review Letters 77, 4 (1996).\n3. Grimme S, Antony J, Ehrlich S, Krieg H. A consistent and accurate ab initio parametrization\nof density functional dispersion correction (DFT-D) for the 94 elements H-Pu. The Journal\nof chemical physics 132, 154104 (2010).\n4. Monkhorst HJ, Pack JD. Special points for Brillouin-zone integrations. Physical Review B\n13, 5 (1976).\n5. Thompson AP, et al. LAMMPS - a flexible simulation tool for particle-based materials\nmodeling at the atomic, meso, and continuum scales. Computer Physics Communications\n271, 108171 (2022).\n6. Aktulga HM, Fogarty JC, Pandit SA, Grama AY. Parallel reactive molecular dynamics:\nNumerical methods and algorithmic techniques. parallel computing 38, 245-259 (2012).\n7. Chenoweth K, van Duin ACT, Goddard WA. ReaxFF Reactive Force Field for Molecular\nDynamics Simulations of Hydrocarbon Oxidation. The Journal of Physical Chemistry A 112,\n1040-1053 (2008).\n8. Rappe AK, Goddard WA, III. Charge equilibration for molecular dynamics simulations. The\nJournal of Physical Chemistry 95, 3358-3363 (1991).\n9. Nakano A. Parallel multilevel preconditioned conjugate-gradient approach to variable-\ncharge molecular dynamics. Computer Physics Communications 104, 59-69 (1997).\n10. Niefind F, et al. Watching (De)Intercalation of 2D Metals in Epitaxial Graphene: Insight into\nthe Role of Defects. Small 20, 2306554 (2024).\n11. Jiang P, Qian X, Yang R. Tutorial: Time-domain thermoreflectance (TDTR) for thermal\nproperty characterization of bulk and thin film materials. J Appl Phys 124, (2018)."}, {"round": 3, "reviewer_comments": {"reviewer_1": "The authors addressed my questions, and I think that the paper can be published"}, "author_response": ""}]

Decision letter