Detailed studies of the processes in low energy H irradiation of Li and Li-compound surfaces

In magnetic confinement fusion energy devices and reactors, the plasma–material interface (PMI) is a dynamic, evolving region of the material that is constantly eroded and re-deposited many times over. In a plasma environment, it is challenging to provide experimental observations of the physical and chemical processes occurring at the surfaces of plasma-facing components (PFCs). Both low and high atomic number candidates for PFCs in fusion devices have benefits and challenges. For instance, low atomic number materials (i.e., Li, Be, B, C) are getters of plasma impurities (e.g., O) but suffer from high sputtering yields from incident plasma particles relative to high atomic number materials. High atomic number materials, such as W, have superior thermomechanical properties, but the core plasma has a very low tolerance for radiative cooling from sputtered, high-Z atoms.^1–7 Recent research using lithium coatings evaporated on a variety of metallic and graphitic surfaces in many tokamaks around the world has provided evidence of the sensitive dependence of plasma behavior on lithiated PFCs.^8–12 Better understanding of the roles of surface chemistry and surface processes at the surfaces of Li PFCs during plasma exposure will be critical for developing a practical operational strategy for tokamak PFCs, e.g., wall conditioning techniques.^2,7

In these devices, the plasma behavior is affected by fuel recycling, which is defined as a ratio of the neutral hydrogen isotope flux from tokamak PFCs into the plasma to the flux of ionized hydrogen isotope flux from the plasma to the PFCs. Lithium is expected to be a low recycling material due to its ability to retain hydrogen from plasma.^13,14 LTX-β is a spherical tokamak designed to systematically study the effect of a low recycling boundary on plasma performance in a tokamak.^15,16 This is done by evaporatively coating the stainless-steel walls of LTX-β with lithium. It had long been theorized that such a low-recycling boundary would improve plasma performance by suppressing high pressure gradients in the edge plasma, thus reducing the edge localized modes (ELMs) and disruptions, which lower the damage of the PFCs, pollution of plasma, and tritium accumulation.¹⁷ Low recycling regimes reduce the influx of cold neutral hydrogen into the plasma, thereby keeping the plasma edge hot.¹⁸ Such an increase in the plasma edge temperature and a consequent reduction in temperature gradients from the core to the edge are expected to enable access to an operating regime with improved particle and energy confinement times,¹⁰ also reducing turbulence driven by temperature gradients. Specifically, reduction in or the elimination of temperature gradients from the core of the plasma to the edge can improve confinement because ion and electron temperature gradients are one of the main mechanisms through which heat escapes a magnetically confined plasma.¹⁹ The flat temperature profiles have recently been demonstrated with LTX, the predecessor of LTX-β.⁹ In addition, in another tokamak, NSTX, Li evaporation decreased the H-mode access power threshold, increased the stored energy, and allowed longer plasma discharges when compared to plasma discharges with no Li conditioning.²⁰ These improvements in plasma performance were correlated with a reduction in plasma impurities and a reduction in fuel recycling and attributed to the formation of Li–O–D complexes.²¹ Other fusion devices that are investigating the effect of Li PFCs on plasma performance include, but are not limited to, EAST,²² HIDRA,²³ ADITYA-Upgrade,²⁴ and RFX-mod.²⁵ 

A novel feature of the LTX-β upgrade is the Sample Exposure Probe (SEP),¹² which includes a vacuum suitcase, which can transport samples under high vacuum to a surface analysis station to conduct the x-ray photoelectron spectroscopy (XPS) analysis of a PFC witness sample on the low field side midplane of LTX-β. The SEP-enabled analysis demonstrated that lithium oxide (Li₂O) grows on top of evaporatively deposited lithium coatings in LTX-β.^11,12 Using the SEP and XPS, various compounds that form on evaporated lithium thin film coatings were successfully identified. In LTX-β, these coatings quickly oxidize to form Li₂O and LiOH. This oxidation occurs due to the presence of water in the background gas of the vacuum vessel, as well as owing to the presence of oxygen from oxide components and oxidized metal PFCs. On the other hand, lithium hydride (LiH) can also be formed at a lithium surface irradiated by hydrogen and its isotopes.¹³ 

It is difficult to study the Li and Li-compound dynamics using theoretical calculations because of Li polarizing features when interacting with other elements. These features are most transparently represented by the quantity called electronegativity, i.e., the chemical property of an element that defines its tendency to attract electrons. The Pauling electronegativity²⁶ of Li is exceptionally low (0.98), one of the lowest in nature, and, therefore, much lower than the other elements in materials of interest in NSTX and LTX-β. The primary impurity for NSTX-U plasma is C, followed by others such as B, O, Fe, Ni, and Li.²⁷ In LTX-beta, only three impurities exist at relatively low levels, Li, C, and O. The difference is primarily due to the near complete coverage of LTX-beta PFCs with Li.^28,11,29 The presence of electrically polarized atoms requires the treatment of the long-range Coulomb and multipole forces, in addition to the short-range potentials usually used in molecular dynamics. This implies a large radius of interactions and, therefore, many “neighboring” atoms that simultaneously interact. This significantly increases computational power requirements. On the other hand, the electronegativities of oxygen are exceptionally high (3.44), while the hydrogen, tungsten, and carbon electronegativities are in the middle (2.2, 2.36, and 2.55, respectively). The low electronegativity of lithium in comparison with these other elements causes its positive polarization in mixtures forming fusion materials, while the other element(s) in these materials becomes negatively charged (the whole system stays neutral). Typically, Li and O strongly polarize each other, creating an ionic solid, while the sign and strength of the polarization of the other atoms in the materials, including hydrogens that impact surfaces from the plasma, depend on the prevailing concentrations as well as on the bulk or surface distributions of Li, O, and possibly other elements.

