Influence of temperature on the displacement threshold energy in graphene

The atomic structure of nanomaterials is often studied using transmission electron microscopy. In addition to image formation, the energetic electrons impinging on the sample may also cause damage. In a good conductor such as graphene, the damage is limited to the knock-on process caused by elastic electron-nucleus scattering. This process is determined by the kinetic energy an atom needs to be sputtered, i.e. its displacement threshold energy Ed. This is typically assumed to have a fixed value for all electron impacts on equivalent atoms within a crystal. Here we show using density functional tight-binding simulations that the displacement threshold energy is affected by thermal perturbations of atoms from their equilibrium positions. This effect can be accounted for in the estimation of the displacement cross section by replacing the constant threshold energy value with a distribution. Our refined model better describes previous precision measurements of graphene knock-on damage, and should be considered also for other low-dimensional materials.

monitored a temperature-dependent behaviour of E d in TiO 2 and discussed the effects of temperature on the primary knock-on atom displacement threshold energies and the defect formation probability. Robinson's simulation results show a dramatic increase in the E d values of O atoms, from 18 ± 3 to 53 ± 5 eV, at temperatures of 300 and 750 K respectively.
In this article, we prove that the displacement threshold energy E d can not be assumed to be constant for each atom at finite temperatures and show that it deviates from the E d value at 0 K. Importantly, we show that in contrast to the results by Robinson 27 , the displacement threshold in our simulations needs to be represented as a distribution around the zero-temperature value. The width of this distribution depends on the temperature. We explain how this is due to the kinematics of the atomic lattice, using density-functional tight-binding (DFTB) simulations of graphene. We further derive an equation for calculating the maximum transferred energy resulting from the relativistic scattering of an electron from a vibrating atomic nucleus and compare it to previous approximations, and provide an improved theoretical model to predict displacement cross sections taking into account the realistic distribution of E d at a given temperature T. Our results show that temperature needs to be taken into consideration when describing knock-on irradiation effects, not only via the atomic velocities or a simple change in the E d value, but also due to the stochastic spread of E d caused by thermal displacements of the atoms.

Knock-on Displacements and their Cross Section
Due to the great difference in mass between an electron and a nucleus (a factor of ca. 22000 for carbon), the amount of transferred energy is limited. The maximum energy is transferred in a head-on collision between the electron and nucleus when the electron backscatters. If the transferred energy to the atom is large enough to produce a vacancy in the lattice which does not spontaneously recombine with the displaced atom, then the atom is considered displaced (knocked out). This minimum energy transfer needed to knock out an atom is called the displacement threshold energy E d .
In a TEM experiment with electron energies below 100 keV it is necessary to consider the vibration of the atoms in the out-of-plane direction (due to the typical experimental geometry). Resulting from the summing of the initial momenta, if an atom is hit by an electron while it happens to move parallel to the incoming electron beam it can acquire a higher transferred energy  E n than if it were at rest, as originally suggested by Brown and Augustinyak 32 and later elaborated and quantified by Meyer et al. 15 and Susi et al. 16 . To describe this situation, we need to consider a relativistic scattering process between an electron and a nucleus. The electron (mass m, energy E e , momentum p e ), a relativistic projectile, collides with a moving non-relativistic target, the nucleus (mass M, energy E n , momentum p n ). Taking into account energy and momentum conservation we can express the momentum of the electron after collision as n n e , which delivers an algebraically solvable equation for the maximum energy that an electron can transfer to a nucleus moving with velocity v parallel to the incident beam: n e e e 2 2

