# Mesoscopic modeling of field evaporation on atom probe tomography

Constantinos Hatzoglou, Benjamin Klaes, Fabien Delaroche, Gérald da Costa, Brian Geiser, Markus Kühbach, Peter Wells, François Vurpillot

► To cite this version:

Constantinos Hatzoglou, Benjamin Klaes, Fabien Delaroche, Gérald da Costa, Brian Geiser, et al. Mesoscopic modeling of field evaporation on atom probe tomography. Journal of Physics D: Applied Physics, 2023, 56 (37), pp.375301. 10.1088/1361-6463/acd649. hal-04344607

HAL Id: hal-04344607

https://normandie-univ.hal.science/hal-04344607

Submitted on 29 May 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Mesoscopic modeling of field evaporation on atom probe tomography

Constantinos Hatzoglou, Benjamin Klaes, Fabien Delaroche, Gérald da Costa, Brian Geiser, Markus Kühbach, Peter Wells, François Vurpillot

► To cite this version:

Constantinos Hatzoglou, Benjamin Klaes, Fabien Delaroche, Gérald da Costa, Brian Geiser, et al. Mesoscopic modeling of field evaporation on atom probe tomography. Journal of Physics D: Applied Physics, 2023, 56 (37), pp.375301. 10.1088/1361-6463/acd649. hal-04344607

HAL Id: hal-04344607

https://normandie-univ.hal.science/hal-04344607v1

Submitted on 29 May 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
PAPER • OPEN ACCESS

# Mesoscopic modeling of field evaporation on atom probe tomography

To cite this article: Constantinos Hatzoglou et al 2023 *J. Phys. D: Appl. Phys.* **56** 375301

