Molecular dynamics simulations of the sputtering of boron and boron oxide surfaces

Abstract

We employ classical molecular dynamics (MD) simulations to study processes governed by particle–surface interactions. The interatomic potential energy functions are described by a neural network potential (NNP) trained on an extensive set of density functional theory (DFT) calculations in a semi-iterative fashion. Potential construction and simulation set up follow the Behler–Parrinello approach. As a specific example we investigate sputtering, reflection, and adsorption phenomena occurring on boron and boron oxide surfaces under the impact of deuterium atoms, systems that reflect recent developments in materials science. Besides the frequent use of boron surfaces as an oxygen-gathering material in technical applications, boron-based compounds will be used in future fusion devices. Understanding their interaction with energetic plasma particles is essential, yet their stability and sputtering behavior at the atomic level have remained largely unexplored. From our simulations, we analyzed sputtering yields, adsorption, and reflection events on both boron and boron oxide surfaces. The production simulations included about 750 atoms and covered a range of incident energies and impact angles. Increasing the deuterium impact energy generally leads to an increase in sputtering yields but with distinct energy-dependent trends. The sputtering yield of boron from boron oxide surfaces remains significantly lower than from surfaces of pristine boron. We use an analytical approach to estimate an effective surface binding energy. The work presented here, especially the various energy- and angle dependent rates, can be used to create a parametric model of the B and B₂O₃ surfaces under the impact of hot particles.

---

1 Introduction

The interactions between surfaces and their environment are crucial for understanding the degradation of materials and devices in various applications, ranging from semiconductor manufacturing to nuclear fusion.^1–3 In thermonuclear devices, particularly in ITER, the interaction between the plasma particles and plasma-facing materials (PFM) is important in determining the efficiency, lifetime, and safety of the fusion process.⁴ Among the components most exposed to the plasma environment are the divertor and the first wall of the blanket. These materials must withstand extreme heat loads, high particle fluxes, and continuous irradiation by energetic plasma particles such as deuterium (D), tritium (T), helium (He), and many other ions present in the plasma.⁵ During the fusion operations, PFMs can degrade the structural integrity and lead to the contamination of the plasma with impurities. To reduce these effects, materials with low atomic numbers (low-Z), like boron and beryllium, are important due to their low sputtering yields, low hydrogen retention at higher temperatures, and their smaller loss of energy when they enter the plasma.

Boron and its compounds, particularly boron trioxide (B₂O₃), have attracted attention as potential candidates for use in PFM protective layers. Boronization – depositing boron onto PFM surfaces – is routinely applied in several experiments to reduce impurity influx and improve plasma confinement, taking advantage of both the low nuclear weight and the oxygen affinity of boron.^6–10 Sputtering of boron hydrates under the impact of hydrogen isotopes is called chemical sputtering. Such molecule formation processes may cause larger degradation of boron-containing surfaces than purely physical sputtering. The detailed atomistic mechanisms governing the sputtering and surface evolution of boron-based materials under deuterium bombardment are not yet fully understood.

A large number of theoretical and computational studies have been devoted to the sputtering of materials under various conditions.¹¹ Fitting experimental sputtering yields using empirical or semi-empirical formulas, such as the Eckstein¹² or Yamamura¹³ models, allow for the interpretation and extrapolation of experimental data. Quantum chemical calculations, particularly DFT and the various AIMD schemes, are the most accurate methods among computational techniques, allowing analysis of sputtering events, surface binding energies, energy transfer, and mechanisms of sputtering during particle impacts, although they remain computationally expensive and impractical for large-scale and/or long-time simulations.

Classical MD simulations with simple analytical potentials are fast and can model large systems.¹⁴ However, their accuracy depends strongly on the quality of the interatomic potentials. Recent progress in machine learning, especially neural network potentials (NNPs) trained on DFT datasets, now allow simulations that retain DFT-level accuracy with near-classical MD efficiency.^15,16 These NNP-based MD simulations are becoming essential in sputtering research, bridging the gap between computational cost and predictive accuracy for complex surface and plasma–wall interaction phenomena. Simulations with the Monte-Carlo based Binary Collision Approximation (BCA) method offer statistical predictions over larger ensembles and longer timescales but often require fitting parameters and miss atomistic detail.¹⁷