2
Setting v = 0 would recover the result of the static nucleus approximation 6 . In experimental studies, the measurable parameter is the displacement cross section. Meyer et al. 15 were the first to quantitatively show that under 80 keV electron irradiation, the defect-free graphene lattice remains undisturbed and that knock-on damage begins a few keV above this energy. To theoretically predict the displacement cross section, they approximated the phonon population of the material using a Debye model. The three-dimensional velocity of the atom was considered instead of the better justified one-dimensional velocity, which led to an overestimate of the out-of-plane mean square velocity by a factor of three. This model was improved upon by Susi et al. 16 replacing the three-dimensional Debye model by the out-of-plane phonon density of states calculated with density-functional theory (DFT).
While these models both take into account the fact that the target atoms are vibrating, they assume that the displacement threshold energy itself is constant for every electron impact. In this work we present a model that includes a temperature dependency of the displacement threshold: where P Ed T (E) is a normal distribution of probabilities for the displacement threshold E d at a certain temperature T integrated over all possible energies E, P vz T (v) is the probability distribution of velocities of the target atoms in the out-of-plane direction introduced in refs 15,16 , and σ d is the displacement cross section derived by Seitz and Kohler 33 , which has been used to calculate cross sections for carbon nanostructures 15,16,34 . www.nature.com/scientificreports www.nature.com/scientificreports/

Methods
We performed density-functional tight-binding (DFTB)-based MD to study the displacement threshold energy of graphene at finite temperatures. Atomic Simulation Environment 35 was used as the front-end for the MD to set up, perform, visualize, and analyze the results with the non self-consistent and non spin-polarized DFTB+ calculator 36-38 (version 1.3). We employed the "matsci-0-3" parameter set in our calculations 39 . Although our absolute value of the displacement threshold energy is lower than that previously obtained with DFT or experimentally estimated 16 , a similar methodology has been useful to study the dynamics of carbon atoms in graphene under electron irradiation 7,40 . We point out that although more accurate DFTB-based methods are available (including self-consistent methods and spin-polarized calculations), we found that they also fail to reproduce the correct threshold but are more than an order of magnitude more demanding computationally and not necessary to establish a qualitative understanding of the variations in the threshold value depending on the temperature of the structure. The initial structures were generated as a 15 × 15 pristine graphene supercell structure (450 atoms) thermalized to 73, 139, 216, 307, 383 and 434 K. Displacement thresholds were estimated by picking one C atom at a time and giving it initial momentum in the z-direction (out-of-plane) until it got ejected (when its distance from the lattice reached at least 5 Å), similar to earlier works 16,34,41,42 . We repeated this process for up to 450 atoms in each thermalized structure. The time step used in the calculations was 0.3 fs and the E d was estimated to within 0.1 eV for each displaced atom.

