Phonon-Mediated Thermal Transport in Nanocrystalline Silicon Using Machine-Learning Interatomic Potentials

Source: arXiv:2607.06470 · Published 2026-07-07 · By Houssem Rezgui, Catalina Coll Benejam, Miguel Pruneda, Clivia M. Sotomayor Torres

TL;DR

This paper addresses the challenge of accurately modeling phonon-mediated thermal transport in structurally complex nanocrystalline silicon, where grain boundaries (GBs) strongly impede heat dissipation. Classical empirical interatomic potentials have known limitations in capturing vibrational properties and phonon scattering at interfaces, leading to less reliable predictions. The authors develop a unified machine-learning interatomic potentials (MLIPs) framework combining Gaussian approximation potential (GAP) and multi-atomic cluster expansion (MACE) models with lattice-dynamical calculations and non-equilibrium molecular dynamics (NEMD). This approach derives harmonic and anharmonic force constants to compute phonon dispersions, lifetimes, and lattice thermal conductivity consistently. They apply it to bulk and nanocrystalline silicon bicrystals with varying GB misorientation, benchmarking against classical potentials (Stillinger–Weber and Tersoff) and experimental data. The results demonstrate that MLIPs provide significantly improved accuracy for bulk vibrational properties and anharmonic phonon scattering. In nanocrystalline silicon, NEMD reveals that thermal boundary resistance at GBs depends sensitively on atomic-scale interfacial structure and roughness, with MLIPs predicting much stronger phonon scattering and lower thermal conductivity than classical potentials. The work establishes MLIPs as a predictive methodology for nanoscale thermal transport in disordered low-dimensional materials.

Key findings

• MLIPs (GAP and MACE) reproduce bulk silicon phonon density of states and dispersion relations closely matching experimental data, unlike Stillinger–Weber (SW) and Tersoff potentials which overestimate optical phonon frequencies and distort branch curvature (Fig. 1).
• Phonon lifetimes predicted by MLIPs at 300 K agree well with physical expectations and experimental trends, whereas SW and Tersoff overestimate lifetimes leading to inflated mean free paths (Fig. 2).
• Temperature-dependent lattice thermal conductivity obtained using MLIPs accurately matches experimental values and captures its decreasing trend with increasing temperature, while empirical potentials systematically overestimate conductivity, especially at low temperatures (Fig. 3a).
• In symmetric tilt grain-boundary bicrystals (3189 atoms), phonon density of states in the GB core shows significant broadening and redistribution compared to near-GB and pristine bulk, increasing with misorientation angle (Fig. 5).
• Phonon lifetimes within the GB core drop to a few picoseconds across frequencies with only weak dependence on misorientation angle beyond 10°, indicating phonon scattering is dominated by local atomic disorder inside the interface (Fig. 6).
• Mean free paths of phonons within the GB core are drastically reduced to tens of nanometers even for low-frequency acoustic modes, with enhanced suppression at higher misorientation (Fig. 7).
• NEMD simulations quantify thermal boundary resistance showing that MLIPs predict substantially higher interface scattering and lower interfacial thermal conductivity compared to classical potentials.
• Cumulative thermal conductivity in the GB core saturates at shorter mean free paths and is nearly temperature independent, in contrast to bulk, revealing interface-dominated heat transport regimes (Fig. 9).

Methodology — deep read

The authors start with the premise that accurate modeling of phonon transport in nanostructured silicon requires faithful capture of harmonic and anharmonic vibrational properties, especially near grain boundaries where local atomic disorder disrupts bulk phonon modes.

1. Threat Model & Assumptions: The study focuses on physical phonon transport mechanisms and does not consider adversarial threats. The aim is to provide predictive physical models for thermal conductivity and phonon scattering at grain boundaries, assuming no extrinsic perturbations beyond the modeled atomic structures. The MLIPs are trained to replicate DFT-calculated energies and forces, treating DFT as ground truth.

