Logo Beilstein Institut

Influence of water contamination on the sputtering of silicon with low-energy argon ions investigated by molecular dynamics simulations

  1. 1,2 ORCID Logo ,
  2. 3 and
  3. 1 ORCID Logo
1Advanced Instrumentation for Nano-Analytics (AINA), Materials Research and Technology Department (MRT), Luxembourg Institute of Science and Technology (LIST), 4422 Belvaux, Luxembourg
2University of Luxembourg, 4365 Esch-sur-Alzette, Luxembourg
3Thermo Fisher Scientific, Hillsboro, OR, 97124, USA
  1. Corresponding author email
Associate Editor: E. Meyer
Beilstein J. Nanotechnol. 2022, 13, 986–1003. https://doi.org/10.3762/bjnano.13.86
Received 17 Mar 2022, Accepted 06 Sep 2022, Published 21 Sep 2022
A non-peer-reviewed version of this article has been posted as a preprint https://doi.org/10.3762/bxiv.2022.14.v1
Full Research Paper
cc by logo
PDF
Album
Supp Info
Cite

Abstract

Focused ion beams (FIB) are a common tool in nanotechnology for surface analysis, sample preparation for electron microscopy and atom probe tomography, surface patterning, nanolithography, nanomachining, and nanoprinting. For many of these applications, a precise control of ion-beam-induced processes is essential. The effect of contaminations on these processes has not been thoroughly explored but can often be substantial, especially for ultralow impact energies in the sub-keV range. In this paper we investigate by molecular dynamics (MD) simulations how one of the most commonly found residual contaminations in vacuum chambers (i.e., water adsorbed on a silicon surface) influences sputtering by 100 eV argon ions. The incidence angle was changed from normal incidence to close to grazing incidence. For the simulation conditions used in this work, the adsorption of water favours the formation of defects in silicon by mixing hydrogen and oxygen atoms into the substrate. The sputtering yield of silicon is not significantly changed by the contamination, but the fraction of hydrogen and oxygen atoms that is sputtered largely depends on the incidence angle. This fraction is the largest for incidence angles between 70 and 80° defined with respect to the sample surface. Overall, it changes from 25% to 65%.

Keywords: angle dependency; argon ions; contamination; focused ion beams; ion bombardment; low energy; molecular dynamics; silicon; simulations; water

Introduction

Focused ion beams (FIB) play an increasingly important role in materials research areas such as nanoanalysis (e.g., secondary ion mass spectrometry (SIMS) [1-3] and sample preparation for transmission electron microscopy (TEM) [4], atom probe tomography (APT) [5], and ion beam analysis used for life sciences applications [6,7]), surface patterning [8], nanolithography [9], nanomachining [10,11], and nanoprinting at room [12] and cryogenic temperatures [13]. The development of nanotechnology relies on lower ion beam energies and smaller spot sizes to reduce the thickness of the layer damaged by ion beams, and to increase the lateral resolution for precise machining and sample characterization. For most of these applications, the quality of the sample surface and its cleanliness are essential and, therefore, highly controlled. Depending on the application, the ion beam energy is in the range of 10 to 30 keV when small spot sizes are required (i.e., spot sizes in the nanometre range) and at a few keV or even in the sub-keV range when low surface damage or minimized atomic mixing is required. One example is low-energy depth profiling SIMS to resolve thin films in multilayered samples [14].

Another example is TEM sample preparation, where the achievement of the highest lateral resolutions in the subsequent TEM analysis requires lamellae thicknesses between 10 and 20 nm, but goes along with a typical amorphous layer of 2 to 4 nm formed during the sample preparation by FIB milling [15-17]. Such an amorphous layer represents a substantial part of the thickness of the sample and information coming from this part does not correspond to the initial sample structure. Minimizing the thickness of this amorphous layer during FIB milling is essential because most samples analysed in high-precision instruments are prepared using this method. This can be best achieved using low-beam energies, ideally in the sub-keV range [18], since low-energy ion beams (under 500 eV) produce a thinner amorphous layer due to their lower penetration depth. Investigations performed with low-energy argon ions [15,19,20] have shown that the current model describing the sputter yields and the sputtering processes (such as sputtering threshold and the amorphization process) does not fit with experimental data, leading to discrepancies that cannot be ignored.

Contaminants on the sample surface can also play a critical role during the milling process: water can be found on nearly every sample surface [21]. For example, in a vacuum chamber at a pressure of 10−8 mbar, there are still 106 to 109 molecules per cm3. This contaminant has a strong impact on the ion beam process by modifying sputtering processes. Even under an ultrahigh vacuum of 10−12 mbar there are still 104 molecules·cm−3 remaining in the experimental chamber, thus making water by far the most common contaminant. These assumptions can be confirmed by SIMS experiments [22]. It is therefore safe to assume that most experimental chambers able to hold vacuum between 10−6 to 10−10 mbar will contain water molecules that will interact with the sample surface. Water has been widely studied [23-26] and is well described by many simulation tools available to model samples in such extreme cases.

Nowadays, very little is known about the influence of contaminations on the amorphization process under ion irradiation. Thanks to molecular dynamics (MD) simulations, a wide range of materials properties and process parameters can be reproduced [27-29]. In this paper, we are using MD simulations to study the sputtering of a surface with water contamination by sub-500 eV ions. The information of interest is the chemical reactions occurring at and below the sample surface, as well as the mixing of the contaminant layer into the sample. The ReaxFF reactive force field is used in this work [30]. ReaxFF force fields are specifically tuned for a set of atomic interactions. They are developed from quantum calculations and are adapted for MD simulations, providing faster calculations than pure quantum electrodynamics (QED)/density functional theory (DFT) and more information than classical MD/binary collision approximation (BCA) simulations. ReaxFF simulations can cover a broad range of applications: from DNA molecule bombardments with heavy atoms [31] to graphene layer deposition on copper surfaces [32]. Numerous simulations using the ReaxFF potential have been conducted in the past years, several describing a given sample with contaminant interactions (such as water or organic particles) [33-35].

To understand how the water layer influences the amorphization process, we designed metrics aiming to characterize the degree of amorphization in the sample. Previous works have tried to give a broader picture of the phenomenon in solids or in crystals [36], by using a very small sample size and performing an exhaustive study on bond orders and bond length average (Gaussian core model), as well as analysing the granularity of the system in liquid/solid phases. However, these methods were designed using a limited sample size of up to a few hundreds of atoms. Other studies have tried to elaborate a local metric for condensed-phase environments [37], giving a good descriptor for phase transitions. While they provide good results for a specific crystalline structure, they are tedious to adapt to other structures. In this paper, we developed a simple metric based on bond length variations, which can be easily adapted to different crystalline structures to characterize sample amorphization and also to study sputtering surfaces exposed to residual water molecules.

Computational Methods

Force fields

