# Ab-Initio Simulation of Field Evaporation

Jiayuwen Qi  Dept. of Materials Science and Engineering, The Ohio State
University, Columbus, OH, USA  Christian Oberdorfer  Dept. of Materials
Science and Engineering, The Ohio State University, Columbus, OH, USA
Emmanuelle A. Marquis  Dept. of Materials Science and Engineering, University
of Michigan, Ann Arbor, MI, USA  Wolfgang Windl<EMAIL_ADDRESS>Dept. of
Materials Science and Engineering, The Ohio State University, Columbus, OH,
USA

###### Abstract

A new simulation approach of field evaporation is presented. The model
combines classical electrostatics with molecular dynamics (MD) simulations.
Unlike previous atomic-level simulation approaches, our method does not rely
on an evaporation criterion based on thermal activation theory, instead,
electric-field-induced forces on atoms are explicitly calculated and added to
the interatomic forces. Atoms then simply move according to the laws of
classical molecular dynamics and are “evaporated” when the external force
overcomes interatomic bonding. This approach thus makes no ad-hoc assumptions
concerning evaporation fields and criteria, which makes the simulation fully
physics-based and “ab-initio” apart from the interatomic potential. As proof
of principle, we perform simulations to determine material dependent critical
voltages which allow assessing the evaporation fields and the corresponding
steady-state tip shapes in different metals. We also extract critical
evaporation fields in elemental metals and sublimation energies in a high
entropy alloy to have a more direct comparison with tabulated values. In
contrast to previous approaches, we show that our method is able to
successfully reproduce the enhanced zone lines observed in experimental field
desorption patterns. We also demonstrate the need for careful selection of the
interatomic potential by a comparative study for the example of Cu-Ni alloys.

††preprint: APS/123-QED

##  I Introduction

Field evaporation is a field-assisted atom removal process in high-electric
field nanoscience with applications in liquid ion sources or materials
characterization. During field evaporation, surface atoms of a specimen are
removed from their original lattice sites and ionized by virtue of a large
electric field [1]. Field evaporation is the central process in atom probe
tomography (APT), where a high voltage is applied to a needle-shaped specimen,
creating an electric field strong enough to evaporate atoms from the apex.
These atoms are then accelerated towards a position-sensitive detector. Their
impact time and positions are collected for reconstruction of the data, where
the time of flight allows determining the specific atom identity while the
combination of arrival time with the $(x,y)$ detector coordinates are used for
reconstruction.

Artifacts are often present and interfere with the fidelity of the
reconstruction [2], especially in heterogeneous materials and complex
microstructures [3]. However, since the field evaporation process is fully
destructive, it is fundamentally impossible to validate the truthfulness of
APT reconstruction experimentally. The only viable way to quantify fidelity
and uncertainty of reconstruction is via a fully physical computational
forward model of field evaporation where every atom is traceable.
Consequently, simulation of field evaporation has become an essential tool to
detect and interpret artifacts and improve the quality of reconstructed data
[4]. Physical modelling and simulation of field evaporation has implications
beyond the field of atom probe tomography, and can contribute to understanding
the electrochemical behavior of surfaces during corrosion [5, 6] and
electrocatalysis [7].

A number of simulation approaches have been developed since the advent of
field ion microscopy in the 1950s spanning from continuum to atomistic scales
[8]. Initially, simulations were performed to unravel the order of evaporation
and the evolution of the surface morphology [9, 10]. These simulations were
based purely on geometry without including field effects. Later continuum
models using the level-set method [11, 12, 13] were developed to simulate the
evolution of the tip shape in field evaporation [14], which provides an
alternative way of guiding APT reconstruction by considering the aberration of
the tip shape from the hemispherical ideal assumption. These continuum models
are computationally efficient and able to model realistic sample geometries.
However, because the atomic structure is not included, they are not able to
capture the atomic behavior or crystallographic features of the sample such as
terraces and facets which are closely related to the trajectories of the
evaporated atoms and formation of poles and zone lines in the desorption
patterns.

At the other end of the simulation spectrum, first-principles studies of field
evaporation within density-functional theory (DFT) can be used to accurately
calculate with quantum mechanical approaches desorption energies without and
with applied electric field and thus determine the critical evaporation field
where the evaporation barrier goes to zero [15, 16, 17, 18, 19, 20, 21, 22].
However, the size of the structural models used in DFT calculations is limited
to at most a few hundred atoms and very short time scales, which makes it
impossible to examine the full physical process of field evaporation for a
realistic sample with realistic electric field. It also restricts these
calculations to simplified, slab-like structures.

To fill the gap between the DFT and continuum methods, Vurpillot et al.
pionnered simulations of a sample with a crystallographic structure and
realistic geometry and field distribution [23]. The atomic structure of an
elemental sample tip was modeled by a 3D rigid cubic grid, while the
surrounding electrostatic potential was determined within classical
electrostatics from Laplace’s equation. The general procedure had three steps:
(1) removal of the surface atom with the highest electric field; (2)
recalculation of the electric field distribution for the updated geometry; and
(3) calculation of the trajectory of the evaporated atom by integration of
Newton’s equation of motion. Building on this approach, Geiser et al.
introduced grids with different sizes, making the model compatible with any
cubic structure emitter [24]. Oberdorfer et al. replaced the rigid lattice
structure with adaptive unstructured meshes based on Voronoi cells, where the
mesh at the tip is determined by the atom positions with increasing coarseness
away from the tip [25]. Beyond Vurpillot’s approach whose use of the Laplace
equation restricted modeling materials with infinite permittivity, Oberdorfer
solved Poisson’s equation to compute the electrostatic potential and electric
field. This model now can be applied to any structure with any surface
morphology [26].

Several variations to this general approach have been examined. These include
solving the electrostatic field from the Robin equation [27], which has not
been widely applied. A kinetic Monte Carlo (KMC) approach to add randomness to
the choice of the evaporation event found little effect at experimental
temperatures [28]. More impactful was the coupling of molecular dynamics (MD)
with the finite element method in different simulation approaches [29, 30, 31,
32]. There, instead of a “static” tip structure, a “dynamic” tip structure is
thermally equilibrated for several ps after an atom is removed from the tip.
Through the addition of surface dynamics, several artefacts could be
reproduced and explained. For example, athermal surface migration of solute
atoms with much stronger bonds [31], and rearrangement of the electrostatic
field following the field evaporation of atoms [33].

