Machine-Learned Interatomic Potential for Predictive Simulation of MoS2 Epitaxy

Emir Bilgili; Nicholas Taormina; Richard Hennig; Simon R. Phillpot; Youping Chen

arxiv: 2512.15952 · v3 · pith:KA45DUUDnew · submitted 2025-12-17 · ❄️ cond-mat.mtrl-sci

Machine-Learned Interatomic Potential for Predictive Simulation of MoS2 Epitaxy

Emir Bilgili , Nicholas Taormina , Richard Hennig , Simon R. Phillpot , Youping Chen This is my paper

Pith reviewed 2026-05-16 21:08 UTC · model grok-4.3

classification ❄️ cond-mat.mtrl-sci

keywords MoS2machine-learned interatomic potentialUF3epitaxial growthmolecular dynamicsdefect energies2D materialsvan der Waals gaps

0 comments

The pith

A machine-learned interatomic potential for MoS2 reproduces defect energies and simulates layered epitaxial growth matching experiments.

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The paper develops a machine-learned interatomic potential using the UF3 framework for multilayer MoS2. This potential matches density functional theory results on lattice constants, interlayer binding energies, phonon spectra, elastic properties, and especially defect and edge formation energies with high correlation. Non-equilibrium molecular dynamics simulations with the potential then show homoepitaxial growth that forms layers separated by van der Waals gaps and triangular domains bounded by zigzag edges, consistent with experimental observations. The potential runs only about twice as slow as the fastest empirical models, opening the door to large-scale predictive simulations of 2D material synthesis.

Core claim

The UF3 machine-learned interatomic potential for MoS2 reproduces DFT lattice constants, binding energies, phonons, and elastic tensors across phases while capturing defect formation energies with R squared of 0.91 and relative zigzag versus armchair edge energies within 5 percent of DFT. Non-equilibrium molecular dynamics using this potential demonstrates layered homoepitaxial growth that produces van der Waals gaps between successive epilayers and triangular domains bounded by zigzag edges.

What carries the argument

The ultra-fast force field (UF3) machine-learned interatomic potential trained on DFT data, which models interatomic forces to enable large-scale non-equilibrium molecular dynamics of MoS2.

If this is right

Large-scale atomistic simulations of MoS2 growth become feasible at modest computational cost.
Growth morphologies under varied temperature or flux conditions can now be predicted before experiment.
Domain shapes and edge structures during synthesis can be modeled from accurate formation energies.
Multilayer stacking with controlled van der Waals gaps can be studied atomistically.

Where Pith is reading between the lines

These are editorial extensions of the paper, not claims the author makes directly.

The same training strategy could be applied to related materials such as WS2 or MoSe2 to study their epitaxial behavior.
Zigzag edge preference during growth may affect the electronic transport or optical properties of the resulting films.
Adding more explicit edge and defect trajectories to the training set could further improve accuracy for complex growth dynamics.

Load-bearing premise

The DFT training configurations cover the atomic environments that appear during dynamic non-equilibrium epitaxial growth, including evolving edges and defects.

What would settle it

A molecular dynamics run under the same growth conditions that fails to produce clear van der Waals gaps between layers or that forms domains bounded by armchair edges instead of zigzag edges would falsify the claim.

Figures

Figures reproduced from arXiv: 2512.15952 by Emir Bilgili, Nicholas Taormina, Richard Hennig, Simon R. Phillpot, Youping Chen.

**Figure 4.** Figure 4: Learned three-body (3B) interaction potentials for unique trios in the UF3 model. Each surface is plotted as a function of the two radial legs 𝑟𝑖𝑗, 𝑟𝑖𝑘, and the bond angle 𝜃𝑖𝑗𝑘. The 3B surfaces exhibit bean-shaped topologies with sharp repulsive regions (penalties > 0.4 eV) away from the ideal trigonal-prismatic geometry and neutral or mildly attractive pockets near equilibrium. These topologies reflect th… view at source ↗