Classical molecular dynamics (CMD) is not capable of accounting for these charge changes in the atoms and the effects of atomic coordinate-dependent polarizations since these are the quantum mechanical effects. Our earlier results based on quantum-classical molecular dynamics (QCMD) simulations,^30–32 which were validated by XPS measurements,^31–34 demonstrated that the oxygen content in the surface layers of carbon (graphitic) PFCs was the main driver for deuterium uptake in lithiated carbon surfaces. Future fusion reactors will likely not use graphitic PFCs, but rather high-Z metal substrates possibly coated by Li. Consequently, understanding the dynamic surface chemistry of Li films and coatings converted into Li₂O, LiH, and LiOH, which are then subjected to hydrogen irradiation by plasma will be of use for high-duty cycle fusion reactor operation. However, QCMD calculations approximate a solution for the quantum mechanical eigenvalue problem for thousands of electrons at each time step requiring hundreds and thousands of times more computing time than CMD, resulting in a formidable computational problem when thousands of time steps, thousands of impacts trajectories, and impact energies larger than a few eV are involved.

To use the faster and simpler CMD method with pre-parameterized short-range bond-order potentials, a semi-empirical method like the Electronegativity Equalization Method (EEM)^35,36 is applied at each step (or after a small number of steps) for the calculation of the dynamics of atomic charges. The Reactive Force Field (ReaxFF) potentials^37–40 have these capabilities and can provide results for large surface slabs and a range of impact energies with about two orders of magnitude shorter computational times than QCMD.

In the present report, we describe our simulation results for physical and chemical processes at Li₂O and LiOH surfaces induced by low-energy H, D, or T impacts. As seen in previous density functional theory (DFT) calculations⁴¹ and demonstrated with CMD using ReaxFF calculations and comparisons with experimental results,⁴² we used the same Reax potentials for impacts of various hydrogen isotopes, adapting only to their mass. Amorphous oxide surfaces were created by the amorphization of either Li₂O or LiOH molecular crystal slabs, thus creating mixtures of atoms with the atomic ratio and mass density corresponding to these crystals. The used CMD method and the development of the needed ReaxFF potentials are described in Sec. II. We compare our current results for Li₂O and LiOH surfaces with earlier results for LiH and Li surfaces in Sec. III. This comparison provides insight into the effects of the various mixtures of these chemically reactive elements on the surface dynamic response when irradiated by low-energy hydrogen atoms. Finally, in Sec. IV, we give our conclusions.

As discussed in Introduction, a distinct consideration of the Li₂O and LiOH surfaces arises from the electronegativities of lithium and oxygen; i.e., they easily charge in contact with other elements, partially giving or receiving the charges in their electron cloud to other atoms in their vicinity, polarizing the mixed material. This causes long-range, Coulomb-type interactions between the surface constituents, which depend on the instantaneous coordinates of the atoms. In addition, the study of chemical sputtering, retention, reflection, and other processes for an impact energy range of 1–100 eV is computationally intensive, not only because the calculation has to be repeated for various energies and many trajectories but also because particles incident at higher energies require a larger computational cell, shown in Fig. 1, i.e., a larger number of atoms. This is so because the maximal penetration depth of an impact particle increases with its initial kinetic energy, especially for small incident angles θ in Fig. 1. The cell in the z-direction must be deep enough to avoid artificial reflections from the lower vacuum boundary. The surface, 2D boundary conditions are provided at the sides of the computational cell, i.e., in the X–Y directions. For these reasons, the combination of CMD and EEM was chosen for these studies. The latter is a semi-empirical method that applies a set of pre-calibrated material- and coordinate-dependent parameters to estimate the change of charges of each atom as their coordinates change during the collision cascade of an impact atom on the surface. The EEM is implemented in the Large Scale Atomic/Molecular Massively Parallel Simulator (LAMMPS)⁴³ with the Reactive Force Field (ReaxFF) method and potentials developed by van Duin and co-workers.^37–40 The EEM slows computation by an order of magnitude, but this is still orders of magnitude faster than QCMD, even when an approximation to density functional theory (DFT) is used for the quantum mechanical part. ReaxFF is a bond order potential with proven capability to treat correctly chemical processes in the material of interest. Since we are considering the impacts of H, D, and T atoms in a low energy (1–100 eV) range, chemical sputtering is the main process of erosion. Consequently, we need a detailed description of the hydrogenated, lithiated, and oxidized surface chemistry. The ReaxFF method, which uses a bond order/bond distance relation, is well suited for this requirement. We have previously verified this approach by comparison with some relevant QCMD results, obtaining good qualitative and even quantitative agreements.^14,31,34