Although the current state of the art in modeling has increased realism by
adding atomic structure, realistic field, and atomic motion, one of the most
crucial approximations, i.e. the convenient but unphysical separation of
dynamics from the process of evaporation, has yet to be addressed. To trigger
an evaporation event, an ad-hoc criterion is usually introduced, such as
comparing the local field at the site of an atom to a tabulated zero-barrier
evaporation field (ZBEF) [25]. However, Yao et al. showed that ZBEF is site-
dependent [19]. Since ZBEF values are typically related to an activation
energy calculated from sublimation energy, workfunction, and ionization
energy, Parviainen et al. introduced a first site-dependence by determining
location-dependent binding energy with classical MD, while keeping ionization
energy and work function as global constants [30]. Besides the difficulty of
obtaining site-dependent activation energies, evaporation models also ignore
the initial detachment velocity of an evaporated atom. It is usually assumed
that the initial detachment velocity is zero, thus the atom trajectory always
follows the local electric field. Furthermore, it is assumed that an initial
velocity component not parallel to the field would be small and therefore have
negligible effects on the detector image except for some minor blurring [34].
In contrast, we show that the component of the electrostatic force not
parallel to the interatomic force on the evaporating atom results in
noticeable lateral acceleration, and gives rise to detector patterns with
enhanced zone lines as we will discuss below.

To address these issues, we present a new field evaporation simulation
approach – “TAPSim-MD”. It combines the finite element APT simulation package
TAPSim [25] with the popular molecular dynamics simulation package LAMMPS
(Large-scale Atomic/Molecular Massively Parallel Simulator) [35] and replaces
criterion-based evaporation with dynamic, force-based evaporation. The
simulation is then executed in the form of an MD run where the atoms move
under the influence of net forces consisting of the sum of the usual
interatomic forces and the electric-field-induced forces. This eliminates the
ad-hoc assumptions of an evaporation criterion, the lack of location
dependence of the zero-barrier evaporation field, and the zero launch velocity
that were part of all previous approaches. Here and in subsequent publications
we will show that this approach is more realistic and can provide new
understanding of composition dependence of evaporation field, local-
environment effects, detailed evaporation mechanisms, their effects on the
observed detector images, and artifacts in reconstruction.

##  II Methodology

The electrostatic force calculation follows Oberdorfer’s method [25], where
the key characteristic was the meshing between Voronoi cells in the vacuum
surrounding the tip and the atomic tip structure, where each atom was chosen
as the center of a Voronoi cell, while the mesh becomes coarser away from the
tip for computational efficiency. Inaccuracies from meshing can be minimized
due to the adaptive distribution of Voronoi cells in addition to the emitter
cells. After the mesh is established, the electrostatic potential is solved
from Poisson’s equation on the mesh following the approach described in [26].
The electric field is calculated from the electrostatic potential through the
numerical differentiation scheme for the irregular grid using the method based
on interpolation between the generator points outlined in [25]. The
polarization charge on each atom is then calculated from Gauss’s law (even
below the surface in case the material is not a metal and has a finite
dielectric constant), and the electrostatic force on each atom is obtained
from the product of the charge and the electric field. In metals that we
examine here, the dielectric constant is infinite, which we approximate by a
very large value of $10^{6}$ to keep the numerics stable. As a side note, care
needs to be taken for materials with strongly field-dependent dielectric
constant such as SrTiO3, where the dielectric constant in the range of
evaporation fields is more than one order of magnitude smaller than in the
field-free case [36].

Oberdorfer’s method uniquely determines realistic field-induced force
contributions, which when added to the interatomic forces, allows an “ab-
initio” simulation of field evaporation without specified evaporation
criterion. Thus, evaporation is not postulated through a criterion, but
observed as a process happening naturally during an MD simulation, where an
atom is pulled of the surface by the field-induced force. This process can be
detected through vanishing of interatomic forces or by defining a critical
velocity for an atom. The sound velocity of the studied material is a
convenient value. This thus eliminates the “selection of surface atoms for
desorption” described in [25]. Calculations of trajectories after evaporation
and recording of the data then follow again the procedure outlined in [25].

While this approach may sound straightforward, marrying two major codes,
especially from the two fundamentally different domains of finite-element
modeling and atomic simulation required a large number of handshake steps and
streamlining, especially since computational performance is essential to allow
for tips that are at least close to realistic dimensions in size.

In our MD simulations, we use the microcanonical ensemble (NVE), a Langevin
thermostat with a damp parameter of 2.50, and the default timestep of 1 fs for
metal units. In order to fix the position of the tip base, we set for the
atoms in the bottom layer of the virtual tip all velocity components to zero,
as well as the $z$-components of their forces. In addition, the “fix recenter”
command is used to constrain the center-of-mass position of the bottom layer
in the $xy$-plane.

Regarding the choice of interatomic potentials, we use embedded atom model
(EAM) potentials in all but one simulations in this paper to balance the
accuracy and speed of computation, since the simulations are computationally
rather intensive (see Sec. V). EAM potentials are known to work very well for
close-packed structures and can also work for bcc when the cutoff is long
enough. When choosing a potential, we first pick a series of viable candidate
potentials based on if their fitting database includes quantities relevant to
field evaporation and run test simulations. Then we make decisions based on
the assessment of the simulations which should be stable and reasonable. As an
example for potential choice, we discuss the process for Cu-Ni alloys in Sec.
IV.

Based on this, we ended up with potentials by Mishin and co-workers for the
fcc-metals Al [37], Cu [38], and Ni [37] as well as for hcp Co [39]. For Au,
our potential of choice is by Foiles, Baskes, and Daw [40]. For Ti, which is
more challenging becuase of its high-temperature phase transformation from hcp
to bcc combined with unique mechanical properties, we chose a very recent
potential by Mendelev et al. [41]. While a viable EAM potential was found for
bcc W [42], a comparable potential could not be found for bcc Fe, which then
was the only case where we had to resort to a MEAM potential [43], whose
performance due to the included three-body terms is however markedly slower.
For the high entropy alloy (HEA), there was only one potential available from
Farkas and Caro [44]. For validation, we also calculated elemental processes
with it as was shown in Fig. 6 and found that they compared favorably to the
expectation from Tsong [45].

##  III Validation of the model

To validate the new capabilities of the proposed approach, we first
demonstrate the electrostatic potential and electric field obtained by
Poisson’s equation solver. We then compare steady-state tip shapes for
different crystal structures of metals and discuss values of critical voltages
that for the first time “TAPSim-MD” is able to quantitatively evaluate, which
allows assessment of evaporation fields across different simulations. We show
that the derived average critical evaporation fields correlate strongly with
tabulated ZBEFs that are used in classical reconstruction. Alloying effects
are investigated in terms of sublimation energies for a high entropy alloy in
comparison to its constituent elemental metals. We also show that comparison
between the correctly modelled interactomic forces and electric-field-induced-
forces explains the origin of enhanced zones lines in field desorption
patterns that for the first time successfully reproduce experimental
observations.

###  III.1 Electrostatic potential and electric field

The electrostatic potential and electric field distributions are examined in a
W tip with the size of 10 nm in radius and 27 nm in height. Figure 1 shows the
calculated (a) electric field and (b) electrostatic potential. The potential
drops rapidly from tip’s surface to the vacuum space by about 75% within 44 nm
along the central axis, which results in a strong electric field around the
tip. The strongest field with the warmest color is around the top surface of
the tip, which is about twice of the field around the side surface of the
shaft.