2. Data: Bulk silicon data comprising 216-atom 3×3×3 supercells, relaxed and equilibrated via MD at 300 K, serve as the main reference dataset with energies and forces extracted. Nanocrystalline structures with symmetric tilt grain boundaries were built in 60×30×30 Å3 boxes containing ~3189 atoms, with four different misorientation angles (10°, 20°, 30°, 40°). Atomic configurations and forces from these simulations were collected, with separate datasets extracted for near-GB and GB core regions.

3. Architecture/Algorithm: Two MLIPs were constructed— a Gaussian Approximation Potential (GAP) using combined two-body and many-body Smooth Overlap of Atomic Positions (SOAP) descriptors, and a MACE model employing equivariant message-passing neural networks preserving physical symmetries. Both potentials were trained on the same DFT reference dataset to ensure comparability. The GAP model uses nmax=7, lmax=6, ζ=2 with 800 sparse points, trained on energies and forces simultaneously. MACE was trained with 20% validation split for up to 200 epochs on GPUs.

4. Training Regime: The GAP model was retrained with the QUIP framework using a combined SOAP and two-body descriptor. MACE used standard train routines optimizing energy and force losses jointly. Both models relied on DFT-based datasets including crystalline and amorphous configurations to cover diverse atomic environments found in bulk and grain boundary regions.

5. Evaluation Protocol: After training, harmonic and anharmonic force constants were computed using finite displacement methods via Phonopy and Phono3py on supercells, with forces derived from the MLIPs. These were used to calculate phonon density of states, dispersions, lifetimes, mean free paths, group velocities, and lattice thermal conductivity. MLIPs were benchmarked against classical Stillinger–Weber and Tersoff potentials and compared to experimental measurements. Non-equilibrium molecular dynamics (NEMD) simulations quantified thermal boundary resistance at GBs, directly revealing phonon scattering effects from interfacial structural disorder.

6. Reproducibility: The paper references use of open-source tools (LAMMPS, QUIP/QUIPY, Phonopy/Phono3py, ASE, Atomsk). Training and evaluation details including hyperparameters and protocols are described with additional information in Supplementary Material. Code and trained model availability are not explicitly stated.

Concrete example: For bulk silicon, the authors equilibrated a 216-atom supercell at 300 K, extracted forces, and trained the GAP and MACE potentials. Using these, they computed second- and third-order force constants with 3×3×3 displacement supercells and calculated phonon dispersion and lifetimes. These predicted phonon lifetimes matched closely with experimental trends, unlike classical potentials which overestimated lifetimes. This accurate description of anharmonic scattering yielded lattice thermal conductivity in excellent agreement with measurements over 100–500 K. Extending to bicrystal models with explicit grain boundaries, MLIPs informed NEMD simulations demonstrated that atomic-scale GB disorder strongly reduces phonon lifetimes, mean free paths, and thermal conductivity inside the interfacial region, consistent with the anticipated physics of enhanced phonon scattering by crystalline mismatch and roughness.

Technical innovations

• Integration of Gaussian Approximation Potentials (GAP) and equivariant neural network-based MACE models to jointly model harmonic and anharmonic phonon properties in bulk and nanocrystalline silicon with high fidelity.
• Use of MLIP-derived second- and third-order force constants within a unified Phonopy/Phono3py lattice-dynamical workflow enabling internally consistent computation of phonon dispersions, lifetimes, and thermal conductivity.
• Systematic application of non-equilibrium molecular dynamics (NEMD) simulations informed by MLIPs to directly quantify thermal boundary resistance and phonon scattering at grain boundaries with variable misorientation and roughness.
• Demonstration that classical empirical potentials (Stillinger–Weber, Tersoff) underestimate anharmonic phonon scattering and overpredict thermal conductivity, especially near interfaces, whereas MLIPs match first-principles and experimental data closely.

Datasets

• Bulk silicon supercells — 216 atoms — DFT-generated reference dataset including crystalline and amorphous configurations
• Symmetric tilt grain-boundary bicrystals — ~3189 atoms per sample — Atomistically constructed and relaxed nanocrystalline silicon models

Baselines vs proposed