The ReaxFF force field differs from other established force fields as it aims to bridge quantum mechanics (QM) and classical MD. Quantum mechanics algorithms are limited to small-sized samples (up to a few hundreds of atoms) due to the difficulty to solve Schrödinger equations for many-body systems whereas classical MD makes the assumption of fixed bonds and usually do not include complex chemical reactions. ReaxFF captures the dynamic nature of bond formation and breaking, based on bond-order calculation from interatomic distances. ReaxFF also includes nonbonded terms such as van der Waals and Coulombic interactions in the energy calculation, resulting in a sum of partial energies [38]. The ReaxFF force field was developed by van Duin et al. [34] for the Si–C–O–H system, and for all simulations described in this paper, ReaxFFSiOCH(2019) is used. To reproduce the interaction between argon ions and silicon, hydrogen, and oxygen atoms, we used DFT to simulate the potential energy between each of these pairs: Ar–Ar, Ar–Si, Ar–H, and Ar–O. Once extracted, the potential energy was fitted using the Morse potential which is described in Equation 1:

[2190-4286-13-86-i1]
(1)

where r represents the interatomic distance in distance units, r0 the interatomic distance in the equilibrium state (also in distance units), D0 the depth well in energy units (defined per atomic interaction and related to the molecule dissociation energy), and α the control parameter on the well thickness, proportional to the curvature of the potential at its origin, in inverted length units [39].

The DFT data has been obtained using the VASP [40-43] package. The Perdew–Burke–Ernzerhof generalized gradient approximation (GGA-PPE) [44] is selected together with the projector augmented wave (PAW) [45,46] method for the pseudopotentials describing core electrons of the different elements. Van der Waals contributions are included using the Tkatchenko–Scheffler method [47]. A plane wave cut-off of 500 eV is used and k-points are automatically generated. The self-consistent calculation is performed using a Davidson block iteration scheme with a threshold of 10−6 eV/atom. The force threshold for the conjugate gradient geometry optimisation is set to 10−4 eV/Å. The interatomic interactions are calculated by changing the distance between Ar and H and between O and Si atoms by 0.1–12.0 Å in 0.1 Å steps in a simulation box of 25 × 10 × 10 Å3 by taking the energy for a couple of Ar–X interactions by the energy of the two atoms taken individually.

The fit of the Morse potential can be seen in Figure 1 for the Ar–Si bond, and the remaining data of the fits can be found in Supporting Information File 1. The Ar–X interactions have a potential well in the order of a few meV, which makes it rather difficult for the fit to precisely match the well. However, this is not problematic as the potential well is in the order of magnitude of thermal vibrations and its impact on the outcome of the simulations is small. The fits match the repulsive part of the potential for the Ar–X couples well, except for high energies. However, for energies up to 100 eV discrepancies are small so that the interpolation to the ZBL potential at short interatomic distances can be omitted.

[2190-4286-13-86-1]

Figure 1: Fit of the DFT data for the argon–silicon potential in (a) the 1–5 Å region, with the potential well, and (b) the entire potential.

The ReaxFF potential is used to model the damage cascade and the collisions between silicon, hydrogen, and oxygen particles. An hybrid pair style [48] was used in the MD simulations to match both the ReaxFF potential (for Si–Si, Si–O, Si–H, O–O, O–H, and H–H bonds) and the Morse potential (for Ar–Si, Ar–O, Ar–H, and Ar–Ar interactions). The charge equilibration for the ReaxFF potential was performed using the algorithms developed by Rappé et al. [49-51], namely QEq, allowing charge distributions to be calculated according to geometry. Due to the dynamic nature of bonds and geometrical configurations, the charge distribution QEq is performed every 10 time steps during the simulation of ion bombardment, allowing to cut down computation costs while maintaining a good fidelity to the charge distribution. Furthermore, we implemented a variable time step, which allowed us to control the variation in either position (with a limit from 0.01 to 0.1 Å) or energy (up to 0.4 eV) in which the time step would be reduced in order to take these variations into account. As a result, a maximal displacement of one interatomic distance observed in the simulation could not occur in a time frame smaller than approximately 40 time steps, which reassured us regarding the 10 time step limit for QEq. During water deposition on top of the silicon sample, the charge equilibration was performed at each step to fully describe the deposition of the contaminant on top of the sample.

Molecular dynamics simulations

Molecular dynamics simulations of argon bombardment on a silicon surface were carried using the LAMMPS code [52]. A pristine silicon sample containing 5248 atoms was created and used as a reference sample. On top of this sample was added a water layer to study the effect of contaminations. Periodic boundary conditions were applied in x- and y-directions, with some vacuum above and below the sample in the z-direction. The different axes are defined in Figure 2.

[2190-4286-13-86-2]

Figure 2: Visualisations of the samples during equilibration steps: (a) represents the silicon sample after Stillinger–Weber equilibration, showing dimers on the surface with a VDW representation; (b) represents the water layer after equilibration (top view), while (c) shows the side view, both in a VDW representation. (d) Final sample once equilibrated with both the water layer and the silicon sample using ReaxFF and energy minimized with water molecules in a VDW representation, and silicon atoms in a CPK representation. The direction of the z axis is defined by the blue arrow, the y axis by the green arrow and the x axis by the red arrow.

The simulations for sample preparation were run in an NVT ensemble using the Langevin thermostat with a damping parameter set to 100 time steps and the temperature set to 300 K [53]. The pristine sample was equilibrated using the Stillinger–Weber potential to obtain the dimers formed on the Si(100) surface [21,54], with a time step of 0.1 fs for 1 ns. Some snapshots of the dimers on the surface can be found in the Figure 2.

The equilibration was pursued with the ReaxFF potential for another 500 ps in the same aforementioned conditions. The largest displacements observed in the sample were vibrational motions due to a temperature of 300 K (between 0.5 to 1 fs), and the time step was chosen to limit the maximum displacement between steps to one tenth of that motion, resulting in a 0.1 fs time step for the simulations. Furthermore, as ReaxFF tends to reach a better stability at a shorter time step, a variable time step was defined between 0.05 and 0.1 fs, with a maximal allowed energy variation per atom between each step of 5 eV. This process allowed to limit the maximum displacement of atoms per time step while the damage cascade happened and to increase it after the energy dissipated into the system. During the ion irradiation simulations, the bottom slab of the sample (one lattice high) was set in an NVT ensemble using the same methodology in order to reproduce a heat-bath effect [55]. The remaining atoms were kept in the NVE ensemble in order to avoid any impact of a thermostat on the collision cascade.

Once properly equilibrated, the simulation of argon bombardment of pristine Si(100) was run for 50 ps, which was enough for the modelling of energy dissipation and sputtering, and then thermalized to 300 K. The last script was set to end as soon as the computed mean temperature reached 300 K. The whole process was repeated until 500 ion impacts had been simulated. Compared to experiments, the time between two ion impacts needed to be reduced by several orders of magnitude due to computational constraints (i.e., a few hundred picoseconds in simulations compared to µs in experiments) leading to an overall higher irradiation density yet keeping the sample at room temperature. To reproduce an ion beam, the initial positions of the argon ions were randomized in the surface plane. Each simulation set consisted of 500 irradiation events at a set incidence angle and at an impact energy of 100 eV. For a surface area of 18.9 nm2, this corresponds to a fluence of 2.6 × 1015 atoms·cm−2. The simulations were carried out for incidence angles of 0, 20, 30, 45, 65, 72, 75, 83, and 85° with respect to the surface normal. These angles were selected to cover a good sample in the 0–85° range, including angles of specific interest, such as 72° and 83° which represent a channelling angle for Si(100) [56] and a grazing incidence angle, respectively, frequently used in TEM lamellae preparation [57]. A visualization of a set of simulations can be seen in Figure 3 and the global parameters of the simulations can be seen in the Table 1. Typical values for fluence in experiments are between 1015 and 1017 ions/cm2.