We will discuss the method in more detail with the example of the most challenging Li₂O surface target. The method, similar in all details to the one described here, was also applied to the present study of the LiOH surface and in the past to the Li and LiH surfaces.^44–46 The initial Li₂O slab has a fluorite crystal structure with a lattice constant of 4.69 Å. The crystalline supercell of 10 000 li₂O units, after 3D energy minimization (adaptation to the ReaxFF potential), thermalization to 300 K, and relaxation to the surface boundary conditions and 2D periodicity, occupied a box of dimensions 4.7 × 4.7 × 11.7 nm³ and mass density 1.93 g/cm³, which is within 3.5% of the experimental value of 2 g/cm³. The amorphization of the energy-optimized crystalline supercell was done using the usual procedure by repeated fast heating to 10 000 K and slow annealing to 300 K with 3D periodicity, followed by relaxation to the surface boundary conditions and 2D periodicity in the x- and y-directions. The depth of this box was sufficient for slowing down and the thermalization of the vast majority of the retained impact atoms. After amorphization, we note that the LiOH slab occupied a box of dimensions 5.6 × 5.6 × 12.2 nm³, with a total of 29 735 atoms, where there was approximately one-third of this number of each of the Li, O, and H atoms. The notation a:M and c:M will be used for amorphous and crystalline slabs of M molecules, respectively.

The methods for the target slab creation and irradiation with H are like those described in Refs. 44–46, so only the necessary details are described here. When an atom of initial energy E and incident angle θ impinges on the slab, the choice of the MD time step Δt depends on the impact velocity, which is a function of the impact energy and mass of the impact atom. A time step that is too long might miss the interactions on the surface and tends to distort energy conservation. On the other hand, a time step that is too short would require a huge number of steps to describe the desired process. The best compromise of Δt = 0.25 fs for all energies 1 ≤ E ≤ 100 eV was found and applied in the calculations reported here. The number of time steps for a single trajectory was chosen for each impact energy and angle based on the requirement that the total time is sufficient for the full realization of the considered processes, reflection and sputtering. Since the total number of steps N was between 4 000 and 10 000, this calculation captures phenomena at the ps timescale (1–2.5 ps). Uniformly scanning the surface area of 3.5 × 3.5 nm² by N = 1 022 repeated impacts was sufficient to obtain a maximal standard error (MSE)³⁰ less than 20% for most calculated quantities (probabilities of reflection, retention, sputtering) when varying the parameters of the system (impact energy and angle, impact particle mass, and surface structure). The error bars in all figures herein are shown with MSE uncertainties. The approximate MSE values for all data are calculated as $p / n$⁠, where n is the size of a desired outcome sample. An outcome probability is calculated by p = n/N_s, where N_s is the number of successful trajectories among the N impacts. A successful trajectory is one that followed an impact particle that either traveled all N steps in the computation box or went over the box border and was removed. Therefore, N_s≤ N. The total length of the computation box was limited to about 15.5 nm in the z-direction, of which about 12 nm was the initial range of the surface atoms, the rest being the domain of motion of reflected and sputtered particles in vacuum. In many cases, the MSEs were within the size of the symbols used to plot the data and may not be visible in the figures. In some cases, when p was small, e.g., for sputtering at the lowest impact energies E, the uncertainty margins were larger than 20%. Each impact was done independently on the same initial surface, changing N times only the point of origin for the incident particle over the plane, approximately 1 nm above the surface. In the process, each atom, i.e., impact atom or surface atoms, was labeled during calculations. Hence, the origin of ejected species (sputtered or reflected) was clearly distinguished. These calculations qualify well for the so-called “embarrassing parallelization,” where each computer node in a supercomputer is used for one trajectory. The reflection, retention, and sputtering probabilities were always obtained “on the same footing.” The total computing time for obtaining and analyzing results in Sec. III was close to 7 × 10⁶ CPU hours.

Although the ratio of the total number of Li and O atoms in the a:Li₂O slab in Fig. 1 is two, the distribution of Li and O atoms is not entirely homogeneous. Figure 2(a) shows that the number of Li to O atoms in the first layer exceeds the expected ratio of 2 by more than three times, reaching two and oscillating around two after a depth of 2 Å. This is a product of nonideal amorphization with open top surface boundary, as well as the fact that the topmost atomic layer of the c:Li₂O slab in Fig. 1 is formed of only Li atoms. We note that, due to polarization, each Li and O atom carries a charge of about +0.5 and −1, respectively. Therefore, inhomogeneities in the number of atoms in the layers result in corresponding inhomogeneities in the charge of the layers as shown in Fig. 2(b). Consequently, in addition to the covalent potential field that each impact atom experiences while penetrating and traveling through the surface, the incident H atom is generally also charged while penetrating through the a:Li₂O surface and, thus, subject to the variable Coulomb charge as shown in Fig. 2(b). However, in the case of a:LiOH, the distribution of the number of Li + H atoms to the number of O atoms is oscillating also about 2 [Fig. 2(c)] but without preferential accumulation of atoms in the top players. Consequently, the top layer charges in Fig. 2(d) oscillate with significantly smaller amplitudes. The surface inhomogeneities presented in Fig. 2 can influence the probabilities of surface processes, as will be discussed in Sec. III.

Here, x_i,QM is the energy calculated by the QM method, x_i,ReaxFF is the energy calculated by ReaxFF parameters, and σ is the accuracy parameter (inverse of training weight) for the i^th data. We used the standalone ReaxFF code to train the force field parameters with the single-parameter search optimization scheme. In the single-parameter search optimization scheme, one force field parameter at a time is varied by a small amount to minimize the overall error [Eq. (1)] in the training data. This method is further described by van Duin et al.,⁴⁹ and the standalone ReaxFF code can be made available upon request.⁵⁰ 