Migunov et al. experimentally determined the electric field and electrostatic
potential around a charged tip from a numerical analysis of charge density
measurements on the tip surface from electron holography [46]. The
distribution of the electrostatic potential and electric field around the atom
probe needle is reprinted in Fig. 1c. The apex radius of the tip used in the
experiment was about 15.6 nm, which is comparable to the virtual tip used in
our simulation. The experimental results showed that the potential drops from
tip’s surface to the vacuum space by about 75% within 25 nm along the central
axis, and the strongest field is found around the top surface of the tip which
is about 4 times of the field around the side surface of the shaft.
Considering the considerable uncertainties in charge density measurement
paired with the strong approximations underlying the numerical analysis, which
include perfect cylindrical symmetry of the sample and ellipsoidal shape of
the tip, the agreement between our simulations and experiment is good.

Figure 1:  Contour plots of the (a) electric field distribution and (b)
electrostatic potential distribution for the center of a 10 nm-radius virtual
W tip. The potential in (b) between two nodes of the mesh is calculated by
linear interpolation. The electric field in (a) is computed from the gradient
of the potential. (c) is the distribution of the electrostatic potential and
electric field around the atom probe needle from electron holography
measurements. The white lines represent electric field lines. The colors
correspond to equipotential contours. Figure (c) is reprinted from [46].

###  III.2 Steady-state tip shape

In APT experiments, the data collected until a steady-state tip shape is
reached are usually discarded. Once the steady-state shape is reached,
evaporation of the surface atoms happens typically layer by layer, resulting
in a periodically returning tip morphology. Steady-state shapes of virtual
tips are also observed in our simulations. Figure 2 includes snapshots of the
steady-state shapes of tips of the three most common crystal structures for
metals in $\langle 001\rangle$ (or $\langle 0001\rangle$) orientation, bcc W,
fcc Al, and hcp Ti. The size of the virtual tips is 4 nm in radius. The
surface atoms are colored by the magnitude of the electric-field-induced
forces (external forces). As expected, the field is concentrated on the edges
of the atomic terraces. Therefore, the evaporation of surface atoms always
tends to progress from the edges towards the center, which causes faceting.
Moreover, the prominent lattice planes revealed by the top view of the tips in
Fig. 2 follow the relative prominence rule that states that planes with the
largest interplanar spacings are the most prominent [47]. For bcc, the most
prominent planes are $\\{110\\}$, and the second prominent planes are
$\\{002\\}$; for fcc, the most prominent planes are $\\{111\\}$, and the
second prominent planes are $\\{002\\}$; for hcp, the most prominent planes
are $\\{0002\\}$, and the second prominent planes are $\\{1\overline{1}01\\}$
(the most prominent planes in general should be $\\{01\overline{1}0\\}$,
however they are the shaft surface of the cylinder, therefore are not
considered here). In Fig. 2, we see indeed $\\{110\\}$ and $\\{002\\}$ for W,
$\\{111\\}$ and $\\{002\\}$ for Al and $\\{0002\\}$ and
$\\{1\overline{1}01\\}$ for Ti virtual tip as prominent planes respectively.

Figure 2: Steady-state tip shapes of different metals. Virtual tip’s size is 4
nm in radius. Surface atoms are colored by the magnitude of the electric-
field-induced forces. The top view of tips are displayed in the second row.

###  III.3 Material-dependent critical voltage

The applied voltage in previous simulation approaches does not have a
quantitative meaning, since codes like TAPSim [25] evaporate atoms based on a
local field-strength criterion such as the relative strength of the local
field vs. the tabulated evaporation field of the respective chemical species
without considering the local bonding environment. Thus, the voltage is
usually set to some reasonable value to create the electric field
distribution, irrespective of the material examined.

In our approach, evaporation results from the competition between interatomic
forces and field-induced forces, which gives a physical meaning and material
dependence to the applied voltage. In order to trigger a single field
evaporation event at a time, the applied voltage is controlled to increase
first rapidly by 10 V per simulation step then slowly by 0.1 V per simulation
step until the first atom evaporates. After that, it is reduced successively
by 0.1 V per step until no more evaporation events are detected, and ramped up
by 0.1 V per step again until the next evaporation event happens. Figure 3
shows that after some initially high values which we call “over voltage”
originating from the initially chosen tip shape as discussed in Sec. III.3.2,
the evaporation voltage oscillates around a constant value which defines the
“critical voltage”. Its value slightly decreases with decreasing tip’s height
(keeping the field constant). Since the interatomic forces are determined from
carefully fitted interatomic potentials and have realistic values, different
materials with different bond strengths will require different voltage values.
Also, the same atomic species may require different voltage values when
presents in different alloys or compounds where the atomic bonding strength
changes and depends on their coordination. Thus, the critical voltage has a
physical meaning and can be benchmarked against experimental observations.

####  III.3.1 Critical voltage results

The results in Fig. 3 are from virtual tips in the shape of a hemisphere on
top of a cylinder to mimic the idealized shape assumed for classical
reconstruction. The size of the virtual tips is 4 nm in radius and 12 nm in
total height. The $\langle 001\rangle$ direction of the tip ($\langle
0001\rangle$ direction for HCP structure) is oriented along the $z$ axis.
Simulated elemental metals as shown in Fig. 3a include the fcc metals Al, Cu,
Ni, and Au; hcp metals Co and Ti; and bcc metals W and Fe.

Figure 3: (a) Critical evaporation voltages of different metals. The
hemispherical tip size is 4 nm in radius and 12 nm in total height. (b) An Au
tip deforms during the simulation which leads to fracture. Surface atoms are
colored by the magnitude of the electric-field-induced forces.

The sequence of critical voltages from high to low is W ($\sim 295.0$ V), Ni
($\sim 241.5$ V), Co ($\sim 215.3$ V), Fe ($\sim 208.9$ V), Cu ($\sim 184.4$
V), Ti ($\sim 170.7$ V), and Al ($\sim 137.6$ V), which are obtained
respectively as the averages of their voltages within the third quarter of the
corresponding time ranges. The sequence of the tabulated ZBEF in Tsong’s paper
[45] from high to low is W (5.7 V/Å), Co (3.6 V/Å), Ni (3.5 V/Å), Fe (3.5
V/Å), Cu (3.0 V/Å), Ti (2.5 V/Å), Al (1.9 V/Å). Despite the reversed sequence
of Ni and Co, the simulation results are otherwise in qualitative agreement
with the theoretical predictions, which is a first indication that the use of
interatomic forces and field-induced forces is able to capture the basic
features of field evaporation. The simulation of Au went to failure during the
simulation as shown in Fig. 3b, due to Au’s high ductility and peculiar
deformation mechanism, which also makes it very difficult to stabilize the
pure gold specimen in experiments [48]. While it is not easy to convert the
voltage value from the simulation to a corresponding ZBEF due to the geometry
dependence of electric fields and the very different, much smaller geometry of
the virtual setup in comparison to experiments, we will show in Sec. III.4
that an approximate virtual evaporation field can be determined from the local
forces and benchmarked against tabulated values in a more direct way.