• Stillinger–Weber potential: phonon lifetime overestimated by ~2× in acoustic range vs GAP and MACE which align with experiment
• Tersoff potential: phonon dispersion optical branch frequencies overestimated by up to 15% vs MLIPs matching experiment (Fig. 1)
• Thermal conductivity at 300 K: SW and Tersoff overestimate by ~30-50% compared to experiment; MLIPs within 10% of measured values (Fig. 3a)
• Cumulative thermal conductivity as a function of phonon mean free path: MLIPs match Li et al. first-principles results and FDTR experiments within experimental uncertainty (Fig. 3b)
• Grain boundary core phonon lifetimes: Reduced by about an order of magnitude relative to bulk for all misorientation angles using MLIPs (Fig. 6)
• Thermal boundary resistance from NEMD: MLIPs predict 2–3× higher resistance than classical potentials, reflecting stronger phonon scattering

Figures from the paper

Figures are reproduced from the source paper for academic discussion. Original copyright: the paper authors. See arXiv:2607.06470.

Fig 2: Phonon lifetimes of bulk silicon as a function of phonon frequency predicted by

Fig 3: Bulk silicon thermal transport properties: (a) thermal conductivity compared with

Fig 4: Atomistic structures of symmetric tilt grain boundaries in silicon for misorientation

Fig 5: Phonon density of states of near-GB and grain-boundary core regions for different

Fig 6: Phonon lifetimes of the grain-boundary core for different misorientation angles: (a)

Fig 6 (page 12).

Fig 7: Phonon mean free path as a function of frequency within the grain-boundary core of

Fig 10: Schematic of the non-equilibrium molecular dynamics setup used to compute

Limitations

• MLIP training relies on DFT reference data that may not fully sample all relevant grain boundary atomic configurations or temperature ranges.
• Finite system sizes (up to ~3000 atoms) in NEMD limit exploration of larger-scale microstructural effects and distribution of GB types encountered experimentally.
• Lack of explicit evaluation of model performance under extreme roughness or amorphization beyond symmetric tilt grain boundaries.
• No adversarial testing or uncertainty quantification was reported to assess robustness of MLIP predictions under perturbations.
• Reproducibility constrained by no clear public release of trained MLIP weights or full training datasets.
• No direct exploration of electron-phonon coupling effects which can be important for thermal transport in devices.

Open questions / follow-ons

• How do MLIPs perform in modeling phonon transport across more complex, non-symmetric grain boundaries with variable atomic scale roughness or impurity segregation?
• Can MLIPs be extended to accurately model thermal transport in polycrystalline silicon with distributions of GBs and defects at larger length scales?
• What is the role of electron-phonon interactions and carrier scattering at grain boundaries in combined thermal and electrical transport?
• How does temperature-dependent structural dynamics beyond harmonic and anharmonic phonons influence interfacial thermal resistance at nanoscale GBs?

Why it matters for bot defense

Though not directly related to bot-defense or CAPTCHAs, this work demonstrates how advanced machine-learning potentials can provide more physically accurate models in complex systems than classical fixed-form potentials. For bot-defense engineers, the methodological rigor exemplified here—careful benchmarking against experiments, unified computational pipelines, and explicit modeling of interfaces—highlights the importance of ensuring simulation fidelity when deploying ML models in security-critical applications. Captchas or bot detection systems relying on ML could similarly benefit from physics-informed modeling, systematic validation, and attention to boundary conditions or adversarial environments. Furthermore, the study underscores the value of combining machine learning with domain-specific knowledge to improve predictive capabilities — a useful principle when considering defenses against sophisticated automated attacks.

Cite

@article{arxiv2607_06470,
  title={ Phonon-Mediated Thermal Transport in Nanocrystalline Silicon Using Machine-Learning Interatomic Potentials },
  author={ Houssem Rezgui and Catalina Coll Benejam and Miguel Pruneda and Clivia M. Sotomayor Torres },
  journal={arXiv preprint arXiv:2607.06470},
  year={ 2026 },
  url={https://arxiv.org/abs/2607.06470}
}

Read the full paper

• arXiv:2607.06470 (abstract)
• arXiv:2607.06470 (PDF)