Figure 3 shows the comparison between DFT and ReaxFF energies. Figures 3(a)–3(c) show the equations of state (EOS) for Li₂O_Fm3m (antifluorite), LiOH_C2me (orthorhombic), and LiOH_p4nm (tetragonal) lattices, while Fig. 3(d) shows the heats of formation for these lattices with Li_BCC, O₂ (gas), and H₂ (gas) as zero-energy references. Figures 3(e) and 3(f) show the pair potential between a Li atom and the O atom of OH₂ and OH, also referred to as a Li–OH₂/Li–OH bond scan. Figures 3(g)–3(i) show the atom pair potential between Li–H,⁴⁶ O–H, and Li–O atom pairs, respectively. The comparison is shown for energies greater than 100 eV as the close-range interactions are crucial for the simulation of a high-energy H/D/T atom impact on the amorphous Li₂O slab. We train the strain-energy relation for Li₂O and LiOH crystals as per the training methodology of Shin et al.⁵¹ The stiffness constants for the strain-energy relationship were obtained from the Materials Project Database.⁵² The detailed strain-energy relations for training Li–O parameters are given in Sec. S1 in the supplementary material. The target data for training these pair potentials were defined as a combination of DFT and Ziegler-Biersack-Littmark (ZBL) energies where ZBL energies were used for distances below 0.6 Å, as explained in Sec. S2 in the supplementary material. For the equation of states, heats of formation energies, Li–OH pair potential, and Li–OH₂ pair potential, the DFT data were obtained from previously published work.³⁹ However, for gas-phase DFT (Li–H, Li–O, and O–H pair potentials), we used Amsterdam density functional (ADF) code within the Amsterdam Modelling Suite.⁵³ We performed spin polarized calculations with generalized gradient approximation (GGA) based Perdew–Burke–Ernzerhof (PBE) exchange-correlation functionals in combination with the TZP basis set (triple Z with 2 polarization functions).

To stress the amorphous nature of our surface slabs, we consistently use the name of a surface compound with prefix “a:.” Thus, amorphous Li₂O is denoted by a:Li₂O, amorphous LiOH is written as a:LiOH, etc. A variation in an isotope of the hydrogen atom, which impacts a surface, means the same isotope variation in hydrogen in the surface, if contained in the surface compound. For example, H + a:LiOH, D + a:LiOD, and T + a:LiOT are the only allowed combinations.

The probabilities of reflection (P_R) of the incident H/D/T particle were calculated for each impact energy and angle of incidence as a ratio of the number of reflected trajectories (N_R) and the total number of successful trajectories (N_S), with P_R= N_R/N_S. Then, the probabilities of retention per H/D/T impact (P_ret) are in all cases defined by particle conservation, i.e., P_ret= 1-P_R. Since reflection in most cases has a smaller probability and larger variability than retention, we will mainly report the reflection probabilities.

The reflection probability of H from amorphous Li₂O (a:Li₂O) (Fig. 4) is small at 1 eV, resulting in retention close to 0.90 for all incident angles θ. However, for 1 eV at the a:LiOH surface, reflection is significantly higher for H impact at θ = 45° and smaller for D impact at θ = 0°, resulting in retention of about 0.75 and 0.98, respectively. Irrespective of θ, and unlike a:LiOH, a:LiH,³¹ and a:Li²⁹ surfaces, reflection probability from a:Li₂O has a sharp peak value at about 5 eV, and at larger impact energies, it steadily decreases, resulting in increased retention of H. For θ = 0° and 45°, the retention probability reaches 0.95 for a 100 eV impact. For near-grazing incidence at θ = 85°, retention stays between 0.4 and 0.6 for E ≥ 10 eV. The P_R of H with θ = 0° from a:Li, a:LiOH, and a:LiH slowly decreases with energy, staying almost parallel to each other, and after reaching the maximums at about E = 25 eV, ensuing retentions of about 0.75, 0.82, and 0.91, respectively. These all stay above the P_R values from a:Li₂O at θ = 0 for E larger than 40 eV. However, at θ = 85°, the reflection probabilities of incident H from a:LiH, a:LiOH, and a:Li exceed those from a:Li₂O for E larger than 10 eV, with the difference increasing with increasing impact energy.

Concerning the isotopic effects of reflection (Fig. 5), at low energies, the reflection probability decreases with isotope mass. That trend reverses at projectile energies of 25 eV. At 100 eV impact energy, the reflection probability of tritium of 0.1 is about three times larger than that of hydrogen, 0.037, resulting in lower retention of tritium than of hydrogen and deuterium. This is not the case with reflection from LiOH, which decreases with the mass of the impacting projectile, as shown by the D-curve in Fig. 4.