####  III.3.2 Effect of tip shape on critical voltage

As can be seen in Fig. 3a, there is a considerable over-voltage stage during
the initial shape evolution that lasts a few ns, whose value is about 20% to
25% higher than the obtained critical voltage. We find that this over-voltage
is closely related to the initial hemispherical tip shape and causes
inefficiency in the simulation. Figure 4 is an example of a 4 nm-radius
$\langle 100\rangle$ oriented W tip starting with a hemispherical shape. The
figure includes snapshots of the tip at 0.11 ns, 2.64 ns, 2.74 ns and 13.01
ns. Since the tip surface is hemispherical, there is no pronounced
concentration of electric field on the surface in the beginning, as indicated
by the homogeneous coloring of the surface atoms observed up to 2.64 ns. In
order to create higher local electric fields and trigger evaporation events,
the voltage needs to increase. In this example, due to the stable
hemispherical shape, although the voltage goes up considerably by about 34%,
only about 1250 atoms are evaporated during the first 2.53 ns. However, during
the following 0.1 ns, about 600 atoms are evaporated very quickly, which
results in a very large $\\{002\\}$ and some smaller $\\{110\\}$ facets. The
$\\{002\\}$ terrace has a large edge where the electric field can concentrate,
and the applied voltage drops promptly by about 17% within 0.1 ns approaching
the critical voltage. The tip then evolves to a semi-stable shape as shown at
13.01 ns. Since this avalanche of atoms cannot be analyzed in a meaningful
way, the initial computation wastes considerable computer time. The effect is
exacerbated as the number of atoms ripped off in the avalanche is too large
and results in a too-flat tip surface with a much smaller curvature compared
to the steady-state shape, therefore more steps are required before the final
tip shape is reached. In this example, the tip is not able to evolve to its
steady-state shape shown in Fig. 2 before the whole tip is evaporated.
Therefore, a starting hemispherical shape, while perhaps intuitive since close
to the faceted steady-state shape, is highly inefficient.

Figure 4: Evolution of a 4 nm-radius $\langle 001\rangle$ hemispherical W tip.
Surface atoms are colored by the magnitude of the electric-field-induced
forces.

In order to overcome this inefficiency, we found a cylinder tip with a flat
end cap to be a much better initial structure. We demonstrate this for a
$\langle 001\rangle$ oriented Al tip with 4 nm radius. We compare the
evolution of evaporation voltage between a cylindrical and hemispherical tip
in Fig. 5. In contrast to the voltage applied to the hemispherical tip, the
voltage applied to the cylindrical tip quickly reaches its critical voltage
within about 2 ns, thanks to the concentration of electric field on the rim of
the cylinder. Starting from the edge, the tip gradually evolves to its steady-
state shape. There is no “over-voltage” stage or avalanche-evaporation of
atoms in this process and the steady-state shape is reached within short
simulation time. Therefore we use cylinders as start shapes of virtual tips in
all other simulations. Different critical voltages obtained in the two cases
result from different tip heights when the voltages reach the respective
plateaus. In addition, since a cylinder is distinctly different from the final
steady-state shape, the consistent evolution of the tip shape to its steady
state in independent simulations further validates the present simulation
method.

Figure 5: Applied voltage vs. time in hemispherical and cylindrical Al virtual
tips. The tip structure of each case is (a) $\langle 001\rangle$ oriented
hemispherical tip and (b) $\langle 001\rangle$ oriented cylindrical tip. The
size of both tips is 4 nm in radius. Surface atoms are colored by the
magnitude of the electric-field-induced forces.

###  III.4 Critical evaporation field and average sublimation energy

We derive critical evaporation fields of different elemental metals from
simulation results to have a direct comparison to tabulated ZBEF data. To
examine alloying effects in an extreme case, simulations are conducted in an
fcc Fe18.8Ni21.9Cr18.9Co21.8Cu18.6 high entropy alloy, which is later referred
to as the FeNiCrCoCu HEA or HEA for simplicity. The tip was cut from the
periodic replica of a 32-atom special quasi-random structure (SQS) cubic cell
with 6 Fe, 6 Cr, 6 Cu, 7 Ni and 7 Co. We compare the obtained sublimation
energies to the tabulated data, as well as the sublimation energies and
critical voltages for the different elements in their elemental form and in
the HEA, where a “mixed” critical voltage of the HEA is found.

Although the local electric field is site-dependent for surface atoms, in
order to compare to the tabulated ZBEF in Tsong’s paper [45], the critical
evaporation field is derived as an averaged value over all evaporation events.
The local electric field of an evaporated atom is chosen as its critical
evaporation field at the time when its interatomic force goes to zero. As
shown in Fig. 6a, typical tabulated ZBEF values [45] are about 3.5 times
higher than the critical electric field obtained in our simulations. Since the
aberration is systematic, it demonstrates that while the reduced scales of the
MD setup cause some scaling in potential and fields, the process of
evaporation is described with high fidelity with linear dependence between
simulated and calculated values.

A similar analysis has also been done for zero-field sublimation energies of
surface atoms of the FeNiCrCoCu HEA cylinder tip. Simulations were performed
for the HEA tip of 5 nm in radius as well as for elemental metals of each
alloy component using the same interatomic potential. Sublimation energies are
computed for $\langle 001\rangle$ surface atoms. In Fig. 6b, circles are
results obtained for elemental metals of HEA components, while diamonds are
obtained in the HEA alloy. Symbols of the same element are in the same color.
Sublimation energies of elemental metals are very close to each other in the
case of Co, Fe, and Ni, whereas Cu bonds are much weaker than the others.
Compared to the homogeneous bonding environment in the elemental metals, there
is a distribution of sublimation energies for each of the components in the
HEA as shown in Fig. 6b, which indicates a mixing effect caused by different
bonding contributions from different neighbor atoms. Average sublimation
energies of each component in the HEA have a similar sequence compared to the
elemental metals, however, the discrepancy of bond strength among components
is mitigated. The standard deviation of sublimation energies of the elemental
metals is 0.50 eV, while that of the HEA atoms is only 0.28 eV. The trend of
homogenization of bonds is also reflected in the critical voltage. As shown in
Tab. 1, the critical voltage of the HEA equals the average of the critical
voltages of elemental metals of the components following the rule of mixture,
which is an important result that we will revisit in Sec. IV in the case of
Cu-Ni alloy, where we discuss the assessment of the reliability of interatomic
potentials.