View the [article online](#) for updates and enhancements.

You may also like

- Understanding atom probe's analytical performance for iron oxides using correlation histograms and ab initio calculations
Se-Ho Kim, Shalini Bhatt, Daniel K Schreiber et al.

- Atomistic simulations of field assisted evaporation in atom probe tomography
S Parviainen, F Djurabekova, S P Fitzgerald et al.

- Nanoscale photoconductive switching effect applied to atom probe tomography
L. Zhao, A. Normand, J. Houard et al.

<figure><img src="image_2.png" /></figure>

The
Electrochemical
Society

Advancing solid state &
electrochemical science & technology

<figure><img src="image_3.png" /></figure>

DISCOVER
how sustainability
intersects with
electrochemistry & solid
state science research

<figure><img src="image_4.png" /></figure>
PAPER • OPEN ACCESS

# Mesoscopic modeling of field evaporation on atom probe tomography

To cite this article: Constantinos Hatzoglou et al 2023 *J. Phys. D: Appl. Phys.* **56** 375301

View the [article online](#) for updates and enhancements.

You may also like

- Understanding atom probe's analytical performance for iron oxides using correlation histograms and ab initio calculations
Se-Ho Kim, Shalini Bhatt, Daniel K Schreiber et al.

- Atomistic simulations of field assisted evaporation in atom probe tomography
S Parviainen, F Djurabekova, S P Fitzgerald et al.

- Nanoscale photoconductive switching effect applied to atom probe tomography
L. Zhao, A. Normand, J. Houard et al.

<figure><img src="image_2.png" /></figure>

The
Electrochemical
Society

Advancing solid state &
electrochemical science & technology

<figure><img src="image_3.png" /></figure>

DISCOVER
how sustainability
intersects with
electrochemistry & solid
state science research

<figure><img src="image_4.png" /></figure>
# Mesoscopic modeling of field evaporation on atom probe tomography

Constantinos Hatzoglou<sup>1,2,*</sup> ◐, Benjamin Klaes<sup>1</sup>, Fabien Delaroche<sup>1</sup>, Gérald Da Costa<sup>1</sup>, Brian Geiser<sup>3</sup>, Markus Kühbach<sup>4,5</sup>, Peter B Wells<sup>6</sup> and François Vurpillot<sup>1</sup>

E-mail: constantinos.hatzoglou@ntnu.no

Received 14 December 2022, revised 30 March 2023

Accepted for publication 17 May 2023

Published 13 June 2023

<figure><img src="image_5.png" /></figure>

CrossMark

## Abstract

Reconstructions in atom probe tomography (APT) are biased by image distortions arising from dynamic changes of the specimen geometry that controls image projection. Despite the strong efforts to build realistic models for understanding and reproducing image artifacts, the current models are too slow or not adapted to be routinely used in image correction approaches. To understand the APT imaging process for real size samples submitted to realistic experimental conditions of electric field and temperature, we propose an alternative simulation tool based on a coarse-grained model of the sample surface. The surface electric field on a meshed surface is calculated by using continuous models describing field evaporation. The dynamic evolution of the sample surface and the image projection are predicted using materials properties. We show that the interplay between temperature and electric field is an important ingredient in predicting the ion projection, in pure metals and in more complex materials. This fast approach accurately reproduces the well-known local magnification and trajectory overlaps effects in the evaporation of small particles. By combining prior knowledge about the sample structure and properties, the model could be used to improve the reconstruction approaches for complex sample geometries.

Keywords: field evaporation, simulation, atom probe tomography

(Some figures may appear in colour only in the online journal)

## 1. Introduction

Atom probe tomography (APT) is a powerful technique to analyze materials in three dimensions (3D) at the atomic scale

[1–3]. The main physical principle is the field evaporation (FE) of atoms which are located at the material surface. FE is the transition, in the presence of a sufficiently strong electric field, of a material directly from the solid phase to the plasma phase without passing through sublimation. For most materials, the critical field giving rise to FE without any thermal activation is in the range 10–60 V nm⁻<fnref n="1" />. This intense electric field is easily achieved at the surface of a material using the classical electrostatic point effect. The material is prepared as a sharp needle with an end radius R in the range of 10–100 nm and is submitted to a high voltage V (a few kV) to

<figure><img src="image_6.png" /></figure>

<footnote n="1"><sup>1</sup> Normandie Université, UNIROUEN, INSA Rouen, CNRS, Groupe de Physique des Matériaux, 76000 Rouen, France</footnote>

<footnote n="2"><sup>2</sup> Department of Materials Science and Engineering, NTNU, Norwegian University of Science and Technology, Trondheim 7491, Norway</footnote>

<footnote n="3"><sup>3</sup> CAMECA Instruments Inc., 5470 Nobel Dr, Madison, WI 53711, United States of America</footnote>

<footnote n="4"><sup>4</sup> Max-Planck-Institut für Eisenforschung GmbH, Max-Planck-Straße 1, 40237 Düsseldorf, Germany</footnote>

<footnote n="5"><sup>5</sup> Consortium FAIRmat c/o Humboldt-Universität zu Berlin, Zum Großen Windkanal 2, 12489 Berlin, Germany</footnote>

<footnote n="6"><sup>6</sup> Intel Corporation, Hillsboro, OR 97124, United States of America</footnote>

<footnote>* Author to whom any correspondence should be addressed.</footnote>

<footnote>Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.</footnote>
# Mesoscopic modeling of field evaporation on atom probe tomography

Constantinos Hatzoglou<sup>1,2,*</sup> ◐, Benjamin Klaes<sup>1</sup>, Fabien Delaroche<sup>1</sup>, Gérald Da Costa<sup>1</sup>, Brian Geiser<sup>3</sup>, Markus Kühbach<sup>4,5</sup>, Peter B Wells<sup>6</sup> and François Vurpillot<sup>1</sup>

E-mail: constantinos.hatzoglou@ntnu.no

Received 14 December 2022, revised 30 March 2023

Accepted for publication 17 May 2023

Published 13 June 2023

<figure><img src="image_5.png" /></figure>

CrossMark

## Abstract

Reconstructions in atom probe tomography (APT) are biased by image distortions arising from dynamic changes of the specimen geometry that controls image projection. Despite the strong efforts to build realistic models for understanding and reproducing image artifacts, the current models are too slow or not adapted to be routinely used in image correction approaches. To understand the APT imaging process for real size samples submitted to realistic experimental conditions of electric field and temperature, we propose an alternative simulation tool based on a coarse-grained model of the sample surface. The surface electric field on a meshed surface is calculated by using continuous models describing field evaporation. The dynamic evolution of the sample surface and the image projection are predicted using materials properties. We show that the interplay between temperature and electric field is an important ingredient in predicting the ion projection, in pure metals and in more complex materials. This fast approach accurately reproduces the well-known local magnification and trajectory overlaps effects in the evaporation of small particles. By combining prior knowledge about the sample structure and properties, the model could be used to improve the reconstruction approaches for complex sample geometries.

Keywords: field evaporation, simulation, atom probe tomography

(Some figures may appear in colour only in the online journal)

## 1. Introduction

Atom probe tomography (APT) is a powerful technique to analyze materials in three dimensions (3D) at the atomic scale

[1–3]. The main physical principle is the field evaporation (FE) of atoms which are located at the material surface. FE is the transition, in the presence of a sufficiently strong electric field, of a material directly from the solid phase to the plasma phase without passing through sublimation. For most materials, the critical field giving rise to FE without any thermal activation is in the range 10–60 V nm⁻<fnref n="1" />. This intense electric field is easily achieved at the surface of a material using the classical electrostatic point effect. The material is prepared as a sharp needle with an end radius R in the range of 10–100 nm and is submitted to a high voltage V (a few kV) to

<figure><img src="image_6.png" /></figure>

<footnote n="1"><sup>1</sup> Normandie Université, UNIROUEN, INSA Rouen, CNRS, Groupe de Physique des Matériaux, 76000 Rouen, France</footnote>

<footnote n="2"><sup>2</sup> Department of Materials Science and Engineering, NTNU, Norwegian University of Science and Technology, Trondheim 7491, Norway</footnote>

<footnote n="3"><sup>3</sup> CAMECA Instruments Inc., 5470 Nobel Dr, Madison, WI 53711, United States of America</footnote>

<footnote n="4"><sup>4</sup> Max-Planck-Institut für Eisenforschung GmbH, Max-Planck-Straße 1, 40237 Düsseldorf, Germany</footnote>

<footnote n="5"><sup>5</sup> Consortium FAIRmat c/o Humboldt-Universität zu Berlin, Zum Großen Windkanal 2, 12489 Berlin, Germany</footnote>

<footnote n="6"><sup>6</sup> Intel Corporation, Hillsboro, OR 97124, United States of America</footnote>

<footnote>* Author to whom any correspondence should be addressed.</footnote>

<footnote>Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.</footnote>
achieve the required electric field. When FE is carefully con-
trolled (by tuning V), atoms are liberated one by one from
the surface, ionized, and projected in the diverging electro-
static field existing close to the tip surface. The sample acts
itself as an electrostatic lens, projecting ions in a determin-
istic way. By collecting ions on a position sensitive detector
placed in front of the sample, the initial positions of all atoms
may be in turn be retrieved. In a single-phase material, FE
is smooth and regular, and the ion projection is straightfor-
ward. APT tomographic reconstructions of these materials
achieve excellent metrological performances [4]. However, in
materials composed of several phases of different chemical
composition, image reconstruction is often plagued by severe
distortions.

The major cause of distortion is the complex, dynamic evolution of samples by FE. This complexity is linked to the interplay between the structural and chemical composition of the materials, and the intensity of the electric field on the sample surface that is governed by the 3D geometry of the sample surface, i.e. the atomic arrangement. We may note that the temperature of the sample surface can also be non-homogenous giving rise to additional distortions. Appropriate reverse projection algorithms must be developed with a correct understanding of the dynamic evolution of the FE at the sample surface. As a result, the evaporation process has to be understood, in order to become fully quantify and thus to enable a possibility to systematic accounting for and correcting of as many of the artifacts inherent in the current reconstruction protocols as possible [5].

A common approach is to use numerical models, which
have proven to be useful to understand some known issues.
These simple approaches using analytical sample shape are
clearly insufficient [6]. Atomic scale approaches to FE were
developed in the early 2000s [7, 8]. These kinds of models
were sufficient to give an overview of the experimental short-
comings, such as improving the understanding and providing
correction methods for trajectory overlaps during evaporation
of small precipitates. For instance, in [8], the artifacts that
occur during FE of nano-oxides embedded in a metallic matrix
were quantified and corrected with this model. More recently,
the artifacts observed in the APT analysis of nanometric voids
induced by irradiation in metallic materials were successfully
explored. Simulated APT experiments demonstrate the local
density variations near voids are controlled by the unique ring
structures as the volume of voids gets exposed and the differ-
ent evaporation fields of the surrounding atoms [9].

Atomic scale models for FE model have shown their value
for understanding imaging artifacts. The basic principle is the
description of the emitter geometry at the atomic scale via an
adaptively, refined mesh. This discretization is used to solve
the electrostatic Laplace equation to obtain the 3D distribu-
tion of the voltage and the electric field. The atom exhibit-
ing the highest field is selected to be field evaporated [10].
The mesh is also useful to compute trajectories of desorbed
ions towards a virtual detector. Since then, a lot of improve-
ments have been introduced by different groups in the APT

community to model emitter structures quite comparable with
real samples [11–13]. Main modifications of FE models were
mesh refinements which enabled computations on more real-
istic specimen shape structure and external boundary condi-
tions for the electrostatic field and chamber [14]. Recently, a
new model was developed by Rolland et al [15] disregarding
the mesh structure. This model is based on a direct evaluation
of the charge distribution over the conducting emitter surface
as a solution of the electrostatic Robin’s equation [16, 17].
In practice, this is equivalent to solving an N-body problem
(where N is the number of atoms at the emitter surface). An
approximation for handling long-range interactions was intro-
duced by Rolland et al [15], developed by Barnes and Hut
[18], in order to reduce sequentially the time complexity for
solving the Robin’s equation from O(N²) to O(Nlog(N)).
Recently, FE simulations are no longer only used to understand
the evaporation process, but also used to perform tomographic
reconstruction [19].

The ultimate dream of using FE models to calculate reverse
projection, atom by atom was proposed by several authors [19,
20]. However, atomic models are still far away from giving
perfect knowledge of the FE process. For instance, zone line
artifacts, linked partly to rolling-up of atoms in the first steps
of flight are not fully treated by current models [21]. The prob-
lem is that atomistic models, so far developed [14, 15], do not
take into account some complex atomic scale phenomena such
as molecular dissociation [22], field induced surface diffusion
[23], or rolling-up motion [21]. To achieve a better descrip-
tion, modeling at a sub-atomic scale is required, which can
be done either by density functional theory (DFT [24, 25]),
or by molecular dynamics simulations [26–29]. Attempting
to incorporate sub-atomic physics in current models drastic-
ally increases computational complexity, resulting in a reduc-
tion of model sizes that can be physically simulated. On the
other hand, some necessary approximations are also used at
the macroscopic scale, so that atom probe chamber details
and limitations of the position-sensitive detector limitations
are often ignored in models. We note several limitations in
APT experimental conditions, such as detection efficiency
and detector spatial precision (detection efficiency between
30%–70%, multi-hit capability, field-of-view ±20–30°, spa-
tial resolution ~0.2 mm). As a result, the exact comparison
of atomic scale FE model and experimental data is clearly not
feasible.

However, it does not mean that a simulation tool cannot
be used to improve the reconstruction algorithm significantly.
Most of the image projection is governed by the sample shape,
thus an accurate understanding of this shape would enable
reconstructions with less distortion. The value of incorporating
atomic details is unclear, since a lot of subtle physical artifacts
may disturb ion trajectories in a non-deterministic way at this
scale. For this reason, but also for the sake of saving computing
time, FE models at the mesoscopic scale have been recently
developed. Mesoscopic models assume an evolving sample
surface with a defined curvature (tensor) field whose evolution
is simulated as a function of time. Numerically analyzable in
achieve the required electric field. When FE is carefully con-
trolled (by tuning V), atoms are liberated one by one from
the surface, ionized, and projected in the diverging electro-
static field existing close to the tip surface. The sample acts
itself as an electrostatic lens, projecting ions in a determin-
istic way. By collecting ions on a position sensitive detector
placed in front of the sample, the initial positions of all atoms
may be in turn be retrieved. In a single-phase material, FE
is smooth and regular, and the ion projection is straightfor-
ward. APT tomographic reconstructions of these materials
achieve excellent metrological performances [4]. However, in
materials composed of several phases of different chemical
composition, image reconstruction is often plagued by severe
distortions.

The major cause of distortion is the complex, dynamic evolution of samples by FE. This complexity is linked to the interplay between the structural and chemical composition of the materials, and the intensity of the electric field on the sample surface that is governed by the 3D geometry of the sample surface, i.e. the atomic arrangement. We may note that the temperature of the sample surface can also be non-homogenous giving rise to additional distortions. Appropriate reverse projection algorithms must be developed with a correct understanding of the dynamic evolution of the FE at the sample surface. As a result, the evaporation process has to be understood, in order to become fully quantify and thus to enable a possibility to systematic accounting for and correcting of as many of the artifacts inherent in the current reconstruction protocols as possible [5].

A common approach is to use numerical models, which
have proven to be useful to understand some known issues.
These simple approaches using analytical sample shape are
clearly insufficient [6]. Atomic scale approaches to FE were
developed in the early 2000s [7, 8]. These kinds of models
were sufficient to give an overview of the experimental short-
comings, such as improving the understanding and providing
correction methods for trajectory overlaps during evaporation
of small precipitates. For instance, in [8], the artifacts that
occur during FE of nano-oxides embedded in a metallic matrix
were quantified and corrected with this model. More recently,
the artifacts observed in the APT analysis of nanometric voids
induced by irradiation in metallic materials were successfully
explored. Simulated APT experiments demonstrate the local
density variations near voids are controlled by the unique ring
structures as the volume of voids gets exposed and the differ-
ent evaporation fields of the surrounding atoms [9].

Atomic scale models for FE model have shown their value
for understanding imaging artifacts. The basic principle is the
description of the emitter geometry at the atomic scale via an
adaptively, refined mesh. This discretization is used to solve
the electrostatic Laplace equation to obtain the 3D distribu-
tion of the voltage and the electric field. The atom exhibit-
ing the highest field is selected to be field evaporated [10].
The mesh is also useful to compute trajectories of desorbed
ions towards a virtual detector. Since then, a lot of improve-
ments have been introduced by different groups in the APT

community to model emitter structures quite comparable with
real samples [11–13]. Main modifications of FE models were
mesh refinements which enabled computations on more real-
istic specimen shape structure and external boundary condi-
tions for the electrostatic field and chamber [14]. Recently, a
new model was developed by Rolland et al [15] disregarding
the mesh structure. This model is based on a direct evaluation
of the charge distribution over the conducting emitter surface
as a solution of the electrostatic Robin’s equation [16, 17].
In practice, this is equivalent to solving an N-body problem
(where N is the number of atoms at the emitter surface). An
approximation for handling long-range interactions was intro-
duced by Rolland et al [15], developed by Barnes and Hut
[18], in order to reduce sequentially the time complexity for
solving the Robin’s equation from O(N²) to O(Nlog(N)).
Recently, FE simulations are no longer only used to understand
the evaporation process, but also used to perform tomographic
reconstruction [19].

The ultimate dream of using FE models to calculate reverse
projection, atom by atom was proposed by several authors [19,
20]. However, atomic models are still far away from giving
perfect knowledge of the FE process. For instance, zone line
artifacts, linked partly to rolling-up of atoms in the first steps
of flight are not fully treated by current models [21]. The prob-
lem is that atomistic models, so far developed [14, 15], do not
take into account some complex atomic scale phenomena such
as molecular dissociation [22], field induced surface diffusion
[23], or rolling-up motion [21]. To achieve a better descrip-
tion, modeling at a sub-atomic scale is required, which can
be done either by density functional theory (DFT [24, 25]),
or by molecular dynamics simulations [26–29]. Attempting
to incorporate sub-atomic physics in current models drastic-
ally increases computational complexity, resulting in a reduc-
tion of model sizes that can be physically simulated. On the
other hand, some necessary approximations are also used at
the macroscopic scale, so that atom probe chamber details
and limitations of the position-sensitive detector limitations
are often ignored in models. We note several limitations in
APT experimental conditions, such as detection efficiency
and detector spatial precision (detection efficiency between
30%–70%, multi-hit capability, field-of-view ±20–30°, spa-
tial resolution ~0.2 mm). As a result, the exact comparison
of atomic scale FE model and experimental data is clearly not
feasible.

However, it does not mean that a simulation tool cannot
be used to improve the reconstruction algorithm significantly.
Most of the image projection is governed by the sample shape,
thus an accurate understanding of this shape would enable
reconstructions with less distortion. The value of incorporating
atomic details is unclear, since a lot of subtle physical artifacts
may disturb ion trajectories in a non-deterministic way at this
scale. For this reason, but also for the sake of saving computing
time, FE models at the mesoscopic scale have been recently
developed. Mesoscopic models assume an evolving sample
surface with a defined curvature (tensor) field whose evolution
is simulated as a function of time. Numerically analyzable in
seconds or minutes on a laptop compared to days and months on cluster computers makes mesoscopic models useful tools for research on new methods of reconstruction [19].

Mesoscopic FE models were first developed by Haley and co-workers using level set approaches. First attempts using 2D approximations were found promising [30], and 3D extrapolation is in progress [31–33]. We propose in this paper an alternative approach to simulate FE at a mesoscopic scale (nm to µm scale). The model discretizes the evolving sample surface as a triangulated surface mesh. The evaporation sequence is simulated by estimating the evaporation rate of all the mesh primitives (vertices and facets of the triangles supporting the mesh). Note that the physical process of the FE is injected using a continuous model of FE, that is parameterized to handle the local nature of the material (FE constants, as a function of the specimen structure and composition), but also the local temperatures. Results of this modeling are presented for microstructures commonly encountered performing APT analysis. Results from this simulation are compared to experimental ones, demonstrating the validity and the relevance of this FE model at the mesoscopic scale.

## 2. Basic principles of the FE modeling

### 2.1. Three-dimensional geometry of the field emitter

The first step of the modeling process is to define the shape and the dimensions of the field emitter in 3D. As mentioned in the introduction, an APT sample is prepared as a thin needle with an initial radius of curvature at the apex in the range 10–100 nm. This allows, through application of an electrical potential (in the range of kV), the generation of a sufficiently intense electric field (few tens of V nm⁻¹) to ionize and then to evaporate the atoms located at the apex surface. Whereas the definition of the sample was previously done with the coordinates of each atom of the emitter [13, 15], the emitter surface in this approach is geometrically triangulated. The mesh size defines the minimum scale of features to be modeled as compared to an atomically refined approach. With this model, the sample information (and computing) is reduced to vertices of triangles of the surface mesh which discretizes the field emitter surface. For comparison, atomistic models need to store the information of atomic positions of the entire volume, to finally extract only surface atom positions [15]. The mesoscopic approach presents a true advantage in terms of memory allocated for calculations. At the atomic scale, the allocated memory is in the range of GB (at least) [8, 13, 15] while here it is, at best, one hundred times smaller (10–50 MB). For the initial emitter surface (i.e. before evaporation), the coordinates of the vertices are defined, such as after meshing the surface, with a triangulation, each mesh has the same and the expected area ($S_{at}$) [34]. Figure 1 represents such an emitter with an initial radius of curvature ($R_c$) of 20 nm, a shank angle ($\gamma$) of 10°, a length ($L$) of 100 nm and a mesh element area ($S_{at}$) of 4 N m² per triangle.

To characterize a surface and its evolution during evaporation (presented in the next section), the radius of curvature is

<figure><img src="image_7.png" /><figcaption>Figure 1. Triangulated APT sample ($R_c$ = 20 nm, $L$ = 100 nm, $g$ = 10° and $S_{at}$ = 4 nm²).</figcaption></figure>

commonly used. Indeed, the correlation between the electric field $F$ and the radius of curvature $R_c$ is often considered using the standard expression [35]:

$$F = \frac{V}{k_f R_c} \qquad (1)$$

with $k_f$ the field factor (linked to the specimen morphology [4, 6]) and $V$ the applied potential to the emitter. The radius of curvature is the reciprocal of the curvature, which is determined here using the approach of S. Rusinkiewicz [36]. Curvature quantification is done generally by calculating the mean curvature (average of the principal curvatures) or the Gaussian curvature (product of the principal curvatures). In this study, the mean curvature of each mesh is computed (as done by [36]) to estimate its associated radius of curvature (figure 2(a)). The values thus obtained (figure 2(b)) agree with the simulation input parameters (figure 1). Indeed, the radius of curvature distribution is a bimodal distribution with one mode centered to the set apex radius of curvature (i.e. 20 nm in figure 1) and the other mode with higher radius of curvature (>40 nm) corresponding to the shank (due to the shank angle, this last curvature radius continuously increases as a function of distance to the tip apex).

### 2.2. FE

The FE process of this simulated system is based on the local FE rate at the surface of a material submitted to an electric field $F$. Assuming that the atoms evaporate from a surface at a
seconds or minutes on a laptop compared to days and months on cluster computers makes mesoscopic models useful tools for research on new methods of reconstruction [19].

Mesoscopic FE models were first developed by Haley and co-workers using level set approaches. First attempts using 2D approximations were found promising [30], and 3D extrapolation is in progress [31–33]. We propose in this paper an alternative approach to simulate FE at a mesoscopic scale (nm to µm scale). The model discretizes the evolving sample surface as a triangulated surface mesh. The evaporation sequence is simulated by estimating the evaporation rate of all the mesh primitives (vertices and facets of the triangles supporting the mesh). Note that the physical process of the FE is injected using a continuous model of FE, that is parameterized to handle the local nature of the material (FE constants, as a function of the specimen structure and composition), but also the local temperatures. Results of this modeling are presented for microstructures commonly encountered performing APT analysis. Results from this simulation are compared to experimental ones, demonstrating the validity and the relevance of this FE model at the mesoscopic scale.

## 2. Basic principles of the FE modeling

### 2.1. Three-dimensional geometry of the field emitter

The first step of the modeling process is to define the shape and the dimensions of the field emitter in 3D. As mentioned in the introduction, an APT sample is prepared as a thin needle with an initial radius of curvature at the apex in the range 10–100 nm. This allows, through application of an electrical potential (in the range of kV), the generation of a sufficiently intense electric field (few tens of V nm⁻¹) to ionize and then to evaporate the atoms located at the apex surface. Whereas the definition of the sample was previously done with the coordinates of each atom of the emitter [13, 15], the emitter surface in this approach is geometrically triangulated. The mesh size defines the minimum scale of features to be modeled as compared to an atomically refined approach. With this model, the sample information (and computing) is reduced to vertices of triangles of the surface mesh which discretizes the field emitter surface. For comparison, atomistic models need to store the information of atomic positions of the entire volume, to finally extract only surface atom positions [15]. The mesoscopic approach presents a true advantage in terms of memory allocated for calculations. At the atomic scale, the allocated memory is in the range of GB (at least) [8, 13, 15] while here it is, at best, one hundred times smaller (10–50 MB). For the initial emitter surface (i.e. before evaporation), the coordinates of the vertices are defined, such as after meshing the surface, with a triangulation, each mesh has the same and the expected area (S<sub>at</sub>) [34]. Figure 1 represents such an emitter with an initial radius of curvature (R<sub>c</sub>) of 20 nm, a shank angle (γ) of 10°, a length (L) of 100 nm and a mesh element area (S<sub>at</sub>) of 4 N m² per triangle.

To characterize a surface and its evolution during evaporation (presented in the next section), the radius of curvature is

<figure><img src="image_7.png" /><figcaption>Figure 1. Triangulated APT sample (R<sub>c</sub> = 20 nm, L = 100 nm, g = 10° and S<sub>at</sub> = 4 nm²).</figcaption></figure>

commonly used. Indeed, the correlation between the electric field F and the radius of curvature R<sub>c</sub> is often considered using the standard expression [35]:

$$F = \frac{V}{k_f R_c} \qquad (1)$$

with k<sub>f</sub> the field factor (linked to the specimen morphology [4, 6]) and V the applied potential to the emitter. The radius of curvature is the reciprocal of the curvature, which is determined here using the approach of S. Rusinkiewicz [36]. Curvature quantification is done generally by calculating the mean curvature (average of the principal curvatures) or the Gaussian curvature (product of the principal curvatures). In this study, the mean curvature of each mesh is computed (as done by [36]) to estimate its associated radius of curvature (figure 2(a)). The values thus obtained (figure 2(b)) agree with the simulation input parameters (figure 1). Indeed, the radius of curvature distribution is a bimodal distribution with one mode centered to the set apex radius of curvature (i.e. 20 nm in figure 1) and the other mode with higher radius of curvature (>40 nm) corresponding to the shank (due to the shank angle, this last curvature radius continuously increases as a function of distance to the tip apex).

### 2.2. FE

The FE process of this simulated system is based on the local FE rate at the surface of a material submitted to an electric field F. Assuming that the atoms evaporate from a surface at a
<figure><img src="image_8.png" /><figcaption>Figure 2. (a) Radius of curvature (nm) associated with each mesh that composed the surface defined in figure 1. (b) Radius of curvature distribution (%) from the tip in figure 2(a).</figcaption></figure>

defined rate, the surface center will move downward along its normal direction. From the evaporation rate equation, assigned to each primitive of the surface mesh (figure 1), the volume variation can be estimated. Using geometrical considerations, it is then possible to estimate the displacement vector associated with each vertex. A new surface is thus obtained. This operation is iterated until a specified end condition is reached (commonly an a priori defined evaporation depth).

Experimentally, the evaporation process is triggered using pulses (electric or laser). Regardless of the pulsing mode, the evaporation rate (K) for a single atom is described empirically by an Arrhenius-type law. Extrapolating to our triangulated system, the evaporation rate of a mesh j (Kj) is then expressed as:

$$K_j \propto \exp\left(-\frac{Q(F_j)}{k_B T_j}\right) \qquad (2)$$

with the activation energy barrier Q (since this process is usually considered thermally assisted), Q (Fj) is the energy barrier considering an external electric field (Fj), kB the Boltzmann constant and Tj is the temperature. A linear approximation can be used to describe the dependence of the energy barrier and the electric field [37]:

$$Q(F_j) = C_j \left(1 - \frac{F_j}{F_{EV,j}}\right) \qquad (3)$$

with Cj as an energetic constant and FEVj the evaporation field. The evaporation field is the value of the field at which the energy barrier is reduced to zero. Note that more refined equations defining Q (F) can be injected in the model without major modifications of the algorithm [1, 38].

In the case of pulsed evaporation, the expressions (2) and (3) are time dependent, since both F or T may be time varying. If the evaporation process is triggered by electrical pulses, then the base temperature of the sample can be considered as a constant. Therefore, the temperature of each mesh (Tj) is set as a constant, defined by the user according to its APT analysis conditions (figure 4(a) with a constant temperature equal to 50 K). Most of the evaporation occurs at the top of electrical pulses. We can then model the evaporation rate per pulse by:

$$\phi_F \sim \tau_F \exp\left(-\frac{C}{k_B T_j}\left(1 - \frac{F_j}{F_{EV,j}}\right)\right). \qquad (4)$$

With τF being a fraction of the voltage pulse applied to the sample [1], and Fj being the peak electric field induced by VDC + Vp.

However, if the analysis is triggered by laser pulses, the temperature can no longer be considered as a constant [39, 40]. Indeed, the interaction of the laser with a nanometric needle induces local enhancement of the temperature, which in turn triggers different local probabilities for the FE. During the pulse, each triangle is submitted to its own temporal evolution of temperature, with a temperature maximum Tmax reached in less than a nanosecond in the general case, followed by a temperature decay that can vary slightly on the surface. Assuming a simple shape of transient temperature pulse of width τT (close to the decay time constant of a few ns), the atom flux φT can be approximated by [41]:

$$\phi_T \sim \left( \frac{k_B T_{\max}}{C \left(1 - \frac{F_j}{F_{EV,j}}\right)} \right) \times \tau_T \exp \left( - \frac{C}{k_B T_{\max}} \left( 1 - \frac{F_j}{F_{EV,j}} \right) \right). \qquad (5)$$
<figure><img src="image_8.png" /><figcaption>**Figure 2.** (a) Radius of curvature (nm) associated with each mesh mesh that composed the surface defined in figure 1. (b) Radius of curvature distribution (%) from the tip in figure 2(a).</figcaption></figure>

defined rate, the surface center will move downward along its normal direction. From the evaporation rate equation, assigned to each primitive of the surface mesh (figure 1), the volume variation can be estimated. Using geometrical considerations, it is then possible to estimate the displacement vector associated with each vertex. A new surface is thus obtained. This operation is iterated until a specified end condition is reached (commonly an *a priori* defined evaporation depth).

Experimentally, the evaporation process is triggered using pulses (electric or laser). Regardless of the pulsing mode, the evaporation rate (K) for a single atom is described empirically by an Arrhenius-type law. Extrapolating to our triangulated system, the evaporation rate of a mesh j (Kj) is then expressed as:

$$K_j \propto \exp\left(-\frac{Q(F_j)}{k_B T_j}\right) \quad (2)$$

with the activation energy barrier Q (since this process is usually considered thermally assisted), Q(Fj) is the energy barrier considering an external electric field (Fj), kB the Boltzmann constant and Tj is the temperature. A linear approximation can be used to describe the dependence of the energy barrier and the electric field [37]:

$$Q(F_j) = C_j \left(1 - \frac{F_j}{F_{EV,j}}\right) \quad (3)$$

with Cj as an energetic constant and FEVj the evaporation field. The evaporation field is the value of the field at which the energy barrier is reduced to zero. Note that more refined equations defining Q(F) can be injected in the model without major modifications of the algorithm [1, 38].

In the case of pulsed evaporation, the expressions (2) and (3) are time dependent, since both F or T may be time varying. If the evaporation process is triggered by electrical pulses, then the base temperature of the sample can be considered as a constant. Therefore, the temperature of each mesh (Tj) is set as a constant, defined by the user according to its APT analysis conditions (figure 4(a) with a constant temperature equal to 50 K). Most of the evaporation occurs at the top of electrical pulses. We can then model the evaporation rate per pulse by:

$$\phi_F \sim \tau_F \exp\left(-\frac{C}{k_B T_j}\left(1 - \frac{F_j}{F_{EV,j}}\right)\right). \quad (4)$$

With τF being a fraction of the voltage pulse applied to the sample [1], and Fj being the peak electric field induced by VDC + Vp.

However, if the analysis is triggered by laser pulses, the temperature can no longer be considered as a constant [39, 40]. Indeed, the interaction of the laser with a nanometric needle induces local enhancement of the temperature, which in turn triggers different local probabilities for the FE. During the pulse, each triangle is submitted to its own temporal evolution of temperature, with a temperature maximum Tmax reached in less than a nanosecond in the general case, followed by a temperature decay that can vary slightly on the surface. Assuming a simple shape of transient temperature pulse of width τT (close to the decay time constant of a few ns), the atom flux φT can be approximated by [41]:

$$\phi_T \sim \left( \frac{k_B T_{\max}}{C \left(1 - \frac{F_j}{F_{EV,j}}\right)} \right) \times \tau_T \exp \left( - \frac{C}{k_B T_{\max}} \left( 1 - \frac{F_j}{F_{EV,j}} \right) \right). \quad (5)$$
<figure><img src="image_9.png" /><figcaption>**Figure 3.** Schematic description of a vertex *i* and its connected triangle fan (CMᵢ). It is also reported for a mesh *j* that is a CMᵢ, its evaporation rate (Kⱼ) and its associated normal direction (nⱼ).</figcaption></figure>

Since the prefactor of the flux equation is linearly dependent of T_max, or F, the exponential term is the main source of flux variation as a function of F or T. In a first approach, the use of a constant prefactor τ'_{T} = (k_B T_{max} / (C(1 - F_{EV,j}))) × τ_{T} of about 1 ns is a fair approximation.

The displacement vector of a defined vertex *i* ($\vec{D}_i$) results then from the evaporation rate of all the mesh ($K_j$) that are connected to it ($j \in CM_i$ in figure 3), considering their outward normal direction ($\vec{n}_j$), since the electric field at the surface ($F_j$) is collinear to the latter:

$$ \vec{D}_i = -\epsilon \times \left( \sum_{j \in CM_i} K_j \vec{n}_j \right). \qquad (6) $$

A constant parameter ($\epsilon$) set by the user is added in the expression determining the vector displacement. Thus, it is considered that at each pulse (i.e. here each evaporation simulation iteration), the sample is evaporated by a certain amount. In order to be closest to the experiment, this quantity must be small ($\epsilon < 0.1$ nm). This quantity is in fact the time step of the numerical scheme. Since the maximum evaporation rate is fixed, the displacement in depth is proportional to the time step between two iterations. The stability of the scheme, of first order in this simple approach, is strongly dependent on this quantity. Systematic tests were performed on different geometries to evaluate the propagation of errors of this approach for nanometer-size meshes with $\epsilon < 0.1$ nm found experimentally as a good compromise between scheme stability and efficiency of the method.

Thus, at each iteration, it is necessary to determine four parameters for each triangle mesh j, in order to estimate the FE process of the sample, these are: Tⱼ, Cⱼ, Fⱼ and FEVⱼ. In fact, only three parameters are here considered, since the energetic constant (Cⱼ) will be fixed as a constant for each mesh. This constant, as documented in [3, 42] was only determined for very few materials, with values in the range 0.5–3 eV. It has to

be set by the user in the range of tenth to few eV to correspond to experimental measurements [37, 43–46].

## 2.3. Temperature and phase

In laser assisted atom probe, the maximum temperature is asymmetric on the tip apex, inducing then schematically an illuminated so-called laser side and a shadow side. The temperature distribution along the apex due to the laser interaction can be accurately reproduced with advanced models, such as proposed by Houard et al [47] (but also exhibiting in return a huge time computing). In this study, however, an analytical model is preferred, with a polynomial description of the distribution of peak temperature on the sample (expression (7)). This expression can reproduce the asymmetric distribution of the temperature (figure 4(b)),

$$ T_{(X_i, Y_i)} = T_0 \left( 1 + \alpha(X_i, Y_i)^\beta \right). \qquad (7) $$

With Xᵢ and Yᵢ the coordinates of the vertex *i*, α and β are parameters defined by the user to introduce different temperature evolution trend and T₀ the base temperature also the temperature in the shadow side of the sample. This last expression provides the temperature that is assigned to each vertex. The temperature that is associated to each triangle (Tⱼ) is the average temperature of its vertices.

The different phases in a material are distinguished by their evaporation fields. In the case of a single-phase material, the evaporation field of all the mesh (FEVⱼ) can be then considered as a constant (in the range of tens V nm⁻¹). However, in multiphase materials, before each evaporation iteration, the evaporation field associated to each vertex needs to be defined. To evaluate these evaporation fields, at each evaporation step, the position of each vertex is compared to the position inside a 3D volume that defines the phases of interest. This volume contains the information of shapes, sizes, (particles, layers, ...) and elemental nature (evaporation field) that are used as input data.

## 2.4. Electric field

It remains now the most important point, which is the calculation of the electric field for each mesh (Fⱼ). In this study, the electric field is obtained from the estimation of the charge distribution on the emitter surface. This model, as is the one developed by Rolland et al [15], is based on the direct evaluation of the charge distribution over the emitter conducting surface, as a solution of the electrostatic Robin’s equation [16, 17]. As demonstrated by Rolland and collaborators, and extrapolated in this modeling approach, the charge qⱼ (a positive non-integral partial charge) on a defined mesh j having an area (S_at,j estimated at each evaporation iteration using the Heron formula) can be written as:
<figure><img src="image_9.png" /><figcaption>**Figure 3.** Schematic description of a vertex *i* and its connected triangle fan (*CM*<sub>*i*</sub>). It is also reported for a mesh *j* that is a *CM*<sub>*i*</sub>, its evaporation rate (*K*<sub>*j*</sub>) and its associated normal direction (*n*<sub>*j*</sub>).</figcaption></figure>

Since the prefactor of the flux equation is linearly dependent of $T_{\text{max}}$, or $F$, the exponential term is the main source of flux variation as a function of $F$ or $T$. In a first approach, the use of a constant prefactor $\tau_T' = (\frac{k_B T_{\text{max}}}{C(1-\frac{F_{\text{EV,j}}}{F_{\text{EV}}})}) \times \tau_T$ of about 1 ns is a fair approximation.

The displacement vector of a defined vertex $i$ ($\vec{D}_i$) results then from the evaporation rate of all the mesh ($K_j$) that are connected to it ($j \in CM_i$ in figure 3), considering their outward normal direction ($\vec{n}_j$), since the electric field at the surface ($F_j$) is collinear to the latter:

$$ \vec{D}_i = -\epsilon \times \left( \sum_{j \in CM_i} K_j \vec{n}_j \right). \qquad (6) $$

A constant parameter ($\epsilon$) set by the user is added in the expression determining the vector displacement. Thus, it is considered that at each pulse (i.e. here each evaporation simulation iteration), the sample is evaporated by a certain amount. In order to be closest to the experiment, this quantity must be small ($\epsilon < 0.1$ nm). This quantity is in fact the time step of the numerical scheme. Since the maximum evaporation rate is fixed, the displacement in depth is proportional to the time step between two iterations. The stability of the scheme, of first order in this simple approach, is strongly dependent on this quantity. Systematic tests were performed on different geometries to evaluate the propagation of errors of this approach for nanometer-size meshes with $\epsilon < 0.1$ nm found experimentally as a good compromise between scheme stability and efficiency of the method.

Thus, at each iteration, it is necessary to determine four parameters for each triangle mesh *j*, in order to estimate the FE process of the sample, these are: $T_j$, $C_j$, $F_j$ and $F_{\text{EV},j}$. In fact, only three parameters are here considered, since the energetic constant ($C_j$) will be fixed as a constant for each mesh. This constant, as documented in [3, 42] was only determined for very few materials, with values in the range 0.5–3 eV. It has to

be set by the user in the range of tenth to few eV to correspond to experimental measurements [37, 43–46].

## 2.3. Temperature and phase

In laser assisted atom probe, the maximum temperature is asymmetric on the tip apex, inducing then schematically an illuminated so-called laser side and a shadow side. The temperature distribution along the apex due to the laser interaction can be accurately reproduced with advanced models, such as proposed by Houard et al [47] (but also exhibiting in return a huge time computing). In this study, however, an analytical model is preferred, with a polynomial description of the distribution of peak temperature on the sample (expression (7)). This expression can reproduce the asymmetric distribution of the temperature (figure 4(b)),

$$ T_{(X_i, Y_i)} = T_0 \left( 1 + \alpha(X_i, Y_i)^\beta \right). \qquad (7) $$

With $X_i$ and $Y_i$ the coordinates of the vertex $i$, $\alpha$ and $\beta$ are parameters defined by the user to introduce different temperature evolution trend and $T_0$ the base temperature also the temperature in the shadow side of the sample. This last expression provides the temperature that is assigned to each vertex. The temperature that is associated to each triangle ($T_j$) is the average temperature of its vertices.

The different phases in a material are distinguished by their evaporation fields. In the case of a single-phase material, the evaporation field of all the mesh ($F_{\text{EV},j}$) can be then considered as a constant (in the range of tens V nm<sup>-1</sup>). However, in multiphase materials, before each evaporation iteration, the evaporation field associated to each vertex needs to be defined. To evaluate these evaporation fields, at each evaporation step, the position of each vertex is compared to the position inside a 3D volume that defines the phases of interest. This volume contains the information of shapes, sizes, (particles, layers, ...) and elemental nature (evaporation field) that are used as input data.

## 2.4. Electric field

It remains now the most important point, which is the calculation of the electric field for each mesh ($F_j$). In this study, the electric field is obtained from the estimation of the charge distribution on the emitter surface. This model, as is the one developed by Rolland et al [15], is based on the direct evaluation of the charge distribution over the emitter conducting surface, as a solution of the electrostatic Robin’s equation [16, 17]. As demonstrated by Rolland and collaborators, and extrapolated in this modeling approach, the charge $q_j$ (a positive non-integral partial charge) on a defined mesh *j* having an area ($S_{\text{at},j}$ estimated at each evaporation iteration using the Heron formula) can be written as:
<figure><img src="image_10.png" /><figcaption>Figure 4. Temperature associated to each mesh that composed the surface defined in figure 1, if the temperature is set as (a) constant ($T = 50$ K) and as (b) polynomial evolution along X direction.</figcaption></figure>

$$ \frac{q_j}{S_{at,j}} = \frac{1}{2\pi} \sum_{k=1, k \neq j}^{N} q_k \frac{\overrightarrow{n_j} \cdot \overrightarrow{r_{j,k}}}{\| \overrightarrow{r_{j,k}} \|^3} \quad (8) $$

With $\overrightarrow{n_j}$ is the outward normal direction to the mesh j (as in expression (6)), $\overrightarrow{r_{j,k}}$ the vector between the mass center of the mesh j and k, and N the number of meshes at the emitter surface. The mass center of each mesh is easily estimated from its vertex coordinates. It is then possible to construct a sequence $q_{j,n}$ for each mesh j to find the charge distribution over the field emitter surface:

$$ \frac{q_{j,n+1}}{S_{at,j}} = \frac{1}{2\pi} \sum_{k=1, k \neq j}^{N} q_{k,n} \frac{\overrightarrow{n_j} \cdot \overrightarrow{r_{j,k}}}{\| \overrightarrow{r_{j,k}} \|^3} \quad (9) $$

With $q_{j,0}$ being an arbitrary positive number. It has been proven that this sequence tends toward the charge distribution at the equilibrium [48]. Moreover, this sequence has the advantage to converge very rapidly to the equilibrium (figure 5), even when starting from a random distribution of charge. Indeed, in the present case (figure 1), less than 15 iterations are required to reach the equilibrium (figure 5).

As expected, the highest charges are located at the apex of the emitter which is also in those regions where the electric field is the highest (electric field amounts to $q_j/S_{at,j}\epsilon_0$ or $\sigma/\epsilon_0$ with $\sigma = q_j/S_{at,j}$ in accordance with basic electrostatic laws). We assumed that the equilibrium charge distribution is reached through equation (9), if the average relative deviation of the charge distribution between two iterations ($\langle \delta q_{j,n+1} \rangle$) is lower than $10^{-5}\%$ with

$$ \langle \delta q_{j,n+1} \rangle = \frac{1}{N} \sum_{j=1}^{N} \left| \frac{q_{j,n+1} - q_{j,n}}{q_{j,n+1}} \right| \quad (10) $$

<figure><img src="image_11.png" /><figcaption>Figure 5. Charge density distribution ($q/S_{at}$ in a.u.) at the emitter surface.</figcaption></figure>

It must be pointed out that the total charge carried by the conductor remains the same throughout the iterative process. To ensure the same behavior for the expression (9), $q_{j,n}/S_{at,j}$ will be computed, at each iteration, such that:

$$ \sum_{j}^{N} q_{j,n+1} = \sum_{j}^{N} q_{j,n} \quad (11) $$
<figure><img src="image_10.png" /><figcaption>Figure 4. Temperature associated to each mesh that composed the surface defined in figure 1, if the temperature is set as (a) constant ($T = 50$ K) and as (b) polynomial evolution along X direction.</figcaption></figure>

$$ \frac{q_j}{S_{at,j}} = \frac{1}{2\pi} \sum_{k=1, k \neq j}^{N} q_k \frac{\overrightarrow{n_j} \cdot \overrightarrow{r_{j,k}}}{\| \overrightarrow{r_{j,k}} \|^3} \quad (8) $$

With $\overrightarrow{n_j}$ is the outward normal direction to the mesh j (as in expression (6)), $\overrightarrow{r_{j,k}}$ the vector between the mass center of the mesh j and k, and N the number of meshes at the emitter surface. The mass center of each mesh is easily estimated from its vertex coordinates. It is then possible to construct a sequence $q_{j,n}$ for each mesh j to find the charge distribution over the field emitter surface:

$$ \frac{q_{j,n+1}}{S_{at,j}} = \frac{1}{2\pi} \sum_{k=1, k \neq j}^{N} q_{k,n} \frac{\overrightarrow{n_j} \cdot \overrightarrow{r_{j,k}}}{\| \overrightarrow{r_{j,k}} \|^3} \quad (9) $$

With $q_{j,0}$ being an arbitrary positive number. It has been proven that this sequence tends toward the charge distribution at the equilibrium [48]. Moreover, this sequence has the advantage to converge very rapidly to the equilibrium (figure 5), even when starting from a random distribution of charge. Indeed, in the present case (figure 1), less than 15 iterations are required to reach the equilibrium (figure 5).

As expected, the highest charges are located at the apex of the emitter which is also in those regions where the electric field is the highest (electric field amounts to $q_j/S_{at,j}\epsilon_0$ or $\sigma/\epsilon_0$ with $\sigma = q_j/S_{at,j}$ in accordance with basic electrostatic laws). We assumed that the equilibrium charge distribution is reached through equation (9), if the average relative deviation of the charge distribution between two iterations ($\langle\delta q_{j,n+1}\rangle$) is lower than $10^{-5}\%$ with

$$ \langle \delta q_{j,n+1} \rangle = \frac{1}{N} \sum_{j=1}^{N} \left| \frac{q_{j,n+1} - q_{j,n}}{q_{j,n+1}} \right| \quad (10) $$

<figure><img src="image_11.png" /><figcaption>Figure 5. Charge density distribution ($q/S_{at}$ in a.u.) at the emitter surface.</figcaption></figure>

It must be pointed out that the total charge carried by the conductor remains the same throughout the iterative process. To ensure the same behavior for the expression (9), $q_{j,n}/S_{at,j}$ will be computed, at each iteration, such that:

$$ \sum_{j}^{N} q_{j,n+1} = \sum_{j}^{N} q_{j,n} \quad (11) $$
<figure><img src="image_12.png" /><figcaption>Figure 6. (a) Velocity $V$ (nm s$^{-1}$) and velocity deviation ($\Delta V$) for an ion as a function of its distance to sample apex (nm). (b) ICF ($\xi$) distribution calculated with equation (11) (%) according to the compression coefficient $\lambda$.</figcaption></figure>

## 2.5. Trajectory computation

Up to now, the model describes the evolution of the emitter morphology (and the electrical field which derives from the charge distribution), which is an important piece of data to understand APT analyses. By solving Newton’s equation of motion for ions from vertices up to a planar detector situated a few centimeters in front of the sample for each step of the FE process, the dynamic and evolving projection APT law can be deduced. A numeric integration of this equation is conventionally done with an embedded fourth order Runge Kutta (RK4) method (at least for the last updated models [12, 13, 15]) or even a Leap Frog (LF) method for the oldest model [8]). For each method, there are pros and cons. In order to get the best compromise between speed and accuracy, the Verlet method was chosen in this model [49]. The Verlet method combines a good spatial accuracy (comparable to the RK4 method) and a low computation time (comparable to the LF method). Moreover, as it is conventionally done in FE modeling an adaptive time step calculation in solving the motion equation is introduced, since the emitter is at nanometric scale (≈100 nm) and the detector is located at metric scale from the emitter (≈0.1 m). Such adaptive time scale (dt) is determined using the velocity ($V$) deviation between two iterations:

$$dt_{(n+1)} = t_0 \cdot \Delta V = t_0 \cdot \frac{V_{(n)}}{V_{(n)} - V_{(n-1)}}. \quad (12)$$

With $t_0$ the initial time calculation (≤10$^{-16}$ for an accurate computation of the trajectory). The velocity deviation provides a more adaptive time step calculation than the velocity itself. Indeed, up to a certain distance from the emitter apex (few times the radius of curvature, figure 6(a)), the ion velocity increases drastically and then the velocity stabilizes (very low

increase), while the velocity deviation ($\Delta V$) evolves continuously up to the detector. With the adaptive approach, it requires less than 100 iterations to calculate the trajectory up to a detector located at 10 cm from the emitter apex. In this way, most of the iterations are used to compute the trajectory close to the sample apex, where they are most prone to lead to trajectory estimation imprecision related to surface topography and the rest of the iterations are used for trajectory computation far away for the sample, where they are almost straight lines [1, 2].

Experimentally, it is well known that due to the sample morphology and the analysis chamber setup (in particular an electrode, on which the electrical pulses are applied), the ion trajectories are not straight lines, but become compressed [50]. This compression is defined and quantified by a parameter, the so-called, image compression factor (ICF or $\xi$). The parameter is expressed for each meshes $j$, considering small angle approximations:

$$\xi_j = \frac{L}{R_{c,j}} \times \sqrt{\frac{S_{at,j}}{S_{D,j}}} \quad (13)$$

with $L$ the flight length (distance between the sample and the detector), $R_{c,j}$ the radius of curvature of the mesh $j$ and $S_{D,j}$ the mesh area on the detector. The experimental value of $\xi$ is in the range of 1.3–1.6. Assuming an emitter of reduced dimensions (less than a micron in depth), the calculated ICF with the model is naturally reduced as compared to experimental values (blue histogram on the figure 6(b)). To achieve the exact experimental ICF value, it would be necessary to model the full geometry of the sample (needle on its coupon or capillary, local electrode atom probe (LEAP) sample holder, and local electrode ahead). We observed that increasing the length of the simulated sample, increases the average ICF,
<figure><img src="image_12.png" /><figcaption>Figure 6. (a) Velocity $V$ (nm s$^{-1}$) and velocity deviation ($\Delta V$) for an ion as a function of its distance to sample apex (nm). (b) ICF ($\xi$) distribution calculated with equation (11) (%) according to the compression coefficient $\lambda$.</figcaption></figure>

## 2.5. Trajectory computation

Up to now, the model describes the evolution of the emitter morphology (and the electrical field which derives from the charge distribution), which is an important piece of data to understand APT analyses. By solving Newton’s equation of motion for ions from vertices up to a planar detector situated a few centimeters in front of the sample for each step of the FE process, the dynamic and evolving projection APT law can be deduced. A numeric integration of this equation is conventionally done with an embedded fourth order Runge Kutta (RK4) method (at least for the last updated models [12, 13, 15]) or even a Leap Frog (LF) method for the oldest model [8]). For each method, there are pros and cons. In order to get the best compromise between speed and accuracy, the Verlet method was chosen in this model [49]. The Verlet method combines a good spatial accuracy (comparable to the RK4 method) and a low computation time (comparable to the LF method). Moreover, as it is conventionally done in FE modeling an adaptive time step calculation in solving the motion equation is introduced, since the emitter is at nanometric scale (≈100 nm) and the detector is located at metric scale from the emitter (≈0.1 m). Such adaptive time scale (dt) is determined using the velocity (V) deviation between two iterations:

$$dt_{(n+1)} = t_0 \cdot \Delta V = t_0 \cdot \frac{V_{(n)}}{V_{(n)} - V_{(n-1)}}. \quad (12)$$

With $t_0$ the initial time calculation (≤10$^{-16}$ for an accurate computation of the trajectory). The velocity deviation provides a more adaptive time step calculation than the velocity itself. Indeed, up to a certain distance from the emitter apex (few times the radius of curvature, figure 6(a)), the ion velocity increases drastically and then the velocity stabilizes (very low

increase), while the velocity deviation ($\Delta V$) evolves continuously up to the detector. With the adaptive approach, it requires less than 100 iterations to calculate the trajectory up to a detector located at 10 cm from the emitter apex. In this way, most of the iterations are used to compute the trajectory close to the sample apex, where they are most prone to lead to trajectory estimation imprecision related to surface topography and the rest of the iterations are used for trajectory computation far away for the sample, where they are almost straight lines [1, 2].

Experimentally, it is well known that due to the sample morphology and the analysis chamber setup (in particular an electrode, on which the electrical pulses are applied), the ion trajectories are not straight lines, but become compressed [50]. This compression is defined and quantified by a parameter, the so-called, image compression factor (ICF or $\xi$). The parameter is expressed for each meshes $j$, considering small angle approximations:

$$\xi_j = \frac{L}{R_{c,j}} \times \sqrt{\frac{S_{at,j}}{S_{D,j}}} \quad (13)$$

with $L$ the flight length (distance between the sample and the detector), $R_{c,j}$ the radius of curvature of the mesh $j$ and $S_{D,j}$ the mesh area on the detector. The experimental value of $\xi$ is in the range of 1.3–1.6. Assuming an emitter of reduced dimensions (less than a micron in depth), the calculated ICF with the model is naturally reduced as compared to experimental values (blue histogram on the figure 6(b)). To achieve the exact experimental ICF value, it would be necessary to model the full geometry of the sample (needle on its coupon or capillary, local electrode atom probe (LEAP) sample holder, and local electrode ahead). We observed that increasing the length of the simulated sample, increases the average ICF,
which corroborates this. Obviously, increasing the length of the tip, also increases the number of meshes and therefore the time computing. In order to maintain a reasonable computation time (sample length around hundreds of nanometers) and to have an ICF equivalent to the experimental one, we performed a post analytical compression of the trajectories, using an ad hoc coefficient $\lambda$. After a certain distance from the sample (in the range of few $\mu$m corresponding to the distance to an electrode and where the trajectories are straight lines), the trajectory calculation is stopped, and the position on the detector is obtained by a projection, using the velocity at this step, considering a compression of its polar angle, by $\theta = \lambda\theta'$ (with $\theta$ and $\theta'$ respectively the polar angle of the velocity where the calculation is stopped and on the detector), while its azimuthal angle remains unchanged $\varphi = \varphi'$ (with $\varphi$ and $\varphi'$ respectively the azimuthal angle of the velocity at the same step as polar angle), as experimentally observed [1]. Using this approach with a coefficient $\lambda$ equal to 1.5 (red histogram on the figure 6(b)), the average ICF (i.e. 1.53) is then comparable to experiment value. This coefficient is adjustable by the user in order to best adjust to the real analysis.

## 2.6. Algorithmic optimization

One of the crucial points of any simulation model for practical application are computing time and memory consumption. In this section, we show some solutions added to the model to address this. We are not benchmarking the performance of them, as done in [51], but just give a brief overview of the improvements, in term of computing time (table 1).

As it was mentioned previously, the electrical field is estimated by solving the Robin's equation with an $O(N^2)$ complexity (with $N$ the number of meshes). Thus, the first introduced optimization, is a parallelized computation (table 1). Observe that in equation (9), step $n+1$ is uniquely dependent on fixed values at step $n$, with demand for an exchange of pieces of information for each point of charge, parallelization is strongly efficient with a scaling close to the number of cores (Using C# parallel loop). The second improvement, as Rolland et al [15] at atomic scale, is a long range interaction approximation, developed by Barnes and Hut [18], in order to drastically reduce the computing time required for solving Robin's equation solving by reducing the complexity to $O(\text{Log}(N))$. The time computing is directly linked to the number of meshes $N$, so we also introduced logically a reduction of this number (table 1). Regarding the electrical field distribution (figure 5), the evaporation will mainly occur at the sample apex. During the whole process of FE, the model consists in moving the meshes as a function of the evaporation rate (the latter of which is field dependent). Obviously, a refined mesh at the bottom of sample is unnecessary as the electric field on these regions is much lower than the electric field on the apex. An average charge on these meshes is also sufficient to describe the long-range interaction for ion trajectory calculations. As a result, we introduced an evolving mesh area that depends on the sample

<table><caption>**Table 1.** Time computing (second) to perform ten field evaporation iterations ($\approx 0.3-0.6$ nm depth erosion in this case) as a function of the number of cores used (the computation time is estimated using a quad core processor). All the sample have a radius of curvature of 30 nm, a length of 100 nm and the meshes at the apex have an area of 2 N m$^2$. The computation is done considering if the meshes area is constant (i.e. $\alpha = 0$ and $\beta = 0$) or evolving according to the equation (7) (which defines the coefficients $\alpha$ and $\beta$).</caption><thead><tr><th></th><th>Single core</th><th>All cores</th></tr></thead><tbody><tr><td>$\alpha = 0$ and $\beta = 0$</td><td>400</td><td>139</td></tr><tr><td>$\alpha = 10$ and $\beta = 2$</td><td>50</td><td>18</td></tr><tr><td>$\alpha = 50$ and $\beta = 2$</td><td>12</td><td>5</td></tr></tbody></table>

depth, using a polynomial law (as for the temperature (7) but this time with a depth dependency) (figure 7(a)). With this approach, a small mesh area at the apex is defined, so that an accurate estimation of the electrical field is achieved, where the evaporation process prevails, and a large mesh area is imposed on the shank, where the evaporation rate is very low, regarding the charge distribution (figure 7(b)). This approach significantly reduces the number of meshes and therefore the computation load.

The input of all these optimizations makes it possible to simulate the evaporation over a depth of 100 nm on a realistic description of an APT sample keeping a fine scale description of the sample apex (with a mesh area at the apex of 2 N m$^2$) with computation time at the hour scale, on a standard computer. The same result on atomistic models takes days [8, 15] or weeks [12, 13].

## 2.7. Comparison with experimental results

The advantage of atomic models [8, 12, 15] is the possibility to directly compare the results with the experimental ones. To perform such comparison with our model, we must transform meshes into 'atoms'. Between two iterations, we associate to each mesh an evaporated volume ($V_{ev}$), which corresponds then to a number of evaporated atoms $N_{ev}$ ($V_{ev} = N_{ev}V_{at}/Q$ with $V_{at}$ the atomic volume and $Q$ the detection efficiency). Then, we randomly place this number of atoms $N_{ev}$ in the regions where the meshes impact the detector. Some examples of this transformation will be shown in the next section.

# 3. Results

## 3.1. Single phase

Figure 8 represents results of the simulation of FE of a single-phase emitter. The initial morphology (before evaporation) was defined with an initial apex radius of curvature of 20 nm, a length of 100 nm and a shank angle of 10°. The mesh elements area at the apex were set to 4 N m$^2$ with the evolving mesh method (see section 2.6). From this model, 20 nm has been evaporated to test the simulation. More than considering
which corroborates this. Obviously, increasing the length of the tip, also increases the number of meshes and therefore the time computing. In order to maintain a reasonable computation time (sample length around hundreds of nanometers) and to have an ICF equivalent to the experimental one, we performed a post analytical compression of the trajectories, using an ad hoc coefficient $\lambda$. After a certain distance from the sample (in the range of few $\mu$m corresponding to the distance to an electrode and where the trajectories are straight lines), the trajectory calculation is stopped, and the position on the detector is obtained by a projection, using the velocity at this step, considering a compression of its polar angle, by $\theta = \lambda\theta'$ (with $\theta$ and $\theta'$ respectively the polar angle of the velocity where the calculation is stopped and on the detector), while its azimuthal angle remains unchanged $\varphi = \varphi'$ (with $\varphi$ and $\varphi'$ respectively the azimuthal angle of the velocity at the same step as polar angle), as experimentally observed [1]. Using this approach with a coefficient $\lambda$ equal to 1.5 (red histogram on the figure 6(b)), the average ICF (i.e. 1.53) is then comparable to experiment value. This coefficient is adjustable by the user in order to best adjust to the real analysis.

## 2.6. Algorithmic optimization

One of the crucial points of any simulation model for practical application are computing time and memory consumption. In this section, we show some solutions added to the model to address this. We are not benchmarking the performance of them, as done in [51], but just give a brief overview of the improvements, in term of computing time (table 1).

As it was mentioned previously, the electrical field is estimated by solving the Robin's equation with an $O(N^2)$ complexity (with $N$ the number of meshes). Thus, the first introduced optimization, is a parallelized computation (table 1). Observe that in equation (9), step $n+1$ is uniquely dependent on fixed values at step $n$, with demand for an exchange of pieces of information for each point of charge, parallelization is strongly efficient with a scaling close to the number of cores (Using C# parallel loop). The second improvement, as Rolland et al [15] at atomic scale, is a long range interaction approximation, developed by Barnes and Hut [18], in order to drastically reduce the computing time required for solving Robin's equation solving by reducing the complexity to $O(\text{Log}(N))$. The time computing is directly linked to the number of meshes $N$, so we also introduced logically a reduction of this number (table 1). Regarding the electrical field distribution (figure 5), the evaporation will mainly occur at the sample apex. During the whole process of FE, the model consists in moving the meshes as a function of the evaporation rate (the latter of which is field dependent). Obviously, a refined mesh at the bottom of sample is unnecessary as the electric field on these regions is much lower than the electric field on the apex. An average charge on these meshes is also sufficient to describe the long-range interaction for ion trajectory calculations. As a result, we introduced an evolving mesh area that depends on the sample

<table><caption>**Table 1.** Time computing (second) to perform ten field evaporation iterations ($\approx 0.3-0.6$ nm depth erosion in this case) as a function of the number of cores used (the computation time is estimated using a quad core processor). All the sample have a radius of curvature of 30 nm, a length of 100 nm and the meshes at the apex have an area of 2 N m$^2$. The computation is done considering if the meshes area is constant (i.e. $\alpha = 0$ and $\beta = 0$) or evolving according to the equation (7) (which defines the coefficients $\alpha$ and $\beta$).</caption><thead><tr><th></th><th>Single core</th><th>All cores</th></tr></thead><tbody><tr><td>$\alpha = 0$ and $\beta = 0$</td><td>400</td><td>139</td></tr><tr><td>$\alpha = 10$ and $\beta = 2$</td><td>50</td><td>18</td></tr><tr><td>$\alpha = 50$ and $\beta = 2$</td><td>12</td><td>5</td></tr></tbody></table>

depth, using a polynomial law (as for the temperature (7) but this time with a depth dependency) (figure 7(a)). With this approach, a small mesh area at the apex is defined, so that an accurate estimation of the electrical field is achieved, where the evaporation process prevails, and a large mesh area is imposed on the shank, where the evaporation rate is very low, regarding the charge distribution (figure 7(b)). This approach significantly reduces the number of meshes and therefore the computation load.

The input of all these optimizations makes it possible to simulate the evaporation over a depth of 100 nm on a realistic description of an APT sample keeping a fine scale description of the sample apex (with a mesh area at the apex of 2 N m$^2$) with computation time at the hour scale, on a standard computer. The same result on atomistic models takes days [8, 15] or weeks [12, 13].

## 2.7. Comparison with experimental results

The advantage of atomic models [8, 12, 15] is the possibility to directly compare the results with the experimental ones. To perform such comparison with our model, we must transform meshes into 'atoms'. Between two iterations, we associate to each mesh an evaporated volume ($V_{ev}$), which corresponds then to a number of evaporated atoms $N_{ev}$ ($V_{ev} = N_{ev}V_{at}/Q$ with $V_{at}$ the atomic volume and $Q$ the detection efficiency). Then, we randomly place this number of atoms $N_{ev}$ in the regions where the meshes impact the detector. Some examples of this transformation will be shown in the next section.

# 3. Results

## 3.1. Single phase

Figure 8 represents results of the simulation of FE of a single-phase emitter. The initial morphology (before evaporation) was defined with an initial apex radius of curvature of 20 nm, a length of 100 nm and a shank angle of 10°. The mesh elements area at the apex were set to 4 N m$^2$ with the evolving mesh method (see section 2.6). From this model, 20 nm has been evaporated to test the simulation. More than considering
<figure><img src="image_13.png" /><figcaption>**Figure 7.** (a) Triangulated APT sample ($R_c$ = 20 nm, $L$ = 100 nm, $g$ = 10° and with evolving meshing $S_{at,0}$ = 4 N m²) and (b) the associated charge density distribution ($q/S_{at}$ in a.u.) at the emitter surface.</figcaption></figure>

only a constant temperature (figure 8(a), $T$ = 50 K), also non-homogeneous temperature has been tested to check the ability to reproduce laser assisted FE (figure 8(b)). A gradient of temperature was imposed along the X axis of the tip following a polynomial law (expression (7)) with $\beta$ = 2). The peak temperature on the side of the tip was set to 200 K ($\alpha$ = 3).

At a constant temperature ($T$ = 50, figure 8(a)), model predictions agree with those from classical experimental observations. During the evaporation process, the tip apex maintains its quasi-hemispherical shape, with an increase of the radius of curvature (initially 20 nm and at the end $27 \pm 2$ nm, figure 8(c)), due to a non-null shank angle as a function of depth erosion. The meshes, that impact the detector ($R$ = 4 cm, at 11 cm from the tip, with $\lambda$ equal to 1.6), are coming from a small spherical region from the apex ($R$ = 10 nm) and a number of impacts distribute centrally around it.

The introduction of a gradient of temperature on the apex temperature, in order to reproduce the laser assisted evaporation, modifies the evaporation sequence (figure 8(b)) and thus reproduces some experimental observations. The evaporation rate is obviously locally increased by the local high temperature (expression (2)), thus the tip shape gradually becomes asymmetric, as observed experimentally by [52]. The impact distribution on the detector is no more centered (highest density of impact is shifted to the direction of the laser side). By back projecting the detector onto the tip surface, it is interesting to note that the meshes impacting the detector are located on the tip on a deformed elliptical region shifted to the shadow side of the sample. This kind of image projection deformations are experimentally ignored. Considering the dimension of this last region ($X \in [-15,5$ nm]), the variation of radius of curvature on either side of the sample, which therefore impacts the detector, is equal to 5 nm (from 22.5 to 27.5 nm,

figure 8(c)) and the electric field variation is equal to $5$ V nm$^{-1}$ (from 35 to 30 V nm$^{-1}$, figure 8(d)). In experiment, this surface electric field variation generates variations in the post ionization of the emitted ions. In turn, this results in an evolution of the charge state ratio from the laser to the shadow side [53–55]. It is important to note that this asymmetry (and the effects resulting from it) become even more pronounced when the temperature gradient increases (figures 8(c) and (d)).

We may note that the empirical relationship $F \sim V/R$ is not verified through this simple simulation which may be surprising. However, this empirical expression is generally given for a sample of cylindrical symmetry, which is not the case here. Our model could be then used to understand composition measurement fluctuations, induced by the laser assisted evaporation [56] or even to estimate the temperature rise at the sample apex induced by its interaction with the laser [57], introducing an accurate evolution [41, 47].

## 3.2. Multiphase FE: planar layers

