-
Notifications
You must be signed in to change notification settings - Fork 16
Benchmarks
This page will serve to collect the results of post-release benchmarks of RustBCA. These include comparisons to other codes, comparisons to (semi-)empirical formulas, and comparisons to experimental data where available. Benchmarks can be requested using the [benchmark] flag when creating a new issue.
Improvement of titanium hydrogenation by low energy ion irradiation, published in the International Journal of Hydrogen Energy, includes a comparison of experimental data (ERDA) to SRIM. The two plots are significantly different from each other, from which the authors conclude the following:
"The experimental data indicates that hydrogen depth profile is slightly broader than the SRIM simulation, where the experimental value for the ion projected range is 81 nm. This result suggests that a radiation induced diffusion of hydrogen into the titanium bulk is taking place during implantation; which could be related to the high mobility of hydrogen atoms in titanium, increased by a thermal mechanism during implantation, even though low energy implantation is used during the experiment." (emphasis mine).
However, using their digitized data, a different conclusion can be reached when compared to RustBCA run with default values for hydrogen and titanium:
5 keV H on Ti implantation distributions. SRIM results obtained by the authors do not agree with the experimental data - RustBCA results agree very well.
Source https://doi.org/10.1016/j.ijhydene.2015.01.166
A critical review and meta-analysis of xenon-on-carbon sputter yield data, a review paper just published by J. E. Polk, provides a great set of data of Xe on C sputtering yields.
Here's the comparison of RustBCA with the experimental data from Polk's paper for three different Xenon accumulations, 0.0 Xe/C, 0.05 Xe/C, and 0.1 Xe/C. This paper by Kenmotsu suggests using ACAT that the low sputtering threshold seen in experiments in comparison to simulations of pure carbon is due to the accumulation of xenon, which significantly reduces the sputtering threshold. Looking at the above data, it's clear that the low sputtering threshold can only be explained by nonzero xenon accumulation, since the original data come from experiments with different conditions and forms of carbon and no other free parameter has as significant an effect on the sputtering threshold as the presence of relatively high-Z xenon in the target.
I ran these simulations using Oen-Robinson electronic stopping, per the recommendations in the ACAT paper. I used the following properties for carbon:
carbon = {
'symbol': 'C',
'name': 'carbon',
'Z': 6,
'm': 12.011,
'n': 1.1331E29,
'Es': 7.37,
'Ec': 1.0,
'Eb': 0.0,
}
And the properties for xenon in scripts/materials.py
. Since I only wanted to measure sputtering yields, I overwrote Ec for both species with the Es of carbon; this has no effect on the calculated sputtering yields, since anything below Es will not sputter anyway.
I ran all three xenon concentrations using both isotropic and planar surface binding models; isotropic would be more appropriate for an atomically rough surface, whereas planar would be more appropriate for an atomically smooth surface. In this case, I think isotropic surface binding better reproduces the data.
Best agreement is for a xenon accumulation of 10% Xe/C, but agreement is still acceptable for 5% Xe/C. In the ACAT paper, they found best agreement with 14%. Although there is good agreement between RustBCA and experiment here, there is one caveat and a cautionary tale. The caveat is that, when viewed in log-log form, it can be seen that RustBCA still overestimates the sputtering threshold by about 4-6 eV:
This may be due to a lower actual surface binding energy than the literature value, multi-body effects, or some other phenomenon not captured in RustBCA.
The cautionary tale lies in the importance of accumulated xenon; the sputtering yield curves are very different for a pure carbon target. When dealing with ions with a mass ratio far from one, the accumulation of those ions must be accounted for when calculating sputtering yields. The authors of the ACAT paper note that a similar effect was not observed for a non-carbon target (Mo), so this will be material dependent.
This is a straightforward comparison to experimental values originally used to benchmark TRIM.SP (the precursor to TRIDYN). For this comparison, the ergonomic python functions were used, with two modifications. First, the electronic stopping mode was changed from LOW_ENERGY_NONLOCAL (Lindhard-Scharff) to LOW_ENERGY_EQUIPARTITION (50% Lindhard-Scharff + 50% Oen-Robinson). Second, the bulk binding energy was set to 3.0 eV and the surface binding energy was set to 8.68 eV to match the TRIM.SP inputs from 1. Agreement w.r.t. energy is very good above 200eV, but both RustBCA and TRIM.SP calculate significantly lower near-threshold yields (see 1 for the TRIM.SP calculations) than observed in experiment. RustBCA calculates slightly lower sputtering yields (~10%) at non-normal incidence compared to TRIM.SP results from 2. However, they do not report any BCA parameters in that paper, so it is possible that the original results could be recovered by reducing Eb, by using isotropic surface binding, using Lindhard-Scharff, etc.
Figure 1: W on W self-sputtering w.r.t. energy. These values had reported error bars of approximately 5-10%. Note the difference in behavior at near-threshold energies; TRIM.SP results were consistent with RustBCA, not the experiment.
Figure 2: W on W self-sputtering w.r.t. angle. These values had no reported error bars. Two sets of data, one including some sort of tungsten accumulation on the surface, were included for 350 eV.
Comparison of Attractive-Repulsive Potential Simulations to Experimental H on Ni Reflection Coefficients
See examples/benchmark_eam.py. EAM and Experimental data from https://doi.org/10.1016/0022-3115(84)90433-1.
This shows that previous assumptions about the energy regime of validity of BCA codes may not be accurate for all quantities. In particular, reflection coefficients from TRIM were inaccurate below 10s of eV; this seems to have been implicitly extended to all BCA codes in the literature. However, RustBCA has qualitatively good results for reflection at energies significantly lower than 10s of eV.
Using the default settings, RustBCA performs well when compared to both experimental data at higher energies and an EAM model at much lower energies. Modification of the default settings, such as using a realistic interaction potential, or changing the weak collision order, is likely to improve the results further.
See this example for the script to run the Morse potential simulations for the above figure.
By introducing the Morse potential for H-Ni (D=0.4219 eV, r0=2.782 A, alpha=1.4198 1/A, parameters from Olander 1971) and increasing the surface binding energy for H to the literature value for H on Ni surfaces of 2.5 eV (from the default assumption for H of 1.5 eV), and choosing a nonlocal electronic stopping model with an appropriate correction applied (ck=1.09), the results at low energy can be significantly improved (note that the latter two changes were also applied for the Kr-C case in the figure above). However, since the high-energy (that is, small-r) behavior of the Morse potential does not approach the correct value (that is, a screened coulomb-like potential) the accuracy decreases significantly above about 200 eV. It is likely that the accuracy can be further improved with a better interaction potential.
Note that error bars that reflect the small number of EAM trajectories from Baskes (1984), 100, have been included in this updated plot. See examples/test_morse.py for the script used to generate this plot.
This potential combines with sigmoidal interpolation the Kr-C potential for short-range forces and the Morse potential for the attractive portion to reproduce reflection coefficients from 0.1 eV to 10 keV; the interpolation parameters are k=7.0 1/A (interpolation characteristic length) and x0=0.75 A (interpolation transition point). Note that I've reduced the surface binding energy to 1.5 eV, the default value - in the case of an attractive-repulsive potential, it's not actually clear what value the surface binding energy should take - since it should include the contribution of all surface atoms except for the one most recently collided with. Olander (1971) gives a value of ~2.4 eV for the surface binding potential. In the literature, Es is often used as a fitting parameter for sputtering data, so using it as one here is not uncalled for, but it would be nice to have a theoretically sound value to make the values more convincing.
This work was presented at APS-DPP 2023 under the title "Modeling of hydrogen plasma-material interactions with reactive metals."
This paper shows a plot with data extracted from an AC∀T paper that they use to suggest that BCA codes do not produce accurate sputtering yields compared to MD codes:
"We see that our results are considerably lower than the ones obtained using BCA... BCA is heavily influenced by the parameterization used for the particular case, and at low energies BCA has been seen to overestimate the results... "
However, if they had used another BCA, such as SDTrimSP or RustBCA, they would have had markedly better agreement:
I used the literature value for tungsten's bulk binding energy, 3.0 eV, the values in materials.py for everything else, and 100,000 ion trajectories using the Python bindings on the development branch (at time of writing). To recreate the RustBCA results:
from libRustBCA.pybca import *
from scripts.materials import *
import numpy as np
tungsten['Eb'] = 3.0
energies = np.logspace(2, 3, 100)
y = [sputtering_yield(argon, tungsten, energy, 0.0, 100000) for energy in energies]
plt.plot(energies, y)
plt.show()
MD paper: https://doi.org/10.1016/j.commatsci.2022.111876
AC∀T paper: https://doi.org/10.15748/jasse.3.165