Although a homogenization effect of bonds is observed in the HEA, the copper
atoms still have the weakest bonds among the alloy components, and their
average sublimation energy has the largest deviation $-8.3$% from the average
of all atoms in the sample $-4.95$ eV. In contrast, the deviations of Ni, Co,
Fe and Cr are 4.8%, 2.6%, 0.59% and 0.42% respectively. The much weaker
bonding of copper atoms makes them the easiest ones to be evaporated in the
HEA, and their preferential evaporation causes a rough tip surface with
dimples and evaporation of clusters (separation at the “weakest link”, (see
the inset in Fig. 6b). With our carefully controlled voltage-driven
evaporation, we never observe cluster evaporation in elemental metals. This
finding explains and is in agreement with experimental observations on a
similar, compositionally homogeneous equiatomic FeNiCrCoMn HEA, where in
voltage mode an unusually high fraction of 30%-40% multiple vs. single hits is
observed [49], which is more than an order of magnitude larger than what is
typically observed for elemental metal tips [50].

As a final note, we do not observe any recognizable zone lines or poles in the
desorption map of the HEA. This agrees with the general expectation, since for
example in highly alloyed steels, zone lines and poles start to wash out and
disappear with increasing solute concentration [51, 52, 53]. However, at this
point, we cannot make any comparison to experiments so far, as no desorption
map of HEA has been reported in the literature.

Figure 6: (a) Average critical evaporation fields of different elemental
metals from simulations of virtual tips with 4 nm in radius vs. tabulated ZBEF
values from Tsong [45]. The error bars represent the standard deviation of the
observed evaporation fields. (b) Averaged sublimation energies of Co, Ni, Fe,
Cu, and Cr metals in elemental form (circles) and as components of an HEA
(diamonds) vs. tabulated sublimation energies from Tsong [45]. Symbols in the
same color represent the same species. The error bars represent the standard
deviation of sublimation energies of components in the HEA. The virtual tips
were 5 nm in radius and all simulations used the same interatomic potential.
The inset is the tip surface during evaporation. An evaporating cluster is
indicated by a red circle. Table 1:  Comparison of elemental vs. alloy
critical voltages of components of the Fe18.8Ni21.9Cr18.9Co21.8Cu18.6 high
entropy alloy.

Component |  Ni |  Cr |  Co |  Fe |  Cu
---|---|---|---|---|---
Concentration (%) |  21.9 |  18.9 |  21.8 |  18.8 |  18.6
Critical voltages of elemental metals (V) |  342.1 |  317.0 |  311.4 |  284.0 |  242.7

Weighted average of critical voltages (V) |  301.2 |  Critical voltage of HEA (V) |  298.4
---|---|---|---

###  III.5 Field desorption pattern

By considering the important role of interatomic interactions during field
evaporation, we are able to reproduce a key feature of field desorption
patterns observed in experiments – enhanced zone lines. As a demonstration, we
make a comparison between Ni desorption patterns obtained by experiment,
classical finite element APT simulation, and the present TAPSim-MD approach
respectively.

In experiments, zone lines with high intensity of hit events can be observed
for many metals, such as for the case of a $\langle 111\rangle$ oriented Ni
field desorption pattern shown in Fig. 7a. This enhanced-zone-line feature of
field desorption patterns was first reported by Waugh et al. about 45 years
ago in many materials [54], and a “roll-up” model was then proposed to explain
this feature [55]. However, this conjecture has not been validated yet, and no
simulation has been previously able to reproduce this feature. Field
desorption patterns simulated by a traditional electrostatic approach [56, 24,
25, 27] always display broad and depleted stripes around zone lines. Without
considering the local bonding environment of the chemical species, the
patterns only depend on the crystal structures of the materials. Figure 7b is
the desorption pattern for the $\langle$111$\rangle$ oriented face centered
cubic structure, which is found by the classical TAPSim method to look
identical for all fcc-materials [25]. In contrast, results obtained here by
the TAPSim-MD approach find for the first time zone-line enhancement in
agreement with the experiments. Figure 7c is the simulated field desorption
pattern of an $\langle$111$\rangle$ oriented Ni tip by TAPSim-MD. The virtual
tip is 8 nm in radius. The desorption pattern was obtained by evaporating
about 50 layers of atoms. Regions with a warmer color have a higher intensity
of atom hit events. The field desorption pattern has a three-fold symmetry,
and the intensity of atom hit events is higher around the three
$\langle$110$\rangle$ zone lines, which agrees with the experimental result in
Fig. 7a. Our simulations demonstrate clearly that the electric field is
perpendicular to the envelope of the surface. On the other hand, interatomic
forces are anisotropic, and their overall directions are governed by the
crystallographic structure. Because of that, the electrostatic desorption
force is not always antiparallel to the adhesive interatomic force at the
moment of desorption, making the atoms around certain zone lines leave the
surface with a non-zero lateral velocity, which accumulates them along certain
zone lines. The full quantitative description of the exact mechanism of the
enhanced zone lines exceeds the scope of the current publication and will be
described in [57].

Figure 7: Field desorption patterns obtained for an $\langle$111$\rangle$
oriented Ni tip in (a) experiment, the figure is reprinted from [54]; (b)
TAPSim simulation for an fcc $\langle$111$\rangle$ oriented structure (tip
radius = 25 nm), the figure is reprinted from [58]; (c) proposed TAPSim-MD
simulation (tip radius = 8 nm).

##  IV Influence of interatomic potentials on simulations

As mentioned before, the selection of the interatomic potential is crucial to
the MD simulation. To demonstrate its impact, we compare the results from two
common EAM potentials [59, 60] for the Cu-Ni binary alloy system and assess
the physical sensibility of their results.

The virtual tips for these simulations are cylinders with 4 nm in radius and
12 nm in height. Nine virtual tips were generated with randomly arranged Ni
impurities with concentrations of 0%, 5%, 15%, 30%, 50%, 70%, 85%, 95%, 100%.
Evaporation from these tips was studied with the interatomic potential by
Foiles et al. [59] in comparison to the potential by Onat et al. [60]. The
important difference between the potentials is that Foiles uses two elemental
EAM potentials where the cross term $\phi_{AB}$ of interactions between two
species is defined as the geometric average of the elemental pair potentials,
whereas Onat fits a dedicated cross term for the alloy potential.

The critical voltages in the simulations with the potentials were determined
in the same way described before. For an isomorphous system such as Cu-Ni,
properties follow in general the mixing rule as we have just shown for the
case of an HEA in Sec. III.4, Fig. 6b. As shown in Fig. 8, the Foiles
potential produces critical voltages that have a pronounced quadratic
dependence on the composition, while the linear relationship that we expect
according to our HEA simulations is reproduced by the Onat potential.
Considering the nature of the cross potential as a geometric average in the
potential from Foiles, the quadratic relation should simply originate from the
number of bonds in a binary solid solution, which quadratically depends on the
composition [61]. However, the potential from Onat shows a linear dependence
consistent with general expectations, and our results for the miscible
elements in the FeNiCrCoCu high entropy alloy also display a “mixed” critical
voltage. Compared to the dedicated cross term fitted in Onat’s potential, the
definition of the cross term as the geometric average in Foiles’s is too
simplified to accurately describe interactions between Cu and Ni atoms, which
is the reason why later alloy approaches in EAM do not use it anymore [62, 63,
64, 65]. Since the ZBEF of a species is directly connected to its sublimation
energy which in turn is closely related to the environment as was already
recognized early in e.g. the “Müller formula” [19], a correct description of
interatomic interactions is significantly important to obtaining reliable
sublimation energies.

Figure 8: Critical voltages of Cu1-xNix alloys obtained by using interatomic
potentials from [59] and [60]. Error bars represent the mean squared errors of
the critical voltages.

##  V Computational performance

With regard to the computational performance, the combination of a traditional
finite element solver with LAMMPS results in a computationally demanding
simulation, whose performance also limits the size of the virtual sample. So
far, the largest tip used in the simulation is 15 nm in radius with 2,005,877
atoms. The computation speed of the MD part largely depends on the type of the
chosen interatomic potentials. We use EAM potentials in most of the cases to
balance the efficiency and accuracy of the simulation. The computational
effort increases significantly with the number of atoms in the virtual tip
using the same interatomic potential on the Owens HPC system at the Ohio
Supercomputer Center [66], making the use of highly parallel supercomputers
indispensable. To talk about scaling more specifically, since the number of
atoms in the top surface of the tip scales with the square of the tip radius,
if we want to evaporate the same number of atom layers for different tip
sizes, the simulation time required for the tip with a larger radius will be
even longer. As for our examples, in the case of $\langle 110\rangle$ oriented
elemental tungsten with the EAM potential from Marinica et al. [42], it takes
approximately 19.97 CPU hours to evaporate a layer of atoms (604 atoms) in Tip
A (28,204 atoms, 4 nm in radius), 331.94 CPU hours to evaporate a layer of
atoms (1,110 atoms) in Tip B (103,363 atoms, 5 nm in radius), and 3767.18 CPU
hours to evaporate a layer of atoms (4,432 atoms) in Tip C (484,213 atoms, 10
nm in radius), suggesting order-$N$ scaling of CPU time with the number of
atoms in a layer, or alternatively, $N^{2}$ scaling with the tip radius. On
our system, the dependence of CPU hours on evaporated number of atoms can be
described with $R^{2}=0.998$ by $t_{\rm CPU}[{\rm CPU-h}]=N-673$.

##  VI Conclusion

We have introduced TAPSim-MD as a fully physical simulation tool for field
desorption. By integrating LAMMPS with the field and surface-charge solver
TAPSim, the new TAPSim-MD approach is able to evaporate atoms in a virtual tip
by performing classical MD simulations driven by both electric-field-induced
forces and interatomic forces. This approach was validated from several
aspects. We have first shown that the electric field distribution obtained by
solving Poisson’s equation compares favorably with results from electron
holography measurements [46]. Second, we demonstrate that the material-
dependent critical voltages and resulting steady-state shapes of virtual tips,
as well as critical evaporation fields found in elemental metals and
sublimation energies obtained in a high entropy alloy are in line with
previous predictions [45] and general expectations. Finally, field desorption
patterns displaying enhanced zone lines are successfully reproduced for the
first time, in contrast to all previous forward-simulation work [8]. By
eliminating the ad-hoc assumptions adopted by previous approaches concerning
evaporation fields and related criteria, our model is able to capture key
characteristics of field evaporation in a fully physics-based “ab-initio” way,
provided a sensible interatomic potential is chosen as illustrated for Cu-Ni
alloys. Therefore, sufficient benchmarking and validation needs to be done to
balance prediction fidelity and computational efficiency.

The development and validation of TAPSim-MD described in this paper as well as
its general nature which is not restricted to a specific material system such
as metals enables applications to in principle any arbitrary material as long
as a reasonable potential can be identified, as well as examining questions
that could not be answered before. Obvious examples include a better
understanding of the influence of inhomogeneous bonding on the evaporation
sequence and related artifacts as well as field desorption maps. In addition,
the full knowledge of the initial atom positions and the resulting detector
hit maps enables uncertainty quantification to benchmark classical
reconstruction methods.

## Acknowledgements

This work is sponsored by AFOSR (PM Dr. Ali Sayir) under Award No.
FA9550-14-1-0249 and FA9550-19-1-0378. Computing resources are supplied and
maintained by the Ohio Supercomputer Center [66] under Grant No. PAA0010.
Results of the simulations were visualized with the programs OVITO [67] and
ParaView [68]. We thank Dr. David Seidman for helpful discussions concerning
the experimental situation of field evaporation of gold.

## References

  * Miller and Forbes [2014] M. K. Miller and R. G. Forbes, Field evaporation and related topics, in  _Atom-Probe Tomography: The Local Electrode Atom Probe_, edited by M. K. Miller and R. G. Forbes (Springer US, Boston, MA, 2014) pp. 111–187.
  * Arslan  _et al._ [2008] I. Arslan, E. A. Marquis, M. Homer, M. A. Hekmaty, and N. C. Bartelt, Towards better 3-d reconstructions by combining electron tomography and atom-probe tomography, Ultramicroscopy 108, 1579 (2008).
  * Marquis  _et al._ [2011] E. A. Marquis, B. P. Geiser, T. J. Prosa, and D. J. Larson, Evolution of tip shape during field evaporation of complex multilayer structures, Journal of Microscopy 241, 225 (2011).
  * Vurpillot  _et al._ [2013] F. Vurpillot, B. Gault, B. P. Geiser, and D. J. Larson, Reconstructing atom probe data: A review, Ultramicroscopy 132, 19 (2013).
  * Scully  _et al._ [2020] J. R. Scully, S. B. Inman, A. Y. Gerard, C. D. Taylor, W. Windl, D. K. Schreiber, P. Lu, J. E. Saal, and G. S. Frankel, Controlling the corrosion resistance of multi-principal element alloys, Scripta Materialia 188, 96 (2020).
  * Frankel  _et al._ [2021] G. S. Frankel, J. D. Vienna, J. Lian, X. Guo, S. Gin, S. H. Kim, J. Du, J. V. Ryan, J. Wang, W. Windl, C. D. Taylor, and J. R. Scully, Recent advances in corrosion science applicable to disposal of high-level nuclear waste, Chemical Reviews 121, 12327 (2021).
  * Novaes  _et al._ [2021] L. F. T. Novaes, J. Liu, Y. Shen, L. Lu, J. M. Meinhardt, and S. Lin, Electrocatalysis as an enabling technology for organic synthesis, Chem. Soc. Rev. 50, 7941 (2021).
  * Vurpillot and Oberdorfer [2015] F. Vurpillot and C. Oberdorfer, Modeling atom probe tomography: A review, Ultramicroscopy 159, 202 (2015).
  * Moore [1962] A. Moore, The structure of atomically smooth spherical surfaces, Journal of Physics and Chemistry of Solids 23, 907 (1962).
  * Suvorov and Razinkova [1977] A. Suvorov and T. Razinkova, Computer simulation of tungsten field evaporation by field-ion microscope, Surface Science 65, 354 (1977).
  * Haley  _et al._ [2013] D. Haley, M. P. Moody, and G. D. Smith, Level set methods for modelling field evaporation in atom probe, Microscopy and Microanalysis 19, 1709–1717 (2013).
  * Xu  _et al._ [2015] Z. Xu, D. Li, W. Xu, A. Devaraj, R. Colby, S. Thevuthasan, B. P. Geiser, and D. J. Larson, Simulation of heterogeneous atom probe tip shapes evolution during field evaporation using a level set method and different evaporation models, Computer Physics Communications 189, 106 (2015).
  * Fletcher  _et al._ [2019] C. Fletcher, M. P. Moody, and D. Haley, Fast modelling of field evaporation in atom probe tomography using level set methods, Journal of Physics D: Applied Physics 52, 435305 (2019).
  * Fletcher  _et al._ [2020] C. Fletcher, M. P. Moody, and D. Haley, Towards model-driven reconstruction in atom probe tomography, Journal of Physics D: Applied Physics 53, 475303 (2020).
  * Sánchez  _et al._ [2004] C. G. Sánchez, A. Y. lozovoi, and A. alavi, Field-evaporation from first-principles, Molecular Physics 102, 1045 (2004).
  * Ono  _et al._ [2005] T. Ono, T. Sasaki, J. Otsuka, and K. Hirose, First-principles study on field evaporation of surface atoms from w(011) and mo(011) surfaces, Surface Science 577, 42 (2005).
  * Peralta  _et al._ [2013] J. Peralta, S. R. Broderick, and K. Rajan, Mapping energetics of atom probe evaporation events through first principles calculations, Ultramicroscopy 132, 143 (2013).
  * Karahka and Kreuzer [2013] M. Karahka and H. Kreuzer, Field evaporation of oxides: A theoretical study, Ultramicroscopy 132, 54 (2013), iFES 2012.
  * Yao  _et al._ [2015] L. Yao, T. Withrow, O. D. Restrepo, W. Windl, and E. A. Marquis, Effects of the local structure dependence of evaporation fields on field evaporation behavior, Applied Physics Letters 107, 241602 (2015).
  * Cui and Ringer [2018] X.-Y. Cui and S. P. Ringer, On the nexus between atom probe microscopy and density functional theory simulations, Materials Characterization 146, 347 (2018).
  * Ashton  _et al._ [2020] M. Ashton, A. Mishra, J. Neugebauer, and C. Freysoldt, Ab initio description of bond breaking in large electric fields, Physical Review Letters 124, 176801 (2020).
  * Ohnuma [2021] T. Ohnuma, First-Principles Calculation of the Evaporation Field and Roll-up Effect of M (M = Fe, Cu, Si, and Mn) on the Fe (001) and Fe Step Structure, Microscopy and Microanalysis , 1 (2021), publisher: Cambridge University Press.
  * Vurpillot  _et al._ [1999a] F. Vurpillot, A. Bostel, A. Menand, and D. Blavette, Trajectories of field emitted ions in 3d atom-probe, The European Physical Journal - Applied Physics 6, 217 (1999a).
  * Geiser  _et al._ [2009] B. P. Geiser, D. J. Larson, S. Gerstl, D. Reinhard, T. F. Kelly, T. J. Prosa, and D. Olson, A system for simulation of tip evolution under field evaporation, Microscopy and Microanalysis 15, 302 (2009).
  * Oberdorfer  _et al._ [2013] C. Oberdorfer, S. M. Eich, and G. Schmitz, A full-scale simulation approach for atom probe tomography, Ultramicroscopy 128, 55 (2013).
  * Oberdorfer and Schmitz [2011] C. Oberdorfer and G. Schmitz, On the field evaporation behavior of dielectric materials in three-dimensional atom probe: A numeric simulation, Microscopy and Microanalysis 17, 15 (2011).
  * Rolland  _et al._ [2015] N. Rolland, F. Vurpillot, S. Duguay, and D. Blavette, A meshless algorithm to model field evaporation in atom probe tomography, Microscopy and Microanalysis 21, 1649 (2015).
  * Gruber  _et al._ [2011] M. Gruber, F. Vurpillot, A. Bostel, and B. Deconihout, Field evaporation: A kinetic monte carlo approach on the influence of temperature, Surface Science 605, 2025 (2011).
  * Djurabekova  _et al._ [2011] F. Djurabekova, S. Parviainen, A. Pohjonen, and K. Nordlund, Atomistic modeling of metal surfaces under electric fields: Direct coupling of electric fields to a molecular dynamics algorithm, Physical Review E 83, 026704 (2011).
  * Parviainen  _et al._ [2015] S. Parviainen, F. Djurabekova, S. P. Fitzgerald, A. Ruzibaev, and K. Nordlund, Atomistic simulations of field assisted evaporation in atom probe tomography, Journal of Physics D: Applied Physics 49, 045302 (2015).
  * Oberdorfer  _et al._ [2018] C. Oberdorfer, T. Withrow, L. J. Yu, K. Fisher, E. A. Marquis, and W. Windl, Influence of surface relaxation on solute atoms positioning within atom probe tomography reconstructions, Materials Characterization 146, 324 (2018).
  * Veske  _et al._ [2018] M. Veske, A. Kyritsakis, K. Eimre, V. Zadin, A. Aabloo, and F. Djurabekova, Dynamic coupling of a finite element solver to large-scale atomistic simulations, Journal of Computational Physics 367, 279 (2018).
  * Katnagallu  _et al._ [2018] S. Katnagallu, M. Dagan, S. Parviainen, A. Nematollahi, B. Grabowski, P. A. J. Bagot, N. Rolland, J. Neugebauer, D. Raabe, F. Vurpillot, M. P. Moody, and B. Gault, Impact of local electrostatic field rearrangement on field ionization, Journal of Physics D: Applied Physics 51, 105601 (2018).
  * Vurpillot  _et al._ [2018] F. Vurpillot, S. Parviainen, F. Djurabekova, D. Zanuttini, and B. Gervais, Simulation tools for atom probe tomography: A path for diagnosis and treatment of image degradation, Materials Characterization 146, 336 (2018).
  * Plimpton [1995] S. Plimpton, Fast parallel algorithms for short-range molecular dynamics, Journal of Computational Physics 117, 1 (1995).
  * van der Berg  _et al._ [1995] R. A. van der Berg, P. W. M. Blom, J. F. M. Cillessen, and R. M. Wolf, Field dependent permittivity in metal‐semiconducting srtio3 schottky diodes, Applied Physics Letters 66, 697 (1995).
  * Mishin  _et al._ [1999] Y. Mishin, D. Farkas, M. J. Mehl, and D. A. Papaconstantopoulos, Interatomic potentials for monoatomic metals from experimental data and ab initio calculations, Phys. Rev. B 59, 3393 (1999).
  * Mishin  _et al._ [2001] Y. Mishin, M. J. Mehl, D. A. Papaconstantopoulos, A. F. Voter, and J. D. Kress, Structural stability and lattice defects in copper: Ab initio, tight-binding, and embedded-atom calculations, Phys. Rev. B 63, 224106 (2001).
  * Pun and Mishin [2012] G. P. P. Pun and Y. Mishin, Embedded-atom potential for hcp and fcc cobalt, Phys. Rev. B 86, 134116 (2012).
  * Foiles  _et al._ [1986] S. M. Foiles, M. I. Baskes, and M. S. Daw, Embedded-atom-method functions for the fcc metals cu, ag, au, ni, pd, pt, and their alloys, Phys. Rev. B 33, 7983 (1986).
  * Mendelev  _et al._ [2016] M. I. Mendelev, T. L. Underwood, and G. J. Ackland, Development of an interatomic potential for the simulation of defects, plasticity, and phase transformations in titanium, The Journal of Chemical Physics 145, 154102 (2016).
  * Marinica  _et al._ [2013] M.-C. Marinica, L. Ventelon, M. R. Gilbert, L. Proville, S. L. Dudarev, J. Marian, G. Bencteux, and F. Willaime, Interatomic potentials for modelling radiation defects and dislocations in tungsten, Journal of Physics: Condensed Matter 25, 395502 (2013).
  * Baskes [1992] M. I. Baskes, Modified embedded-atom potentials for cubic materials and impurities, Phys. Rev. B 46, 2727 (1992).
  * Farkas and Caro [2018] D. Farkas and A. Caro, Model interatomic potentials and lattice strain in a high-entropy alloy, Journal of Materials Research 33, 3218 (2018).
  * Tsong [1978] T. T. Tsong, Field ion image formation, Surface Science 70, 211 (1978).
  * Migunov  _et al._ [2015] V. Migunov, A. London, M. Farle, and R. E. Dunin-Borkowski, Model-independent measurement of the charge density distribution along an fe atom probe needle using off-axis electron holography without mean inner potential effects, Journal of Applied Physics 117, 134301 (2015).
  * Miller [2000a] M. K. Miller, Field ion microscopy, in  _Atom Probe Tomography: Analysis at the Atomic Level_ (Springer US, Boston, MA, 2000) Chap. 3, pp. 45–83.
  * Ast and Seidman [1968] D. G. Ast and D. N. Seidman, The field ion microscopy of gold, Applied Physics Letters 13, 348 (1968).
  * Muniandy  _et al._ [2021] Y. Muniandy, M. He, M. Eizadjou, E. P. George, J. J. Kruzic, S. P. Ringer, and B. Gludovatz, Compositional variations in equiatomic crmnfeconi high-entropy alloys, Materials Characterization 180, 111437 (2021).
  * Gault  _et al._ [2012] B. Gault, M. P. Moody, J. M. Cairney, and S. P. Ringer, Experimental protocols in atom probe tomography, in  _Atom Probe Microscopy_, edited by B. Gault, M. P. Moody, J. M. Cairney, and S. P. Ringer (Springer New York, New York, NY, 2012) p. 150.
  * Miller and Hetherington [1991] M. Miller and M. Hetherington, Local magnification effects in the atom probe, Surface Science 246, 442 (1991).
  * Miller [2000b] M. Miller, The development of atom probe field-ion microscopy, Materials Characterization 44, 11 (2000b).
  * Yao [2016] L. Yao, A filtering method to reveal crystalline patterns from atom probe microscopy desorption maps, MethodsX 3, 268 (2016).
  * Waugh  _et al._ [1975] A. R. Waugh, E. D. Boyes, and M. J. Southon, Field-desorption microscopy and the atom probe, Nature 253, 342 (1975).
  * Waugh  _et al._ [1976] A. R. Waugh, E. D. Boyes, and M. J. Southon, Investigations of field evaporation with a field-desorption microscope, Surface Science 61, 109 (1976).
  * Vurpillot  _et al._ [1999b] Vurpillot, Bostel, and Blavette, The shape of field emitters and the ion trajectories in three‐dimensional atom probes, Journal of Microscopy 196 (1999b).
  * Qi  _et al._ [2022] J. Qi, C. Oberdorfer, E. A. Marquis, and W. Windl, Origin of enhanced zone lines in field desorption maps (2022).
  * Oberdorfer  _et al._ [2015] C. Oberdorfer, S. M. Eich, M. Lütkemeyer, and G. Schmitz, Applications of a versatile modelling approach to 3d atom probe simulations, Ultramicroscopy 159, 184 (2015).
  * Foiles [1985] S. M. Foiles, Calculation of the surface segregation of ni-cu alloys with the use of the embedded-atom method, Phys. Rev. B 32, 7685 (1985).
  * Onat and Durukanoğlu [2013] B. Onat and S. Durukanoğlu, An optimized interatomic potential for cu–ni alloys with the embedded-atom method, Journal of Physics: Condensed Matter 26, 035404 (2013).
  * Oberdorfer and Windl [2019] C. Oberdorfer and W. Windl, Bond-order bond energy model for alloys, Acta Materialia 179, 406 (2019).
  * Ward  _et al._ [2012] L. Ward, A. Agrawal, K. M. Flores, and W. Windl, Rapid production of accurate embedded-atom method potentials for metal alloys (2012).
  * Roth  _et al._ [2017] P. C. Roth, H. Shan, D. Riegner, N. Antolin, S. Sreepathi, L. Oliker, S. Williams, S. Moore, and W. Windl, Performance analysis and optimization of the rampage metal alloy potential generation software, in  _Proceedings of the 4th ACM SIGPLAN International Workshop on Software Engineering for Parallel Systems_, SEPS 2017 (Association for Computing Machinery, New York, NY, USA, 2017) p. 11–20.
  * Hegde  _et al._ [2022] A. Hegde, E. Weiss, W. Windl, H. Najm, and C. Safta, Bayesian calibration of interatomic potentials for binary alloys, Computational Materials Science 214, 111660 (2022).
  * Johnson [1989] R. Johnson, Alloy models with the embedded-atom method, Physical Review B 39, 12554 (1989).
  * Center [1987] O. S. Center, Ohio supercomputer center (1987).
  * Stukowski [2009] A. Stukowski, Visualization and analysis of atomistic simulation data with OVITO–the open visualization tool, Modelling and Simulation in Materials Science and Engineering 18, 015012 (2009).
  * Ayachit [2015] U. Ayachit,  _The ParaView Guide: A Parallel Visualization Application_ (Kitware, Inc., Clifton Park, NY, USA, 2015).