Table 1: References of the simulation parameters.

  Surface (nm2) No. of impacts Energy (eV) Fluence (ions/cm2)
simulations 18.9 500 100 2.6 × 1015
[2190-4286-13-86-3]

Figure 3: Representations of (a) clean and (b) contaminated samples. The silicon atoms are shown in yellow, oxygen atoms in red, hydrogen atoms in white, and argon atoms in pink. Both samples were bombarded at 45° and 100 eV and are viewed from the (110) direction after 500 irradiation events.

To prepare a contaminated sample, a slab of water was equilibrated first and then deposited onto the Si sample. The water molecules were equilibrated using the ReaxFFOH(2017) [24] potential in a fully periodic space using the NPT ensemble. Since the adsorption of water molecules on a surface is a random process, the orientation of the molecules could freely evolve. Snapshots of water equilibration can be found in Figure 2. Next, the potential was switched to the ReaxFF SiOCH potential [34] used for the bombardment simulations, followed by another energy minimization. The change to ReaxFFSiOCH(2019) did not induce any significant change during the potential change, the main visible one being a slight variation of the bond length O–H in the water molecules that did not directly interact with silicon particles (i.e., water molecules on top of the contaminant layer). Thereafter, the water layer was deposited on top of the pristine silicon sample in the NVE ensemble for up to 1 μs in order to reproduce the gradual contamination of the sample by residual gas molecules, by giving a vertical energy of 0.5 eV to the water slab. After water deposition, a series of irradiation simulations was carried out using the same conditions than those for the pristine sample.

For the estimation of water contamination, we based our simulations on instrument chambers with a vacuum level between 10−8 and 10−6 mbar and a sample area of 18.9 nm2. According to [58,59] water forms about 10% of the residual gas in the experimental chamber. By using the Hertz–Knudsen formula (Hertz–Knudsen formula, p. 159 in [60]) we can determine, for a specific pressure, the temperature and atomic mass of a molecule and the number of collisions per second and per square centimetre:

[2190-4286-13-86-i2]
(2)

where Za is the number of collisions per second, P is the pressure in Torr, T is the temperature in Kelvin, and m is the atomic mass of the molecule of interest in amu. With P = 7.5 × 10−8 Torr, T = 300 K, and m = 18.02 amu we obtain a value for the collisions per square centimetre of Za = 3.58 × 1013 collisions·s−1·cm−2. Once scaled to the sample size, we obtain Za = 6.61 collisions·s−1·surface−1, from which about 10% would be water collisions. Hence, for typical vacuum conditions a monolayer of water is formed within minutes.

Results and Discussion

Sample characteristics and contaminant effects

At first, bond length distributions were calculated to observe the effect of the water contaminant on these distributions. The nominal Si–Si bond in a diamond lattice is 2.358 Å. Figure 4 shows the distributions for both samples before ion bombardment. For the Si–Si bonds, the most probable bond length is at the correct length, but a significant number of bonds are in the 2.1 Å region. This can be explained by dimers appearing on top and bottom surfaces which have a nominal bond length (for Si dimers) of 2.1 Å.

[2190-4286-13-86-4]

Figure 4: Bond length distributions for the pristine sample (blue) and the contaminated sample (pink) regardless of the particle type with a cut-off at 3.2 Å (at which distance the distribution for the 2nd neighbours starts).

The distribution for the contaminated sample shows two additional regions: the region between 0 and 2 Å and the region between 2.6 and 3.3 Å. The smallest bond length around 1 Å is due to O–H bonds [61] in water molecules. The second peak observed at approx. 1.4 to 1.7 Å is associated to both Si–H (1.56 Å) [62] and Si–O (1.58 Å) [63] bonds. The Si–H bonds are located near the surface between Si dimers and some of the H atoms of water after the reaction between silicon and oxygen, leading to the splitting of water molecules. Studies have suggested that oxygen dissociates on a single Si–Si dimer, allowing for the recombination of dissociated hydrogens with Si surface atoms [64]. Some counts are also measured between 2.6 and 3.2 Å, which come from the amorphization of regular Si–Si bonds due to argon collisions. Since most of the Si–O and Si–H bonds are shorter than these counts, it is safe to assume that the local disorder induced by reactions between water molecules and the sample surface tends to disorganize the silicon lattices. These bonds are longer than those of the crystalline structure, and their distribution follows a Gaussian distribution. Both distributions are depicted in Figure 5, showing the evolution of the counts with respect to the number of bombardments (i.e., with increasing fluence).

[2190-4286-13-86-5]

Figure 5: Evolution of bond lengths in the contaminated sample under bombardment. (a) The first graph shows bonds between 2.6 to 3.2 Å describing the amorphous layer formed by Si–Si bonds, and (b) the second graph shows bonds between 1 and 2 Å describing O–X and H–X bonds.

Furthermore, Figure 5 shows how 500 impacts modified bond length distributions: the counts related to amorphous bonds (between 2.6 and 3.2 Å) have clearly augmented, with a reduction in the O–H bond counts (for the contaminated sample). Hence, most of the water adsorbed on the sample surface is either fractioned and mixed in the sample with repeated argon atom impacts, creating more Si–H and Si–O bonds, or sputtered from the sample surface.

Description of the amorphization

The two samples at several stages of the sputtering process are shown in Figure 6 to illustrate the evolution of the amorphization process. The penetration depth of argon atoms is limited to the first few lattices inside the silicon, for both the contaminated and clean samples. For the clean sample, the amorphization depth linearly evolves until the implanted argon atoms reach a saturation concentration due to the diffusion and desorption of excess atoms. In the contaminated sample, the observations are similar for argon and for the contaminants, and hydrogen is implanted deeper than oxygen. Both atoms also tend to be implanted deeper than argon, proportionally to their respective masses, but some occurrences of deep argon implantation can happen. The amorphization can be divided into two phases: in an initial phase, most silicon atoms stay in their lattice position, even after receiving some energy from the incident argon ion, and only a small number of defects is formed. Water molecules are usually fragmented by an argon ion on the sample surface. The hydrogen and oxygen atoms are lighter than argon and are pushed deeper into the sample, globally increasing the damage per collision. The argon ion is either backscattered or implanted. Previously implanted argon atoms can also diffuse to the sample surface when receiving enough energy. The different processes gradually move deeper into the sample with increasing fluence. The second phase is characterized by several impacts which happen in close surroundings of a given region, and most atoms are displaced from their initial equilibrium position. This phase describes an amorphous layer. The depth of the amorphous layer varies as a function of the incidence angle and is maximum at around 70° with respect to the surface normal. In the amorphous region, contaminants and argon atoms reach a saturation concentration. Once the argon saturation is reached, the thickness of the amorphous layer reaches a maximum, which is defined by the range of the argon atoms. At the same time, material is removed by sputtering from the surface.

[2190-4286-13-86-6]