Many studies in the sputtering field have focused on fitting interatomic potentials and on investigating the erosion behavior of boron and its compounds. Domínguez-Gutiérrez and Krstić¹⁸ investigated the erosion of boron and boron carbide in contact with deuterium plasma using reactive force fields. They calculated sputtering yields under low-energy particles (up to 30 eV) and found a significant effect of chemical sputtering in boron and boron carbides. Deng and Du¹⁹ developed and validated composition-dependent empirical potentials specifically for boron oxide in multicomponent glasses, enabling large-scale MD studies of B₂O₃ structural dynamics and rare atomistic configurations. However, accurate interatomic potentials for B/B₂O₃-plasma particle are still lacking.

In this paper MD simulations of sputtering processes from boron and boron oxide under deuterium bombardment are performed at various incident energies and impact angles. A specifically created neural network potential drives the simulations. The data and findings presented below contribute to an understanding of surface erosion processes, enhancing the development of materials for future fusion applications. The simulation trajectories provide detailed information of incident particle reflection, retention, and adsorption ratios as functions of impact energy and angle, together with effective surface binding energies (E_SB), allowing insight into the general field of surface interactions under particle bombardment.

2 Methods

2.1 Behler–Parrinello method

In the Behler–Parrinello approach the total energy of the system is expressed as the sum of its atomic energies E_i.²⁰ The atomic energies are derived from element-specific, atom-centered high-dimensional neural networks (HDNN). The Cartesian coordinates of the neighbors of each atom are transformed into weighted radial and angular symmetry functions, to ensure invariance against rotation, translation and exchange within the same atomic species. A cut-off radius of 7 Å restricts the environment of each atom. The atomic coordinates were transformed through symmetry functions, using parameters provided in Tables S1 and S2 of the SI, and served as input for the feedforward neural network. The network had two hidden layers with 25 nodes each, using soft-plus activation functions and offsets to predict the material energies. Once the parameters of the NNPs are trained (see Section 2.3), the MD simulations are conducted using a modified LAMMPS code.^21,22

2.2 Generation of training data

Density functional theory calculations generate the energies and forces used to train and test the NN parameters. The Behler–Parinello formalism thus acts as a bridge to produce DFT-quality forces in the simulations which are much faster than it would be possible with direct DFT. Initial training and test sets are generated from systems where the atomic positions are allowed to slightly deviate from their DFT-optimized equilibrium structures. Here, after a first training, new data were generated from actual sputtering simulations. This yielded trajectories with various incident angles of deuterium on small supercells. During these simulations, the lowest atomic layer was fixed, while other atoms were allowed to move unrestricted. The collected structures from these runs were recomputed at the DFT level. In all DFT calculations we employed the PW91 functional²³ within the generalized gradient approximation (GGA) framework as implemented in the Vienna Ab Initio Simulation Package (VASP).^24,25 A large supercell and a Gamma-centered k-point mesh of 3 × 3 × 3 were used with a plane wave basis set and a cut-off energy of 350 eV. The convergence threshold of the forces was set to 10⁻⁴ eV Å⁻¹. The simulations were spin-unpolarized. This methodology is consistent with previous studies on the interaction of plasma particles with various targets.

2.3 Network training

The training data were then used to train Neural Network Potentials for the boron–deuterium (B(100)–D)²⁶ and boron-oxide–deuterium (B₂O₃(001)–D)²⁷ systems. The final reference datasets include a mix of surface and bulk structures, with 4861 configurations for the pure boron systems and 5695 for the boron oxide systems, containing a total of 655 433 and 605 454 atoms, respectively. Of these, 90% of the structures were used to train the final NNP, while the remaining 10% were used for testing. After 100 training steps, the model achieved root-mean-square errors (RMSE) of 6.72 and 10.1 meV per atom for the predicted energies, and 0.43 and 0.66 eV Å⁻¹ for the atomic forces in the boron and boron oxide test sets, respectively.