**Figure 11.** Figure 11: Non-equilibrium MoS₂ homoepitaxial growth simulation at 1450 K using the UF3 MLIP. (a–c) Snapshots of the atomic positions of an evolving monolayer triangular domain during epitaxial growth, exposing S-terminated zigzag edges, consistent with experimental observations. The substrate consists of 104,976 atoms and is approximately 17.0 × 14.7 nm² in-plane and 8 nm thick; the growth rate was ~0.006 ML/ns. Vi… view at source ↗

read the original abstract

A machine-learned interatomic potential (MLIP) for multilayer MoS2 was developed using the ultra-fast force field (UF3) framework. The UF3 MLIP reproduces key properties in strong agreement with DFT including lattice constants, interlayer binding energies, and phase-stability. Furthermore, the potential reasonably captures the phonon spectra and the highly anisotropic elastic tensor across monolayer (1H) and bulk (2H, 3R) MoS2 phases. Critically, defect and edge formation energies are captured with high fidelity, exhibiting a strong correlation with DFT (R^2 = 0.91) across ten defective monolayers and reproducing the relative difference between the free energies of zigzag and armchair edges within 5% of DFT. Non-equilibrium molecular dynamics simulations reveal layered homoepitaxial growth consistent with experimental observations, demonstrating the formation of van der Waals gaps between successive epilayers and triangular domains bounded by zigzag edges. The robust UF3 MLIP, which is only ~2X slower than the fastest empirical potentials, enables large-scale atomistic simulations of MoS2 epitaxial growth.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Desk Editor's Note private letter to a colleague

This UF3 MLIP for MoS2 matches DFT on defects and edges with solid numbers and runs large-scale growth MD that produces the expected vdW gaps and zigzag domains, but the dynamic environments in epitaxy are not directly validated.

read the letter

The paper gives a working UF3 machine-learned potential for MoS2 that matches DFT closely on several key properties and then uses it for non-equilibrium simulations of epitaxial growth. It applies the UF3 framework to multilayer phases of MoS2. The validation covers lattice constants, interlayer binding, phonon spectra, and the anisotropic elastic tensor for 1H, 2H, and 3R phases. Defect formation energies across ten defective monolayers correlate with DFT at R squared of 0.91. The relative free energies of zigzag versus armchair edges come within 5 percent of DFT values. Those numbers look reliable. The new part is running molecular dynamics under deposition conditions to watch homoepitaxial layering. The simulations produce van der Waals gaps between layers and triangular domains with zigzag edges, which lines up with what experiments see. The potential runs only about twice as slow as the fastest empirical potentials, so it opens up larger systems. The soft spot is whether the training set covers the atomic environments that appear during actual growth. The validation uses equilibrium defect and edge structures. Growth involves moving adatoms, fluctuating edges, and stacking faults under continuous flux. If the potential gets the forces or barriers wrong in those transient states, the layer-by-layer ordering and domain shapes could be artifacts rather than true predictions. The paper does not report direct checks on adatom diffusion barriers or growth rates against experiment or higher-level calculations. This work is for people who need atomistic models of 2D material processing at scales beyond DFT. A reader working on TMD epitaxy or defect engineering in MoS2 would get a usable potential and some demonstration runs. It is not a new method but a careful domain application. The evidence is quantitative where it counts, and the growth claim has some independent grounding from the defect validation. I would send it to peer review. The manuscript is ready for referees to check the training details and the extrapolation assumptions.

Referee Report

1 major / 2 minor

Summary. The manuscript develops a UF3 machine-learned interatomic potential for multilayer MoS2. It reports quantitative agreement with DFT for lattice constants, interlayer binding energies, phonon spectra, anisotropic elastic tensors in 1H/2H/3R phases, and defect/edge formation energies (R²=0.91 across ten defective monolayers; zigzag-armchair edge energy difference within 5%). Non-equilibrium molecular dynamics simulations with this potential are used to model homoepitaxial growth, producing van der Waals gaps between epilayers and triangular domains bounded by zigzag edges, stated to be consistent with experimental observations. The potential is noted to be only ~2X slower than empirical potentials.