The distributions of the scattered energies of the reflected atoms, averaged over the successful trajectories, show a linear dependence on the impact energy for both a:Li₂O [Fig. 6(a)] and a:LiOH [Fig. 6(b)] surfaces. This is also a case with a:LiH⁴⁶ and a:Li⁴⁴ surfaces. The slopes of these dependencies increase with impact angle. The slopes decrease with increasing isotopic mass, as seen for the a:Li₂O surface at θ = 45° and a:LiOH surface at θ = 0°. At almost grazing incidence (85°), the reflection energies are close to elastic scattering, with the ratio of E_av/E for impact H close to 95% (a:Li2O and a:LiH) and 85% (a:LiOH). This is also indicated by scattering angles of reflected particles, averaging over all successful trajectories at a given energy, shown for a:Li2O in Fig. 7. Having a very small component of velocity toward the surface in the case of almost grazing incidence means that impacting H will not penetrate and interact with atoms in the surface (see Fig. 9, for example), so its motion will be mainly elastic. This observation also corresponds to the Rutherford scattering formula for charged particles, which indicates specular scattering angles for elastic scattering. The average reflection angles show specular behavior at θ = 85° and E > 5 eV, while in the other cases, the average reflection angles primarily fall in the range of 45°–60°. Like the case of average reflection energies, the average reflection angles increase when the isotopic mass of the impact atom decreases.

The distributions of the deepest reach of impact particles inside the solid before reflection are of interest for estimating the depth of a lithium compound film for a considered impact energy below which reflection and sputtering mechanisms depend on the reflection from the substrate.⁴⁴ These are shown in Fig. 8 for a characteristic example of 50 eV impact energy and 45° impact angle for all three isotopes of H.

The reflection depths for T, D, and H at the 95%-mark [Fig. 8(a)] are 0.7, 1, and 1.35 nm, respectively. These are counted from the top layer in Fig. 1 (“0” at reflection depths in Fig. 8), the position of which is taken as an average over a bumpy top surface. The reflection depth values are about a factor of two smaller than the H(D) reflection depth in an amorphous LiH (LiD) surface⁴⁶ and five to six times smaller than the reflection depth of H in a:Li.⁴⁴ This indicates a larger energy loss of H atoms and isotopes in Li₂O, probably due to the larger mass density of Li₂O (about 2 g/cm³) but also due to the larger chemical and polarization forces, which need to be studied in the future. Since the closer values of masses of the impact and surface atoms mean higher kinetic energy transfer in their collisions, the trend of a decrease in the reflection depth with increasing isotopic mass is due to the more efficient kinetic energy transfer to the surface atoms.

A good manifestation of the reflection depths is their average values over all reflected trajectories for a given impact energy and angle. As shown in Fig. 9, the change in the mean reflection depth (MRD) for H impacts at a:Li₂O for an incident angle of 0° up to an impact energy of 75 eV is about −0.09 Å/eV, while for a 45° angle of incidence, this value is −0.07 Å/eV, calculated as the average slope of the curves. For the a:LiOH surface, the slope of MRD as a function of E is −0.18 Å/eV, about a factor of 2 larger than the corresponding quantity for a:Li₂O. However, for almost grazing incidence of H, the slope is negligible over the entire energy range; i.e., the average reflection depth does not change with impact energy and stays close to zero for both a:LiOH and Li₂O surfaces. Note that the average values, in all cases, capture less than 50% of the reflected atom distribution due to the asymmetry in the distributions of the curves, as seen in Fig. 8(b).

We now discuss the distinct behavior of the reflection probabilities as functions of the impact energy in Figs. 4 and 5. We believe that a significant role in the appearance of the peaked structure at 5 eV is played by the inhomogeneities in the content of Li and O atoms in the top surface layers, as shown in Fig. 2(a). While the ratio of the number of Li to O atoms is approximately 2 in the bulk, the first Angstrom of the surface that faces vacuum has a ratio of about 7, which reduces to 2 at 2 Å depth. Since Li atoms are positively charged (∼0.5 a.u. of charge), and oxygen atoms are charged negatively (∼−1 a.u.), the whole system is neutral. Consequently, the first layer is positively charged, and the third layer becomes negatively charged. The impact atom becomes slightly negatively charged due to interaction with the surface (−0.18 a.u. at the top layers; −0.3 a.u. in bulk). For a 5 eV impact energy, none of the H, D, nor T atoms that reflect off the surface penetrate the solid deeper than 2 Å (Fig. 9), and it seems the charged top surface layers play a role in the peak of the reflection. The role of the surface layer decreases with increasing impact energy as the H atoms penetrate deeper into the bulk, and so the reflection probability decreases (Figs. 4 and 5). In the case of a:LiOH, both Li and H in the surface charge positively, while O charges negatively. However, the number of Li + H atoms at the top layers of the a:LiOH surface is not significantly different than the number of O atoms, as shown in Fig. 2(c). Consequently, there are no significant total charges at the top layers [Fig. 2(d)] and the reflection probability of H does not show a peaked trend at lower energies as is a case of a:Li₂O.

The question arises why the reflection probability of H from a:Li₂O is smaller than that from a:LiOH, as seen in Fig. 4. Since the reflection and retention probabilities complement to 1, the question can be reformulated to “Why is probability of hydrogen to be absorbed, i.e., retained at impact energies above 25 eV is bigger in a:Li2O than in a:LiOH?” Retention of neutral hydrogen atoms from a:Li2O and a:LiOH is a complex process, which depends on several factors. The average distance of atoms in a:Li2O is smaller (0.205 nm) than that in a:LiOH (0.234 nm), which gives the possibility to H atoms to move easier through a:LiOH without interacting, having a larger average reflection depth (Fig. 9) and larger retained penetration depth, as will be shown in Sec. III C. More chances for H to interact means also increased the possibility for bonding with surface atoms when it slows down by collision cascade to almost thermal energies. Since impact H in the case of a:Li₂O gets negatively charged even in bulk, it can possibly bond only to positively charged Li atom complexes. However, H atoms impacting to a:LiOH are dominantly positively charged in bulk and can only bond to twice less O atoms in a:LiOH then there are Li atoms in a:Li₂O.