Fig. 1 shows the accuracy of the NNP with respect to the DFT data together with the distribution of data points for the B–D dataset; the corresponding plots for the B₂O₃–D dataset, are presented in Fig. S1 (SI). Considering the ranges of both energies and forces, the figures confirm that indeed the DFT data is carried over to the NNP model with no visible loss of accuracy. The expression to create the symmetry-adapted input for the NNP are given in the SI.

Fig. 1 Correlation between DFT-calculated and NNP-predicted atomic energies (left) and forces (right) for the boron–deuterium dataset, accompanied by histograms illustrating the distribution of their values. The structures shown below and at the bottom scatter points provide representative examples corresponding to the distributions presented above.

Further validation of the trained NNPs was done by calculating the surface atom removal energy and the adsorption curves of deuterium atoms on boron (Fig. S2, SI) and boron oxide (Fig. S3, SI) surfaces. These calculations also served as cross-validation tests of the trained NNPs, allowing us to estimate the RMSE of models on these structures. The obtained RMSE values for energies and forces are given in Tables S3 and S4 (SI). The results show that the RMSE values are lower than those from the training dataset tests, which confirms that the NNPs can be reliably used for sputtering simulations.

2.4 MD simulation protocol

Molecular dynamics can simulate sputtering in two ways: cumulative (many impacts in one simulation) or non-cumulative (one impact in one simulation). The cumulative case resembles real processes, while the non-cumulative case gives clearer statistics by isolating single impacts without accumulated damage. Employing the trained NNPs, we performed non-cumulative MD simulations of bombardment by deuterium (D) projectiles with energies of ranging from 20 to 200 eV at impact angles of 0°, 30°, and 60°. 2000 simulations for each set of parameters were performed and the resulting trajectories were used in the statistical analysis of the events. The initial simulation boxes were created from B₇₆₈ and (B₂O₃)₁₅₀ surfaces, oriented along the z direction. The box dimensions (in x, y, z) were 19.60 × 17.05 × 32.1 Å for the boron system and 21.80 × 18.80 × 36.8 Å for the boron oxide system (Fig. 2). Initially, the atomic positions were relaxed at 300 K in the NVT ensemble for 0.1 ps. The obtained surfaces were then used for sputtering simulations in the NVE ensemble. The D atoms as projectiles were deposited as described in Section 2.2. The timestep was dynamically adjusted between 10⁻⁴ and 1 fs depending on the velocity of the projectile, with the step size being resampled every five steps such that no atom moved more than 0.05 Å between these five steps. The simulation proceeds until one of several halting conditions is met.

Fig. 2 Simulation cells used in the non-cumulative MD simulations. Left: B(100); right: B₂O₃(001).

3 Results

3.1 Sputtering yields of boron

In Fig. 3a the sputtering yields of boron from a pure boron surface under deuterium bombardment are shown together with the available literature data. The experimental data from Hechtl et al.²⁸ and computational data from Yamamura,¹³ along with the Eckstein curve²⁹ are included. Our values for different incident angles are indicated by colored triangles. It is seen that at perpendicular impact the yields increase as the impact energy increases. The statistical error of the sputtering yield was estimated using a bootstrapping procedure,³⁰ where random resampling with replacement was applied to the set of MD trajectories. The variance of the resampled averages was then used to calculate the uncertainty of the sputtering yield.

Fig. 3 Sputtering yields of boron from a boron surface under bombardment by deuterium (a) and boron atoms (b). Triangles represent our NNP simulations. 0° means perpendicular impact.

At lower impact energies our sputtering yields are smaller than the experimental and computational benchmark points. The large error bars for our data points are caused by the underlying Poisson statistics due to infrequent events. Increasing impact angles causes higher sputtering yields. At 60° the highest sputtering yields are found. They reach a plateau above the benchmark points, while for 30° a small maximum at 120 eV occurs with a decline for higher energies.