The first layer case that we studied with our model is a horizontal layer (figure 9(a)), which is also the most studied case with FE simulations [30–32, 58–63]. A key advantage of our model is its ability to model samples with dimensions comparable with experimental ones. Thus, we simulate FE for a sample with a radius of curvature of 50 nm, a length of 200 nm and a shank angle of 5°, on which we introduce a layer with a thickness of 30 nm. The mesh elements area is set to 4 N m² at the apex, using the evolving mesh method along the shank. The simulations assume constant sample temperature of 50 K. Two evaporation field values for the layer ($F_{Ev,B}$) have been studied: 0.7 and 1.4 times the matrix one ($F_{Ev,A}$). We will focus on describing the first case ($F_{Ev,B}$ = 0.7 × $F_{Ev,A}$) as the second
<figure><img src="image_13.png" /><figcaption>**Figure 7.** (a) Triangulated APT sample ($R_c$ = 20 nm, $L$ = 100 nm, $g$ = 10° and with evolving meshing $S_{at,0}$ = 4 N m²) and (b) the associated charge density distribution ($q/S_{at}$ in a.u.) at the emitter surface.</figcaption></figure>

only a constant temperature (figure 8(a), $T$ = 50 K), also non-homogeneous temperature has been tested to check the ability to reproduce laser assisted FE (figure 8(b)). A gradient of temperature was imposed along the X axis of the tip following a polynomial law (expression (7)) with $\beta$ = 2). The peak temperature on the side of the tip was set to 200 K ($\alpha$ = 3).

At a constant temperature ($T$ = 50, figure 8(a)), model predictions agree with those from classical experimental observations. During the evaporation process, the tip apex maintains its quasi-hemispherical shape, with an increase of the radius of curvature (initially 20 nm and at the end $27 \pm 2$ nm, figure 8(c)), due to a non-null shank angle as a function of depth erosion. The meshes, that impact the detector ($R$ = 4 cm, at 11 cm from the tip, with $\lambda$ equal to 1.6), are coming from a small spherical region from the apex ($R$ = 10 nm) and a number of impacts distribute centrally around it.

The introduction of a gradient of temperature on the apex temperature, in order to reproduce the laser assisted evaporation, modifies the evaporation sequence (figure 8(b)) and thus reproduces some experimental observations. The evaporation rate is obviously locally increased by the local high temperature (expression (2)), thus the tip shape gradually becomes asymmetric, as observed experimentally by [52]. The impact distribution on the detector is no more centered (highest density of impact is shifted to the direction of the laser side). By back projecting the detector onto the tip surface, it is interesting to note that the meshes impacting the detector are located on the tip on a deformed elliptical region shifted to the shadow side of the sample. This kind of image projection deformations are experimentally ignored. Considering the dimension of this last region ($X \in [-15,5$ nm]), the variation of radius of curvature on either side of the sample, which therefore impacts the detector, is equal to 5 nm (from 22.5 to 27.5 nm,

figure 8(c)) and the electric field variation is equal to $5$ V nm$^{-1}$ (from 35 to 30 V nm$^{-1}$, figure 8(d)). In experiment, this surface electric field variation generates variations in the post ionization of the emitted ions. In turn, this results in an evolution of the charge state ratio from the laser to the shadow side [53–55]. It is important to note that this asymmetry (and the effects resulting from it) become even more pronounced when the temperature gradient increases (figures 8(c) and (d)).

We may note that the empirical relationship $F \sim V/R$ is not verified through this simple simulation which may be surprising. However, this empirical expression is generally given for a sample of cylindrical symmetry, which is not the case here. Our model could be then used to understand composition measurement fluctuations, induced by the laser assisted evaporation [56] or even to estimate the temperature rise at the sample apex induced by its interaction with the laser [57], introducing an accurate evolution [41, 47].

## 3.2. Multiphase FE: planar layers

The first layer case that we studied with our model is a horizontal layer (figure 9(a)), which is also the most studied case with FE simulations [30–32, 58–63]. A key advantage of our model is its ability to model samples with dimensions comparable with experimental ones. Thus, we simulate FE for a sample with a radius of curvature of 50 nm, a length of 200 nm and a shank angle of 5°, on which we introduce a layer with a thickness of 30 nm. The mesh elements area is set to 4 N m² at the apex, using the evolving mesh method along the shank. The simulations assume constant sample temperature of 50 K. Two evaporation field values for the layer ($F_{Ev,B}$) have been studied: 0.7 and 1.4 times the matrix one ($F_{Ev,A}$). We will focus on describing the first case ($F_{Ev,B}$ = 0.7 × $F_{Ev,A}$) as the second
<figure><img src="image_14.png" /><figcaption>**Figure 8.** Cross section of the sample during evaporation (charge density in color, a.u.), the impact density on the detector (number of impact per grid unit area) and the meshes that impact the detector, according to the applied temperature (*T*): (a) *T* is constant and equal to 50 K and (b) evolving temperature modeled with a polynomial law, corresponding to an evolution from 50 to 200 K along the *X* axis (expression (6) with α = 3 and β = 2). (c) Radius of curvature (nm) and (d) electric field (V nm⁻¹) along the *X* axis, as function of the applied temperature: constant (*T* = 50 K), evolving with β = 2 and α = 2 (*T* = 50–150 K) or α = 3 (*T* = 50–200 K).</figcaption></figure>

case is similar, except for the observation that the evaporation sequence is reversed, as we will be discussed later.

Regarding the evolution of the field at the surface during the evaporation (figure 9(b)), we define two regions of interest (ROIs): the interface A/B (figure 9(c), ROI 1) and the interface B/A (figure 9(d), ROI 2). During the evaporation of the interface A/B, we observe an increase of the field on the A phase, due to its protruding position. This expands the magnification locally and decreases the field of view (figure 10(a)) since only

a limited region of the tip apex is collected on the finite size detector, (here detector radius is *R* = 4 cm at 11 cm from the tip, with λ equal to 1.6). At this step, we also observe a slight increase of the electric field at the interface B/A, induced by a local decrease of the radius of curvature, because the evaporation rate of A phase is lower than the B phase. During the evaporation of interface B/A, the reverse phenomenon is observed. The B phase becomes flatter (high radius and low field), due to its higher evaporation rate compared to the phase A below,
<figure><img src="image_14.png" /><figcaption>**Figure 8.** Cross section of the sample during evaporation (charge density in color, a.u.), the impact density on the detector (number of impact per grid unit area) and the meshes that impact the detector, according to the applied temperature (*T*): (a) *T* is constant and equal to 50 K and (b) evolving temperature modeled with a polynomial law, corresponding to an evolution from 50 to 200 K along the *X* axis (expression (6) with *α* = 3 and *β* = 2). (c) Radius of curvature (nm) and (d) electric field (V nm⁻¹) along the *X* axis, as function of the applied temperature: constant (*T* = 50 K), evolving with *β* = 2 and *α* = 2 (*T* = 50–150 K) or *α* = 3 (*T* = 50–200 K).</figcaption></figure>

case is similar, except for the observation that the evaporation sequence is reversed, as we will be discussed later.

Regarding the evolution of the field at the surface during the evaporation (figure 9(b)), we define two regions of interest (ROIs): the interface A/B (figure 9(c), ROI 1) and the interface B/A (figure 9(d), ROI 2). During the evaporation of the interface A/B, we observe an increase of the field on the A phase, due to its protruding position. This expands the magnification locally and decreases the field of view (figure 10(a)) since only

a limited region of the tip apex is collected on the finite size detector, (here detector radius is *R* = 4 cm at 11 cm from the tip, with λ equal to 1.6). At this step, we also observe a slight increase of the electric field at the interface B/A, induced by a local decrease of the radius of curvature, because the evaporation rate of A phase is lower than the B phase. During the evaporation of interface B/A, the reverse phenomenon is observed. The B phase becomes flatter (high radius and low field), due to its higher evaporation rate compared to the phase A below,
<figure><img src="image_15.png" /><figcaption>**Figure 9.** (a) Initial sample morphology. (b) Evolution of the charge density (according to the phase: A or B, in a.u.), the simulated potential (V), the one if there is no B phase (reported as reference), as the function of the number of evaporated atoms. Morphology, radius of curvature (nm), and electric field (V nm⁻¹) at two different steps defined on the (b): (c) interface A/B defined as ROI 1 (region of interest) and (d) interface B/A defined as ROI 2.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Figure 10.** Perfect 3D reconstruction if the evaporation field of the layer is equal to (a) 0.7 or (b) 1.4 times the matrix one. Reconstructions were performed by extracting tip shape every 3–4 nm.</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>**Figure 9.** (a) Initial sample morphology. (b) Evolution of the charge density (according to the phase: A or B, in a.u.), the simulated potential (V), the one if there is no B phase (reported as reference), as the function of the number of evaporated atoms. Morphology, radius of curvature (nm), and electric field (V nm⁻¹) at two different steps defined on the (b): (c) interface A/B defined as ROI 1 (region of interest) and (d) interface B/A defined as ROI 2.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Figure 10.** Perfect 3D reconstruction if the evaporation field of the layer is equal to (a) 0.7 or (b) 1.4 times the matrix one. Reconstructions were performed by extracting tip shape every 3–4 nm.</figcaption></figure>
inducing then a drop of the magnification and an increase of the field of view (figure 10(a)). The radius of curvature of phase A decreases regarding its high evaporation field.

For the second case (i.e. $F_{Ev,B} = 1.4 \times F_{Ev,A}$), the situation is the exact opposite. The evaporation of the interface A/B corresponds to the evaporation of the previous B/A, and reciprocally for the B/A interface. We obtain then an increase of the field of view at the A/B interface and the opposite at the B/A interface (figure 10(b)).

Note that these results and observations agree with those which have already been obtained in the literature, but they have been obtained (i) on a real size tip and this (ii) with a relatively short computation time (2 h for each evaporation field).

### 3.3. Nanoparticles

The FE process of second-phase particles has been widely studied with atomic modeling [8, 10, 13, 64, 65], for almost 20 years now [7]. APT reconstruction particles are known in the literature to be strongly biased by artifacts named local magnification effects and trajectory overlaps [8, 66, 67]. Models enabled to obtain better quantitative understanding of these phenomena at least in idealized geometries. However, the cumbersome nature of atomic scale simulations often limited the result to a small range of particle size, or FE properties.

Local magnification occurs when different phases with individually different evaporation properties and fields are evaporated at the same time. Those phases (for example A and B) will develop locally different radii of curvature, (respectively $R_{C,A}$ and $R_{C,B}$), governed by their individual evaporation field and evaporation properties, which leads then to different local magnification (respectively $M_A$ and $M_B$). In a simple approximation, the following expression correlates their magnification, according to the local radius of curvature:

$$M_B = M_A \times \frac{R_{C,A}}{R_{C,B}}. \qquad (14)$$

This leads to trajectory aberrations. These are responsible for morphological distortions and in the worst case cause trajectory overlaps (when trajectories from different phases impact the detector at the same place), inducing then chemical biases. Those aberrations lead then to erroneous reconstructed atomic densities as was demonstrated analytically in [5, 64].

We have simulated the FE of a spherical particle, with a radius of 4 nm that is embedded in a sample with a radius of curvature of 20 nm, a length of 100 nm and a shank angle of 10° (figure 1). The meshes area is equal to 1 N m² at the apex, with an evolution along the shank. The simulations operate at a constant temperature equal to 50 K. Different evaporation fields of the particles ($F_{Ev,B}$) have been studied, from 0.4 to 1.4 times the matrix one ($F_{Ev,A}$) (by step of 0.1). After a transformation of meshes into atoms, associating an elemental nature corresponding to the phase from which they originate, a 3D reconstruction is performed, from their positions onto the detector, using standard reconstruction protocol [68, 69]. Concentration profiles were computed: one cutting through the particle along the evaporation direction and another one pointing in a direction perpendicular to the latter. Whatever is the

evaporation field of the particle, the size of the particle (Z), along the evaporation direction is found almost constant and equal to the expected one (i.e. a radius of 4 nm). This observation matches previous observations from atomic modeling [13, 67]. For each profile in the perpendicular direction, we calculated the reduced density profile ($\rho$, atomic density normalized by the matrix one), a measurement of the average core concentration in atoms of the particle, and a calculation of the shape factor (S, square of the dimension along the perpendicular direction (i.e. X) divided by the Z dimension, i.e. $(X/Z)^2$). The shape factor is a descriptor for the morphological distortion [5]. All those measurements are reported in the figure 11.

As expected, no local magnification, no concentration bias, and no morphological distortion occur when the two phases have the same evaporation fields ($F_{Ev,B} = F_{Ev,A}$, $X_B = 100$ at.%, $\rho = 1$, $S = 1$) (figure 9). Otherwise, as the results show, we can distinguish three different principal cases, according to the evaporation field of the particle. The origins and the interpretations of these three behaviors will be based on the simulation results reported in the figure 12.

$$\bullet F_{Ev,B} > F_{Ev,A}.$$

This case is referred in the literature as a high field particle case, the particle has a higher evaporation field than its surrounding matrix. Those particles exhibit no chemical bias in the core ($X_B = 100$ at.%), but a shape dilatation ($S > 1$) surplus atomic density fluctuation, with an higher one specifically located at the interface ($\rho > 1$) and a lower in the core ($\rho < 1$) (figure 11). Those observations have already been made experimentally for MgO particles in a Au matrix [70], nitride dispersion strengthened steel particles [71] or even Al₃(Zr₁-xTiₓ) particles in an α-Al matrix [72]. When this kind of particle merges at the surface, its evaporation rate is lower than the matrix. Consequently, the particle is protruding (with a small radius of curvature and so a high electric field compare to the matrix) and the matrix at the interface is exhibiting a very high radius of curvature. This interface can also have negative radii of curvature, since the concavity can be modified, when the evaporation field difference is very important. This surface topography induces divergent trajectories from the particle, as well as of the matrix at the interface, but to a lower extent. In effect, a size expansion of the particle is observed due to trajectory aberrations of the core (then a low atomic density of impact on the detector) and a trajectory overlap between atoms located at the interface (matrix and particle with then a high atomic density).

$$\bullet F_{Ev,A} > F_{Ev,B} \ge 0.8 \times F_{Ev,A}.$$

This case is referred to in the literature as the low field particle case. Those particles exhibit no chemical bias in the core ($X_B = 100$ at.%), but a shape compression ($S < 1$) and a high atomic density ($\rho > 1$) (figure 11). Those observations have already been done experimentally for an α' particle in FeCr [73, 74]. Their evaporation rate is higher than the matrix. In effect, the matrix at the interface is then becoming protruding and the particle develops a convex interface (negative radii)
inducing then a drop of the magnification and an increase of the field of view (figure 10(a)). The radius of curvature of phase A decreases regarding its high evaporation field.

For the second case (i.e. $F_{Ev,B} = 1.4 \times F_{Ev,A}$), the situation is the exact opposite. The evaporation of the interface A/B corresponds to the evaporation of the previous B/A, and reciprocally for the B/A interface. We obtain then an increase of the field of view at the A/B interface and the opposite at the B/A interface (figure 10(b)).

Note that these results and observations agree with those which have already been obtained in the literature, but they have been obtained (i) on a real size tip and this (ii) with a relatively short computation time (2 h for each evaporation field).

### 3.3. Nanoparticles

The FE process of second-phase particles has been widely studied with atomic modeling [8, 10, 13, 64, 65], for almost 20 years now [7]. APT reconstruction particles are known in the literature to be strongly biased by artifacts named local magnification effects and trajectory overlaps [8, 66, 67]. Models enabled to obtain better quantitative understanding of these phenomena at least in idealized geometries. However, the cumbersome nature of atomic scale simulations often limited the result to a small range of particle size, or FE properties.

Local magnification occurs when different phases with individually different evaporation properties and fields are evaporated at the same time. Those phases (for example A and B) will develop locally different radii of curvature, (respectively $R_{C,A}$ and $R_{C,B}$), governed by their individual evaporation field and evaporation properties, which leads then to different local magnification (respectively $M_A$ and $M_B$). In a simple approximation, the following expression correlates their magnification, according to the local radius of curvature:

$$M_B = M_A \times \frac{R_{C,A}}{R_{C,B}}. \qquad (14)$$

This leads to trajectory aberrations. These are responsible for morphological distortions and in the worst case cause trajectory overlaps (when trajectories from different phases impact the detector at the same place), inducing then chemical biases. Those aberrations lead then to erroneous reconstructed atomic densities as was demonstrated analytically in [5, 64].

We have simulated the FE of a spherical particle, with a radius of 4 nm that is embedded in a sample with a radius of curvature of 20 nm, a length of 100 nm and a shank angle of 10° (figure 1). The meshes area is equal to 1 N m² at the apex, with an evolution along the shank. The simulations operate at a constant temperature equal to 50 K. Different evaporation fields of the particles ($F_{Ev,B}$) have been studied, from 0.4 to 1.4 times the matrix one ($F_{Ev,A}$) (by step of 0.1). After a transformation of meshes into atoms, associating an elemental nature corresponding to the phase from which they originate, a 3D reconstruction is performed, from their positions onto the detector, using standard reconstruction protocol [68, 69]. Concentration profiles were computed: one cutting through the particle along the evaporation direction and another one pointing in a direction perpendicular to the latter. Whatever is the

evaporation field of the particle, the size of the particle (Z), along the evaporation direction is found almost constant and equal to the expected one (i.e. a radius of 4 nm). This observation matches previous observations from atomic modeling [13, 67]. For each profile in the perpendicular direction, we calculated the reduced density profile ($\rho$, atomic density normalized by the matrix one), a measurement of the average core concentration in atoms of the particle, and a calculation of the shape factor (S, square of the dimension along the perpendicular direction (i.e. X) divided by the Z dimension, i.e. $(X/Z)^2$). The shape factor is a descriptor for the morphological distortion [5]. All those measurements are reported in the figure 11.

As expected, no local magnification, no concentration bias, and no morphological distortion occur when the two phases have the same evaporation fields ($F_{Ev,B} = F_{Ev,A}$, $X_B = 100$ at.%, $\rho = 1$, $S = 1$) (figure 9). Otherwise, as the results show, we can distinguish three different principal cases, according to the evaporation field of the particle. The origins and the interpretations of these three behaviors will be based on the simulation results reported in the figure 12.