Significance. If the potential's accuracy extends reliably to the non-equilibrium regimes of growth, the work supplies an efficient tool for large-scale predictive simulations of MoS2 epitaxy that can access morphologies and dynamics beyond direct DFT reach. The quantitative multi-property validation and the reported reproduction of experimentally relevant growth features (vdW gaps, domain shapes) represent a concrete advance for 2D materials modeling.

major comments (1)

[Validation of defect/edge properties and epitaxial growth MD results] The central claim that non-equilibrium MD produces predictive layered homoepitaxial growth (vdW gaps and zigzag-bounded triangular domains) depends on the MLIP's fidelity outside the training distribution. Validation is limited to static equilibrium quantities: defect formation energies (R²=0.91) and edge energy differences (within 5% of DFT). No explicit checks are reported on force accuracy, adatom migration barriers, or forces in transient configurations such as partial edges or deposition-flux environments encountered in the growth trajectories. This leaves open whether the observed ordering and morphologies are robust predictions or artifacts of extrapolation.

minor comments (2)

[Abstract] In the abstract, 'phase-stability' should read 'phase stability'.
[Methods] The manuscript notes free parameters (UF3 hyperparameters and cutoffs) but does not tabulate their final values; adding these would aid reproducibility.

Simulated Author's Rebuttal

1 responses · 0 unresolved

We thank the referee for their constructive and detailed review of our manuscript. We address the major comment point by point below, providing the strongest honest defense of our work while acknowledging limitations where they exist.

read point-by-point responses

Referee: [Validation of defect/edge properties and epitaxial growth MD results] The central claim that non-equilibrium MD produces predictive layered homoepitaxial growth (vdW gaps and zigzag-bounded triangular domains) depends on the MLIP's fidelity outside the training distribution. Validation is limited to static equilibrium quantities: defect formation energies (R²=0.91) and edge energy differences (within 5% of DFT). No explicit checks are reported on force accuracy, adatom migration barriers, or forces in transient configurations such as partial edges or deposition-flux environments encountered in the growth trajectories. This leaves open whether the observed ordering and morphologies are robust predictions or artifacts of extrapolation.

Authors: We thank the referee for raising this important point regarding transferability to non-equilibrium growth conditions. Our validation strategy prioritizes properties most directly relevant to epitaxy: defect formation energies across ten defective monolayers (R²=0.91) and the zigzag-armchair edge energy difference (within 5% of DFT), both of which govern domain shape and layer stacking. These quantities were computed on configurations outside the primary training set, providing evidence of extrapolation capability. The non-equilibrium MD trajectories spontaneously produce van der Waals gaps and zigzag-bounded triangular domains that match experimental morphologies, which would be unlikely if the potential were severely misbehaving in transient states. That said, we did not report explicit force RMSE on deposition-flux or partial-edge configurations, nor adatom migration barriers. We will revise the manuscript to add (i) force-error statistics on a held-out set containing edge and defective structures and (ii) a brief discussion of training-data coverage for growth-relevant environments, thereby strengthening the transferability argument without altering the central claims. revision: partial

Circularity Check

0 steps flagged

No circularity: MLIP trained on DFT data; epitaxial growth is extrapolation from validated potential

full rationale

The paper develops a UF3 MLIP fitted to DFT configurations for MoS2, then validates it on held-out static properties including defect formation energies (R²=0.91) and edge energy differences (within 5% of DFT). Non-equilibrium MD simulations of homoepitaxial growth are run with this potential to observe vdW gaps and zigzag-bounded domains. No self-definitional reductions, fitted inputs renamed as predictions, or load-bearing self-citations appear in the derivation chain. The growth trajectories constitute independent dynamic extrapolation beyond the enumerated training/validation sets, making the central claim self-contained against external benchmarks.