In addition to deuterium bombardment, we investigated the self-sputtering of boron surfaces using also the B–D trained NNP. The results, presented in Fig. 3b, show the sputtering yields of boron when bombarded by boron atoms. The agreement between the present results (blue triangles) and the former experimental²⁸ and computational¹³ data, as well as with Eckstein's curve is much better, probably due to better training data and the system involving only one element. These B self-sputtering simulations were only performed for perpendicular impact.

3.2 Sputtering of boron oxide

As mentioned, boron surfaces often serve the purpose of gathering oxygen, with boron oxides being formed. Simulations of the sputtering of a B₂O₃ surface by deuterium provides insight into the mechanism of sputtering of boron and oxygen atoms from them. The sputtering yields of boron atoms are illustrated in Fig. 4a, with different triangles representing the incident angles. The Eckstein curve for pure boron surfaces is included for comparison as it is also in the oxygen sputtering panel (Fig. 3b). The highest yield is observed at an impact angle of 60°. At an angle of 30°, the boron sputtering yield reaches a maximum at 200 eV. Below 150 eV, the calculated sputtering yields of boron from B₂O₃ surface are seen to be lower than those from pure boron surfaces. Oxygen atoms form strong bonds with boron, creating a robust surface that requires higher energy to remove boron atoms.

Fig. 4 Sputtering yields of boron atoms (a) and oxygen atoms (b). For D colored triangles indicate different incident angles. The Eckstein curve given for reference shows yields for pure boron surfaces only.

We have also analyzed oxygen sputtering from B₂O₃. Fig. 4b shows the yields of oxygen atoms. The yields of oxygen atoms increase monotonously with the energies and angles. Initially, the oxygen sputtering yields are lower than those of boron, but they become comparable at 200 eV. The differences between boron and oxygen can be attributed to the energy transfer factor:²⁹ γ = 4M₁M₂/(M₁ +M₂)² where M₁ and M₂ are the atomic masses of the projectile and target atom, respectively. Boron has a lower mass compared to the oxygen; accordingly, boron sputtering yields become higher when the surface contains both elements.

The results also show the significant effect of the impact angle on sputtering yields. At higher impact angles, the projectiles transfer some of their energy near the surface, increasing the probability of atom sputtering. This results in higher yields compared to normal incidence (0°). Cases with normal (=perpendicular) impact led to deeper penetration (details in Section 3.4) and less efficient sputtering. This is especially true at high energies, above 200 eV, where projectiles can penetrate deeper into the material and generate collision cascades that extend beyond the finite simulation cell. Consequently, part of the energy transfer and defect production is not fully captured, which can lead to underestimated sputtering yields and inaccurate defect statistics. Such effects can be more accurately described by cumulative simulations.³¹

3.3 Reflection, adsorption and retention probabilities of the projectile

While sputtering yields are readily measurable, for example by weight loss of the surface and an obvious sign of degradation, other processes take also place if D atoms interact with the surfaces. We analyzed D reflection probabilities, probabilities of D adsorption on the surface and retention of deuterium in the bulk. This is especially important as it can lead to the formation of boranes and boron hydrides and, via oxygen–hydrogen interactions, to the formation of hydroxyl groups. These chemical modifications can increase surface degradation of PFMs. We find that lower projectile energies, reflection and adsorption phenomena are dominant (Fig. 5), as is evident from the relatively low sputtering yields.

Fig. 5 Deuterium atoms after collision with a boron surface: probabilities of reflection (left), adsorption (middle), and retention (right).

At higher energies, the deuterium atoms possess sufficient kinetic energy to penetrate deeper into the samples, causing more collision cascades and can effectively distort the surface structure. Collision cascades can also create subsurface damage and vacancies, further contributing to the sputtering processes. This is true for different impact angles even though the probability of reflection and adsorption is higher at larger incident angles. At normal incidence (0°), the projectiles penetrate the surface more efficiently.

Differences between the two surfaces manifest themselves mainly in the reflection probabilities. Reflection of deuterium from boron oxide (Fig. S4 in SI) is lower than from pure boron surfaces. This is not unexpected, since with D implanted in B₂O₃, strong oxygen–deuterium interactions can appear. This in turn supplies energy for bond breaking and for chemical reactions, which can be observed indirectly by the increased sputtering yields in B₂O₃ compared to pure boron.