Figure 6: Evolution of the clean (a, b, c) and contaminated (d, e, f) samples under ion irradiation, before the 1st, 250th, and 500th bombardment showing the evolution of the sample surface. The argon particle on top of the simulation is the ion which is bombarded onto the sample.

Amorphization coefficient

To properly describe the amorphization process, we choose to monitor how the bond length distribution changes with the number of impacts. In particular, we calculated the coefficient described in Equation 3:

[2190-4286-13-86-i3]
(3)

To characterize the amorphization process as a function of depth, the simulation box was split into different slabs along the z-direction, each slab having the height of one-unit cell and the initial slab being situated at the bottom of the sample with an undisturbed surface. Then, the number of bonds and bond lengths in each slice are calculated to obtain the µ coefficient associated to each slab. By comparing the different values of µ to the radial distribution functions (RDFs) of the different slabs, the data can be grouped into three different categories: the crystalline slabs, the semi-amorphous slab(s), and the amorphous slabs (see discussion on RDF for more details). The 0.94 < µ < 1 interval describes a region where no defects are present. Changes in µ are due to variations in bond length coming from thermal vibrations around the equilibrium bond length. This region is crystalline. The 0.89 < µ < 0.94 interval describes regions that contain defects but are not completely amorphous (i.e., the crystalline structure is preserved at least in a small part of the slab). The collisions have only caused minor damages. For simplicity, we call this transitional phase a “partially amorphous” or “semi-amorphous” region in the discussion. The µ < 0.89 region designates slabs with almost no remaining crystalline structure. There is a large variance in bond lengths, generally describing a damaged/amorphous area with a lot of disorder.

A full breakdown of a sample bombarded by 0°, 45°, and 75° argon ions and the associated µ coefficients can be found in Supporting Information File 1, where each slab is detailed and compared with the associated RDFs. As per design, the “semi-amorphous” layer is supposed to be one-slab thick, as it is the transition between the amorphized region and the crystalline/intact part of the sample. In some cases, deeper penetration through the channelling of some atoms can cause a few slabs to be in the semi-amorphous regime. Figure 7 shows the three regions in the sample, and Figure 8 represents the evolution of µ in the sample with respect to the depth of the slab as well as with respect to the angle, describing how the incidence angle of the ion influences the amorphization depth. The sample shows a continuous transition from crystalline to amorphous. A large difference is observed between pristine and contaminated samples close to the grazing incidence. For the pristine sample the material is barely amorphized, whereas the presence of the contaminant leads to a full amorphization. In general, the incidence angle changes the depth at which the transition between amorphous and partially amorphized layer happens. Close to normal incidence the demarcation is clear. At larger angles, the transition is progressively shifted closer towards the surface. Especially for angles above 83° there is almost no amorphous layer, the upper part of the sample is at most partially amorphized by the collisions. This can be explained by the fact that at such high angles the argon particles are mainly backscattered and interact very weakly with the sample. In the case of the contaminated sample, most of the damage is done via the fractioning of water molecules on the sample surface. The implantation of argon ions is minimal.

[2190-4286-13-86-7]

Figure 7: Representation of each region in the sample, the amorphous region is the one with the highest disorder, the crystalline region is the one with the lowest disorder, and the partially amorphized region represents the transition between the two regions.

[2190-4286-13-86-8]

Figure 8: Evolution of the μ coefficient with respect to the angle for (a) pristine and (b) contaminated samples. The last data point at 55 Å has a low value in the 0.825–0.850 region due to the non-periodicity along the z axis in the simulation, causing the formation of dimers on the bottom of the sample.

Compared to other existing methods, our solution only uses the nearest neighbours of a given atom, which allows for a faster computation speed when using a large number of atoms. The results show that with increasing disorder the amorphization can be described by the change in bond lengths and bond angles. Therefore, the variation and distribution of angles in the amorphous region evolve to a point where all the angular distribution information becomes irrelevant. A methodology based on the averaged bond order parameters would hence become very time consuming for a larger-scale sample (compared to the methodology by Lechner and Dellago [36]). On the other hand, a method using a local-order metric as described by Martelli et al. [37] can be used for a whole sample. However, in the current study the information on the transition of the crystal structure between different regions cannot be obtained. For these reasons, we developed our methodology (i.e., a bond-centric analysis based on the definition of regions). Furthermore, by using a local-order metric, a high degree of overlapping would make it difficult to observe the precise region where the transition between each previously defined region would happen.

Radial distribution function

The RDF has been calculated for different slabs to support the findings of the previous section. The algorithm used to calculate RDF in finite non-periodic samples can be found in the article written by Kopera and Retsch [65]. With the previously defined slabs, it is possible to calculate the RDF for each region, showing clear differences in the local structural coherency. Figure 9 shows the RDF for each slab in pristine and contaminated samples. A first observation is the difference between RDFs of each slab: the amorphous layer has almost no coherency whereas the crystalline slabs have a strong local order and patterns repeating past the 2nd nearest neighbours (i.e., for crystalline structures the peaks above 5 Å are still well defined while different peaks are not separated for the amorphous structure). The semi-amorphous region shows some coherency after the 2nd neighbour peak but with some characteristics of amorphous RDF.

[2190-4286-13-86-9]

Figure 9: Radial distribution function for (a) clean and (b) contaminated samples for all previously determined regions. The limit fixed to 12 Å represents 2.5 lattices (or more than 10 neighbours).

We can also observe the difference between clean and contaminated samples in Figure 10: For the clean sample, a smaller peak appears between the 2nd and 3rd peaks in the 2.6 to 3.2 Å region. These bonds are due to defect production in silicon after argon irradiation. For the contaminated sample, the same peak can be observed and another peak appears in the 1.5 Å region, which is created by the Si–H and Si–O bonds due to fragmentation of water molecules. From these observations we can determine the influence of the incidence angle on the amorphization process by comparing the RDF of a slab (Figure 11). When comparing the RDF at grazing angles for the clean and contaminated samples, there is less damage directly done by argon ions. However, more damage comes from hydrogen and oxygen atoms mixed into the sample. Due to grazing incidence, most argon ions are backscattered; however, the deposited energy is high enough to split the water molecules and push hydrogen and oxygen atoms deeper into the sample. As hydrogen is very light it can travel further into the sample, channel through the lattices, and preferentially go into interstitial sites. The same reasoning applies to oxygen except that it implants at a lower depth due to a higher mass which makes channelling less likely. The Si–O bonds are also stronger than the Si–H bonds, which leads to different implantation depths for both particles: oxygen particles are caught more easily by silicon atoms and, once formed, it requires a higher energy to break stronger Si–O bonds. Hence, this leads to a shallower implantation of oxygen into the sample. At higher angles, the RDF of an amorphous slab is almost identical to that of the semi-amorphous region, meaning that the disorder created by the collisions is marginal and the variations of the µ coefficient of the slab are minimal. In the next section, the implantation depths of argon, hydrogen, and oxygen are discussed regarding their impact on the amorphization process.

[2190-4286-13-86-10]

Figure 10: Comparison of the radial distribution function for both samples in the amorphous layer. The same cut-off of 12 Å is applied.

[2190-4286-13-86-11]

Figure 11: Comparison of the radial distribution function in the amorphous slab of the contaminated sample with respect to the angle of impact. The angles selected for this study are listed above and aim to cover values between 0° and 85°.