Axiom & Free-Parameter Ledger

1 free parameters · 1 axioms · 0 invented entities

The central claim rests on the assumption that DFT data used for training captures all relevant physics for epitaxy; no new physical entities are postulated.

free parameters (1)

UF3 hyperparameters and cutoff radii
Chosen during training to fit DFT forces and energies; specific values not stated in abstract.

axioms (1)

domain assumption DFT calculations provide accurate reference data for training
Invoked implicitly when claiming agreement with DFT as ground truth.

pith-pipeline@v0.9.0 · 5508 in / 1183 out tokens · 33315 ms · 2026-05-16T21:08:14.162005+00:00 · methodology

discussion (0)

Lean theorems connected to this paper

Citations machine-checked in the Pith Canon. Every link opens the source theorem in the public Lean library.

IndisputableMonolith/Cost/FunctionalEquation.lean washburn_uniqueness_aczel unclear

?

unclear
Relation between the paper passage and the cited Recognition theorem.

UF3 MLIP reproduces key properties... defect and edge formation energies... R²=0.91... Non-equilibrium molecular dynamics simulations reveal layered homoepitaxial growth... van der Waals gaps... triangular domains bounded by zigzag edges.
IndisputableMonolith/Foundation/RealityFromDistinction.lean reality_from_one_distinction unclear

?

unclear
Relation between the paper passage and the cited Recognition theorem.

The UF3 framework represents atomic interactions through cubic B-spline interpolations of two- and three-body terms.

What do these tags mean?

matches: The paper's claim is directly supported by a theorem in the formal canon.
supports: The theorem supports part of the paper's argument, but the paper may add assumptions or extra steps.
extends: The paper goes beyond the formal theorem; the theorem is a base layer rather than the whole result.
uses: The paper appears to rely on the theorem as machinery.
contradicts: The paper's claim conflicts with a theorem or certificate in the canon.
unclear: Pith found a possible connection, but the passage is too broad, indirect, or ambiguous to say the theorem truly supports the claim.

Reference graph

Works this paper leans on

8 extracted references · 8 canonical work pages

[1]

Manzeli, D

Kinetic dependence homoepitaxially grown MoS 2 epilayers. a) The growth rate is ~0.05 ML/ns , the growth temperature is 900K, and the defect density of the triangular domain is ~2.1 1/nm2. b) The growth rate is ~0.02 ML/ns, the growth temperature is 1200K and the defect density of the triangular domain is ~0.3 5 1/nm2. c) The growth rate is ~0.0 06 ML/ns,...

work page doi:10.1038/natrevmats.2017.33 2017
[2]

(5) Dong, J.; Zhang, L.; Dai, X.; Ding, F

https://doi.org/10.1007/s40820-023-01273-5. (5) Dong, J.; Zhang, L.; Dai, X.; Ding, F. The Epitaxy of 2D Materials Growth. Nat. Commun. 2020, 11 (1),

work page doi:10.1007/s40820-023-01273-5 2020
[3]

(6) Aljarb, A.; Fu, J.-H.; Hsu, C.-C.; Chuu, C.-P.; Wan, Y .; Hakami, M.; Naphade, D

https://doi.org/10.1038/s41467-020-19752-3. (6) Aljarb, A.; Fu, J.-H.; Hsu, C.-C.; Chuu, C.-P.; Wan, Y .; Hakami, M.; Naphade, D. R.; Yengel, E.; Lee, C.-J.; Brems, S.; Chen, T.-A.; Li, M.-Y .; Bae, S.-H.; Hsu, W.-T.; Cao, Z.; Albaridy, R.; Lopatin, S.; Chang, W.-H.; Anthopoulos, T. D.; Kim, J.; Li, L.-J.; Tung, V . Ledge-Directed Epitaxy of Continuously ...

work page doi:10.1038/s41467-020-19752-3 2020
[4]

(28) Wei, Z.; Zhang, C.; Kan, Y .; Zhang, Y .; Chen, Y