$$\bullet F_{Ev,B} > F_{Ev,A}.$$

This case is referred in the literature as a high field particle case, the particle has a higher evaporation field than its surrounding matrix. Those particles exhibit no chemical bias in the core ($X_B = 100$ at.%), but a shape dilatation ($S > 1$) surplus atomic density fluctuation, with an higher one specifically located at the interface ($\rho > 1$) and a lower in the core ($\rho < 1$) (figure 11). Those observations have already been made experimentally for MgO particles in a Au matrix [70], nitride dispersion strengthened steel particles [71] or even Al₃(Zr₁-xTiₓ) particles in an α-Al matrix [72]. When this kind of particle merges at the surface, its evaporation rate is lower than the matrix. Consequently, the particle is protruding (with a small radius of curvature and so a high electric field compare to the matrix) and the matrix at the interface is exhibiting a very high radius of curvature. This interface can also have negative radii of curvature, since the concavity can be modified, when the evaporation field difference is very important. This surface topography induces divergent trajectories from the particle, as well as of the matrix at the interface, but to a lower extent. In effect, a size expansion of the particle is observed due to trajectory aberrations of the core (then a low atomic density of impact on the detector) and a trajectory overlap between atoms located at the interface (matrix and particle with then a high atomic density).

$$\bullet F_{Ev,A} > F_{Ev,B} \geq 0.8 \times F_{Ev,A}.$$

This case is referred to in the literature as the low field particle case. Those particles exhibit no chemical bias in the core ($X_B = 100$ at.%), but a shape compression ($S < 1$) and a high atomic density ($\rho > 1$) (figure 11). Those observations have already been done experimentally for an α' particle in FeCr [73, 74]. Their evaporation rate is higher than the matrix. In effect, the matrix at the interface is then becoming protruding and the particle develops a convex interface (negative radii)
<figure><img src="image_17.png" /><figcaption>Figure 11. Reduced atomic density, measured through the particle (with an initial radius of 4 nm) in the direction perpendicular to the evaporation direction, according the distance from the particle core (nm), average core concentration in particle atoms (at.%) and the shape factor, as a function of the evaporation field ratio of the particle with the surrounding matrix (from 0.4 to 1.4 by 0.1). The dash lines in the reduced atomic density map represents the interface between the particle and the matrix, define at half height of the concentration in the particle core and the one in the matrix.</figcaption></figure>

and a ‘flat’ core (high radii), thus a low electric field. This topography induces on the one hand convergent trajectories from the particle and on the other hand divergent trajectories from the matrix at the interface. This divergence is not enough to overlap with the trajectories for ions originating from the particle (so no chemical bias [75]). In this simulated case, it should be noted that the difference in topography within the particle, between the core and the interface, causes a variation in the atomic density, both exhibiting an over atomic density, but a higher one at the interface (especially visible for $F_{Ev,B} = 0.9 \times F_{Ev,A}$). In both cases, the observed high atomic density is thus only induced by the morphological compression.

$$ \bullet F_{Ev,B} < 0.8 \times F_{Ev,A}. $$

This latter case is not very widespread in the literature [1], since it is difficult to probe in experiment. Indeed, there is a chemical composition bias ($X_B < 100$ at.%), an over atomic density ($\rho > 1$) and a morphological variation depending of the evaporation field: compression if $0.5 \times F_{Ev,A} < F_{Ev,B} < 0.8 \times F_{Ev,A}$ ($S < 1$) and a dilation otherwise ($S > 1$). The common point in this evaporation field range is the particle trajectories crossing-over due to an intense convergence, caused by an entire convex particle (negative radii). Up to now, this phenomena has been observed experimentally in oxide dispersion strengthened (ODS) steel particles by [64] and in SiGe layer embedded in SiO$_2$ [76]. In all these cases, the matrix interface divergence (induced by it protruding) is enough to overlap with particle trajectories (resulting in chemical bias [5]). If $F_{Ev,B} > 0.5 \times F_{Ev,A}$, the cross-over occurs far from the tip. Consequently, the defocusing which is induced by the crossover is insufficient to compensate for the initial

focusing. This results in a compressed particle ($S < 1$). This is the case of Cu particles in Fe [77]. As the FE decreases, the cross-over occurs closer to the tip leading to an apparent divergence of the trajectories ($S < 1$) which increases the overlap (so the chemical bias).

Our model allows to reproduce not only experimental observations but also results that were obtained by atomic simulation [1]. With our model, we observe that ion cross-over occurs for particle evaporation field lower than $0.7-0.8 \times F_{Ev,A}$ (figure 11). With an atomic model, Vurpillot observes that this phenomena occurs for particle evaporation field values lower than $0.5-0.7 \times F_{Ev,A}$ [1]. It is obvious that the particle size ratio with the tip size [5, 64] is correlated with the evaporation field at which this crossing over occurs. In this case, our model presents the following advantages: (i) The model reproduces experimental and atomic scale simulation results in reasonable computing hour (less than one hour per case) and (ii) to simulate real sample dimension (with a radius of 50 nm), also in a reasonable time computing but mainly without any memory allocation issue (200 MB compared to 30 GB in atomic model [12]). The investigation of the size effect (of the particles and the tip) will be the subject of a future study.

### 3.4. Vertical layer

Our last studied case is a setup with a vertical layer studied for laser assisted evaporation. Our objective is to see if the model is capable of reproducing observations made by Lee et al [78] as well to confirm their interpretations if possible. They made TEM observations of a vertical layer Si embedded in SiO$_2$ (figure 13(a)) and its opposite (SiO$_2$ embedded in Si, figure 14(a)), after APT analysis.
<figure><img src="image_17.png" /><figcaption>Figure 11. Reduced atomic density, measured through the particle (with an initial radius of 4 nm) in the direction perpendicular to the evaporation direction, according the distance from the particle core (nm), average core concentration in particle atoms (at.%) and the shape factor, as a function of the evaporation field ratio of the particle with the surrounding matrix (from 0.4 to 1.4 by 0.1). The dash lines in the reduced atomic density map represents the interface between the particle and the matrix, define at half height of the concentration in the particle core and the one in the matrix.</figcaption></figure>

and a ‘flat’ core (high radii), thus a low electric field. This topography induces on the one hand convergent trajectories from the particle and on the other hand divergent trajectories from the matrix at the interface. This divergence is not enough to overlap with the trajectories for ions originating from the particle (so no chemical bias [75]). In this simulated case, it should be noted that the difference in topography within the particle, between the core and the interface, causes a variation in the atomic density, both exhibiting an over atomic density, but a higher one at the interface (especially visible for $F_{Ev,B} = 0.9 \times F_{Ev,A}$). In both cases, the observed high atomic density is thus only induced by the morphological compression.

$$ \bullet F_{Ev,B} < 0.8 \times F_{Ev,A}. $$

This latter case is not very widespread in the literature [1], since it is difficult to probe in experiment. Indeed, there is a chemical composition bias ($X_B < 100$ at.%), an over atomic density ($\rho > 1$) and a morphological variation depending of the evaporation field: compression if $0.5 \times F_{Ev,A} < F_{Ev,B} < 0.8 \times F_{Ev,A}$ ($S < 1$) and a dilation otherwise ($S > 1$). The common point in this evaporation field range is the particle trajectories crossing-over due to an intense convergence, caused by an entire convex particle (negative radii). Up to now, this phenomena has been observed experimentally in oxide dispersion strengthened (ODS) steel particles by [64] and in SiGe layer embedded in SiO$_2$ [76]. In all these cases, the matrix interface divergence (induced by it protruding) is enough to overlap with particle trajectories (resulting in chemical bias [5]). If $F_{Ev,B} > 0.5 \times F_{Ev,A}$, the cross-over occurs far from the tip. Consequently, the defocusing which is induced by the crossover is insufficient to compensate for the initial

focusing. This results in a compressed particle ($S < 1$). This is the case of Cu particles in Fe [77]. As the FE decreases, the cross-over occurs closer to the tip leading to an apparent divergence of the trajectories ($S < 1$) which increases the overlap (so the chemical bias).

Our model allows to reproduce not only experimental observations but also results that were obtained by atomic simulation [1]. With our model, we observe that ion cross-over occurs for particle evaporation field lower than $0.7-0.8 \times F_{Ev,A}$ (figure 11). With an atomic model, Vurpillot observes that this phenomena occurs for particle evaporation field values lower than $0.5-0.7 \times F_{Ev,A}$ [1]. It is obvious that the particle size ratio with the tip size [5, 64] is correlated with the evaporation field at which this crossing over occurs. In this case, our model presents the following advantages: (i) The model reproduces experimental and atomic scale simulation results in reasonable computing hour (less than one hour per case) and (ii) to simulate real sample dimension (with a radius of 50 nm), also in a reasonable time computing but mainly without any memory allocation issue (200 MB compared to 30 GB in atomic model [12]). The investigation of the size effect (of the particles and the tip) will be the subject of a future study.

### 3.4. Vertical layer

Our last studied case is a setup with a vertical layer studied for laser assisted evaporation. Our objective is to see if the model is capable of reproducing observations made by Lee et al [78] as well to confirm their interpretations if possible. They made TEM observations of a vertical layer Si embedded in SiO$_2$ (figure 13(a)) and its opposite (SiO$_2$ embedded in Si, figure 14(a)), after APT analysis.
<figure><img src="image_18.png" /><figcaption>**Figure 12.** (from top to bottom) Electric field (V·nm⁻¹), radius of curvature (nm) at the sample surface (in top view), the impact density on the detector (number of impact per grid unit area), the trajectories from the particle and the 3D reconstruction (V = 10 × 2 × 17 nm³, atoms of matrix in blue and those of the particles in orange), according to its evaporation field F<sub>Ev,B</sub>: (a) F<sub>Ev,B</sub> > F<sub>Ev,A</sub>, (b) F<sub>Ev,A</sub> > F<sub>Ev,B</sub> ≥ 0.8 × F<sub>Ev,A</sub> and (c) F<sub>Ev,B</sub> < 0.8 × F<sub>Ev,A</sub> with F<sub>Ev,A</sub> the evaporation field of the matrix.</figcaption></figure>

In the first case, they observe that the SiO₂ layer radius of
curvature on either side is not the same (from 25 nm in the
shadow side to 35 nm in the illuminated side). The authors

attribute this to the laser assisted evaporation, considering that
the laser side develops higher radii [60]. To test this, we simu-
late a vertical layer of Si (thickness of 20 and a sample radius
<figure><img src="image_18.png" /><figcaption>**Figure 12.** (from top to bottom) Electric field (V·nm⁻¹), radius of curvature (nm) at the sample surface (in top view), the impact density on the detector (number of impact per grid unit area), the trajectories from the particle and the 3D reconstruction (V = 10 × 2 × 17 nm³, atoms of matrix in blue and those of the particles in orange), according to its evaporation field F<sub>Ev,B</sub>: (a) F<sub>Ev,B</sub> > F<sub>Ev,A</sub>, (b) F<sub>Ev,A</sub> > F<sub>Ev,B</sub> ≥ 0.8 × F<sub>Ev,A</sub> and (c) F<sub>Ev,B</sub> < 0.8 × F<sub>Ev,A</sub> with F<sub>Ev,A</sub> the evaporation field of the matrix.</figcaption></figure>

In the first case, they observe that the SiO₂ layer radius of
curvature on either side is not the same (from 25 nm in the
shadow side to 35 nm in the illuminated side). The authors

attribute this to the laser assisted evaporation, considering that
the laser side develops higher radii [60]. To test this, we simu-
late a vertical layer of Si (thickness of 20 and a sample radius
<figure><img src="image_19.png" /><figcaption>**Figure 13.** (a) TEM observation of Si layer embedded in SiO2 [78] after APT analysis. Simulated sample morphology, considering constant (*T* = 50 K) or (b) a gradient of *T* (in *X* direction, from 50 to 200 with *α* = 2 and *β* = 3). (d) Radius of curvature at the sample surface (c) in top view. Average radius of curvature (nm) along the *X* axis in the region, noted in the figure (c) as (d) shadow side or as (e) laser side. (g) Electric field (V nm⁻¹) at the sample surface (c) in top view. (h) Average electric field profile along the *X* direction of the figure (g). On the figures (e), (f) and (h), it is also reported the results obtained for another temperature evolution from 50 to 150 K (with *α* = 2 and *β* = 3). Reprinted from [78], Copyright (2014), with permission from Elsevier.</figcaption></figure>

of 25 nm [78]) with an evaporation field of Si 0.7 times lower than the surrounding SiO₂. The simulations have been performed at constant temperature (T = 50 K, figure 12(a)) and with an evolutive temperature (from 50 to 200 K, with α = 2 and β = 3, figure 13(b)). Clearly, we note that the temperature simulation makes it possible to reproduce a surface morphology very similar to that of the experiment, which is not the case for constant temperature (same radii in both side for SiO₂ layers). We measure a variation of 11 nm in the radii of curvature on either side of the tip (at X = -20 nm, the radius is equal to 16 nm, figure 13(e), and at X = 20 nm, the radius is equal to 27 nm, figure 13(f)), thus confirming the conclusions of Lee *et al*. We also reproduce the difference of the evaporation depth for both SiO₂ layers (about 3–4 nm). Our results, however, are not quantitatively similar to the experimental ones. The variation of radii is quite the same (about 10 nm), but the simulated radii are lower (in the shadow side, 16 nm in simulation and 25 nm experimentally). In order to be

quantitative, it will undoubtedly be necessary to increase the size of the tip, even to modify the evaporation field of the Si layer, but more fundamentally to model more accurately the temperature distribution in the tip apex (figures 13(d)–(h)).

For the opposite case (SiO₂ in Si, figure 14), we make similar observations. Indeed, we note that the temperature simulation (figure 14(c)) makes it possible to reproduce a surface morphology very similar to that of the experiment (figure 14(a)), which is not the case for constant temperature (figure 14(b)). However, in this case, our results are far from being quantitatively comparable with experiment. They measure an increase of 10 nm of the SiO₂ layer radius (from 35 to 45 nm, figure 14(a)), whereas in simulation, we measure an increase of 4 nm, in the worst case (i.e. at the highest temperature in the figure 14(d), for X = -6 and X = 6 nm). Moreover, we observe a difference in the Si layer top depth (about 3–4 nm), as experimentally (figure 14(a)), but in the opposite direction (the uppermost layer is the laser side experimentally
<figure><img src="image_19.png" /><figcaption>**Figure 13.** (a) TEM observation of Si layer embedded in SiO2 [78] after APT analysis. Simulated sample morphology, considering constant ($T = 50$ K) or (b) a gradient of $T$ (in $X$ direction, from 50 to 200 with $\alpha = 2$ and $\beta = 3$). (d) Radius of curvature at the sample surface (c) in top view. Average radius of curvature (nm) along the $X$ axis in the region, noted in the figure (c) as (d) shadow side or as (e) laser side. (g) Electric field ($V \cdot nm^{-1}$) at the sample surface (c) in top view. (h) Average electric field profile along the $X$ direction of the figure (g). On the figures (e), (f) and (h), it is also reported the results obtained for another temperature evolution from 50 to 150 K (with $\alpha = 2$ and $\beta = 3$). Reprinted from [78], Copyright (2014), with permission from Elsevier.</figcaption></figure>

of 25 nm [78]) with an evaporation field of Si 0.7 times lower than the surrounding SiO₂. The simulations have been performed at constant temperature ($T = 50$ K, figure 12(a)) and with an evolutive temperature (from 50 to 200 K, with $\alpha = 2$ and $\beta = 3$, figure 13(b)). Clearly, we note that the temperature simulation makes it possible to reproduce a surface morphology very similar to that of the experiment, which is not the case for constant temperature (same radii in both side for SiO₂ layers). We measure a variation of 11 nm in the radii of curvature on either side of the tip (at $X = -20$ nm, the radius is equal to 16 nm, figure 13(e)), and at $X = 20$ nm, the radius is equal to 27 nm, figure 13(f)), thus confirming the conclusions of Lee *et al*. We also reproduce the difference of the evaporation depth for both SiO₂ layers (about 3–4 nm). Our results, however, are not quantitatively similar to the experimental ones. The variation of radii is quite the same (about 10 nm), but the simulated radii are lower (in the shadow side, 16 nm in simulation and 25 nm experimentally). In order to be

quantitative, it will undoubtedly be necessary to increase the size of the tip, even to modify the evaporation field of the Si layer, but more fundamentally to model more accurately the temperature distribution in the tip apex (figures 13(d)–(h)).

For the opposite case (SiO₂ in Si, figure 14), we make similar observations. Indeed, we note that the temperature simulation (figure 14(c)) makes it possible to reproduce a surface morphology very similar to that of the experiment (figure 14(a)), which is not the case for constant temperature (figure 14(b)). However, in this case, our results are far from being quantitatively comparable with experiment. They measure an increase of 10 nm of the SiO₂ layer radius (from 35 to 45 nm, figure 14(a)), whereas in simulation, we measure an increase of 4 nm, in the worst case (i.e. at the highest temperature in the figure 14(d), for $X = -6$ and $X = 6$ nm). Moreover, we observe a difference in the Si layer top depth (about 3–4 nm), as experimentally (figure 14(a)), but in the opposite direction (the uppermost layer is the laser side experimentally
<figure><img src="image_20.png" /><figcaption>**Figure 14.** (a) TEM observation of SiO₂ layer embedded in Si [78] after APT analysis. Simulated sample with constant *T* (*T* = 50 K) or (b) assuming a temperature gradient (in *X* direction, from 50 to 200 with α = 2 and β = 3). (d) Radius of curvature at the sample surface (c) in top view. (e) Average radius of curvature (nm) along the *X* axis in the region squared in the figure (d). (f) Electric field (V nm⁻¹) at the sample surface (c) in top view. (g) Average electric field profile along the *X* direction of the figure (f). On the figures (e) and (g), it is also reported the results obtained for another temperature evolution from 50 to 150 K (with α = 2 and β = 3). Reprinted from [78], Copyright (2014), with permission from Elsevier.</figcaption></figure>

and inversely in simulation). This difference is due to a non-centered SiO₂ layer between the two layers of Si, as observed by [78] before evaporation, in this case. However, it is obvious regarding the simulation results that the initial sample morphology (and even the evaporation field and temperature) must be modified to quantitatively correspond with experimental observation, even if, the main observations are well reproduced with the mesoscopic model.