On the other hand, D and T incident atoms lose more energy in collisions with the heavier Li atoms than H atoms, causing more spread in the reflection depth distributions and a smaller probability of reflection with isotopic mass, as seen in Fig. 5 for impacts at 5 eV. The same behavior of decreasing energy of the reflected atoms with increasing isotopic mass is shown in Fig. 6(a). On the other hand, when the energy of the incident atom increases, the role of the first layer in the reflection decreases, implying a change of the trend with isotopic mass. As seen in Fig. 5, for E ≥ 50 eV, the reflection probability increases with isotopic mass, likely caused by a smaller depth of penetration with increasing isotopic mass. Namely, the peaks in the reflection depth close to the surface are shifted deeper in bulk with the decrease in the isotopic mass, as illustrated in Fig. 8(b).

The sputtering probability per incident atom (P_sp) is calculated for each impact energy and angle of incidence as a ratio of the number of sputtered trajectories (n_sp) and the total number of completed (successful) trajectories (N_S), i.e., P_sp= n_sp/N_S. When the sputtering probability is calculated for Li, O, or H (separately), each value of n_sp is replaced by the number of sputtered Li, O, or H atoms. Li and H are sputtered in the form of atoms and diatomics. When sputtered Li or H appears in the form of molecules, the sputtering probability is calculated by counting all ejected Li or H atoms. For all studied cases, the dominant sputtered particles were H or Li atoms, with the rare appearance of sputtered O atoms.

We note that the sputtered particles from a:Li and a:Li₂O are mainly Li atoms. A comparison of the sputtered atoms of Li from a:Li₂O, a:LiOH, a:LiH, and a:Li, when irradiated by H, is shown in Fig. 10. At 0° angle of incidence, the sputtering probabilities of Li from Li₂O are below 4% and are between those from a:Li (largest) and from a:LiH (smallest). The sputtering of Li from a:LiOH is the lowest, below 0.8%. The sputtering of H₂ molecules from a:LiH reaches 1%, while the sputtering of oxygen from a:LiOH is 0.25% at most. The total sputtering probability from a:LiH and a:LiOH is by far dominated by hydrogen atoms and is reaching 10% and 5% at their respective maximum impact energies of E = 25 eV and E = 75 eV, respectively. At 85° angle of incidence, Li sputtering probabilities from a:Li₂O are below 6% and are close to those of H sputtering from a:LiH and total sputtering from a:LiOH, all being much smaller than Li sputtering from a:Li, which is close to 20% at higher impact energies.

Isotopic effects for the total sputtering probability from a:Li₂O as a function of impact energy of H, D, and T are shown in Fig. 11. As one would expect, the total sputtering probability is lowest for the impact of H, and it increases for D and T impacts. The total sputtering probabilities reach similar values for D and T, consistent with the small difference in the reflection probabilities of D and T, as shown in Fig. 5.

As seen in Fig. 12(a), the average sputtered energies of Li atoms from a:Li₂O, similar to sputtering from a:LiH in Ref. 46, are not strong functions of impact energy and are within a range of a few eV. The exception is sputtering with H at an 85° angle of incidence, which exhibits an almost linear increase in the average sputtered energy of Li atoms with H impact energy and reaches close to 10 eV for 100 eV impacts. This behavior is similar to the sputtering of H atoms from a:LiH.⁴⁶ However, kinetic energies of the sputtered particles from a:LiOH at Fig. 12(b) are dominated by H atoms and reach close to 6 eV for H impact and about 5 eV for D impact at θ = 0° and E = 100 eV. The sputtered energies are even smaller θ = 85°, reaching almost 3 eV.

The depth in the surface from which sputtered particles originate for various impact energies and angles could define the depth to which the erosion of the surface by fusion plasma is taking place and is potentially an important quantity. Similarly, like in the case of the average reflection depth in Fig. 9, we calculate the mean sputtering depth (MSD), averaging the sputtering depth overall sputtered trajectories.

These are presented in Fig. 13 for various impact energies. The sputtering depth is expected to be the largest for impact angle θ = 0°. As shown in the figure, the absolute value |MSD| is smallest for a:Li₂O, being less than 2.5 Å, in the whole range of energies, while it is the largest for a:LiOH, reaching up to 8 Å. For a:Li and a:LiH, maximal |MSD| values are about 5 and 4 Å, at E = 100 eV. As expected, for almost grazing incidence, |MSD| is not larger than 4 Å for a:LiOH. It is notable that the absolute sputtering depth of deuterium at a:LiOH for θ = 0° is smaller than the one of hydrogen within the margins of the maximal standard error in Fig. 13. Like in the case of the mean reflection depth, |MSD| is close to the median of a distribution of sputtering depths, similar to those in Fig. 8(b). Thus, the maximal sputtering depth is about factor <2 larger than the data reported in Fig. 13.