3.4 Energy profile of sputtered boron atoms

The energy distribution of boron atoms which have been sputtered from a boron surface by either D or B particles, or by D from a B₂O₃ surface, is illustrated in Fig. 6 where all outgoing energies and angles are sampled together. The shaded areas under the distributions represent their normalized yields. It turns out that these energy profiles closely follow a lognormal distribution³² which is represented by the solid lines. The positions of the maxima in the final fit to such a distribution approximate the surface binding energy, defined as the amount of energy required to remove an atom from the surface. The effective surface binding energy (E_SB) and indeed the energy distribution of sputtered atoms itself is important because it can be used to predict how far they spread out into the plasma and where they are possibly prone to surface re-deposition.

Fig. 6 Energy distributions of sputtered boron atoms. The shaded areas represent the normalized yield. The solid lines are lognormal distributions fitted to the data.

The effective surface binding energy for boron sputtered from deuterium impact was found to be 5.14 eV; the reference value of this quantity in a BCA simulation was, for example, taken as 5.73 eV.²⁹ Under self-bombardment conditions, the effective E_SB for boron was 5.15 eV, closely matching the value from deuterium bombardments. This consistency suggests that the underlying atomic-scale bonding characteristics of the boron surface remain the dominant factor influencing the effective E_SB, even for different projectile masses and energy transfer dynamics.

Similarly, sputtering of boron oxide by D can be used to extract effective surface binding energies of boron and oxygen, leading to E^B_SB = 4.46 eV and E^O_SB = 3.53 eV, respectively. This effective E_SB for boron is 10% less than a pure boron surface. This difference can be explained in terms of the energy transfer value as well. However, the joined B and O energy spectra (green line in Fig. 6) give a value of 4.36 eV, which means the contribution of oxygen to the total sputtering yield is negligible. Thus, from the viewpoint of sputtering, boron determines the effective surface binding energy of the oxide surface.

3.5 Chemical sputtering

Another important process possibly influencing the effective E_SB values is chemical sputtering. Both boron and boron oxide surfaces can strongly interact with hydrogen isotopes, which can lead to the ejection of hydrides and alter both the sputtering mechanism and surface composition. To investigate this, we analyzed the sputtering of B–B, B–D, B–O, and O–D molecular fragments. We find that the sputtering yields of B–D and O–D molecules each contribute less than 5% to the total sputtering yield. On the other hand, B–O sputtering becomes significant at low impact energies, reaching up to 20% of the boron sputtering yield. Fig. 7 shows energy-resolved sputtering yields of B, O and B–O sputtering events. Oxygen-induced chemical sputtering reaches a maximum around 150 eV, highlighting its role in modifying the sputtering behavior at intermediate energies. The relative importance of chemical sputtering is visible from Fig. 7 as it increases with increasing energy up to 150 eV. One could expect that chemical sputtering would be even more prominent in the scenario of cumulative simulations where adsorbed deuterium molecules from former irradiation events are present.

Fig. 7 Sputtering yields of B, O and B–O from B₂O₃ surfaces under deuterium bombardment.

4 Conclusions

Employing HDNNP-based MD simulations, we have investigated the sputtering of boron and boron oxide under deuterium bombardment, as well as boron self-sputtering on pure boron surfaces. Simulations with various impact angles and energies allowed us to analyze sputtering yields and related quantities. Sputtering from boron oxides shows element-specific behavior. Energy spectra of the sputtered atoms were computed, allowing estimation of the effective surface binding energies. The E_SB values for boron oxide were found to be lower than those of pure boron. However, the sputtering yields of boron from the oxide surface remain comparable to those from the pure boron surface. Additionally, on the oxide surface, chemical sputtering involving the formation of B–O radicals occurs. Along with the computed probability densities, these data are useful for predicting the behavior of boron–oxygen surface layers and assessing the long-term stability of boron-based plasma-facing materials.

Author contributions