Ion implantation

The implantation depth of each chemical element (e.g., hydrogen, oxygen, and argon) and the number of implanted atoms greatly vary with respect to the incidence angle. At grazing incidence, argon is not implanted but rather backscattered. For the contaminated sample, the argon can fragment water molecules even at grazing incidence. Hydrogen and oxygen atoms displaced at the same surface by an incident argon ion are either sputtered or implanted into the sample, causing significantly deeper damage than that observed for the pristine sample, as shown in Figure 12. At low fluences, most atoms are not involved in collision cascades. However, with increasing fluence (i.e., with the number of impacts going up to 500) almost all atoms of the surface will be displaced from their initial position and all water molecules will be eventually destroyed. After getting some energy in an atomic collision, the direction of hydrogen or oxygen atoms allows for channelling in some situations, leading to a deeper penetration into the bulk without causing damage.

[2190-4286-13-86-12]

Figure 12: Comparison between samples after irradiation with an incident beam at 100 eV and 85° for (a) the pristine sample and (b) the contaminated sample, where the silicon atoms are in the compact representation (yellow). The oxygen (red) and hydrogen (white) atoms are represented with their van der Waals interaction sphere.

The biggest difference between pristine and contaminated sample is also observed in these high angles (80° and above), where in the pristine sample most of the interactions are located on the surface, and in the contaminated sample hydrogen and oxygen are causing damage deeper in the sample. Since these angles are commonly used for the final polishing in TEM lamellae preparation, it is interesting to understand the impact of water in this specific case.

To quantify the impact of contaminants on the irradiation process, we looked at the distributions of implanted argon particles and compared them to the distribution of contaminants for each incidence angle. The result of these comparisons can be found in Figure 13. In Figure 14, the total number of retained argon atoms in the sample is plotted as a function of the angle. At normal incidence the Ar implantation depth is the highest and it decreases with the incidence angle. For angles above 80°, argon ions are not implanted at all into the sample; they are either backscattered or ejected from the sample during subsequent bombardments. Like the contaminant distribution, we can see some implantation depth deeper than the depth of the amorphous slab.

[2190-4286-13-86-13]

Figure 13: Distribution of implanted contaminants, (a) oxygen, (b) hydrogen and (c) argon, with respect to the selected angles. The values above the peaks represent the mean implantation depth. Different colours are used for different angles to increase the readability of the plots.

[2190-4286-13-86-14]

Figure 14: Number of argon atoms retained in the sample, regardless of the depth of implantation, with respect to the angle.

Furthermore, the mean implantation depth of hydrogen is twice as deep as that of oxygen. This behaviour can be explained by two factors: (i) Hydrogen atoms are involved twice as often as oxygen in displacements, due to the number of atoms of each species in the water molecule. This is true in the early stage of the sputtering process since the water molecules are intact and, therefore, the collision will either sputter or fragment the molecules. During latter events, hydrogen and oxygen might be pushed into the sample. At higher fluences, where most of the hydrogen atoms are dissociated from oxygen atoms, the probability of hydrogen being involved in a collision cascade is twice as high as that for oxygen, since the ratio between the two species stays approximately constant throughout the process. (ii) Since hydrogen is lighter than oxygen, the particle can easily channel through the silicon lattice and penetrate deeper into the sample. This phenomenon can happen if the hydrogen is impacted at the correct angle, either by initial collisions with argon ions or by recoil atoms during the collision cascade. This will allow it to travel in a straight path through the silicon lattice, reaching and causing damage in deeper regions.

Overall, the contaminant will be displaced into the amorphous region of the sample, causing this region to be thicker than that of the pristine sample. This effect can impact the preparation of samples such as ultrathin lamellae for TEM.

Reaction products

Water molecules are either fragmented or sputtered by ion bombardment so that oxygen and hydrogen atoms are implanted into the sample or sputtered. The implanted or adsorbed species and fragments will interact with the silicon atoms of the sample, forming reaction products, which contributes to the amorphization of the sample surface. These products can range from simple (e.g., Si–O or Si–H pairs) to complex (e.g., Si–H3–O products at higher incidence angles) products and are found near the surface region. Since the fragments are highly displaced in the near-surface, amorphous region and are continuously displaced by further argon impacts, the reaction products continuously change with fluence. In Figure 15, we showed the two most common products (i.e., Si–O and Si–H) and plotted the results with a sliding mean averaged over 50 points to get the trend and the evolution with respect to the incidence angle and the fluence. We observe a clear trend with a progressive reduction in the Si–O and Si–H counts for intermediate angles. Yet, at either low angles (0° – 30°) or high angle (>80°) we observe that the counts remain high and in the same order of magnitude. Hence, the implanted O and H atoms are either replaced or remain trapped in the amorphous layer. Most of the water molecules will be sputtered at higher incidence angles. The implanted species will tend to significantly move inside the amorphous region. This behaviour is due to either collisions displacing large amounts of atoms in the amorphized region, or to relaxation of previous displacement cascades. Due to the strong interaction between silicon and oxygen particles, Si–O products tend to be sputtered as a cluster or the oxygen migrates into the amorphous region.

[2190-4286-13-86-15]

Figure 15: Evolution of the Si–O (a) and Si–H (b) products with respect to the fluence and the angle.

The other relevant reaction products (e.g., Si–OH, Si–H2, and Si–O2) are displayed in Supporting Information File 1. In addition, a graph with the raw counts is displayed to give a better picture of the evolution of the reaction products.

Sputtering yields and sputtered clusters

Sputtering yields are important for experiments as they define how quickly a sample is eroded. Here, we compare the sputtering yields of pristine and contaminated samples averaged over 500 impacts (Figure 16). No significant difference can be observed between silicon sputtering yields coming from pristine and contaminated samples. As expected, both sputtering yields increase with an increasing incidence angle to reach a maximum of approx. 65°. For partial sputtering yields of hydrogen and oxygen, the curves are flatter and reach their maximum values for incidence angles of approx. 70° (i.e., for higher grazing incidence angles than those of silicon). The total sputtering yield of the contaminated sample (i.e., the sum of silicon, oxygen, and hydrogen partial sputtering yields) is much higher than that of the pristine sample, showing that water molecules are quite easily sputtered.

[2190-4286-13-86-16]

Figure 16: Sputtering yields of pristine and contaminated samples with respect to the incidence angle. For the sputtering yields of the contaminated sample, a distinction between partial Si, O, and H sputtering yields are shown in red, blue, and green, respectively, and the total sputtering yield is shown in magenta.

A second parameter of interest related to sputtering is the sputtered fraction of hydrogen and oxygen atoms (Figure 17). At normal incidence, approx. 25% of the hydrogen and oxygen atoms are sputtered, meaning that by far the largest fraction is implanted and only a few atoms stay on the sample surface. This can be explained by the fact that the velocity vector of argon ions is oriented straight into the sample, and this orientation is maintained for the energy transfer to the target atoms. Interestingly, the partial sputtering yields of hydrogen and oxygen atoms are the highest at incidence angles between 70 and 80° (i.e., above angles where the silicon sputtering yields are the highest), leading to the sputtering of approx. 75% of the adsorbed molecules. At these angles, most of the argon ions are backscattered. However, the energy transfer is still high enough to fragment the water molecules and induce the sputtering of hydrogen and oxygen atoms. Hence, if the implantation of surface contaminations is to be minimized in experiments, the incidence angle should be chosen in this range.