The mesoscopic simulation model allows to reproduce the evaporation of layers (but also of particles). The model also reproduces the laser assisted evaporation with a simple approach (polynomial evolution). The next optimization of this model will be to introduce a temperature evolution which can take into account [41, 47] the thermal diffusivity variation between different phases.

**4. Conclusion**

A new model of FE was presented. The approach is based on the classical relationship between electric field (estimated from the Robin equation) and an evaporation rate formalism that is applied this time on a triangulated surface mesh of an emitter. The model is compared to predictions from a previous model at the atomic scale. This model, despite its simple approach, has many advantages: (1) low memory allocation requirement (around 100 MB), (2) low computation time (around 1 h), (3) real sample dimension can be introduced (regarding previous advantages), (4) can simulate various microstructural objects (particles, layer, ...), with results comparable to the experiments and (5) also to reproduce laser assisted evaporation. This last point remains to be optimized
<figure><img src="image_20.png" /><figcaption>**Figure 14.** (a) TEM observation of SiO₂ layer embedded in Si [78] after APT analysis. Simulated sample with constant *T* (*T* = 50 K) or (b) assuming a temperature gradient (in *X* direction, from 50 to 200 with α = 2 and β = 3). (d) Radius of curvature at the sample surface (c) in top view. (e) Average radius of curvature (nm) along the *X* axis in the region squared in the figure (d). (f) Electric field (V nm⁻¹) at the sample surface (c) in top view. (g) Average electric field profile along the *X* direction of the figure (f). On the figures (e) and (g), it is also reported the results obtained for another temperature evolution from 50 to 150 K (with α = 2 and β = 3). Reprinted from [78], Copyright (2014), with permission from Elsevier.</figcaption></figure>

and inversely in simulation). This difference is due to a non-centered SiO₂ layer between the two layers of Si, as observed by [78] before evaporation, in this case. However, it is obvious regarding the simulation results that the initial sample morphology (and even the evaporation field and temperature) must be modified to quantitatively correspond with experimental observation, even if, the main observations are well reproduced with the mesoscopic model.

The mesoscopic simulation model allows to reproduce the evaporation of layers (but also of particles). The model also reproduces the laser assisted evaporation with a simple approach (polynomial evolution). The next optimization of this model will be to introduce a temperature evolution which can take into account [41, 47] the thermal diffusivity variation between different phases.

**4. Conclusion**

A new model of FE was presented. The approach is based on the classical relationship between electric field (estimated from the Robin equation) and an evaporation rate formalism that is applied this time on a triangulated surface mesh of an emitter. The model is compared to predictions from a previous model at the atomic scale. This model, despite its simple approach, has many advantages: (1) low memory allocation requirement (around 100 MB), (2) low computation time (around 1 h), (3) real sample dimension can be introduced (regarding previous advantages), (4) can simulate various microstructural objects (particles, layer, ...), with results comparable to the experiments and (5) also to reproduce laser assisted evaporation. This last point remains to be optimized
in order to quantitatively reproduce the laser assisted experi-
mental results.

The model accurately solves the electric field at the sur-
face, and the ion trajectory from the sample to the detector.
Given its performance in computing time and memory alloca-
tion, this model could be used as a basis for a 3D reconstruction
algorithm [19, 20].

Data availability statement

The data cannot be made publicly available upon publication
because they contain commercially sensitive information. The
data that support the findings of this study are available upon
reasonable request from the authors.

Acknowledgments

The authors are grateful for scientific and technical input
and support from Semiconductor Research Corporation (SRC)
under task ID 2679.001.

ORCID iD

Constantinos Hatzoglou https://orcid.org/0000-0001-7101-0028

References