S. S.: investigation, formal analysis, data curation, writing, visualisation M. P.: conceptualization, funding acquisition, methodology, administration T. M.: resources, writing K. H.: validation, resources, writing.

Conflicts of interest

There are no conflicts of interest.

Data availability

The data supporting this article have been included as part of the SI. See DOI: https://doi.org/10.1039/d5ra05589j.

Acknowledgements

The work has partially been carried out within the framework of the EUROfusion Consortium. MP received funding from the Euratom Research and Training Programme by Grant Agreement No. 101052200-EUROfusion. The views and opinions expressed herein do not necessarily reflect those of the European Commission. The computational results have been obtained using the HPC infrastructure LEO of the University of Innsbruck and the Vienna Scientific Clusters VSC.

References

1. S. Brezinsek, Plasma-surface interaction in the Be/W environment: conclusions drawn from the JET-ILW for ITER, J. Nucl. Mater., 2015, 463, 11–21.
2. R. Neu, Preparing the scientific basis for an all metal ITER, Plasma Phys. Contr. Fusion, 2011, 53, 124040.
3. T. Hino and T. Yamashina, Review on Plasma Facing Materials and Suitable Divertor Configuration of a Fusion Experimental Reactor, Mater. Trans., JIM, 1993, 34(11), 1106–1110.
4. J. Linke, J. Du, T. Loewenhoff, G. Pintsuk, B. Spilker, I. Steudel and M. Wirtz, Challenges for plasma-facing components in nuclear fusion, Matter Radiat. Extremes, 2019, 4, 056201.
5. S. J. Shin, L. B. B. Aji, J. H. Bae, A. M. Engwall, J. A. Hammons, G. V. Taylor, L. R. Sohngen, P. B. Mirkarimi and S. O. Kucheyev, Magnetron sputter deposition of boron carbide in Ne and Ar plasmas, J. Appl. Phys., 2024, 135, 085303.
6. R. A. Langley, Review of wall conditioning and wall materials for fusion research devices, MRS Bull., 1990, 15, 42–46.
7. K. Ohya and N. Saji, Erosion and Carbon Deposition on Boron Exposed to Plasmas, J. Plasma Fusion Res., 1999, 75(4), 416–425.
8. T. Hino, K. Iwamoto, Y. Hirohata, T. Yamashina, A. Sagara, N. Noda, N. Inoue, Y. Kubota, N. Natsir, O. Motojima, T. Matsuda, T. Sogabe, K. Kuroda and M. Yabe, Properties of boron coatings used as plasma facing material of fusion device, Thin Solid Films, 1994, 253, 518–521.
9. D. Zhao, S. Brezinsek, R. Yi, J. Oelmann, L. Cai, F. Wu, G. Sergienko, M. Rasinski, M. Mayer, C. P. Dhard and D. Naujoks, Investigation of boron distribution and material migration on the W7-X divertor by picosecond LIBS, Phys. Scr., 2022, 97(2), 024005.
10. J. Linke, H. Bolt, R. Doerner, H. Grübmeier, Y. Hirooka, H. Hoven, C. Mingam, H. Schulze, M. Seki, E. Wallura, T. Weber and J. Winter, Performance of boron/carbon first wall materials under fusion relevant conditions, J. Nucl. Mater., 1990, 176–177, 856–863.
11. M. Shoji, G. Kawamura, J. Romazanov, A. Kirschner, A. Eksaeva, D. Borodin, S. Masuzaki and S. Brezinsek, Boron transport simulation using the ERO2.0 code for real-time wall conditioning in the large helical device, Nucl. Mater. Energy, 2020, 25, 100853.
12. W. Eckstein and R. Preuss, New fit formulae for the sputtering yield, J. Nucl. Mater., 2003, 320, 209–213.
13. Y. Yamamura, K. Sakaoka and H. Tawara, Computer Simulation and Data Complilation of Sputtering Yield by Hydrogen Isotopes (1H+, 2D+,3T+) and Helium (4He+), Ion Impact from Monatomic Solids at Normal Incidence, NIFS-DATA.
14. K. Nordlund, C. Björkas, T. Ahlgren, A. Lasa and A. E. Sand, Multiscale modelling of plasma-wall interactions in fusion reactor conditions, J. Phys. D Appl. Phys., 2014, 47, 224018.
15. S. Shermukhamedov and M. Probst, Modelling the impact of argon atoms on a tungsten surface, Eur. Phys. J. D, 2022, 76, 169.
16. S. Shermukhamedov, L. Chen and R. Nazmutdinov, Sputtering and reflection from a beryllium surface: effects of hydrogen isotope mass, impact position and surface binding energy, Nucl. Fusion, 2022, 62, 066024.
17. E. Wolfgang, R. Dohmen, A. Mutzke and R. Schneider, SDTrimSP: A Monte-Carlo Code for Calculating Collision Phenomena in Randomized Targets, IPP Reports, 2007, vol. 12/3, p. 245.
18. F. J. Domínguez-Gutiérrez and P. S. Krstić, Chemical sputtering of boronized and oxidized carbon surfaces irradiated by low-energy deuterium atoms, J. Appl. Phys., 2017, 121, 215302.
19. L. Deng and J. Du, Development of boron oxide potentials for computer simulations of multicomponent oxide glasses, J. Am. Ceram. Soc., 2019, 102, 2482–2505.
20. J. Behler and M. Parrinello, Generalized neural-network representation of high-dimensional potential-energy surfaces, Phys. Rev. Lett., 2007, 98, 1–4.
21. A. Singraber, J. Behler and C. Dellago, Library-Based LAMMPS Implementation of High-Dimensional Neural Network Potentials, J. Chem. Theory Comput., 2019, 15, 1827–1840.
22. A. Singraber, T. Morawietz, J. Behler and C. Dellago, Parallel Multistream Training of High-Dimensional Neural Network Potentials, J. Chem. Theory Comput., 2019, 15, 3075–3092.
23. J. P. Perdew, K. Burke and M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 1996, 77, 3865–3868.
24. G. Kresse and J. Hafner, Ab initio molecular dynamics for liquid metals, Phys. Rev. B:Condens. Matter Mater. Phys., 1993, 47, 558–561.
25. G. Kresse and J. Hafner, Ab initio molecular-dynamics simulation of the liquid-metalamorphous- semiconductor transition in germanium, Phys. Rev. B:Condens. Matter Mater. Phys., 1994, 49, 14251–14269.
26. W. Hayami and S. Otani, Effect of surface energy on the growth of boron nanocrystals, J. Phys.:Conf. Ser., 2008, 176, 012017.
27. N. W. Assaf, M. De La Pierre, M. K. Altarawneh, M. W. Radny, Z. T. Jiang and B. Z. Dlugogorski, Structure, Stability, and (Non)Reactivity of the Low-Index Surfaces of Crystalline B2O3-I, J. Phys. Chem. C, 2017, 121, 11346–11352.
28. C. G.-R. E. Hechtl, A. Mazanec, W. Eckstein and J. Roth, Sputtering behavior of boron and boron carbide, J. Nucl. Mater., 1992, 196–198, 713–716.
29. R. Behrisch and W. Eckstein, Sputtering by Particle Bombardment: Experiments and Computer Calculations from Threshold to MeV energies, Springer, Berlin, 2007, vol. 110.
30. S. Gonçalves and H. White, Bootstrap standard error estimates for linear regression, J. Am. Stat. Assoc., 2005, 100, 970–979.
31. A. Lasa, K. Schmid and K. Nordlund, Modelling of W-Be mixed material sputtering under D irradiation, Phys. Scr., 2014, 2014(T159), 014059.
32. S. Shermukhamedov, T. Maihom and M. Probst, Modelling the Impact of Argon Atoms on a WO3 Surface by Molecular Dynamics Simulations, Molecules, 2024, 29(24), 5928.

---

Footnote

† This work is dedicated to Professor Jumras Limtrakul on the occasion of his 72^nd birthday, in recognition of his contributions to the establishment of computational chemistry in Thailand.

---