https://ikee.lib.auth.gr/record/334898/files/Diouslous.pdf. (28) Wei, Z.; Zhang, C.; Kan, Y .; Zhang, Y .; Chen, Y . Developing Machine Learning Potential for Classical Molecular Dynamics Simulation with Superior Phonon Properties. Comput. Mater. Sci. 2022, 202, 111012. https://doi.org/10.1016/j.commatsci.2021.111012. (29) Han, L.; Qin, G.; Jia, B.; Chen,...

work page doi:10.1016/j.commatsci.2021.111012 2022
[5]

(33) Nair, A

https://doi.org/10.1038/s41524-024-01357-9. (33) Nair, A. K.; Da Silva, C. M.; Amon, C. H. Machine-Learning-Derived Thermal Conductivity of Two-Dimensional TiS2/MoS2 van Der Waals Heterostructures. APL Mach. Learn. 2024, 2 (3), 036115. https://doi.org/10.1063/5.0205702. (34) Wan, S.; Tong, R.; Han, B.; Zhang, H. Machine Learning Interatomic Potential for ...

work page doi:10.1038/s41524-024-01357-9 2024
[6]

(39) Tipton, W

https://doi.org/10.48550/arXiv.2511.08330. (39) Tipton, W. W.; Hennig, R. G. A Grand Canonical Genetic Algorithm for the Prediction of Multi-Component Phase Diagrams and Testing of Empirical Potentials. J. Phys. Condens. Matter 2013, 25 (49), 495401. https://doi.org/10.1088/0953-8984/25/49/495401. (40) Revard, B. C.; Tipton, W. W.; Yesypenko, A.; Hennig, ...

work page doi:10.48550/arxiv.2511.08330 2013
[7]

(51) Tree-Structured Parzen Estimator: Understanding Its Algorithm Components and Their Roles for Better Empirical Performance

https://doi.org/10.1088/0953-8984/16/23/005. (51) Tree-Structured Parzen Estimator: Understanding Its Algorithm Components and Their Roles for Better Empirical Performance . OptunaHub. https://hub.optuna.org/samplers/tpe_tutorial/ (accessed 2025-10-06). (52) Plimpton, S. Fast Parallel Algorithms for Short -Range Molecular Dynamics. J. Comput. Phys. 1995, ...

work page doi:10.1088/0953-8984/16/23/005 2025
[8]

(55) George, J.; Hautier, G.; Bartók, A

https://doi.org/10.5281/zenodo.5668319. (55) George, J.; Hautier, G.; Bartók, A. P.; Csányi, G.; Deringer, V . L. Combining Phonon Accuracy with High Transferability in Gaussian Approximation Potential Models. J. Chem. Phys. 2020, 153 (4), 044104. https://doi.org/10.1063/5.0013826. (56) Wei, X.; Wang, Y .; Shen, Y .; Xie, G.; Xiao, H.; Zhong, J.; Zhang, G...

work page doi:10.5281/zenodo.5668319 2020

[1] [1]

Manzeli, D

Kinetic dependence homoepitaxially grown MoS 2 epilayers. a) The growth rate is ~0.05 ML/ns , the growth temperature is 900K, and the defect density of the triangular domain is ~2.1 1/nm2. b) The growth rate is ~0.02 ML/ns, the growth temperature is 1200K and the defect density of the triangular domain is ~0.3 5 1/nm2. c) The growth rate is ~0.0 06 ML/ns,...

work page doi:10.1038/natrevmats.2017.33 2017

[2] [2]

(5) Dong, J.; Zhang, L.; Dai, X.; Ding, F

https://doi.org/10.1007/s40820-023-01273-5. (5) Dong, J.; Zhang, L.; Dai, X.; Ding, F. The Epitaxy of 2D Materials Growth. Nat. Commun. 2020, 11 (1),

work page doi:10.1007/s40820-023-01273-5 2020

[3] [3]