The principal mechanisms of sputtering in the range of impact energies below 100 eV for close to perpendicular impacts of H to the a:Li surface, established in Ref. 44, can be applied also to the present work. A reflected atom, on its way toward the top layer, collides with the atoms of the surface and possibly breaks off a bond, and if the transferred energy is sufficient, a surface particle is sputtered. The “reflected” particle may leave the surface if it did not lose its energy or stay on the surface. Although it is difficult to show a quantitative relation between sputtering and reflection probabilities, we can safely state that in all studied cases, we saw a positive correlation between the two.

It is not necessary to separately discuss in detail the probabilities of retention for various cases of impact particles since the probability of retention is complementary to the probability of reflection, i.e., P_ret= 1-P_R. The information on P_R has already been discussed in Sec. III A.

Information on the sizes needed for the a:Li₂O slab used for MD calculations is contained in the distributions of the retained particles. These distributions could also serve as initial conditions for the calculation of the diffusion of the accumulated hydrogen upon impact. The maximum reflection depth limits, as discussed in Figs. 8 and 9, do not exceed 3 and 4 nm for a:Li₂O and a:LiOH (for H at E = 100 eV and θ = 0°), indicating that the surface slabs >12nm thick, which were used in this work, are more than sufficient for studying the reflection and sputtering under the considered conditions.

The distributions of the retained particles (Fig. 14) were studied in the a:Li₂O and a:LiOH slabs described in Sec. II. While all cases show a maximum in the retained particle distributions at less than 3.5 nm depth, the widths and heights differ significantly between cases. This is best visualized by using the normalized cumulative counts of retained atoms, showing that 90% of incident H atoms at 100 eV and 0° angle of incidence are retained in a:Li₂O at depths less than 7 nm and 90% of incident D and T atoms at 100 eV and 45° angle of incidence are retained at depths less than 4.5 and 4 nm, respectively. This number is bigger for H impact at a:LiOH, amounting to about 9 nm. About 90% of incident H atoms at 45° angle of incidence are retained at depths less than 6 nm in a:Li₂O, while that number is 2 nm for an 85° angle of incidence. The values for a:Li₂O are very similar to those for the retention distributions of a:LiH in Ref. 46. The tails of the distributions extend no more than 12 nm. As a comparison, for an a:Li surface, 90% of incident H atoms at 100 eV and 0° angle of incidence are retained at depths less than 15 nm with the tail of the distribution extending to 22 nm, while 90% of incident D atoms at 100 eV and 0° angle of incidence are retained at depths less than 12.5 nm with the tail of the distribution extending to 20 nm.⁴⁴ These significantly smaller depths reached by the retained atoms in a:Li₂O (and a:LiH⁴⁶) compared to pure Li are partially a consequence of the higher mass density of Li₂O (2 g/cm³) and LiH (0.78 g/cm³) compared to that of Li of 0.5 g/cm³.

In LTX-β and other fusion devices that have lithium-coated plasma–material interfaces, the coatings could quickly oxidize to form Li₂O and LiOH. The oxidation occurs due to the presence of oxygen from oxide components and oxidized metal PFCs, as well as due to the presence of water in the background gas of the vacuum vessel. Our computational model, an amorphous mixture of Li and O obtained by the amorphization of either crystalline Li₂O or crystalline LiOH, show distinctive features upon irradiation by hydrogen and its isotopes compared to Li and LiH. Both Li and O readily polarize each other and can polarize both the incident and the surface hydrogen due to the large difference in their electronegativities. We have included the charging effects into CMD using improved ReaxFF/EEM potentials, pre-calibrated by quantum-mechanical calculations (DFT) to consider the changes in the dynamic polarization of included atoms caused by changes in their coordinates. We investigated a range of collision parameters: impact energies of 1–100 eV at angles of incidence of 0° to −85° with the target temperature at 300 K.

Reflection of H from a:Li₂O has distinctive features, with a prominent peak in the reflection probability at low energies (<5 eV), which then strongly decreases toward higher energies, resulting in the retention probability exceeding 0.95 at 100 eV. This peak is not obtained for the reflection of H from a:LiOH, a:LiH, and a:Li, which show almost flat dependencies at higher impact energies, reaching retention probabilities of 0.90s for a:LiH, 0.80s for a:LiOH, and 0.70s for a:Li with normal H incidence. At lower energies (<25 eV), the reflection probability decreases with the mass of the impacting isotopes of H. However, at higher energies, this trend reverses in the case of a:Li₂O. The specific effects with a:Li₂O are ascribed to inhomogeneities in the topmost surface layer, which is dominated by Li atoms. Such inhomogeneities at the surface are expected and must be treated with care in modeling at lower incident energies. When the impact energy of the incident H increases, the role of these top-surface effects in reflection gradually diminishes. The kinetic energy of the reflected particles increases linearly with impact energy and slightly decreases with the increase in their isotopic mass. However, reflection becomes close to elastic at almost grazing angle of incidence, showing specular angular features. In these cases, retention probabilities of all Li compounds studied here are between 0.3s (a:Li, a:LiOH) and 0.2s (a:LiH), with the exception of a:Li₂O, which has retention in range of 0.6s for E = 100 eV. The mean reflection depth does not exceed 2 nm at any impact angle and decreases with increasing isotopic mass while strongly depending on the incidence angle in the considered energy range.