[2190-4286-13-86-17]

Figure 17: Fraction of the contaminants sputtered after 500 impacts with respect to the angle. The values are represented as the percentage of total particles removed by counting the contaminant particles in the initial sample and deducting the counts in the final steps for each angle.

In the current simulations, the water layer was added before the argon irradiation was started. When renewing the water layer during the irradiation process, as it happens during experiments, the partial sputtering yield of silicon might be lowered by the presence of the surface contamination, at least if the contamination level is above a critical threshold.

In order to study in our simulations the evolution of a given sample under continuous irradiation, including the effect of water contamination, we kept the energy at 100 eV and argon ions were bombarding the sample one after the other, in a dynamic sputtering mode. This differs from the literature where monatomic and cluster bombardment at high energies can be found (e.g., [66-68]). These studies also use single ion bombardment and no dynamic sputtering. Yet, despite these differences, we observed a similar range of damage and amorphization as in [68], which in our study was increased by the presence of water molecules. Furthermore, the predicted sputtering yields in [69] for energies between 200 and 500 eV at normal incidences using an angular distribution of +/−5° are bigger than those in our observations, which at normal incidence give a very low sputtering yield. Energies of approx. 100 eV are closer to the sputtering threshold, and the incidence angle plays a critical role as our data, at slightly off-normal conditions, matches their predicted yields.

A third parameter of interest is the distribution of the sputtered clusters. Overall, several cases are observed: (i) an atom is sputtered alone, (ii) several atoms are sputtered during the same collision cascade independently of one another, and (iii) a cluster of atoms is sputtered. In Figure 18 the evolution of the sputtered clusters is plotted with respect to the incidence angle. The highest sputtering yield is observed for intact water molecules, usually occurring in the first simulations when the water layer is relatively undamaged. At higher incidence angles, the amount of water molecules which is sputtered intact is the highest: since the argon ion is at grazing incidence, the energy transfer is efficient for molecular lift-off and far less for bond breaking, which predominantly happens towards normal incidence. This explains also the increased amount of sputtered O–H clusters at smaller angles where the argon ions more aggressively damage the water layer, increasing the probability of water fragmentation. For silicon, the most common sputtered cluster is Si2. Its sputtering yield follows the same trend of Si. The second most common fragment for silicon is SiO, which is especially stable due to the nature and the strength of the Si–O bonds. Finally, there are low counts for Si–H and more complex clusters. They represent an almost negligible fraction of the total counts of the sputtered particles.

[2190-4286-13-86-18]

Figure 18: Sputtering yields for (a) silicon-related and (b) oxygen-related clusters.

For most simulations there is no intact water molecule remaining at the end of the sputtering events. This is shown in Figure 19 which describes the number of intact clusters with respect to fluence or to the number of impacts, respectively. From this observation we can deduce that most damages occur during the first irradiation simulations, either by sputtering or fragmentation. Once most of the water molecules have been fragmented, the atomic mixing is facilitated by the presence of fragments in the amorphous region. Furthermore, at higher incidence angles we observe an increased rate of water dissociation.

[2190-4286-13-86-19]

Figure 19: Graph showing the probability to dissociate at least one water cluster per impact, with respect to angle and fluence. In some cases, one collision could dissociate several clusters of water.

Conclusion

In this work, we analysed the influence of the incidence angle on 100 eV sputtering processes. We used MD simulations for irradiation with argon, with fluences up to 2.6 × 1015 ions/cm2 and ions of silicon samples with a clean surface and covered with a monolayer of water molecules. Contaminants adsorbed on the sample surface amplify the damage done by the argon ions. A methodology based on the modification of bond lengths has been developed to characterize the irradiation-induced damage. Three regions have been identified in the sample: (i) a totally amorphous layer closer to the surface, (ii) a transitive, partially amorphized layer presenting defects, and (iii) the crystalline region, which remains unaffected underneath by ion irradiation. Closer to grazing incidence (e.g., above 75°), argon implantation becomes sparse; however, fragmentation and sputtering of water molecules reach maximum values. Hence, sputtering yields and molecular fragmentation show that specific angles can be selected to reduce the mixing of contaminants into the sample in ultralow energy sputtering. Compared to the pristine sample, the mixing of hydrogen and oxygen atoms into the sample increases the overall degree of amorphization. Furthermore, the water layer on the silicon surface was prepared before ion irradiation and not renewed during the milling process. However, in experiments, contaminations are continuously adsorbed on the sample surface, leading to a renewal of the contamination layer. Furthermore, at higher impact energy the implantation depth of the ions is increased, leading most likely to a different amorphization process. Both aspects will be explored in future studies.

Supporting Information

Supporting Information File 1: Additional figures and tables with DFT data fitting plots, g(r) per slab and grouping details, reaction products, and ReaxFF potential.
Format: PDF Size: 3.8 MB Download

Acknowledgements

Dr. Chad Rue, Dr. Richard J. Young from Thermo Fisher Scientific and Dr. Tom Wirtz from Luxembourg Institute of Science and Technology are acknowledged and thanked for all the fruitful discussions and their comments during the writing process.

Funding

This project was supported by the National Research Fund (FNR), Luxembourg (C17/MS/11682850/ULOWBEAM).