[1] Lefebvre W, Vurpillot F and Sauvage X 2016 *Atom Probe Tomography: Put Theory into Practice* (Boston, MA:Elsevier) (available at: https://shop.elsevier.com/books/atom-probe-tomography/lefebvre/978-0-12-804647-0)

[2] Gault B, Moody M P, Cairney J M and Ringer S P 2012 *Atom Probe Microscopy* (Berlin: Springer)

[3] Miller M K 2000 *Atom Probe Tomography: Analysis at the Atomic Level* (Berlin: Springer)

[4] Hatzoglou C, Da Costa G and Vurpillot F 2019 Enhanced dynamic reconstruction for atom probe tomography *Ultramicroscopy* **197** 72–82

[5] Hatzoglou C, Radiguet B and Pareige P 2017 Experimental artefacts occurring during atom probe tomography analysis of oxide nanoparticles in metallic matrix: quantification and correction *J. Nucl. Mater.* **492** 279–91

[6] Mesa G, Dobado-Fuentes E and Sáenz J J 1996 Image charge method for electrostatic calculations in field-emission diodes *J. Appl. Phys.* **79** 39–44

[7] Vurpillot F 2001 Étude de la fonction de transfert pointe-image de la sonde atomique tomographique *Thèse de doctorat* Université de Rouen, France

[8] Blavette D, Vurpillot F, Pareige P and Menand A 2001 A model accounting for spatial overlaps in 3D atom-probe microscopy *Ultramicroscopy* **89** 145–53

[9] Wang X *et al* 2020 Interpreting nanovoids in atom probe tomography data for accurate local compositional measurements *Nat. Commun.* **11** 1–11

[10] Vurpillot F, Gaillard A, Da Costa G and Deconihout B 2013 A model to predict image formation in atom probe tomography *Ultramicroscopy* **132** 152–7

[11] Geiser B, Larson D, Gerstl S, Reinhard D, Kelly T, Prosa T and Olson D 2009 A system for simulation of tip evolution under field evaporation *Microsc. Microanal.* **15** 302–3

[12] Oberdorfer C, Eich S M, Lütkemeyer M and Schmitz G 2015 Applications of a versatile modelling approach to 3D atom probe simulations *Ultramicroscopy* **159** 184–94

[13] Oberdorfer C and Schmitz G 2011 On the field evaporation behavior of dielectric materials in three-dimensional atom probe: a numeric simulation *Microsc. Microanal.* **17** 15–25

[14] Vurpillot F and Oberdorfer C 2015 Modeling atom probe tomography: a review *Ultramicroscopy* **159** 202–16

[15] Rolland N, Vurpillot F, Duguay S and Blavette D 2015 A meshless algorithm to model field evaporation in atom probe tomography *Microsc. Microanal.* **21** 1649–56

[16] Robin G 1886 Sur la distribution de l'électricité à la surface des conducteurs fermés des conducteurs ouverts *Ann. Sci. Ec. Norm.* **3** 3–58

[17] Robin G 1887 Distribution de l'électricité sur une surface fermée convexe *C. R. Acad. Sci. SeR.* **104** 1834–6

[18] Barnes J and Hut P 1986 A hierarchical O(N log N) force-calculation algorithm *Nature* **324** 446–9

[19] Beinke D, Oberdorfer C and Schmitz G 2016 Towards an accurate volume reconstruction in atom probe tomography *Ultramicroscopy* **165** 34–41

[20] Ling Y-T, Cools S, Bogdanowicz J, Fleischmann C, Beenhouwer J D, Sijbers J and Vandervorst W 2022 A bottom-up volume reconstruction method for atom probe tomography *Microsc. Microanal.* **28** 1102–15

[21] Waugh A R, Boyes E D and Southon M J 1976 Investigations of field evaporation with a field-desorption microscope *Surf. Sci.* **61** 109–42

[22] Saxey D W 2011 Correlated ion analysis and the interpretation of atom probe mass spectra *Ultramicroscopy* **111** 473–9

[23] Gault B, Danoix F, Hoummada K, Mangelinck D and Leitner H 2012 Impact of directional walk on atom probe microanalysis *Ultramicroscopy* **113** 182–91

[24] Karahka M and Kreuzer H J 2013 Field evaporation of oxides: a theoretical study *Ultramicroscopy* **132** 54–59

[25] Silaeva E P, Arnoldi L, Karahka M L, Deconihout B, Menand A, Kreuzer H J and Vella A 2014 Do dielectric nanostructures turn metallic in high-electric dc fields? *Nano Lett.* **14** 6066–72

[26] Zanuttini D, Blum I, Rigutti L, Vurpillot F, Douady J, Jacquet E, Anglade P-M and Gervais B 2017 Simulation of field-induced molecular dissociation in atom-probe tomography: identification of a neutral emission channel *Phys. Rev. A* **95** 061401

[27] Vurpillot F, Parviainen S, Djurabekova F, Zanuttini D and Gervais B 2018 Simulation tools for atom probe tomography: a path for diagnosis and treatment of image degradation *Mater. Charact* **146** 336–46

[28] Oberdorfer C, Withrow T, Marquis E A and Windl W 2019 Atomistic-simulation based modeling of atom probe tomography *Microsc. Microanal.* **25** 284–5

[29] Qi J, Oberdorfer C, Windl W and Marquis E A 2022 *Ab initio* simulation of field evaporation *Phys. Rev. Mater.* **6** 093602

[30] Xu Z, Li D, Xu W, Devaraj A, Colby R, Thevuthasan S, Geiser B P and Larson D J 2015 Simulation of heterogeneous atom probe tip shapes evolution during field evaporation using a level set method and different evaporation models *Comput. Phys. Commun.* **189** 106–13

[31] Haley D, Moody M P and Smith G D W 2013 Level set methods for modelling field evaporation in atom probe *Microsc. Microanal.* **19** 1709–17

[32] Haley D, Bagot P A J and Moody M P 2018 Extending continuum models for atom probe simulation *Mater. Charact* **146** 299–306

[33] Fletcher C, Moody M P and Haley D 2019 Fast modelling of field evaporation in atom probe tomography using level set methods *J. Phys. Appl. Phys.* **52** 435305

[34] Deserno M 2004 How to generate equidistributed points on the surface of a sphere
in order to quantitatively reproduce the laser assisted experi-
mental results.

The model accurately solves the electric field at the sur-
face, and the ion trajectory from the sample to the detector.
Given its performance in computing time and memory alloca-
tion, this model could be used as a basis for a 3D reconstruction
algorithm [19, 20].

Data availability statement

The data cannot be made publicly available upon publication
because they contain commercially sensitive information. The
data that support the findings of this study are available upon
reasonable request from the authors.

Acknowledgments

The authors are grateful for scientific and technical input
and support from Semiconductor Research Corporation (SRC)
under task ID 2679.001.

ORCID iD

Constantinos Hatzoglou https://orcid.org/0000-0001-7101-0028

References

[1] Lefebvre W, Vurpillot F and Sauvage X 2016 *Atom Probe Tomography: Put Theory into Practice* (Boston, MA:Elsevier) (available at: https://shop.elsevier.com/books/atom-probe-tomography/lefebvre/978-0-12-804647-0)

[2] Gault B, Moody M P, Cairney J M and Ringer S P 2012 *Atom Probe Microscopy* (Berlin: Springer)

[3] Miller M K 2000 *Atom Probe Tomography: Analysis at the Atomic Level* (Berlin: Springer)

[4] Hatzoglou C, Da Costa G and Vurpillot F 2019 Enhanced dynamic reconstruction for atom probe tomography *Ultramicroscopy* **197** 72–82

[5] Hatzoglou C, Radiguet B and Pareige P 2017 Experimental artefacts occurring during atom probe tomography analysis of oxide nanoparticles in metallic matrix: quantification and correction *J. Nucl. Mater.* **492** 279–91

[6] Mesa G, Dobado-Fuentes E and Sáenz J J 1996 Image charge method for electrostatic calculations in field-emission diodes *J. Appl. Phys.* **79** 39–44

[7] Vurpillot F 2001 *Étude de la fonction de transfert pointe-image de la sonde atomique tomographique Thèse de doctorat* Université de Rouen, France

[8] Blavette D, Vurpillot F, Pareige P and Menand A 2001 A model accounting for spatial overlaps in 3D atom-probe microscopy *Ultramicroscopy* **89** 145–53

[9] Wang X *et al* 2020 Interpreting nanovoids in atom probe tomography data for accurate local compositional measurements *Nat. Commun.* **11** 1–11

[10] Vurpillot F, Gaillard A, Da Costa G and Deconihout B 2013 A model to predict image formation in atom probe tomography *Ultramicroscopy* **132** 152–7

[11] Geiser B, Larson D, Gerstl S, Reinhard D, Kelly T, Prosa T and Olson D 2009 A system for simulation of tip evolution under field evaporation *Microsc. Microanal.* **15** 302–3

[12] Oberdorfer C, Eich S M, Lütkemeyer M and Schmitz G 2015 Applications of a versatile modelling approach to 3D atom probe simulations *Ultramicroscopy* **159** 184–94

[13] Oberdorfer C and Schmitz G 2011 On the field evaporation behavior of dielectric materials in three-dimensional atom probe: a numeric simulation *Microsc. Microanal.* **17** 15–25

[14] Vurpillot F and Oberdorfer C 2015 Modeling atom probe tomography: a review *Ultramicroscopy* **159** 202–16

[15] Rolland N, Vurpillot F, Duguay S and Blavette D 2015 A meshless algorithm to model field evaporation in atom probe tomography *Microsc. Microanal.* **21** 1649–56

[16] Robin G 1886 Sur la distribution de l'électricité à la surface des conducteurs fermés des conducteurs ouverts *Ann. Sci. Ec. Norm.* **3** 3–58

[17] Robin G 1887 Distribution de l'électricité sur une surface fermée convexe *C. R. Acad. Sci. SeR.* **104** 1834–6

[18] Barnes J and Hut P 1986 A hierarchical O(N log N) force-calculation algorithm *Nature* **324** 446–9

[19] Beinke D, Oberdorfer C and Schmitz G 2016 Towards an accurate volume reconstruction in atom probe tomography *Ultramicroscopy* **165** 34–41

[20] Ling Y-T, Cools S, Bogdanowicz J, Fleischmann C, Beenhouwer J D, Sijbers J and Vandervorst W 2022 A bottom-up volume reconstruction method for atom probe tomography *Microsc. Microanal.* **28** 1102–15

[21] Waugh A R, Boyes E D and Southon M J 1976 Investigations of field evaporation with a field-desorption microscope *Surf. Sci.* **61** 109–42

[22] Saxey D W 2011 Correlated ion analysis and the interpretation of atom probe mass spectra *Ultramicroscopy* **111** 473–9

[23] Gault B, Danoix F, Hoummada K, Mangelinck D and Leitner H 2012 Impact of directional walk on atom probe microanalysis *Ultramicroscopy* **113** 182–91

[24] Karahka M and Kreuzer H J 2013 Field evaporation of oxides: a theoretical study *Ultramicroscopy* **132** 54–59

[25] Silaeva E P, Arnoldi L, Karahka M L, Deconihout B, Menand A, Kreuzer H J and Vella A 2014 Do dielectric nanostructures turn metallic in high-electric dc fields? *Nano Lett.* **14** 6066–72

[26] Zanuttini D, Blum I, Rigutti L, Vurpillot F, Douady J, Jacquet E, Anglade P-M and Gervais B 2017 Simulation of field-induced molecular dissociation in atom-probe tomography: identification of a neutral emission channel *Phys. Rev. A* **95** 061401

[27] Vurpillot F, Parviainen S, Djurabekova F, Zanuttini D and Gervais B 2018 Simulation tools for atom probe tomography: a path for diagnosis and treatment of image degradation *Mater. Charact* **146** 336–46

[28] Oberdorfer C, Withrow T, Marquis E A and Windl W 2019 Atomistic-simulation based modeling of atom probe tomography *Microsc. Microanal.* **25** 284–5

[29] Qi J, Oberdorfer C, Windl W and Marquis E A 2022 *Ab initio* simulation of field evaporation *Phys. Rev. Mater.* **6** 093602

[30] Xu Z, Li D, Xu W, Devaraj A, Colby R, Thevuthasan S, Geiser B P and Larson D J 2015 Simulation of heterogeneous atom probe tip shapes evolution during field evaporation using a level set method and different evaporation models *Comput. Phys. Commun.* **189** 106–13

[31] Haley D, Moody M P and Smith G D W 2013 Level set methods for modelling field evaporation in atom probe *Microsc. Microanal.* **19** 1709–17

[32] Haley D, Bagot P A J and Moody M P 2018 Extending continuum models for atom probe simulation *Mater. Charact* **146** 299–306

[33] Fletcher C, Moody M P and Haley D 2019 Fast modelling of field evaporation in atom probe tomography using level set methods *J. Phys. Appl. Phys.* **52** 435305

[34] Deserno M 2004 How to generate equidistributed points on the surface of a sphere
[35] Gomer R 1992 *Field Emissions and Field Ionization* (American Inst. of Physics) (available at: https://books.google.no/books/about/Field_Emissions_and_Field_ Ionization.html?id=0iIkSfKA6NUC&redir_esc=y)

[36] Rusinkiewicz S 2004 Estimating curvatures and their derivatives on triangle meshes *Proc. 2nd Int. Symp. on 3D Data Processing, Visualization and Transmission, 2004.* 3DPVT 2004 pp 486–93

[37] Kellogg G L 1984 Measurement of activation energies for field evaporation of tungsten ions as a function of electric field *Phys. Rev. B* **29** 4304–12

[38] Kreuzer H J and Nath K 1987 Field evaporation *Surf. Sci.* **183** 591–608

[39] Arnoldi L, Vella A, Houard J and Deconihout B 2012 Antenna effect in laser assisted atom probe tomography: how the field emitter aspect ratio can enhance atomic scale imaging *Appl. Phys. Lett.* **101** 153101

[40] Wang H, Houard J, Arnoldi L, Hideur A, Silaeva E P, Deconihout B and Vella A 2016 Effect of the laser pulse width on the field evaporation behavior of metals and oxides *Ultramicroscopy* **160** 18–22

[41] Vurpillot F, Houard J, Vella A and Deconihout B 2009 Thermal response of a field emitter subjected to ultra-fast laser illumination *J. Phys. Appl. Phys.* **42** 125502

[42] Miller M K and Forbes R G 2009 Atom probe tomography *Mater. Charact.* **60** 461–9

[43] Peralta J, Broderick S R and Rajan K 2013 Mapping energetics of atom probe evaporation events through first principles calculations *Ultramicroscopy* **132** 143–51

[44] Carrasco T, Peralta J, Loyola C and Broderick S R 2018 Modeling field evaporation degradation of metallic surfaces by first principles calculations: a case study for Al, Au, Ag, and Pd *J. Phys. Conf. Ser.* **1043** 012039

[45] Ohnuma T 2019 Surface diffusion of Fe and Cu on Fe (001) under electric field using first-principles calculations *Microsc. Microanal.* **25** 547–53

[46] Ernst N 1979 Experimental investigation on field evaporation of singly and doubly charged rhodium *Surf. Sci.* **87** 469–82

[47] Houard J, Vella A, Vurpillot F and Deconihout B 2011 Three-dimensional thermal response of a metal subwavelength tip under femtosecond laser illumination *Phys. Rev. B* **84** 033405

[48] Cade R 1972 An existence theorem for Robin's equation *Math. Proc. Camb. Philos. Soc.* **72** 489–98

[49] Verlet L 1967 Computer "experiments" on classical fluids. I. Thermodynamical properties of Lennard-Jones molecules *Phys. Rev.* **159** 98–103

[50] Loi S T, Gault B, Ringer S P, Larson D J and Geiser B P 2013 Electrostatic simulations of a local electrode atom probe: the dependence of tomographic reconstruction parameters on specimen and microscope geometry *Ultramicroscopy* **132** 107–13

[51] Kühbach M, Bajaj P, Zhao H, Çelik M H, Jägle E A and Gault B 2021 On strong-scaling and open-source tools for analyzing atom probe tomography data *npj Comput. Mater.* **7** 1–10

[52] Koelling S, Innocenti N, Schulze A, Gilbert M, Kambham A K and Vandervorst W 2011 *In-situ* observation of non-hemispherical tip shape formation during laser-assisted atom probe tomography *J. Appl. Phys.* **109** 104909

[53] Kingham D R 1982 The post-ionization of field evaporated ions: a theoretical explanation of multiple charge states *Surf. Sci.* **116** 273–301

[54] Haydock R and Kingham D R 1980 Post-ionization of field-evaporated ions *Phys. Rev. Lett.* **44** 1520–3

[55] Amirifar N et al 2015 Quantitative analysis of doped/undoped ZnO nanomaterials using laser assisted atom probe tomography: influence of the analysis parameters *J. Appl. Phys.* **118** 215703

[56] Santhanagopalan D, Schreiber D K, Perea D E, Martens R, Janssen Y, Kalifah P and Meng Y S 2015 Effects of laser energy and wavelength on the analysis of LiFePO4 using laser assisted atom probe tomography *Ultramicroscopy* **148** 57–66

[57] Kumar A, Bogdanowicz J, Demeulemeester J, Bran J, Melkonyan D, Fleischmann C and Vandervorst W 2018 Measurement of the apex temperature of a nanoscale semiconducting field emitter illuminated by a femtosecond pulsed laser *J. Appl. Phys.* **124** 245105

[58] Rolland N, Larson D J, Geiser B P, Duguay S, Vurpillot F and Blavette D 2015 An analytical model accounting for tip shape evolution during atom probe analysis of heterogeneous materials *Ultramicroscopy* **159** 195–201

[59] Rolland N, Vurpillot F, Duguay S and Blavette D 2015 Dynamic evolution and fracture of multilayer field emitters in atom probe tomography: a new interpretation *Eur. Phys. J. Appl. Phys.* **72** 21001

[60] Marquis E A, Geiser B P, Prosa T J and Larson D J 2011 Evolution of tip shape during field evaporation of complex multilayer structures *J. Microsc.* **241** 225–33

[61] Larson D J, Geiser B P, Prosa T J and Kelly T F 2012 On the use of simulated field-evaporated specimen apex shapes in atom probe tomography data reconstruction *Microsc. Microanal.* **18** 953–63

[62] Vurpillot F, Rolland N, Estivill R, Duguay S and Blavette D 2016 Accuracy of analyses of microelectronics nanostructures in atom probe tomography *Semicond. Sci. Technol.* **31** 074002

[63] Rolland N, Vurpillot F, Duguay S, Mazumder B, Speck J S and Blavette D 2017 New atom probe tomography reconstruction algorithm for multilayered samples: beyond the hemispherical constraint *Microsc. Microanal.* **23** 247–54

[64] Hatzoglou C, Radiguet B, Vurpillot F and Pareige P 2018 A chemical composition correction model for nanoclusters observed by APT—application to ODS steel nanoparticles *J. Nucl. Mater.* **505** 240–8

[65] Lawitzki R, Stender P and Schmitz G 2021 Compensating local magnifications in atom probe tomography for accurate analysis of nano-sized precipitates *Microsc. Microanal.* **27** 499–510

[66] Miller M K and Hetherington M G 1991 Local magnification effects in the atom probe *Surf. Sci.* **246** 442–9

[67] Vurpillot F, Bostel A and Blavette D 2000 Trajectory overlaps and local magnification in three-dimensional atom probe *Appl. Phys. Lett.* **76** 3127–9

[68] Bas P, Bostel A, Deconihout B and Blavette D 1995 A general protocol for the reconstruction of 3D atom probe data *Appl. Surf. Sci.* **87** 298–304

[69] Vurpillot F, Gault B, Geiser B P and Larson D J 2013 Reconstructing atom probe data: a review *Ultramicroscopy* **132** 19–30

[70] Devaraj A, Colby R, Vurpillot F and Thevuthasan S 2014 Understanding atom probe tomography of oxide-supported metal nanoparticles by correlation with atomic-resolution electron microscopy and field evaporation simulation *J. Phys. Chem. Lett.* **5** 1361–7

[71] Kalokhtina O 2012 Etude par Sonde Atomique Tomographique de la formation de nano-particules dans les aciers ODS et NDS *Thèse de doctorat* Université de Rouen

[72] Knipling K E, Dunand D C and Seidman D N 2007 Atom probe tomographic studies of precipitation in Al-0.1 Zr-0.1 Ti (at.%) alloys *Microsc. Microanal.* **13** 503–16

[73] Novy S, Pareige P and Pareige C 2009 Atomic scale analysis and phase separation understanding in a thermally aged Fe-20 at.% Cr alloy *J. Nucl. Mater.* **384** 96–102
[35] Gomer R 1992 *Field Emissions and Field Ionization* (American Inst. of Physics) (available at: https://books.google.no/books/about/Field_Emissions_and_Field_Ionization.html?id=0iIkSfKA6NUC&redir_esc=y)

[36] Rusinkiewicz S 2004 Estimating curvatures and their derivatives on triangle meshes *Proc. 2nd Int. Symp. on 3D Data Processing, Visualization and Transmission, 2004.* 3DPVT 2004 pp 486–93

[37] Kellogg G L 1984 Measurement of activation energies for field evaporation of tungsten ions as a function of electric field *Phys. Rev. B* **29** 4304–12

[38] Kreuzer H J and Nath K 1987 Field evaporation *Surf. Sci.* **183** 591–608

[39] Arnoldi L, Vella A, Houard J and Deconihout B 2012 Antenna effect in laser assisted atom probe tomography: how the field emitter aspect ratio can enhance atomic scale imaging *Appl. Phys. Lett.* **101** 153101

[40] Wang H, Houard J, Arnoldi L, Hideur A, Silaeva E P, Deconihout B and Vella A 2016 Effect of the laser pulse width on the field evaporation behavior of metals and oxides *Ultramicroscopy* **160** 18–22

[41] Vurpillot F, Houard J, Vella A and Deconihout B 2009 Thermal response of a field emitter subjected to ultra-fast laser illumination *J. Phys. Appl. Phys.* **42** 125502

[42] Miller M K and Forbes R G 2009 Atom probe tomography *Mater. Charact.* **60** 461–9

[43] Peralta J, Broderick S R and Rajan K 2013 Mapping energetics of atom probe evaporation events through first principles calculations *Ultramicroscopy* **132** 143–51

[44] Carrasco T, Peralta J, Loyola C and Broderick S R 2018 Modeling field evaporation degradation of metallic surfaces by first principles calculations: a case study for Al, Au, Ag, and Pd *J. Phys. Conf. Ser.* **1043** 012039

[45] Ohnuma T 2019 Surface diffusion of Fe and Cu on Fe (001) under electric field using first-principles calculations *Microsc. Microanal.* **25** 547–53

[46] Ernst N 1979 Experimental investigation on field evaporation of singly and doubly charged rhodium *Surf. Sci.* **87** 469–82

[47] Houard J, Vella A, Vurpillot F and Deconihout B 2011 Three-dimensional thermal response of a metal subwavelength tip under femtosecond laser illumination *Phys. Rev. B* **84** 033405

[48] Cade R 1972 An existence theorem for Robin's equation *Math. Proc. Camb. Philos. Soc.* **72** 489–98

[49] Verlet L 1967 Computer "experiments" on classical fluids. I. Thermodynamical properties of Lennard-Jones molecules *Phys. Rev.* **159** 98–103

[50] Loi S T, Gault B, Ringer S P, Larson D J and Geiser B P 2013 Electrostatic simulations of a local electrode atom probe: the dependence of tomographic reconstruction parameters on specimen and microscope geometry *Ultramicroscopy* **132** 107–13

[51] Kühbach M, Bajaj P, Zhao H, Çelik M H, Jägle E A and Gault B 2021 On strong-scaling and open-source tools for analyzing atom probe tomography data *npj Comput. Mater.* **7** 1–10

[52] Koelling S, Innocenti N, Schulze A, Gilbert M, Kambham A K and Vandervorst W 2011 *In-situ* observation of non-hemispherical tip shape formation during laser-assisted atom probe tomography *J. Appl. Phys.* **109** 104909

[53] Kingham D R 1982 The post-ionization of field evaporated ions: a theoretical explanation of multiple charge states *Surf. Sci.* **116** 273–301

[54] Haydock R and Kingham D R 1980 Post-ionization of field-evaporated ions *Phys. Rev. Lett.* **44** 1520–3

[55] Amirifar N et al 2015 Quantitative analysis of doped/undoped ZnO nanomaterials using laser assisted atom probe tomography: influence of the analysis parameters *J. Appl. Phys.* **118** 215703

[56] Santhanagopalan D, Schreiber D K, Perea D E, Martens R, Janssen Y, Kalifah P and Meng Y S 2015 Effects of laser energy and wavelength on the analysis of LiFePO4 using laser assisted atom probe tomography *Ultramicroscopy* **148** 57–66

[57] Kumar A, Bogdanowicz J, Demeulemeester J, Bran J, Melkonyan D, Fleischmann C and Vandervorst W 2018 Measurement of the apex temperature of a nanoscale semiconducting field emitter illuminated by a femtosecond pulsed laser *J. Appl. Phys.* **124** 245105

[58] Rolland N, Larson D J, Geiser B P, Duguay S, Vurpillot F and Blavette D 2015 An analytical model accounting for tip shape evolution during atom probe analysis of heterogeneous materials *Ultramicroscopy* **159** 195–201

[59] Rolland N, Vurpillot F, Duguay S and Blavette D 2015 Dynamic evolution and fracture of multilayer field emitters in atom probe tomography: a new interpretation *Eur. Phys. J. Appl. Phys.* **72** 21001

[60] Marquis E A, Geiser B P, Prosa T J and Larson D J 2011 Evolution of tip shape during field evaporation of complex multilayer structures *J. Microsc.* **241** 225–33

[61] Larson D J, Geiser B P, Prosa T J and Kelly T F 2012 On the use of simulated field-evaporated specimen apex shapes in atom probe tomography data reconstruction *Microsc. Microanal.* **18** 953–63

[62] Vurpillot F, Rolland N, Estivill R, Duguay S and Blavette D 2016 Accuracy of analyses of microelectronics nanostructures in atom probe tomography *Semicond. Sci. Technol.* **31** 074002

[63] Rolland N, Vurpillot F, Duguay S, Mazumder B, Speck J S and Blavette D 2017 New atom probe tomography reconstruction algorithm for multilayered samples: beyond the hemispherical constraint *Microsc. Microanal.* **23** 247–54

[64] Hatzoglou C, Radiguet B, Vurpillot F and Pareige P 2018 A chemical composition correction model for nanoclusters observed by APT—application to ODS steel nanoparticles *J. Nucl. Mater.* **505** 240–8

[65] Lawitzki R, Stender P and Schmitz G 2021 Compensating local magnifications in atom probe tomography for accurate analysis of nano-sized precipitates *Microsc. Microanal.* **27** 499–510

[66] Miller M K and Hetherington M G 1991 Local magnification effects in the atom probe *Surf. Sci.* **246** 442–9

[67] Vurpillot F, Bostel A and Blavette D 2000 Trajectory overlaps and local magnification in three-dimensional atom probe *Appl. Phys. Lett.* **76** 3127–9

[68] Bas P, Bostel A, Deconihout B and Blavette D 1995 A general protocol for the reconstruction of 3D atom probe data *Appl. Surf. Sci.* **87** 298–304

[69] Vurpillot F, Gault B, Geiser B P and Larson D J 2013 Reconstructing atom probe data: a review *Ultramicroscopy* **132** 19–30

[70] Devaraj A, Colby R, Vurpillot F and Thevuthasan S 2014 Understanding atom probe tomography of oxide-supported metal nanoparticles by correlation with atomic-resolution electron microscopy and field evaporation simulation *J. Phys. Chem. Lett.* **5** 1361–7

[71] Kalokhtina O 2012 Etude par Sonde Atomique Tomographique de la formation de nano-particules dans les aciers ODS et NDS *Thèse de doctorat* Université de Rouen

[72] Knipling K E, Dunand D C and Seidman D N 2007 Atom probe tomographic studies of precipitation in Al-0.1 Zr-0.1 Ti (at.%) alloys *Microsc. Microanal.* **13** 503–16

[73] Novy S, Pareige P and Pareige C 2009 Atomic scale analysis and phase separation understanding in a thermally aged Fe-20 at.% Cr alloy *J. Nucl. Mater.* **384** 96–102
[74] Hatzoglou C, Radiguet B, Da Costa G, Pareige P, Roussel M, Hernandez-Mayoral M and Pareige C 2019 Quantification of APT physical limitations on chemical composition of precipitates in Fe–Cr alloys *J. Nucl. Mater.* **522** 64–73

[75] Kuksenko V, Pareige C, Genevois C, Cuvilly F, Roussel M and Pareige P 2011 Effect of neutron-irradiation on the microstructure of a Fe–12 at.% Cr alloy *J. Nucl. Mater.* **415** 61–66

[76] Melkonyan D, Fleischmann C, Arnoldi L, Demeulemeester J, Kumar A, Bogdanowicz J, Vurpillot F and Vandervorst W 2017 Atom probe tomography analysis of SiGe fins

embedded in SiO₂: facts and artefacts *Ultramicroscopy* **179** 100–7

[77] Shu S, Wirth B D, Wells P B, Morgan D D and Odette G R 2018 Multi-technique characterization of the precipitates in thermally aged and neutron irradiated Fe-Cu and Fe-Cu-Mn model alloys: atom probe tomography reconstruction implications *Acta Mater.* **146** 237–52

[78] Lee J H, Lee B H, Kim Y T, Kim J J, Lee S Y, Lee K P and Park C G 2014 Study of vertical Si/SiO₂ interface using laser-assisted atom probe tomography and transmission electron microscopy *Micron* **58** 32–37
[74] Hatzoglou C, Radiguet B, Da Costa G, Pareige P, Roussel M, Hernandez-Mayoral M and Pareige C 2019 Quantification of APT physical limitations on chemical composition of precipitates in Fe–Cr alloys *J. Nucl. Mater.* **522** 64–73

[75] Kuksenko V, Pareige C, Genevois C, Cuvilly F, Roussel M and Pareige P 2011 Effect of neutron-irradiation on the microstructure of a Fe–12 at.% Cr alloy *J. Nucl. Mater.* **415** 61–66

[76] Melkonyan D, Fleischmann C, Arnoldi L, Demeulemeester J, Kumar A, Bogdanowicz J, Vurpillot F and Vandervorst W 2017 Atom probe tomography analysis of SiGe fins

embedded in SiO₂: facts and artefacts *Ultramicroscopy* **179** 100–7

[77] Shu S, Wirth B D, Wells P B, Morgan D D and Odette G R 2018 Multi-technique characterization of the precipitates in thermally aged and neutron irradiated Fe-Cu and Fe-Cu-Mn model alloys: atom probe tomography reconstruction implications *Acta Mater.* **146** 237–52

[78] Lee J H, Lee B H, Kim Y T, Kim J J, Lee S Y, Lee K P and Park C G 2014 Study of vertical Si/SiO₂ interface using laser-assisted atom probe tomography and transmission electron microscopy *Micron* **58** 32–37