The incident hydrogen atoms charge slightly negatively in the polarized mixture of Li and O of a:Li₂O, which has a counterintuitive consequence: not bonding of retained hydrogen to negatively charged oxygen atoms, but rather settling down at the interstitial local-minima positions, formed by the groups of lithium and oxygen atoms. Consequently, no sputtered oxygen atoms are seen. On the other hand, retained incident hydrogen in an amorphous LiOH surface could be both positively and negatively charged, allowing various types of interactions with the surface atoms, including dominant sputtering of H atoms and some of Li and O atoms. The probability of sputtering of Li at normal incidence does not exceed 1% for any of the considered lithium surface compounds (except for a:Li, which reaches a few %). The sputtering of H atoms from a:LiOH and a:LiH achieves a value of 5% and 10% at normal incidence, respectively. However, at θ = 85°, the sputtering probabilities of Li from a:Li₂O, a: LiOH, and a:LiH are larger, up to about 6%, while for a:Li, the sputtering probability is over 20% for the highest impact energies. Unlike the reflection, the sputtering probability increases with isotope mass in the whole energy range. The energies of the sputtered Li from a:Li₂O range between 2–3 eV for E = 10 eV and 5–10 eV for E = 100 eV. Interestingly, the sputtered energy of Li from a:Li₂O increases while the sputtering of H from a:LiOH decreases with the increase in the angle of incidence and incident isotopic mass. Mean sputtering depths by H for E = 100 eV and θ = 0° do not exceed 2 Å for a:Li₂O and 8 Å for a:LiOH (a:Li and a:LiH are between these values).

The distribution of retained incident particles could serve as an initial condition for diffusion evolving at a much longer timescale than covered by this computation. The depths of the retained incident H particles are larger for the a:LiOH surface than for a:LiH and a:Li₂O surfaces but do not exceed 10 nm for 90% of the retained particles, which decreases with impact angle and isotopic mass. These retention depth values are smaller than the retention depths in low-density a:Li (15 nm for 90% of the retained particles). The electronic stopping power due to inelastic electronic transitions could be included in CMD indirectly as an external friction force⁵⁴ at the projectile, which can influence the retention depth. We neglected this effect, assuming that it is small in the range of impact energies of neutral particles studied in this work but certainly needs to be explored in the future.

The presence of oxygen in the lithium coatings, which is practically unavoidable, might significantly influence the processes at the plasma–material interface. However, there are also other distinctions between the lithium surface compounds that cause various polarizations of the incident and surface atoms and different combinations of the long- and short-range forces, which all can influence the probability of reflection, retention, and sputtering reported herein. However, the bonding chemistry of the retained particles depends on the chemical and structural composition of the surface and its evolution as the hydrogen accumulates in it. This subject is not covered by the present work, requires a combination of chemical dynamics based on the classical ReaxFF/EEM force field and quantum mechanics, and will be the subject of our future research.

See the supplementary material for strain-energy relation and use of ZBL potential for training Li2O and LiOH crystals.

This material is based upon work by the U.S. Department of Energy, Office of Science/Fusion Energy Sciences, under Award No. DE-SC0019308. AM and SA were supported through PPPL by DOE contract No. DE-AC02-09CH11466. SD and AvanD were supported by DOE Award No. DE-SC0022013 for the development of the Li–O–H potentials. PK acknowledges Princeton University for the use of the Stellar HPC cluster and to XSEDE for access to the SDSC Expanse HPC through Grant No. TG-DMR110037.

The authors have no conflicts to disclose.

PK did the conceptualization, developed computational method and analysis software, did computation, analysis, validation, and visualization and curation of the data, provided some of the computational resources, wrote original draft. ETO, SA, and AM took part in writing—review and editing, visualization, validation, and data curation. SD and AvanD developed and verified advanced ReaxFF potentials for Li–O–H. BK supervised the work, did review and editing of the manuscript, performed project administration, funding acquisition, and provided resources and validation.

P. S. Krstic: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Methodology (lead); Resources (lead); Software (equal); Validation (equal); Visualization (lead); Writing – original draft (lead). E. T. Ostrowski: Validation (supporting); Visualization (supporting); Writing – review & editing (equal). S. Dwivedi: Conceptualization (supporting); Data curation (supporting); Software (equal); Visualization (equal); Writing – review & editing (supporting). S. Abe: Validation (equal); Visualization (equal); Writing – review & editing (equal). A. Maan: Conceptualization (equal); Investigation (equal); Validation (equal); Writing – review & editing (equal). A. C. T. van Duin: Conceptualization (equal); Software (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). B. E. Koel: Conceptualization (equal); Funding acquisition (equal); Project administration (equal); Resources (equal); Supervision (equal); Writing – review & editing (equal).

P. S. Krstic: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Methodology (lead); Resources (lead); Software (equal); Validation (equal); Visualization (lead); Writing – original draft (lead). E. T. Ostrowski: Validation (supporting); Visualization (supporting); Writing – review & editing (equal). S. Dwivedi: Conceptualization (supporting); Data curation (supporting); Software (equal); Visualization (equal); Writing – review & editing (supporting). S. Abe: Validation (equal); Visualization (equal); Writing – review & editing (equal). A. Maan: Conceptualization (equal); Investigation (equal); Validation (equal); Writing – review & editing (equal). A. C. T. van Duin: Conceptualization (equal); Software (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). B. E. Koel: Conceptualization (equal); Funding acquisition (equal); Project administration (equal); Resources (equal); Supervision (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.