(6) Aljarb, A.; Fu, J.-H.; Hsu, C.-C.; Chuu, C.-P.; Wan, Y .; Hakami, M.; Naphade, D

https://doi.org/10.1038/s41467-020-19752-3. (6) Aljarb, A.; Fu, J.-H.; Hsu, C.-C.; Chuu, C.-P.; Wan, Y .; Hakami, M.; Naphade, D. R.; Yengel, E.; Lee, C.-J.; Brems, S.; Chen, T.-A.; Li, M.-Y .; Bae, S.-H.; Hsu, W.-T.; Cao, Z.; Albaridy, R.; Lopatin, S.; Chang, W.-H.; Anthopoulos, T. D.; Kim, J.; Li, L.-J.; Tung, V . Ledge-Directed Epitaxy of Continuously ...

work page doi:10.1038/s41467-020-19752-3 2020

[4] [4]

(28) Wei, Z.; Zhang, C.; Kan, Y .; Zhang, Y .; Chen, Y

https://ikee.lib.auth.gr/record/334898/files/Diouslous.pdf. (28) Wei, Z.; Zhang, C.; Kan, Y .; Zhang, Y .; Chen, Y . Developing Machine Learning Potential for Classical Molecular Dynamics Simulation with Superior Phonon Properties. Comput. Mater. Sci. 2022, 202, 111012. https://doi.org/10.1016/j.commatsci.2021.111012. (29) Han, L.; Qin, G.; Jia, B.; Chen,...

work page doi:10.1016/j.commatsci.2021.111012 2022

[5] [5]

(33) Nair, A

https://doi.org/10.1038/s41524-024-01357-9. (33) Nair, A. K.; Da Silva, C. M.; Amon, C. H. Machine-Learning-Derived Thermal Conductivity of Two-Dimensional TiS2/MoS2 van Der Waals Heterostructures. APL Mach. Learn. 2024, 2 (3), 036115. https://doi.org/10.1063/5.0205702. (34) Wan, S.; Tong, R.; Han, B.; Zhang, H. Machine Learning Interatomic Potential for ...

work page doi:10.1038/s41524-024-01357-9 2024

[6] [6]

(39) Tipton, W

https://doi.org/10.48550/arXiv.2511.08330. (39) Tipton, W. W.; Hennig, R. G. A Grand Canonical Genetic Algorithm for the Prediction of Multi-Component Phase Diagrams and Testing of Empirical Potentials. J. Phys. Condens. Matter 2013, 25 (49), 495401. https://doi.org/10.1088/0953-8984/25/49/495401. (40) Revard, B. C.; Tipton, W. W.; Yesypenko, A.; Hennig, ...

work page doi:10.48550/arxiv.2511.08330 2013

[7] [7]

(51) Tree-Structured Parzen Estimator: Understanding Its Algorithm Components and Their Roles for Better Empirical Performance

https://doi.org/10.1088/0953-8984/16/23/005. (51) Tree-Structured Parzen Estimator: Understanding Its Algorithm Components and Their Roles for Better Empirical Performance . OptunaHub. https://hub.optuna.org/samplers/tpe_tutorial/ (accessed 2025-10-06). (52) Plimpton, S. Fast Parallel Algorithms for Short -Range Molecular Dynamics. J. Comput. Phys. 1995, ...

work page doi:10.1088/0953-8984/16/23/005 2025

[8] [8]

(55) George, J.; Hautier, G.; Bartók, A

https://doi.org/10.5281/zenodo.5668319. (55) George, J.; Hautier, G.; Bartók, A. P.; Csányi, G.; Deringer, V . L. Combining Phonon Accuracy with High Transferability in Gaussian Approximation Potential Models. J. Chem. Phys. 2020, 153 (4), 044104. https://doi.org/10.1063/5.0013826. (56) Wei, X.; Wang, Y .; Shen, Y .; Xie, G.; Xiao, H.; Zhong, J.; Zhang, G...

work page doi:10.5281/zenodo.5668319 2020