Results
The histograms in Fig. 1 show the spread in the displacement threshold E d at different temperatures. The distribution is narrow at low temperatures and gets wider for higher temperatures. The spread is symmetric around the 0 K value of 20.0 eV, as estimated with the DFTB method, which allows us to fit a normal distribution to the histograms. In Fig. 2, we plot the distribution full width at half maximum (FWHM) that varies between 0.91 and 1.47 eV depending on the temperature, starting at lower values for low temperatures and reaching a constant level at around 300 K.
In order to understand what causes the change in the displacement threshold, we analyzed the trajectories of the simulations to gather information on the geometry and forces that might influence the outcome. We extracted   www.nature.com/scientificreports www.nature.com/scientificreports/ information from the thermalized state (before any displacement occurred) for each atom that was about to be ejected. Our analysis suggests that there is no direct correlation between the ejecting atom's initial in-plane/ out-of-plane velocity, z-position, or forces in-or out-of-plane and the energy E d needed to displace it. However, we do find a linear correlation between both the mean z-coordinate of the displaced atom i and its neighbors against the displacement threshold E d , as shown for 300 K in Fig. 3. As can be seen from the scatter in this data, neither one of these parameters can be used to predict the displacement threshold E d for a given atom due to the complex nature of the dynamical process.
We explored the mechanism further by investigating a non-thermalized graphene sheet. Manually displacing just the ejecting atom itself in the z-direction with respect to its neighbours only made it easier to eject, regardless of the direction of the displacement, and thus cannot explain the observed threshold energy spread. However, if we set the initial positions of both the ejecting atom and its three neighbors some tenths of an Ångström above or below the rest of the atoms, we did find an asymmetric change in the displacement threshold energy. Alternatively, to confirm the connection between the neighbors and the ejecting atom, we set initial momentum values to the C atoms' nearest neighbors and started the displacement simulations. Neighbors with negative initial momenta help the C atom escape the lattice (in the positive direction) with less energy, while those with momenta towards the displacement direction contributed to a higher displacement threshold. Setting an initial momentum value for the ejecting atom itself did not explain the spread in the threshold energies, which explains our lack of success in attempting to correlate the displacement threshold only with the starting parameters of the ejecting atom itself. Hence, the ability of the nearest neighbors to follow and pull back the ejected atom, similar to what was noticed for graphene edges in ref. 42 , determines how easy it is to eject the atom.
Finally, to estimate the reliability of the structures thermalized via molecular dynamics within the DFTB method, we compared the root mean square velocity of the atoms for in-plane and out-of-plane directions at each temperature in Fig. 4. As expected, these follow the Maxwell-Boltzmann statistics. In the same plot, we show velocities calculated from the phonon density of states shown in ref. 16 . Clearly, the classical description within standard MD leads to an underestimation of the velocities at least up to 900 K. This means that the results we have obtained from DFTB can not be expected to give quantitatively accurate results for any given temperature. Instead, they must be considered to provide qualitative information to suggest magnitudes and trends for experiments.
Finally, we compared the experimental displacement cross section data acquired with the Nion UltraSTEM 100, presented in ref. 16 to a two-parameter fit to Eq. 3. We numerically integrated Eq. 3 for different displacement thresholds E d and distribution widths (FWHM), and calculated the mean square error weighted by the uncertainties of the experimental values to find the best fit. We found two nearly equal minima with a relative difference of just 3.4%, one at E d = 21.2 eV and w = 0.5 eV, and the other one at E d = 21.3 eV and w = 1.0 eV (see heat map in Fig. 5). In Fig. 6, we plot the experimental data and the one-parameter fit from ref. 16 along with our new two-parameter fit with the lowest error. Although the differences are relatively small, the new fit clearly improves upon the old values as compared to the experimental results. We point out that although this fit includes no information from the DFTB simulations, the values found for the width of the E d distribution are similar to those found in the simulations (ca. 1 eV).   Fit: E d =21.3, w=1.0 eV Fit: E d =21.14 eV [9] σ (barn) Acceleration voltage (kV) Figure 6. Comparison of fitted cross section values to STEM data from ref. 16 . The black squares represent the cross section values from the experiment including the uncertainties shown with error bars, while the circle and the triangle symbols correspond to the one-parameter and the two-parameter fit, respectively. (2019) 9:12981 | https://doi.org/10.1038/s41598-019-49565-4 www.nature.com/scientificreports www.nature.com/scientificreports/

Discussion
The simulation results presented here demonstrate that the displacement threshold energy should not be assumed to have a fixed value for low-dimensional materials when accurate estimates for the displacement cross section are needed. Simulating graphene sheets at various finite temperatures shows that the thermal perturbation of atoms from their equilibrium positions leads to a symmetric spread of their displacement threshold energies by ca. 1 eV around the 0 K value. The main reason for the spread of the displacement threshold is the variation in the momenta of the atoms neighboring the displaced atom at the moment of the impact. If these are opposite to the displacement direction, the atom will be easier to knock out. On the contrary, if their momenta are in the same direction, the atom will require more energy to knock out. We derived a formula for calculating the maximum transferred energy in the relativistic scattering of an electron at a moving nucleus, and improved the model for evaluating the theoretical displacement cross section by including the effect temperature has on the displacement threshold, providing a two-parameter fit which better describes experimental data. An accurate description of the displacements of carbon atoms from graphene is a step forward in a more precise understanding of knock-on damage also in other materials, although the magnitude is expected to be smaller in bulk materials due to the limited available space around the displaced atom.