References

  1. Wirtz, T.; De Castro, O.; Audinot, J.-N.; Philipp, P. Annu. Rev. Anal. Chem. 2019, 12, 523–543. doi:10.1146/annurev-anchem-061318-115457
    Return to citation in text: [1]
  2. Pillatsch, L.; Östlund, F.; Michler, J. Prog. Cryst. Growth Charact. Mater. 2019, 65, 1–19. doi:10.1016/j.pcrysgrow.2018.10.001
    Return to citation in text: [1]
  3. Audinot, J.-N.; Philipp, P.; De Castro, O.; Biesemeier, A.; Hoang, Q. H.; Wirtz, T. Rep. Prog. Phys. 2021, 84, 105901. doi:10.1088/1361-6633/ac1e32
    Return to citation in text: [1]
  4. Wirth, R. Chem. Geol. 2009, 261, 217–229. doi:10.1016/j.chemgeo.2008.05.019
    Return to citation in text: [1]
  5. Thermoscientific 2018, 1–4. https://assets.thermofisher.com/TFS-Assets/MSD/Application-Notes/atom-probe-tomography-sample-prep-plasma-fib-appnote.pdf.
    Return to citation in text: [1]
  6. McPhail, D. S. J. Mater. Sci. 2006, 41, 873–903. doi:10.1007/s10853-006-6568-x
    Return to citation in text: [1]
  7. Tseng, A. A. Small 2005, 1, 924–939. doi:10.1002/smll.200500113
    Return to citation in text: [1]
  8. Kawasegi, N.; Yamaguchi, M.; Kozu, T.; Morita, N.; Nishimura, K. Precis. Eng. 2017, 50, 337–343. doi:10.1016/j.precisioneng.2017.06.007
    Return to citation in text: [1]
  9. Hasan, R. M. M.; Luo, X. Nanolithography: Status and challenges. In 2017 23rd International Conference on Automation and Computing (ICAC), 2017; pp 1–6. doi:10.23919/iconac.2017.8081979
    Return to citation in text: [1]
  10. Sloyan, K.; Melkonyan, H.; Apostoleris, H.; Dahlem, M. S.; Chiesa, M.; Al Ghaferi, A. Nanotechnology 2021, 32, 472004. doi:10.1088/1361-6528/ac1d75
    Return to citation in text: [1]
  11. Llobet, J.; Krali, E.; Wang, C.; Arbiol, J.; Jones, M. E.; Pérez-Murano, F.; Durrani, Z. A. K. Appl. Phys. Lett. 2015, 107, 223501. doi:10.1063/1.4936757
    Return to citation in text: [1]
  12. Córdoba, R.; Mailly, D.; Rezaev, R. O.; Smirnova, E. I.; Schmidt, O. G.; Fomin, V. M.; Zeitler, U.; Guillamón, I.; Suderow, H.; De Teresa, J. M. Nano Lett. 2019, 19, 8597–8604. doi:10.1021/acs.nanolett.9b03153
    Return to citation in text: [1]
  13. De Teresa, J. M.; Orús, P.; Córdoba, R.; Philipp, P. Micromachines 2019, 10, 799. doi:10.3390/mi10120799
    Return to citation in text: [1]
  14. Philipp, P.; Ngo, Q. K.; Shtein, M.; Kieffer, J.; Wirtz, T. Anal. Chem. (Washington, DC, U. S.) 2013, 85, 381–388. doi:10.1021/ac302939m
    Return to citation in text: [1]
  15. Prokhodtseva, A.; Mulders, J.; Vystavel, T. Microsc. Microanal. 2017, 23 (Suppl. 1), 298–299. doi:10.1017/s1431927617002173
    Return to citation in text: [1] [2]
  16. Lécz, Z.; Budai, J.; Andreev, A.; Ter-Avetisyan, S. Phys. Plasmas 2020, 27, 013105. doi:10.1063/1.5123542
    Return to citation in text: [1]
  17. Azuma, Y.; Kurokawa, A. X-Ray Spectrom. 2019, 48, 345–350. doi:10.1002/xrs.3006
    Return to citation in text: [1]
  18. Barry, P. R.; Philipp, P.; Wirtz, T.; Kieffer, J. J. Mass Spectrom. 2014, 49, 185–194. doi:10.1002/jms.3317
    Return to citation in text: [1]
  19. Mach, J.; Šamořil, T.; Voborný, S.; Kolíbal, M.; Zlámal, J.; Spousta, J.; Dittrichová, L.; Šikola, T. Rev. Sci. Instrum. 2011, 82, 083302. doi:10.1063/1.3622749
    Return to citation in text: [1]
  20. Rzeznik, L.; Fleming, Y.; Wirtz, T.; Philipp, P. Beilstein J. Nanotechnol. 2016, 7, 1113–1128. doi:10.3762/bjnano.7.104
    Return to citation in text: [1]
  21. Yoshigoe, A.; Teraoka, Y. Surf. Sci. 2003, 532–535, 690–697. doi:10.1016/s0039-6028(03)00412-6
    Return to citation in text: [1] [2]
  22. Pal, G.; Yadav, R. C.; Akhter, J.; Das, T.; Sarkar, A.; Mallik, C.; Bhandari, R. K. J. Phys.: Conf. Ser. 2012, 390, 012045. doi:10.1088/1742-6596/390/1/012045
    Return to citation in text: [1]
  23. Forrest, E.; Schulze, R.; Liu, C.; Dombrowski, D. Int. J. Heat Mass Transfer 2015, 91, 311–317. doi:10.1016/j.ijheatmasstransfer.2015.07.112
    Return to citation in text: [1]
  24. Zhang, W.; van Duin, A. C. T. J. Phys. Chem. B 2017, 121, 6021–6032. doi:10.1021/acs.jpcb.7b02548
    Return to citation in text: [1] [2]
  25. Yeon, J.; van Duin, A. C. T.; Kim, S. H. Langmuir 2016, 32, 1018–1026. doi:10.1021/acs.langmuir.5b04062
    Return to citation in text: [1]
  26. Fogarty, J. C.; Aktulga, H. M.; Grama, A. Y.; van Duin, A. C. T.; Pandit, S. A. J. Chem. Phys. 2010, 132, 174704. doi:10.1063/1.3407433
    Return to citation in text: [1]
  27. Smith, R.; Harrison, D. E.; Garrison, B. J. Phys. Rev. B 1989, 40, 93–101. doi:10.1103/physrevb.40.93
    Return to citation in text: [1]
  28. Smith, R.; Beardmore, K.; Gras-Marti, A.; Kirchner, R.; Webb, R. P. Nucl. Instrum. Methods Phys. Res., Sect. B 1995, 102, 211–217. doi:10.1016/0168-583x(95)80143-a
    Return to citation in text: [1]
  29. Stillinger, F. H.; Weber, T. A. Phys. Rev. Lett. 1989, 62, 2144–2147. doi:10.1103/physrevlett.62.2144
    Return to citation in text: [1]
  30. Senftle, T. P.; Hong, S.; Islam, M. M.; Kylasa, S. B.; Zheng, Y.; Shin, Y. K.; Junkermeier, C.; Engel-Herbert, R.; Janik, M. J.; Aktulga, H. M.; Verstraelen, T.; Grama, A.; van Duin, A. C. T. npj Comput. Mater. 2016, 2, 15011. doi:10.1038/npjcompumats.2015.11
    Return to citation in text: [1]
  31. Abolfath, R. M.; van Duin, A. C. T.; Brabec, T. J. Phys. Chem. A 2011, 115, 11045–11049. doi:10.1021/jp204894m
    Return to citation in text: [1]
  32. Kowalik, M.; Hossain, M. J.; Lele, A.; Zhu, W.; Banerjee, R.; Granzier-Nakajima, T.; Terrones, M.; Hudson, E. W.; van Duin, A. C. T. Catalysts 2021, 11, 208. doi:10.3390/catal11020208
    Return to citation in text: [1]
  33. Chenoweth, K.; van Duin, A. C. T.; Goddard, W. A. J. Phys. Chem. A 2008, 112, 1040–1053. doi:10.1021/jp709896w
    Return to citation in text: [1]
  34. Newsome, D. A.; Sengupta, D.; Foroutan, H.; Russo, M. F.; van Duin, A. C. T. J. Phys. Chem. C 2012, 116, 16111–16121. doi:10.1021/jp306391p
    Return to citation in text: [1] [2] [3]
  35. Soria, F. A.; Zhang, W.; Paredes-Olivera, P. A.; van Duin, A. C. T.; Patrito, E. M. J. Phys. Chem. C 2018, 122, 23515–23527. doi:10.1021/acs.jpcc.8b07075
    Return to citation in text: [1]
  36. Lechner, W.; Dellago, C. J. Chem. Phys. 2008, 129, 114707. doi:10.1063/1.2977970
    Return to citation in text: [1] [2]
  37. Martelli, F.; Ko, H.-Y.; Oğuz, E. C.; Car, R. Phys. Rev. B 2018, 97, 064105. doi:10.1103/physrevb.97.064105
    Return to citation in text: [1] [2]
  38. van Duin, A. C. T.; Strachan, A.; Stewman, S.; Zhang, Q.; Xu, X.; Goddard, W. A. J. Phys. Chem. A 2003, 107, 3803–3811. doi:10.1021/jp0276303
    Return to citation in text: [1]
  39. Costa Filho, R. N.; Alencar, G.; Skagerstam, B.-S.; Andrade, J. S., jr. EPL 2013, 101, 10009. doi:10.1209/0295-5075/101/10009
    Return to citation in text: [1]
  40. Kresse, G.; Hafner, J. Phys. Rev. B 1993, 47, 558–561. doi:10.1103/physrevb.47.558
    Return to citation in text: [1]
  41. Kresse, G.; Hafner, J. Phys. Rev. B 1994, 49, 14251–14269. doi:10.1103/physrevb.49.14251
    Return to citation in text: [1]
  42. Kresse, G.; Furthmüller, J. Comput. Mater. Sci. 1996, 6, 15–50. doi:10.1016/0927-0256(96)00008-0
    Return to citation in text: [1]
  43. Kresse, G.; Furthmüller, J. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 11169–11186. doi:10.1103/physrevb.54.11169
    Return to citation in text: [1]
  44. Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865–3868. doi:10.1103/physrevlett.77.3865
    Return to citation in text: [1]
  45. Blöchl, P. E. Phys. Rev. B 1994, 50, 17953–17979. doi:10.1103/physrevb.50.17953
    Return to citation in text: [1]
  46. Kresse, G.; Joubert, D. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 1758–1775. doi:10.1103/physrevb.59.1758
    Return to citation in text: [1]
  47. Tkatchenko, A.; Scheffler, M. Phys. Rev. Lett. 2009, 102, 073005. doi:10.1103/physrevlett.102.073005
    Return to citation in text: [1]
  48. pair _ style hybrid command, LAMMPS documentation. https://docs.lammps.org/pair_hybrid.html (accessed Sept 5, 2022).
    Return to citation in text: [1]
  49. Rappe, A. K.; Goddard, W. A., III. J. Phys. Chem. 1991, 95, 3358–3363. doi:10.1021/j100161a070
    Return to citation in text: [1]
  50. Nakano, A. Comput. Phys. Commun. 1997, 104, 59–69. doi:10.1016/s0010-4655(97)00041-6
    Return to citation in text: [1]
  51. Aktulga, H. M.; Fogarty, J. C.; Pandit, S. A.; Grama, A. Y. Parallel Comput. 2012, 38, 245–259. doi:10.1016/j.parco.2011.08.005
    Return to citation in text: [1]
  52. Plimpton, S. J. Comput. Phys. 1995, 117, 1–19. doi:10.1006/jcph.1995.1039
    Return to citation in text: [1]
  53. Grønbech-Jensen, N.; Farago, O. Mol. Phys. 2013, 111, 983–991. doi:10.1080/00268976.2012.760055
    Return to citation in text: [1]
  54. Kamiyama, E.; Sueoka, K. J. Electrochem. Soc. 2010, 157, H323. doi:10.1149/1.3281184
    Return to citation in text: [1]
  55. Chen, J.; Zhang, G.; Li, B. J. Phys. Soc. Jpn. 2010, 79, 074604. doi:10.1143/jpsj.79.074604
    Return to citation in text: [1]
  56. Simonton, R. B.; Kamenitsa, D. E.; Ray, A. M.; Park, C.; Klein, K. M.; Tasch, A. F. Nucl. Instrum. Methods Phys. Res., Sect. B 1991, 55, 39–44. doi:10.1016/0168-583x(91)96132-5
    Return to citation in text: [1]
  57. Stegmann, H.; Ritz, Y.; Utess, D.; Hübner, R.; Zschech, E. AIP Conf. Proc. 2009, 1173, 275–279. doi:10.1063/1.3251233
    Return to citation in text: [1]
  58. Ren, Y.; Lei, D.; Yao, F.; Wang, Z. AIP Conf. Proc. 2019, 2126, 120017. doi:10.1063/1.5117635
    Return to citation in text: [1]
  59. Abdel-Samad, S.; Abdel-Bary, M.; Kilian, K. Vacuum 2005, 78, 83–89. doi:10.1016/j.vacuum.2005.01.001
    Return to citation in text: [1]
  60. van der Heide, P. Secondary Ion Mass Spectrometry; John Wiley & Sons: Hoboken, NJ, USA, 2014; Vol. 407. doi:10.1002/9781118916780
    Return to citation in text: [1]
  61. Demaison, J.; Herman, M.; Lievin, J. Int. Rev. Phys. Chem. 2007, 26, 391–420. doi:10.1080/01442350701371919
    Return to citation in text: [1]
  62. Tuttle, B.; Van de Walle, C. G. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 12884–12889. doi:10.1103/physrevb.59.12884
    Return to citation in text: [1]
  63. Baur, W. H. Acta Crystallogr., Sect. B: Struct. Crystallogr. Cryst. Chem. 1977, 33, 2615–2619. doi:10.1107/s0567740877009029
    Return to citation in text: [1]
  64. Niu, C.-Y.; Wang, J.-T. J. Chem. Phys. 2013, 139, 194709. doi:10.1063/1.4832340
    Return to citation in text: [1]
  65. Kopera, B. A. F.; Retsch, M. Anal. Chem. (Washington, DC, U. S.) 2018, 90, 13909–13914. doi:10.1021/acs.analchem.8b03157
    Return to citation in text: [1]
  66. Samela, J.; Nordlund, K.; Keinonen, J.; Popok, V. N. Nucl. Instrum. Methods Phys. Res., Sect. B 2007, 255, 253–258. doi:10.1016/j.nimb.2006.11.080
    Return to citation in text: [1]
  67. Allen, L. P.; Insepov, Z.; Fenner, D. B.; Santeufemio, C.; Brooks, W.; Jones, K. S.; Yamada, I. J. Appl. Phys. 2002, 92, 3671–3678. doi:10.1063/1.1506422
    Return to citation in text: [1]
  68. Samela, J.; Nordlund, K.; Keinonen, J.; Popok, V. N.; Campbell, E. E. B. Eur. Phys. J. D 2007, 43, 181–184. doi:10.1140/epjd/e2007-00104-y
    Return to citation in text: [1] [2]
  69. Zalm, P. C. J. Appl. Phys. 1983, 54, 2660–2666. doi:10.1063/1.332340
    Return to citation in text: [1]

© 2022 Defoort-Levkov et al.; licensee Beilstein-Institut.
This is an open access article licensed under the terms of the Beilstein-Institut Open Access License Agreement (https://www.beilstein-journals.org/bjnano/terms), which is identical to the Creative Commons Attribution 4.0 International License (https://creativecommons.org/licenses/by/4.0). The reuse of material under this license requires that the author(s), source and license are credited. Third-party material in this article could be subject to other licenses (typically indicated in the credit line), and in this case, users are required to obtain permission from the license holder to reuse the material.

Other Beilstein-Institut Open Science Activities
Journal Logo
Institut Logo
Preprint Logo
Symposia Logo