# Unveiling the crystallization kinetics in Ge-rich GexTe alloys by large
scale simulations with a machine-learned interatomic potential

Dario Baratella  Omar Abou El Kheir  Marco Bernasconi∗ Department of Materials
Science, University of Milano-Bicocca, Via R. Cozzi 55, I-20125 Milano, Italy

###### Abstract

A machine-learned interatomic potential for Ge-rich GexTe alloys has been
developed aiming at uncovering the kinetics of phase separation and
crystallization in these materials. The results are of interest for the
operation of embedded phase change memories which exploits Ge-enrichment of
GeSbTe alloys to raise the crystallization temperature. The potential is
generated by fitting a large database of energies and forces computed within
Density Functional Theory with the neural network scheme implemented in the
DeePMD-kit package. The potential is highly accurate and suitable to describe
the structural and dynamical properties of the liquid, amorphous and
crystalline phases of the wide range of compositions from pure Ge and
stoichiometric GeTe to the Ge-rich Ge2Te alloy. Large scale molecular dynamics
simulations revealed a crystallization mechanism which depends on temperature.
At 600 K, segregation of most of Ge in excess occurs on the ns time scale
followed by crystallization of nearly stoichiometric GeTe regions. At 500 K,
nucleation of crystalline GeTe occurs before phase separation, followed by a
slow crystal growth due to the concurrent expulsion of Ge in excess.

00footnotetext: Email<EMAIL_ADDRESS>
##  I Introduction

Chalcogenide GeSbTe (GST) phase change alloys are employed in emerging new
electronic non-volatile memories called Phase Change Memories (PCMs) Wuttig
and Yamada (2007); Noé  _et al._ (2017); Zhang  _et al._ (2019); Pirovano  _et
al._ (2004). These devices feature a rapid and reversible transition of the
active material between the crystalline and amorphous phases upon Joule
heating. The two phases correspond to the two states of the memory which can
be discriminated thanks to a large contrast in the electrical resistivity.
Readout of the memory consists of the measurement of the electrical resistance
at low bias, while programming is achieved by applying current pulses to
amorphize the crystal via melting (reset) or to recrystallize the amorphous
phase (set).

PCMs based on the flagship compound Ge2Sb2Te5 (GST225), which is a
pseudobinary alloy on the GeTe-Sb2Te3 tie line, was brought to market in 2017
by Intel and Micron in a 3D cross-bar architecture (3D XPointTM) Choe (2017).
More recently, PCMs have also gained attention for embedded applications, in
particular for the automotive sector Cappelletti  _et al._ (2020); Redaelli
_et al._ (2022). The crystallization temperature Tx of the most commonly used
GST225 (420-440 K) Noé  _et al._ (2016) is, however, too low for applications
in embedded memories that require data retention at higher temperatures. In
fact, a compulsory fabrication step of embedded memories is the soldering
process in which the device is exposed to 530 K for a few minutes after the
microcontroller code was written in the memory Redaelli  _et al._ (2022). To
meet this requirement, various materials alternative to GST225 have been
explored, including InSbTe, GaSbTe, and InGeTe alloys Kim  _et al._ (2007,
2009); Cheng  _et al._ (2011a). Additionally, doping with N or O atoms
Morikawa  _et al._ (2007) and enrichment with Sb and Ge have been proposed
Cheng  _et al._ (2011b). Among these options, Ge-rich GeSbTe alloys emerged
recently as the most promising materials with Tx exceeding 600 K for Ge-rich
alloys on the Ge-Sb2Te3 tie-line Zuliani  _et al._ (2013, 2015); Palumbo  _et
al._ (2017); Navarro  _et al._ (2013a). Several high density embedded memories
for microcontrollers have been since then reported in literature Sousa  _et
al._ (2015); Arnaud  _et al._ (2020); Kiouseloglou  _et al._ (2014); Grossier
_et al._ (2023).

The raise in Tx of Ge-rich GST alloys is ascribed to the segregation of Ge
atoms and phase separation into pure Ge and Ge-poor GST alloy during the
crystallization process Zuliani  _et al._ (2015). The crystallization is
supposed to be slowed down by the mass transport involved in the phase
separation, resulting in longer times for crystal nucleation. Evidences of the
phase separation were reported in both the set process Zuliani  _et al._
(2015); Agati  _et al._ (2019); Luong  _et al._ (2021); Rahier  _et al._
(2023); Privitera  _et al._ (2018, 2020); Prazakova  _et al._ (2022); Díaz
Fattorini  _et al._ (2022); Cecchi  _et al._ (2022) and the forming operation
of the memory (initialization by the first programming pulse) Palumbo  _et
al._ (2017). However, the inhomogeneity due to phase separation could cause a
high cell-to-cell variability. Moreover, a drift in the electrical resistance
with time was reported in both the set and reset states. A resistance drift in
the reset state (amorphous phase) is common to all phase change materials as
it is ascribed to structural relaxation towards a more stable amorphous
structure. The resistance drift in the set state (crystalline), instead, is
peculiar to Ge-rich GST alloys as it is due to recrystallization of a residual
amorphous region close to the bottom electrode in the phase separated system
Laurin  _et al._ (2023). Gaining a deeper understanding of the segregation and
crystallization processes is thus mandatory to mitigate these detrimental
effects.

Several details of the overall process are, however, unclear. Although several
works have revealed the presence of amorphous Ge and of a cubic crystalline
phase of GST Zuliani  _et al._ (2013); Agati  _et al._ (2019); Luong  _et al._
(2021); Rahier  _et al._ (2023); Privitera  _et al._ (2018, 2020); Prazakova
_et al._ (2022); Díaz Fattorini  _et al._ (2022); Cecchi  _et al._ (2022), the
composition of the cubic phase is unclear because analytical tools such as
electron energy loss spectroscopy or energy dispersive x-ray measurements
provide only a composition averaged over different grains including pure Ge
grains. The presence of crystalline GeTe Rahier  _et al._ (2023) and of Sb-
rich alloys Petroni  _et al._ (2022, 2023) was also detected.

In this respect, atomistic simulations can provide useful insights on both the
crystallization and segregation processes. Simulations based on Density
Functional Theory (DFT), for instance, revealed that the amorphous phase of
GST alloys is unstable with respect to Ge segregation for Ge content above 50
$\%$ Sun  _et al._ (2021). A high throughput DFT study also revealed that off-
stoichiometric compositions (off the GeTe-Sb2Te3 tie-line) might crystallize
in the cubic phase at the operation conditions of the memories Abou El Kheir
and Bernasconi (2021). These DFT works, however, have only addressed the
thermodynamics of the phase transformation and not its kinetics because of the
limitations in the time and length scales of DFT simulations.

A route to overcome the limitations of DFT methods is the development of
machine-learned interatomic potentials generated by fitting a database of DFT
energies and forces. In the past, an interatomic potential for GeTe Sosso  _et
al._ (2012a) was generated within the Neural Network (NN) scheme proposed by
Behler and Parrinello Behler and Parrinello (2007), while more recently
interatomic potentials for stoichiometric GST225 Abou El Kheir  _et al._
(2024) were generated with the NN scheme implemented in the DeePMD code Zhang
_et al._ (2018); Zeng  _et al._ (2023); Wang  _et al._ (2018) and with the
Gaussian Approximation Potential Zhou  _et al._ (2023).

In the perspective of developing an interatomic potential for large scale
simulations of Ge-rich GST alloys, we here address the study of the simpler
Ge-rich GexTe binary alloy that shares several properties with the ternary GST
system and for which detailed experimental data are available from time
resolved reflectivity, x-ray diffraction, transmission electron microscopy and
Raman spectroscopy Carria  _et al._ (2012); Gourvest  _et al._ (2009); Raoux
_et al._ (2009a). Integration of Ge-rich GexTe alloys in PCMs was also
reported in Ref. Navarro  _et al._ (2013b). We considered in particular the
composition Ge2Te which is close to the Ge63Te37 alloys studied experimentally
in Ref. Carria  _et al._ (2012).

We developed a NN potential suitable to describe Ge2Te and the products of its
crystallization process, stoichiometric GeTe and pure Ge, by using the DeePMD
code Zhang  _et al._ (2018); Zeng  _et al._ (2023); Wang  _et al._ (2018)
already employed in our recent work on GST225 Abou El Kheir  _et al._ (2024).
We validated the potential by analyzing the structural and dynamical
properties of Ge2Te, GeTe, and Ge in the crystalline, liquid, and amorphous
phases. Then, we exploited the potential to perform large-scale simulations of
the crystallization and phase separation processes.

##  II Methods

We generated the NN potential for Ge-rich GexTe alloys by fitting a database
of DFT energies, forces, and virial tensors of 115000 atomic configurations of
Ge, GeTe and Ge2Te in the amorphous, liquid and crystalline phases within the
framework implemented in the DeePMD code Zhang  _et al._ (2018); Zeng  _et
al._ (2023); Wang  _et al._ (2018). Amorphous and liquid Ge, GeTe and Ge2Te
were modeled in cells of 144, 108 and 100 atoms (composition Ge67Te33),
respectively. The crystalline phases of Ge and GeTe were modeled in cells with
64 and 96 atoms. The atomic configurations were extracted from DFT molecular
dynamics (MD) simulations by using the CP2k code VandeVondele  _et al._
(2005). Kohn-Sham orbitals were expanded in a Triple-Zeta-Valence-plus-
Polarization (TZVP) basis set while the electronic density was expanded in
plane waves up to a kinetic cutoff of 100 Ry. The Brillouin Zone (BZ)
integration was restricted to the $\Gamma$-point in all MD simulations. The
time step was set to 2 fs and configurations were extracted every 100-200 fs.
The number of configurations for each phase and composition is given in Table
1. We recalculated energies, forces, and stresses for each configuration added
to the database with a higher accuracy by increasing the kinetic energy cutoff
for the plane-waves expansion of the electronic density to 400 Ry and by using
a 4x4x4 k-point mesh for the BZ integration. A Fermi-Dirac smearing in the
occupation of Kohn-Sham states was used with an electronic temperature of 300
K.

Phase | Number of configurations
---|---
Ge l & a | 27000
GeTe l & a | 35000
Ge2Te l & a | 46000
Crystalline Ge | 3400
Crystalline GeTe | 3600
Table 1: Number of configurations included in the database for each
composition and phase. l & a stands for liquid and amorphous phases.

In NN schemes for the generation of interatomic potentials, the total energy
of the system is written as the sum of individual atomic energies that depend
on the local environment of each atom. In the DeePMD scheme the local
environment is encoded by local descriptors which are generated by an embedded
neural network. A second neural network (fitting network) is built for the
calculation of energies and forces with the local descriptors as input layer.
We designed the embedded network with three hidden layers of 40, 80, and 160
nodes. The cutoff radius rc was set to 7 Å, which is beyond the third
coordination shell of our systems, while the smoothing radius rs was set to 2
Å Zeng  _et al._ (2023). The maximum number of neighbors was set to 80. In the
embedded network, we also exploited the attention mechanism that was recently
implemented in the DeePMD code Zeng  _et al._ (2023). Finally, the network for
the fitting of energy and forces consists of 3 hidden layers with 320 nodes
each. In the embedding and fitting network, we have used the hyperbolic
tangent as an activation function.

The NN potential was generated in an iterative manner. A first version of the
potential was generated by using a small training database of about 6000
configurations of Ge, GeTe and Ge2Te. Then, the potential is used to perform
MD simulations that provide high-energy configurations to enrich the training
database for a second generation of the potential and at the same time we
enriched the database with additional atomic configurations extracted from
DFT-MD simulations at different densities and exploring a wider temperature
range. We iterated this procedure until the potential reached a satisfactory
level of accuracy.

Structural and dynamical properties obtained from NN simulations were compared
with results from DFT-MD simulations performed with the CP2k code with the
same parameters given above. NN-MD simulations were performed by using the
LAMMPS code as MD driver with the DeePMD plugin Thompson  _et al._ (2022).

To assess whether an atom is crystalline in the simulations of the
crystallization process, we used the Steinhardt order parameter
Q${}^{dot}_{n}$ Mickel  _et al._ (2013). In general, the Q${}^{dot}_{n}$
parameter of order $n$ is defined for each atom $i$ by

$Q^{dot}_{n}=\frac{1}{N_{i}}\sqrt{\sum_{j=1}^{N_{i}}\sum_{m=-n}^{n}q_{nm,i}q^{*}_{nm,j}}$
(1)

$q_{nm,i}=\frac{1}{N_{i}}\sum_{j=1}^{N_{i}}Y_{nm}\left(\mathbf{\hat{r}}_{ij}\right),$
(2)

where $N_{i}$ in the number of neighbors of atom $i$ up to a given cutoff, $j$
is the neighbors index, $Y_{nm}$ is the $n$ order spherical harmonic with
degree $m$, and $\mathbf{\hat{r}}_{ij}$ is the unit vector connecting the two
atoms. For the crystallization of GeTe, we considered the Q${}^{dot}_{4}$
order parameter, but also the Q${}^{dot}_{6}$ has been considered for some
analysis as it will be discussed in the relevant section. The distribution of
Q${}^{dot}_{4}$ in crystalline and amorphous phases of GeTe, shown in Fig. S1
in the Supplementary Information, suggests a threshold of Q${}^{dot}_{4}$ =
0.87 for an atom to be crystalline.

To assess whether a Ge atom is segregated in regions of amorphous Ge (a-Ge),
we calculated the SOAP (Smooth Overlap of Atomic Position) similarity kernel
$k_{j}$ for each atom Bartók  _et al._ (2013). $k_{j}$ quantifies the
similarity of the atomic environment around atom $j$ with a reference atomic
environment which is taken here as the average environment in a-Ge at 600 K.
$k_{j}$ is a number that ranges from 0, when the atomic environment is totally
different from that of the reference system, to 1 when it is identical. In the
SOAP formalism the local atomic density around each atom $j$ is expressed as a
sum of Gaussian functions (here with broadening $\sigma_{at}$= 0.3 Å) centered
on the position of its neighbors up to a given cutoff (9 Å here). Then the
density around atom $j$ is expanded in spherical harmonics and radial basis
functions $g_{b}(|\mathbf{r}|)$ as
$\rho_{j}(\mathbf{r})=\sum_{blm}c_{blm}g_{b}(|\mathbf{r}|)Y_{lm}(\mathbf{\hat{r}})$.
The coefficients of this expansion define the so-called power spectrum matrix
$p(j)_{b_{1}b_{2}l}=\pi\sqrt{8/(2l+1)}\sum_{m}(c_{b_{1}lm})^{*}c_{b_{2}lm}$
whose elements are turned into a vector $\mathbf{p}_{j}$ from which the SOAP
kernel kj is calculated as
$k_{j}=(\mathbf{p}_{j}/|\mathbf{p}_{j}|\cdot\mathbf{p}_{ref}/|\mathbf{p}_{ref}|)^{\xi}$,
where $\mathbf{p}_{ref}$ is the average power spectrum of the atoms in the
atomic environment of the reference system and $\xi$ is an integer set to 2 in
our case. The SOAP similarity kernel has been used in a variety of studies to
analyze the atomic structure De  _et al._ (2016); Musil  _et al._ (2021) and
very recently to discriminate between the crystalline and amorphous/liquid
phases in GST alloys Xu  _et al._ (2022). In this work, we calculated the SOAP
kernel using the DScribe Python package Laakso  _et al._ (2023); Himanen  _et
al._ (2020) within the ASE Python library Larsen  _et al._ (2017).

##  III Results and Discussion

###  III.1 Validation of the neural network potential

The accuracy of the NN potential is assessed by the cumulative distribution of
the error on energies, forces and virial tensors shown in Fig. 1 for Ge2Te,
GeTe and Ge. For 80% of the configurations, the error on the energies, forces
and virial is smaller than 8 meV/atom, 150 meV/Å and 25 meV/atom. In
particular, among our systems, GeTe has the largest errors, while Ge2Te and Ge
have significantly lower errors. Overall, the root mean square error (RMSE) on
the energy is 4.4 meV/atom, on forces is 105 meV/Å and on virial is 14.5
meV/atom. We remark that the typical average error obtained with DeePMD for
highly disordered phases of multicomponent systems like ours (i.e liquid
and/or amorphous phases) are in the range 2-7 meV/atom and 90-145 meV/ÅUrata
(2022); Zhai and Wang (2023); Liu  _et al._ (2023); Huang  _et al._ (2023).

Figure 1: Cumulative fraction of error on energies, forces and virial of the
NN potential in Ge, GeTe and Ge2Te.

The NN potential has then been validated on the properties of the liquid,
amorphous and crystalline phases as described in the separate sections below
for pure Ge, stoichiometric GeTe and Ge2Te.

####  III.1.1 The liquid phase

We computed the structural properties of liquid Ge at 1250 K from NN-MD
simulations in a 2400-atom model and from DFT-MD simulations in a 300-atom
model. In the perspective to study later the amorphous phase, we used for the
liquid the experimental density of the amorphous phase of 0.0438 atom/Å3
Laaziri  _et al._ (1995). The models were equilibrated first at 2000 K for 10
ps and then at 1250 K for 20 ps. The structural properties were evaluated from
the last 10 ps of the NVT run at 1250 K. The pair correlation function, the
distribution of the coordination numbers and the bond angle distribution
function of liquid Ge from NN and DFT simulations are compared in Fig. 2.

Figure 2: Structural properties of liquid Ge at 1250 K from NN-MD (black
lines, 2400 atoms) and DFT-MD (red lines, 300 atoms). The experimental density
of the amorphous phase of 0.0438 atom/Å3 is used. a) Pair correlation
function. b) Distribution of the coordination numbers by assuming a bonding
cutoff of 3 Å. Ge is primarily 4-fold and 5-fold coordinated, with a minority
of 3-fold and 6-fold coordinated atoms. The average coordination number of Ge
is 4.28 with the NN potential and 4.31 with DFT. c) Bond angle distribution
function. The peak around 100∘ is a feature of the Ge atoms in a tetrahedral
environment, while the peak at 60∘ highlights overcoordinated 5-fold and
6-fold atoms.

The structural properties of liquid GeTe (4096-atom model) and of liquid Ge2Te
(2400-atom model) from NN-MD simulations were compared with DFT-MD results for
300-atom models for both compositions. For GeTe, we used the experimental
equilibrium density of the liquid phase at 1150 K of 0.03294 atom/Å3 Weber
_et al._ (2017). Since no experimental data are available for the density of
liquid Ge2Te, and in the perspective to study the amorphous phase, we
generated an amorphous model with DFT molecular dynamics by quenching from
1200 K to 300 K in 100 ps in the NPT ensemble which yielded an equilibrium
density at 300 K of 0.0355 atom/Å3. To this aim, we used the Grimme (D3)
semiempirical potential Grimme  _et al._ (2010) for van der Waals (vdW)
interactions which is needed to avoid the coalescence of nanovoids in the
liquid phase as discussed in Ref. Sosso  _et al._ (2012b). Liquid Ge2Te was
then simulated without vdW corrections at the same density of the amorphous
phase of 0.0355 atoms/Å3. The liquid models were equilibrated first at 2000 K
for 10 ps and then for 20 ps at 1150 K for GeTe and at 1200 K for Ge2Te, which
are well above the melting temperature (990 K for the melting of GeTe Klemm
and Frischmuth (1934) and for the incongruent melting of Ge2Te, the liquidus
temperature of Ge2Te is instead 1150 K Schlieper  _et al._ (1999)). Structural
properties were calculated in the last 10 ps. The pair correlation functions,
the distribution of the coordination numbers and the bond angle distribution
of liquid Ge2Te and GeTe from NN and DFT simulations are compared in Fig. 3.
The bonding cutoff of 3 Å for Ge-Ge, 3.22 Å for Ge-Te and 3 Å for Te-Te have
been used. The average partial coordination numbers are shown in Table 2.

Overall, the NN potential reproduces very well the structural properties of
the liquid phase for all the three compositions, i.e. pure Ge, GeTe and Ge2Te.

Figure 3: Structural properties of liquid Ge2Te at 1200 K (left panels) and liquid GeTe at 1150 K (right panels) from DFT and NN simulations. a) Partial pair correlation functions. b) Bond angle distribution. c) Distribution of coordination numbers.  |  | Ge | Te | Total
---|---|---|---|---
Ge2Te | Ge | 2.54 (2.53) | 1.50 (1.51) | 4.04 (4.04)
| Te | 3.05 (3.07) | 0.07 (0.07) | 3.12 (3.14)
GeTe | Ge | 1.25 (1.31) | 2.73 (2.63) | 3.98 (3.94)
| Te | 2.73 (2.63) | 0.16 (0.19) | 2.89 (2.82)
Table 2: Average partial coordination numbers in liquid Ge2Te at 1200 K and
liquid GeTe at 1150 K. DFT data are reported in parenthesis.

Regarding the dynamical properties, we computed the diffusion coefficient
($D$) in liquid GeTe and Ge2Te from the mean square displacement (MSD) and the
Einstein relation $MSD$=$6Dt$ in 30 ps simulations as given in Table 3.

|  | D (10-5 cm2/s)
---|---|---
|  | Ge | Te
Ge2Te | 1100 K | 9.5 (8.7) | 6.1 (5.9)
GeTe | 1120 K | 8.0 (6.3) | 6.0 (5.0)
Table 3: Diffusion coefficients in liquid Ge2Te at 1100 K and liquid GeTe at
1120 K from NN and DFT (in parenthesis) simulations.

Both in Ge2Te and GeTe, the NN potential slightly overestimates the diffusion
coefficients by about 10%. Part of this discrepancy may come from the use of
different simulation cells in the NN and DFT calculations. In fact, the
following scaling of $D$ with the size $L$ of the cell is expected Dünweg and
Kremer (1993)

$D(L)=D_{\infty}-\frac{2.387k_{B}T}{6\pi\eta L},$ (3)

where $\eta$ is the viscosity. For $\eta$ of about 1 mPa$\cdot$s (see the NN
results at 1000 K of Ref. Sosso  _et al._ (2012b)) the correction to $D$ for
liquid GeTe at 1120 K is about 0.9 10-5 cm2/s for a 300-atom cell ($L$=20.88
Å) and 0.375 10-5 cm2/s for a 4096-atom cell ($L$=49.91 Å).

The diffusion coefficient (total and resolved per species) has then been
computed on a wide temperature range in the liquid and supercooled liquid
phases from NN simulations, 50 ps long at each temperature, as shown in Fig.
4. The diffusion coefficient of Ge2Te is fitted by a Cohen-Grest function
Cohen and Grest (1979)
$\log(D(T))=A-{2B}/{(T-T_{0}+[(T-T_{0})^{2}+4CT]^{1/2})}$ which is suitable
for the supercooled phase of fragile liquids. According to Angell Angell
(1995), a supercooled liquid is classified as fragile if the viscosity $\eta$
remains very low down to temperatures close to the glass transition
temperature Tg, where a steep, super-Arrhenius behavior is observed up to the
high value of $\eta$ expected at Tg. A similar super-Arrhenius behavior is
expected for the diffusion coefficient $D$ given its inverse proportionality
with the viscosity $\eta$ according to the Stokes-Einstein relation. On the
contrary the viscosity of an ideal strong liquid follows an Arrhenius behavior
from the melting down to Tg. The fitting parameters for Ge2Te amount to $A$=
-3.52 (-3.49), $B$= 435 K (413 K), $C$= 0.36 K (0.52 K) and $T_{0}$= 315 K
(321 K) with the values in parenthesis referring to the diffusion coefficient
of Ge atoms only. For GeTe, the diffusion coefficient follows an Arrhenius
behaviour D = D0e${}^{-E_{a}/(k_{B}T)}$ down to 500 K as shown in Fig. 4 and
in agreement with previous NN results in Ref. Sosso  _et al._ (2012b), due to
the breakdown of the Stokes-Einstein relation at low temperatures. The fitting
parameters for GeTe are Ea = 0.34 eV and D0 = 2.47 $\cdot$ 10-3 cm2/s. The
super-Arrhenius behaviour in $D$ for GeTe would appear at lower temperatures,
below 500 K, which implies a higher fragility of GeTe with respect to Ge2Te
considering that Tg is expected to be higher for Ge2Te than for GeTe. The
reported glass transition temperature for pure germanium of about 750 K
Imabayashi and Ishimaru (2017) is indeed higher than the value of 423 K
reported for GeTe Chen  _et al._ (2017). For the sake of comparison the
diffusion coefficient of pure Ge as a function of temperature is given in Fig.
S2 in the Supplementary Information.

Figure 4: Diffusion coefficient ($D$) as a function of temperature from NN
simulations of GeTe and Ge2Te for (upper panels) Ge atoms only, (central
panels) Te atoms only and (bottom panels) the total values averaged over all
atoms. The continuous lines are Cohen-Grest fit (see text) of the data for
Ge2Te and a simple Arrhenius fit for GeTe.

####  III.1.2 The amorphous phase

An amorphous model of Ge was generated by cooling the liquid model from 1250 K
to 300 K in 140 ps in both the DFT and NN simulations at the experimental
density of a-Ge of 0.0438 atoms/Å3 Laaziri  _et al._ (1995). Structural
properties were averaged over the last 10 ps of an equilibration run at 300 K
lasting 20 ps. NN and DFT results for the pair correlation function, the
distribution of the coordination numbers and the bond angle distribution
function are compared in Fig.5.

Figure 5: Structural properties of amorphous Ge at 300 K from NN and DFT
simulations. a) Pair correlation function. b) Distribution of coordination
numbers by assuming a bond cutoff of 3 Å. c) Bond angle distribution function.

Amorphous models of Ge2Te (a-Ge2Te) and GeTe (a-GeTe) were obtained by
quenching the liquid to 300 K in 100 ps starting from 1200 K for Ge2Te and
from 1150 K for GeTe in both NN and DFT simulations. The pair correlation
functions, the distribution of the coordination numbers, and the bond angle
distribution functions are compared in Fig.6, while the average partial
coordination numbers are given in Table 4. The bonding cutoffs are the same
used for the liquid as given above.

The agreement between the NN and DFT results is overall very good for Ge, GeTe
and Ge2Te. We just note a slight misfit in the number of 2-fold Te atoms in
a-Ge2Te and in the position and intensity of the peak of the Ge-Ge correlation
function in a-GeTe, possibly due in part to the variability from model to
model in small 300-atom cells as discussed in a previous work on
stoichiometric GeTe Sosso  _et al._ (2012a).

Figure 6: Structural properties of amorphous Ge2Te (left panels) and GeTe (right panels) at 300 K from DFT and NN simulations. a) Partial pair distribution functions. b) Bond angle distribution functions. c) Distribution of the coordination numbers. d) Distribution of the $q$ order parameter for tetrahedricity, also for pure Ge. |  | Ge | Te | Total
---|---|---|---|---
Ge2Te | Ge | 2.59 (2.58) | 1.54 (1.46) | 4.13 (4.04)
| Te | 3.12 (2.96) | 0.00 (0.00) | 3.12 (2.96)
GeTe | Ge | 0.79 (0.87) | 3.34 (3.18) | 4.13 (4.04)
| Te | 3.34 (3.18) | 0.00 (0.04) | 3.34 (3.21)
Table 4: Average partial coordination numbers in amorphous Ge2Te and GeTe at
300 K from NN and DFT (in parenthesis) simulations.

Turning now to a comparison of the structural properties of Ge2Te and GeTe, we
notice that by increasing the fraction of Ge from GeTe to Ge2Te, we observe
also an increase in the fraction of Ge-Ge bonds that favor a tetrahedral
environment of Ge Deringer  _et al._ (2014). This feature shows up as an
increased fraction of 4-fold coordinated Ge in the distribution of the
coordination numbers (see Fig. 6) and in the shift to higher angles of the
peak in the bond angle distribution. In general, we expect a coexistence of
tetrahedral, pyramidal (3-fold coordinated) and defective octahedral
coordination (octahedral angles but coordination lower than six) for Ge atoms.
The fraction of Ge atoms in tetrahedral configurations can be quantified by
using the $q$-parameter for tetrahedricity introduced in Ref. Errington and
Debenedetti (2001) and defined by
$q_{j}=1-\frac{3}{8}\Sigma_{i<k}(\frac{1}{3}+cos(\theta_{ijk}))^{2}$, where
the sum runs over the couples of atoms bonded to a central atom $j$ and
forming a bonding angle $\theta_{ijk}$. The order parameter evaluates to $q$=1
for the ideal tetrahedral geometry, to $q$=0 for the 6-fold coordinated
octahedral site, to $q$=5/8 for a 4-fold coordinated defective octahedral
site, and $q$=7/8 for a pyramidal geometry. The distribution of the
$q$-parameter for 4-coordinated Ge atoms is shown in Fig. 6d for amorphous
Ge2Te, GeTe and Ge at 300 K. The $q$-parameter features two peaks in a-GeTe
due to the coexistence of tetrahedral and defective octahedral configurations
Mazzarello  _et al._ (2010). As discussed in Ref. Spreafico  _et al._ (2011),
the integration of the $q$-parameter in the range 0.8-1.0 gives a good measure
of the fraction of tetrahedral Ge which is 24% in GeTe, similar to previous
DFT works Akola and Jones (2007); Mazzarello  _et al._ (2010). On the other
hand, a-Ge2Te features a single peak centered at the position corresponding to
the tetrahedral configuration (as in a-Ge). The peak is, however, very broad
which suggests that also in Ge2Te a percentage of Ge atoms are in defective
octahedral configurations. By integrating the $q$-distribution for values
greater than 0.8, one obtains a percentage of tetrahedral Ge of approximately
52% in a-Ge2Te , which is in between the values of 24% and 68% for a-GeTe and
a-Ge. Analogously to the liquid phase, the agreement between the NN and DFT
results for the amorphous phase is overall very good for Ge, GeTe and Ge2Te.

####  III.1.3 The crystalline phase

We calculated the energy of crystalline cubic Ge and of trigonal
($\alpha$-phase) GeTe at different volumes and we fitted the energy as a
function of volume by using the Birch-Murnaghan equation of state as shown in
Fig. S3 in the Supplementary Information. The resulting fitting parameters
from NN and DFT calculations are compared in Table 5 while the equilibrium
structural parameters of $\alpha$-GeTe (space group $R3m$)Goldak  _et al._
(1966) are given in Table 6. Crystalline $\alpha$-GeTe, with two atoms per
unit cell, can be viewed as a distorted rocksalt geometry with an elongation
of the cube diagonal along the [111] direction and an off-center displacement
of the inner Te atom along the [111] direction giving rise to a 3+3
coordination of Ge with three short and stronger bonds and three long and
weaker bonds (see Table 6). In the conventional hexagonal unit cell of the
trigonal phase, the structure can be also seen as an arrangement of GeTe
bilayers along the $c$ direction with shorter intrabilayer bonds and longer
interbilayers bonds. We mention that the trigonal ferroelectric phase
transforms into the cubic paraelectric phase ($\beta$-phase, space group
Fm$\bar{3}$m) above the Curie temperature of 705 K Chattopadhyay  _et al._
(1987). In the cubic phase, the alternation of long and short bonds survives
in a disordered manner along all equivalent $<$111$>$ directions as revealed
by extended x-ray absorption fine structure (EXAFS), x-ray total diffraction
measurements Fons  _et al._ (2010); Matsunaga  _et al._ (2011) and MD
simulations Xu  _et al._ (2018). However, more recent molecular dynamics
simulations Dangić  _et al._ (2022) suggest that the order-disorder character
of the phase transition is weaker than as inferred from EXAFS data. The
$\beta$-phase is the structure the a-GeTe crystallizes into at the operation
conditions of the memory devices.

|  | cubic Ge |  |
---|---|---|---|---
Method | E0 (eV/atom) | V0 (Å3/atom) | B (GPa) | B′
NN | -106.972 | 24.22 | 56.99 | 4.56
DFT | -106.972 | 24.22 | 57.58 | 4.64
|  |  $\alpha$-GeTe |  |
Method | E0 (eV/atom) | V0 (Å3/atom) | B (GPa) | B′
NN | -164.671 | 28.33 | 29.27 | 8.42
DFT | -164.671 | 28.34 | 28.88 | 8.55
Table 5: Fitting parameters of the Birch–Murnaghan equation of state of cubic Ge and trigonal GeTe. V0 is the equilibrium volume, E0 the equilibrium energy, B0 the bulk modulus and B′ its derivative respect to pressure. Method |  $a$ (Å) | $\alpha$ |  $V$ (Å3) | $x$ |  $d_{short}$ (Å) |  $d_{long}$ (Å)
---|---|---|---|---|---|---
NN | 4.42 | 56.95° | 57.00 | 0.2340 | 2.86 | 3.29
DFT | 4.33 | 58.14° | 55.00 | 0.2358 | 2.85 | 3.21
Table 6: Equilibrium structural parameters of the $\alpha$-phase of GeTe in
the rhombohedral setting from NN and DFT calculations. $a$ is the lattice
parameter, $\alpha$ is the angle of the trigonal cell, $V$ is the unit cell
volume (two atoms), $x$ assigns the position of Ge atom at $(x,x,x)$ and Te
atom at $(-x,-x,-x)$ in crystallographic units. $d_{short}$ and $d_{long}$ are
the lengths of the short and long Ge-Te bonds.

In summary, although the RMSE for the energies and forces are not very small
(4.4 meV/atom and 105 meV/Å), albeit similar to other NN potentials in
literature for disordered multi-component materials, the validation of the
potential over the properties of liquid, amorphous and crystalline phases is
excellent. Overall, we judge that our potential is sufficiently accurate to
address the study of the crystallization process as discussed in the next
section.

###  III.2 Simulation of the crystallization process

We exploited the NN potential discussed above to perform simulations of the
crystallization process in both GeTe and Ge2Te. We first simulated
stoichiometric GeTe to compare the results with previous NN simulations in
Ref. Sosso  _et al._ (2013).

####  III.2.1 GeTe

We generated a 32768-atom model of a-GeTe at the experimental amorphous
density of 0.0333 atom/Å3 Ghezzi  _et al._ (2011). The model was equilibrated
first at 2000 K for 10 ps and at 1150 K for 40 ps, and it was then cooled to
either 600 or 500 K in 60 ps. Finally, we performed a NVT simulation at the
two target temperatures. The number of crystalline atoms as a function of time
at 600 and 500 K are compared in Fig. 7a with the results from previous NN
simulations Sosso  _et al._ (2013) at the same temperatures and with 4096
atoms. At 600 K, our data are very similar to previous results Sosso  _et al._
(2013) with a fraction of crystallized atoms of almost 75 % after 2 ns. Some
small crystalline nuclei of cubic GeTe form and begin to grow after the first
100 ps, as shown in the snapshot of Fig. 7b from the simulation at 600 K. At
500 K, the number of crystalline atoms increases faster with time in our
simulation than in the previous one from Ref. Sosso  _et al._ (2013), because
our system is eight times larger than that of Ref. Sosso  _et al._ (2013) and
then the nucleation time (time needed for the formation of an overcritical
crystalline nucleus) is shorter. We can thus conclude that our potential
reproduces the crystallization kinetics in a manner very similar to the NN
potential of Ref. Sosso  _et al._ (2013) which was previously used to address
several details of the crystallization process in GeTe Sosso  _et al._ (2013);
Gabardi  _et al._ (2017, 2019); Perego  _et al._ (2023); Sosso  _et al._
(2015); Acharya  _et al._ (2024).

Figure 7: a) Fraction of crystalline atoms (defined by Q${}^{dot}_{4}$) as a
function of time in a 32768-atom model of GeTe in simulations at 500 and 600
K, compared with previous results (4096 atoms) from Ref. Sosso  _et al._
(2013). b) Snapshot of the crystalline atoms in the 32768-atom model of GeTe
after 100 ps at 600 K.

Regarding the mechanism of crystal nucleation and growth, we mention that in a
very recent work Wintersteller  _et al._ (2024) it has been proposed that
crystallization is triggered by a pre-ordering of the face-centered-cubic
(fcc) sublattice of the Te atoms that precedes the ordering of the square
rings of Ge-Te bonds. This picture was inferred from the analysis of the
temperature/time dependence of the Te-Te coordination numbers (first neighbors
on the anionic fcc sublattice) obtained from the fitting of extended x-ray
absorption fine structure (EXAFS) spectrum at the Te K-edge. To check this
picture, we computed the Q${}^{dot}_{6}$ order parameter for the Te-Te
coordination involving only Te atoms and Te-Te contacts at distances up to 4.4
Å which is slightly above the nearest neighbor distance on the Te fcc
sublattice. The Q${}^{dot}_{6}$ for Te atoms only would thus detect the
ordering of the Te fcc sublattice. On the contrary, the Q${}^{dot}_{4}$
parameters involving pairs of atoms up to 3.6 Å measures the ordering of the
simple cubic lattice made of the two fcc sublattices of Ge and Te, i.e. it
detects the ordering of the square rings made of Ge-Te bonds. A pre-ordering
of the Te fcc sublattice should then result into an increase in the number of
crystalline Te atoms assigned by the Q${}^{dot}_{6}$ order parameter before
cristallinity would be detected by the Q${}^{dot}_{4}$ order parameter which
is tuned on the ordering of Ge-Te bonds. The evolution in time of the
crystalline Te atoms assigned by the Q${}^{dot}_{6}$ and Q${}^{dot}_{4}$
parameters is shown in Fig. S4 in the Supplementary Information. No delay is
observed in the onset of the growth of Q${}^{dot}_{4}$ with respect to
Q${}^{dot}_{6}$ which means that there is no pre-ordering of the Te fcc
sublattice. We remark, however, that the experimental analysis in Ref.
Wintersteller  _et al._ (2024) was performed at 440 K on a transformation
occurring on the time scale of hours, while our analysis addresses the
nucleation and growth occurring on the time scale of hundreds of ps at 500-600
K. Therefore, we can not exclude that pre-ordering of the Te fcc sublattice
could take place at the experimental conditions mentioned above, albeit
according to our results this process does not seem to occur at the operation
conditions of PCMs. For our purpose, we considered this validation sufficient
to move to the study of the crystallization in Ge2Te which is discussed
hereafter.

####  III.2.2 Ge2Te

As mentioned in the Introduction, crystallization of the Ge63Te37 alloy
deposited by magnetron sputtering was studied by several means in Ref. Carria
_et al._ (2012). There, it was found that the crystallization process consists
of three separate steps: first, Ge segregates by forming a-Ge regions followed
by the crystallization of GeTe and later of pure Ge. These evidences resulted
from time resolved measurements of optical reflectivity, x-ray diffraction and
Raman spectra during an annealing ramp (0.1 K/s) Carria  _et al._ (2012). The
transformation starts at about 600 K, GeTe crystallization was observed at
about 630 K, while crystalline Ge appears at an even later time at about 650
K. The crystallized sample was also re-amorphized by laser melting and the
resulting amorphous sample crystallizes at the same temperature of the asdep
one. In another work Raoux  _et al._ (2009b), a longer crystallization time of
about one order of magnitude was reported instead for asdep Ge62Te38 sample
than for the melt quenched one. This difference was ascribed to the presence
of an amorphous-crystal interface in the laser melted sample which favors
crystal growth without the need of crystal nucleation Raoux  _et al._ (2009b).
A previous Raman study confirms a multistep crystallization process for
Ge76Te24 and Ge62Te38 Gourvest  _et al._ (2009). GeTe crystallites were
detected at 470 K in Ge62Te38 whereas pure Ge crystallizes at a later time at
650 K Gourvest  _et al._ (2009). A similar behaviour was reported from time-
resolved diffraction measurements in Ref. Raoux  _et al._ (2009a) for Ge62Te38
and in Ref. Navarro  _et al._ (2013b) for Ge62Sb38 and Ge69Sb31. However, the
possibility of the formation of large a-Ge region due to Ge segregation prior
to crystal nucleation of GeTe, was not discussed explicitly in Refs. Gourvest
_et al._ (2009); Raoux  _et al._ (2009a); Navarro  _et al._ (2013b). A
simultaneous crystallization of GeTe and Ge was reported instead at high
temperature of 620 K for Ge70Te30 in Ref. Raoux  _et al._ (2009a) under an
annealing ramp of 1 K/s. The reason behind the discrepancies among different
reports is unclear. Different degree of homogeneity of the alloy under
different preparation conditions or different levels of surface oxidation may
affect the crystallization temperature. Indeed, oxidation was shown to
sizeably lower Tx in GeTe Berthier  _et al._ (2017) and in Ge-rich GST alloys
Agati  _et al._ (2020).

On these premises, we started our analysis by first simulating Ge2Te at 600 K
in a 19200-atom model at the theoretical amorphous density of 0.0355 atom/Å3
(see Sec. III.A.1). By assuming that Ge would segregate as a-Ge at the
experimental density of 0.0438 atom/Å3, a phase separation of Ge2Te at
constant volume and at the average density of 0.0355 atom/Å3 would lead to
GeTe regions at the density of 0.0326 atom/Å3 which is sufficiently close to
the experimental density of a-GeTe of 0.0333 atom/Å3 Ghezzi  _et al._ (2011).
Therefore all simulations were performed at constant volume which is also
closer to the operation conditions of PCMs. The 19200-atom model was
equilibrated first at 2000 K for 10 ps and at 1200 K for 40 ps and then
quenched to 600 K in 60 ps. At this temperature we performed a NVT simulation
lasting 14 ns.

At 600 K, Ge atoms immediately start to segregate from Ge2Te. Phase separation
occurs because the sum of the free energies of GeTe and Ge is lower than that
of Ge2Te in the amorphous of supercooled liquid phases. Indeed, we verified
that the reaction energy of the transformation Ge2Te $\rightarrow$ GeTe + Ge
at 600 K is $\Delta$E = 7 $\pm$ 2 meV/atom, as obtained from 6 independent
simulations of supercooled liquid models of these compositions at 600 K.
Although we did not attempt to estimate the reaction entropy, the results
above confirm that Ge2Te is unstable against phase separation. To quantify the
fraction of segregated Ge atoms, we computed the SOAP similarity kernel
$k_{j}$ for a-Ge (SOAP $k_{j}$) as described in Sec. II. SOAP $k_{j}$ for a Ge
atom is close to one for configurations similar to the average local
configuration in a-Ge and it decreases by increasing the dissimilarity with
a-Ge. The distribution of $k_{j}$ at the beginning (t = 0 ns) and at the end
(t = 14 ns) of the simulation is compared in Fig. 8a to the distribution
$k_{j}$ of a-Ge at the same temperature. We could qualitatively assess that a
Ge atom is segregated when its SOAP $k_{j}$ is higher than 0.92, as the
distribution of $k_{j}$ in a-Ge ranges from 0.92 to 1. In the simulation of
Ge2Te, the $k_{j}$ distribution initially features a single broad peak at
around 0.9, while after 14 ns two peaks are present at 0.6 and at 0.95 which
highlights the occurrence of a phase separation. As shown in Fig. 8b the
number of segregated Ge atoms increases initially very fast and then it levels
off at about 40 $\%$ after 10 ns. Complete phase separation into Ge and GeTe
would correspond to a fraction of 50 $\%$ of segregated Ge atoms. The
segregation of Ge is very clear in the snapshots reported in Fig. 8c where we
show only Ge atoms in an environment close to a-Ge ($k_{j}$ $>$ 0.92) at the
beginning of the simulation and after 14 ns.

Segregation of Ge leads to the formation of Ge-rich regions with Ge88Te12
composition and average density of 0.042 atom/Å3, and of Ge-poor regions with
Ge53Te47 composition and average density of 0.0327 atom/Å3 which is very close
to the experimental density of the amorphous phase of stoichiometric GeTe
(0.0333 atom/Å3). The composition of the Ge-poor region was computing by
averaging over about 4000 atoms in a region at least 10 Å far from the region
of segregated Ge (see Fig. 8). Similarly, the composition of the Ge-rich
region was computed by averaging over about 3000 atoms in the lower part of
the a-Ge-like region (see Fig. 8).

We also performed a second simulation starting from an independent 19200-atom
model at 600 K with very similar results. The snapshot of the segregated atoms
at the beginning and the end (14 ns) of this second simulation is shown in
Fig. S5 in the Supplementary Information. In the process of formation of the
a-Ge regions, Ge atoms is the most diffusing specie as one would envisage from
the diffusion coefficients in Fig. 4. The overall average MSD over 10 ns for
each species is indeed very close to the value given by the Einstein relation.
At 600 K the square root of the average MSD after 10 ns is about 8.5 nm for Ge
atoms and about 5.9 nm for Te atoms. We remark that the edge of the simulation
cell is 8.1 nm.

Figure 8: a) Distribution of the SOAP $k_{j}$ at the beginning (top, t = 0 ns
) and after 14 ns (center) in the simulation of Ge2Te (19200 atoms) at 600 K
compared to the distribution of an independent model of a-Ge at the same
temperature (bottom). b) Number of segregated Ge atoms ($k_{j}$ $>$ 0.92) as a
function of time in the 19200-atom model of Ge2Te at 600 K (red) and 500 K
(blue). c) Snapshots of Ge2Te at the beginning (left) and after 14 ns (right)
at 600 K. Only atoms with an environment close to that of a-Ge are shown (SOAP
$k_{j}$ $>$ 0.92).

Although most of Ge atoms in excess is segregated, we do not see the formation
of crystalline nuclei in the Ge53Te47 region in the simulations lasting 14 ns.
We can conclude that at this temperature the nucleation time in a 19200-atom
cell is longer than 14 ns. Raising the temperature from 600 K to 630 K for a 4
ns does not lead to further segregation of Ge nor to crystal nucleation. In
the attempt to reduce the crystal nucleation time, we simulated a larger
30000-atom model prepared with the same protocol of the smaller ones. We chose
an orthorhombic cell with edges a=b= 101.82 Å and c= 81.45 Å. Phase separation
occurs in a similar manner as shown in Fig. S5 in the Supplementary
Information. In the larger model, however, an overcritical crystalline nucleus
of GeTe appears at 600 K after 6 ns and it grows until the crystal percolates
through the simulation box as shown in Fig. 9.

Once formed, the overcritical crystalline nucleus grows, but at a lower rate
compared to stoichiometric GeTe at the same temperature. At 600 K, we
calculated the crystal growth velocity $v_{g}$=$dR/dt$, where $R$ is the
radius of the crystalline nucleus given in turn by
$R\left(t\right)=\left({3N(t)}/4\pi\rho\right)^{\frac{1}{3}}$, where $N$ is
the number of atoms in the nucleus and $\rho$ is the density of the
crystalline phase (0.0351 atoms/Å3). This assumption is valid only in the
early stage of crystallization when the nuclei do not interact with each other
or with their periodic image. As an example, the evolution of $R(t)$ as a
function of time is shown Fig. S6 in the Supplementary Information for a
single crystalline nucleus in GeTe and Ge2Te (30000-atom cell) at 600 K. In
GeTe, we computed $v_{g}$ by averaging over three/(four) crystalline nuclei at
600 K/(500 K). The resulting vg are reported in Table 7. The crystal growth
velocity is lower in Ge2Te because a further segregation of Ge must take place
during the crystal growth in region with average composition of Ge53Te47. The
average composition of the crystallized region is in fact Ge50Te50.

Figure 9:  a) Crystalline atoms in the segregated large model (30000 atoms) of Ge2Te at t = 6 ns and at t= 10 ns at 600 K. b) Fraction of crystalline atoms (defined by Q${}^{dot}_{4}$, black line) and number of segregated Ge atoms ($k_{j}>$ 0.92, red line) as a function of time. |  $v_{\rm g}$ (m/s)
---|---
| 500 K | 600 K
GeTe (homogeneous) | 1.2 | 4.8
Ge2Te (phase separated, large) |  | 3.1
Ge2Te (phase separated, small) | 0.8 |
Ge2Te (homogeneous, small) | 0.2 |
Ge53Te47 (homogeneous) | 0.8 |
Table 7: Crystal growth velocities vg of overcritical nuclei at 500 K or 600 K
in the different models. Homogeneous GeTe; phase separated Ge2Te (large
30000-atom cell); Ge2Te (small 19200-atom cell) phase separated at 600 K,
quenched at 500 K to promote crystal nucleation; homogeneous Ge2Te before Ge
segregation (small 19200-atom cell); homogeneous Ge53Te47.

We also attempted to reduce the crystal nucleation time in the 19200-atom
smaller segregated model by quenching to 500 K. Indeed, at this temperature,
we see the formation of three crystalline nuclei after some hundreds of ps in
the Ge53Te47 region, as shown in Fig. S7 in the Supplementary Information.
This occurs at positions where the local Ge-Ge coordination is lower than the
average in a-GeTe (0.79). The strong shortening of the nucleation time at 500
K, once the system was mostly phase separated at 600 K, might be due to the
smaller size of the critical nucleus at lower temperatures. Since the local Ge
content must be lower than the average composition in the Ge-poor region (at
incomplete phase separation) to promote crystal nucleation, a smaller region
with lower Ge content is needed if the critical nucleus is smaller, which
would reduce the time needed for the sizable fluctuation in the local Ge
content to occur, and consequently also the nucleation time is reduced. The
crystal growth velocity vg at 500 K for this system (phase separated at 600 K
and then quenched at 500 K) is shown in Table 7. We remark that the a-Ge
regions remain amorphous after crystallization of GeTe at 500 and 600 K. The
same crystallization behavior of the system phase separated at 600 K and then
quenched at 500 K, was observed in the second simulation with an independent
model mentioned above. We also repeated a simulation at 500 K for a
homogeneous system (19200 atom) with composition Ge53Te47 generated by
quenching from the melt and we similarly observed the formation of 5
crystalline nuclei after 0.5 ns. The crystal growth velocity is the same of
that obtained for the phase separated system, as given in Table 7.

Finally, to assess if the system might undergo phase separation also at lower
temperatures, we quenched the liquid model from 1200 K to 500 K in 70 ps and
we then performed a NVT simulation at this temperature for 14 ns. The
percentage of segregated Ge atoms saturates around 6 % after about 5 ns, as
shown in Fig. 8b and no crystallization was observed. We can conceive that the
system does not segregate in Ge2Te at 500 K due to the lower atomic mobility
of Ge atoms at this temperature. The diffusion coefficient as a function of
temperature, shown in Fig. 4, features a super-Arrhenius behavior with a drop
from $1.2\cdot 10^{-5}$ cm2/s at 600 K to $2\cdot 10^{-6}$ cm2/s at 500 K.
However, the slow down with temperature of the kinetics of Ge segregation is
not simply due to the decrease of the growth velocity of the a-Ge regions.
Would this be the case, we should see an increase in the growth time by about
a factor six at 500 K with respect to 600 K (actually cube root of six if we
consider an Avrami theory with just the growth of pre-formed nuclei). In fact,
while at 600 K we see the formation of several overcritical nuclei of a-Ge
like in the time span of few hundreds ps, at 500 K no overcritical nuclei of
a-Ge form in 15 ns. Therefore nucleation of a-Ge region seems to control the
kinetic of phase separation.

However, when we repeated the simulation at 500 K for the second independent
model (19200 atoms) generated in the same manner of the first, a surprise
came. After 12 ns, we observed the formation of a single overcritical nucleus
of crystalline GeTe that grew very slowly by expelling the excess of Ge.
Surprisingly, crystal nucleation occurred without prior phase segregation in a
region where the local content of Ge-Ge bonds is lower than the average in the
Ge2Te homogeneous model and also lower than the average of homogeneous GeTe. A
snapshot of the distribution of Ge atoms with a fraction of Ge-Ge bonds lower
than the average in a-GeTe (0.79) is shown in Fig. S8 in the Supplementary
Information. The crystal growth velocity is low as shown in Table 7 because of
the need of expel a large amount of Ge. The evolution in time of the number of
crystalline atoms and of segregated Ge atoms in the simulation at 500 K of the
second 19200-atom model are compared in Fig. 10. A similar behaviour was
observed for a larger 30000-atom model quenched from the melt at 500 K in
which crystal nucleation occurs before phase separation, albeit after about 4
ns (Fig. 10).

Figure 10: Evolution in time of the number of crystalline atoms and of
segregated Ge atoms at 500 K in the simulation of (upper panel) the second
19200-atom model and (lower panel) the larger 30000-atom model. In both models
crystal nucleation occurs before phase separation.

The dependence of the nucleation time on the different simulation conditions
can be seen in Fig. 11 where the internal energy is reported as a function of
time. A sudden drop of the internal energy corresponds to the onset of
crystallization in GeTe and in Ge2Te at 600 K. At 500 K the crystal growth is
slow and the decrease in energy due to crystallization superimposes to the
effect due to equilibration of the supercooled liquid which is slow as well
for the off-stoichiometric system at this temperature.

We have also studied the crystallization at 500 K and 600 K of an amorphous
30000-atom model generated by quenching from 1200 K to 300 K in 100 ps and
then annealed abruptly at the two target temperatures. The results are similar
to those of the liquid model supercooled from the melt discussed above. The
overheated amorphous phase crystallizes concurrently with the phase separation
at 500 K, while phase segregation precedes crystal nucleation at 600 K. The
evolution in time of the fraction of segregated atoms and of the fraction of
crystalline atoms for this last model at 500 and 600 K are shown in Fig. S9 in
the Supplementary Information.

Figure 11: Evolution in time of the internal energy per atom for a)
stoichiometric GeTe at 600 K (32768-atom cell), b) Ge2Te (30000-atom cell) at
600 K, c) Ge2Te in the small (19200-atom, red line) and large (30000-atom,
black line) models at 500 K. Vertical lines in panel c) indicate the onset of
crystallization in the two models.

Therefore, we have seen two different regimes at different temperatures. At
high temperature (600 K), the nucleation time of a-Ge is very short and phase
separation with segregation of most of Ge in excess occurs in a few ns. This
yields a GeTe-like region with composition Ge53Te47 where the nucleation time
for crystalline GeTe is a few ns in a 30000-atom cell, to be compared with the
value of a few hundreds of ps in stoichiometric GeTe at the same temperature
and for a similar cell size (32768 atoms). At 600 K, phase separation occurs
first and GeTe crystallization follows in a two steps mechanism. At lower
temperatures (500 K), phase separation does not occur on the time scale of our
simulations because the overcritical nuclei of a-Ge do not form, albeit Ge is
still very mobile. Still, nucleation of crystalline GeTe is possible at 500 K
without phase separation, but with a low crystal growth velocity because Ge in
excess must be expelled during growth. At 500 K, crystallization is a single
step process with Ge segregation occurring simultaneously with the growth of
GeTe crystallites.

We attempt now to make a contact with the experimental results of Refs Carria
_et al._ (2012); Gourvest  _et al._ (2009); Raoux  _et al._ (2009a). The
single step crystallization that we see at 500 K seems consistent with the
crystallization seen in Refs. Gourvest  _et al._ (2009); Raoux  _et al._
(2009a) at 470 K. We must remark, however, that our models are generated by
quenching from the melt which corresponds to the operation conditions of the
memory devices, while in Refs. Gourvest  _et al._ (2009); Raoux  _et al._
(2009a) crystallization was studied experimentally in the asdep samples. On
the other hand, the two steps crystallization we see at 600 K, seems
consistent with the two steps process (first segregation of Ge and then
crystallization of GeTe) seen in Ref Carria  _et al._ (2012) at 600-630 K. Our
simulations can not, however, explain why a single step crystallization at a
lower temperature about 500 K was not seen in Ref. Carria  _et al._ (2012). We
can not exclude that our NN potential could overestimate the nucleation rate
of GeTe, as it seems to occur Lee  _et al._ (2020) for the NN potential for
stoichiometric GeTe that we developed previously Sosso  _et al._ (2012a).

##  IV Conclusions

In summary, we have generated a machine-learned potential for the Ge-rich
GexTe binary alloy by using the NN method implemented in the DeePMD package.
We assessed that the NN potential can accurately reproduce the structural
properties of the amorphous, liquid, and crystalline phases of Ge, GeTe, and
Ge2Te alloys. Then, this potential was used to study the crystallization
mechanism of Ge2Te by MD simulations of 19200-atom and 30000-atom cells
lasting up to 40 ns. In stoichiometric GeTe, the system crystallizes in a few
ns at both 600 and 500 K, consistently with previous simulations in literature
Sosso  _et al._ (2013). In Ge2Te instead, we see a different behaviour at
different temperatures. At 500 K, we saw in some models the formation of a
crystalline nucleus of GeTe before phase separation while in an other model
(19200-atom cell) crystal nucleation does not occur on our simulation time (14
ns). Nucleation starts in a point where the fraction of Ge-Ge bonds is lower
than the average of stoichiometric a-GeTe. Crystal growth is, however, slow as
all Ge in excess must be expelled during growth. At this temperature,
overcritical a-Ge nuclei do not form on the time scale of 15 ns if not
preempted by crystal nucleation. On the contrary, at 600 K phase separation
takes place first with the formation of large and separated regions with
average compositions Ge88Te12 and Ge53Te47. The slow down with temperature of
the segregation process seems to be due to the increase of the time needed for
the formation of overcritical nuclei of a-Ge. Crystallization then occurs at
600 K in the Ge-poor region with a nucleation time still much longer than that
of stoichiometric GeTe.

Once a supercritical nucleus is formed in Ge53Te47 regions, crystal growth
proceeds with a further expulsion of Ge in excess (because of a larger driving
force for Ge segregation in the crystalline phase) which, however, implies a
still lower crystal growth velocity with respect to stoichiometric GeTe. We
remark that the a-Ge regions remain amorphous after crystallization of GeTe at
500 and 600 K on the time scale of our simulations.

Overall, our simulations represent a first step towards the atomistic modeling
of the crystallization in Ge-rich GeSbTe alloys which are of interest for
applications in embedded memories. Firstly, we have shown that it is possible
to devise a reliable NN potential suitable to describe the full range of
composition from Ge2Te to GeTe and pure Ge. The database that we generated for
this purpose will contribute to the training set for the development of a NN
potential for the Ge-rich GeSbTe alloys which is in progress. Secondly, we
have shown that phase separation with segregation of most of Ge in excess
occurs on the relatively short time scale of a few ns in Ge2Te at high
temperature (600 K). The complex process of Ge segregation and crystallization
in non-stoichiometric alloys can thus be tackled by MD simulations.

## References

  * Wuttig and Yamada (2007) M. Wuttig and N. Yamada, Nat. Mater. 6, 824 (2007).
  * Noé  _et al._ (2017) P. Noé, C. Vallée, F. Hippert, F. Fillot, and J.-Y. Raty, Semicond. Sci. Technol. 33, 013002 (2017).
  * Zhang  _et al._ (2019) W. Zhang, R. Mazzarello, M. Wuttig, and E. Ma, Nat. Rev. Mater. 4, 150 (2019).
  * Pirovano  _et al._ (2004) A. Pirovano, A. Lacaita, A. Benvenuti, F. Pellizzer, and R. Bez, IEEE Trans. Electron Devices 51, 452–459 (2004).
  * Choe (2017) J. Choe, “Intel 3D XPoint Memory Die Removed from Intel OptaneTM PCM (Phase Change Memory),” (2017).
  * Cappelletti  _et al._ (2020) P. Cappelletti, R. Annunziata, F. Arnaud, F. Disegni, A. Maurelli, and P. Zuliani, J. Phys. D 53, 193002 (2020).
  * Redaelli  _et al._ (2022) A. Redaelli, E. Petroni, and R. Annunziata, Mater. Sci. Semicond. Process. 137, 106184 (2022).
  * Noé  _et al._ (2016) P. Noé, C. Sabbione, N. Bernier, N. Castellani, F. Fillot, and F. Hippert, Acta Mater. 110, 142 (2016).
  * Kim  _et al._ (2007) Y. Kim, U. Hwang, Y. J. Cho, H. Park, M.-H. Cho, P.-S. Cho, and J.-H. Lee, Appl. Phys. Lett. 90, 021908 (2007).
  * Kim  _et al._ (2009) E. T. Kim, J. Y. Lee, and Y. T. Kim, Phys. Status Solidi RRL 3, 103 (2009).
  * Cheng  _et al._ (2011a) H.-Y. Cheng, S. Raoux, and J. L. Jordan-Sweet, Appl. Phys. Lett. 98, 121911 (2011a).
  * Morikawa  _et al._ (2007) T. Morikawa, K. Kurotsuchi, M. Kinoshita, N. Matsuzaki, Y. Matsui, Y. Fujisaki, S. Hanzawa, A. Kotabe, M. Terao, H. Moriya, _et al._ , in  _2007 IEEE International Electron Devices Meeting_ (IEEE, 2007).
  * Cheng  _et al._ (2011b) H. Y. Cheng, T. H. Hsu, S. Raoux, J. Wu, P. Y. Du, M. Breitwisch, Y. Zhu, E. K. Lai, E. Joseph, S. Mittal, R. Cheek, A. Schrott, S. C. Lai, H. L. Lung, and C. Lam, in  _2011 International Electron Devices Meeting_ (2011) p. 3.4.1.
  * Zuliani  _et al._ (2013) P. Zuliani, E. Varesi, E. Palumbo, M. Borghi, I. Tortorelli, D. Erbetta, G. D. Libera, N. Pessina, A. Gandolfo, C. Prelini, L. Ravazzi, and R. Annunziata, IEEE Transactions on Electron Devices 60, 4020 (2013).
  * Zuliani  _et al._ (2015) P. Zuliani, E. Palumbo, M. Borghi, G. Dalla Libera, and R. Annunziata, Solid State Electron. 111, 27 (2015).
  * Palumbo  _et al._ (2017) E. Palumbo, P. Zuliani, M. Borghi, and R. Annunziata, Solid State Electron. 133, 38 (2017).
  * Navarro  _et al._ (2013a) G. Navarro, M. Coué, A. Kiouseloglou, P. Noé, F. Fillot, V. Delaye, A. Persico, A. Roule, M. Bernard, C. Sabbione, _et al._ , in  _2013 IEEE International Electron Devices Meeting_ (IEEE, 2013) p. 21.
  * Sousa  _et al._ (2015) V. Sousa, G. Navarro, N. Castellani, M. Coue, O. Cueto, C. Sabbione, P. Noe, L. Perniola, S. Blonkowski, P. Zuliani, _et al._ , in  _2015 Symposium on VLSI Technology (VLSI Technology)_ (2015) p. T98.
  * Arnaud  _et al._ (2020) F. Arnaud, P. Ferreira, F. Piazza, A. Gandolfo, P. Zuliani, P. Mattavelli, E. Gomiero, G. Samanni, J. Jasse, C. Jahan, J. P. Reynard, R. Berthelon, O. Weber, A. Villaret, B. Dumont, J. C. Grenier, R. Ranica, C. Gallon, C. Boccaccio, A. Souhaite, L. Desvoivres, D. Ristoiu, L. Favennec, V. Caubet, S. Delmedico, N. Cherault, R. Beneyton, S. Chouteau, P. O. Sassoulas, L. Clement, P. Boivin, D. Turgis, F. Disegni, J. L. Ogier, X. Federspiel, O. Kermarrec, M. Molgg, A. Viscuso, R. Annunziata, A. Maurelli, P. Cappelletti, and E. Ciantar, in  _2020 IEEE International Electron Devices Meeting (IEDM)_ (2020).
  * Kiouseloglou  _et al._ (2014) A. Kiouseloglou, G. Navarro, V. Sousa, A. Persico, A. Roule, A. Cabrini, G. Torelli, S. Maitrejean, G. Reimbold, B. De Salvo, F. Clermidy, and L. Perniola, IEEE Transactions on Electron Devices 61 (2014).
  * Grossier  _et al._ (2023) N. Grossier, F. Disegni, A. Ventre, A. Barcella, R. Mariani, V. Marino, S. Mazzara, A. Scavuzzo, M. Bansal, B. Soni, A. Anand, S. Banzal, D. Joshi, R. Narwal, M. Niranjani, K. Trivedi, P. Ferreira, R. Ranica, L. Vullo, A. Cathelin, A. Maurelli, S. Pezzini, and M. Peri, in  _2023 IEEE Symposium on VLSI Technology and Circuits (VLSI Technology and Circuits)_ (2023).
  * Agati  _et al._ (2019) M. Agati, M. Vallet, S. Joulié, D. Benoit, and A. Claverie, J. Mater. Chem. C 7, 8720 (2019).
  * Luong  _et al._ (2021) M. A. Luong, M. Agati, N. Ratel Ramond, J. Grisolia, Y. Le Friec, D. Benoit, and A. Claverie, Phys. Status Solidi RRL 15, 2000471 (2021).
  * Rahier  _et al._ (2023) E. Rahier, M.-A. Luong, S. Ran, N. Ratel-Ramond, S. Saha, C. Mocuta, D. Benoit, Y. Le-Friec, and A. Claverie, Phys. Status Solidi RRL 17, 2200450 (2023).
  * Privitera  _et al._ (2018) S. Privitera, V. Sousa, C. Bongiorno, G. Navarro, C. Sabbione, E. Carria, and E. Rimini, J. Phys. D 51, 145103 (2018).
  * Privitera  _et al._ (2020) S. Privitera, I. López García, C. Bongiorno, V. Sousa, M. Cyrille, G. Navarro, C. Sabbione, E. Carria, and E. Rimini, J. Appl. Phys. 128, 155105 (2020).
  * Prazakova  _et al._ (2022) L. Prazakova, E. Nolot, E. Martinez, D. Rouchon, F. Fillot, N. Bernier, R. Elizalde, M. Bernard, and G. Navarro, Materialia 21, 101345 (2022).
  * Díaz Fattorini  _et al._ (2022) A. Díaz Fattorini, C. Chèze, I. Lø’pez García, C. Petrucci, M. Bertelli, F. Righi Riva, S. Prili, S. M. S. Privitera, M. Buscema, S. Sciuto, Antonella amd Di Franco, G. D’Arrigo, M. Longo, S. De Simone, V. Mussi, E. Placidi, M.-C. Cyrille, N.-P. Tran, R. Calarco, and F. Arciprete, Nanomaterials 12, 1340 (2022).
  * Cecchi  _et al._ (2022) S. Cecchi, I. Lopez Garcia, A. M. Mio, E. Zallo, O. Abou El Kheir, R. Calarco, M. Bernasconi, G. Nicotra, and S. M. Privitera, Nanomaterials 12, 631 (2022).
  * Laurin  _et al._ (2023) L. Laurin, M. Baldo, E. Petroni, G. Samanni, L. Turconi, A. Motta, M. Borghi, A. Serafini, D. Codegoni, M. Scuderi, S. Ran, A. Claverie, D. Ielmini, R. Annunziata, and A. Redaelli, in  _2023 IEEE International Reliability Physics Symposium (IRPS)_ (2023) pp. 1–7.
  * Petroni  _et al._ (2022) E. Petroni, A. Serafini, D. Codegoni, P. Targa, L. Mariani, M. Scuderi, G. Nicotra, and A. Redaelli, Frontiers in Physics 10 (2022).
  * Petroni  _et al._ (2023) E. Petroni, M. Patelmo, A. Serafini, D. Codegoni, L. Laurin, M. Baldo, and A. Redaelli, Phys. Status Solidi RRL 17, 2200458 (2023).
  * Sun  _et al._ (2021) L. Sun, Y.-X. Zhou, X.-D. Wang, Y.-H. Chen, V. L. Deringer, R. Mazzarello, and W. Zhang, npj Comput. Mater. 7, 1 (2021).
  * Abou El Kheir and Bernasconi (2021) O. Abou El Kheir and M. Bernasconi, Nanomaterials 11, 2382 (2021).
  * Sosso  _et al._ (2012a) G. C. Sosso, G. Miceli, S. Caravati, J. Behler, and M. Bernasconi, Phys. Rev. B 85, 174103 (2012a).
  * Behler and Parrinello (2007) J. Behler and M. Parrinello, Phys. Rev. Lett. 98, 146401 (2007).
  * Abou El Kheir  _et al._ (2024) O. Abou El Kheir, L. Bonati, M. Parrinello, and M. Bernasconi, npj Comput. Mater. 10, 33 (2024).
  * Zhang  _et al._ (2018) L. Zhang, J. Han, H. Wang, R. Car, and W. E, Phys. Rev. Lett. 120, 143001 (2018).
  * Zeng  _et al._ (2023) J. Zeng, D. Zhang, D. Lu, P. Mo, _et al._ , J. Chem. Phys. 159, 054801 (2023).
  * Wang  _et al._ (2018) H. Wang, L. Zhang, J. Han, and E. Weinan, Comput. Phys. Commun. 228, 178 (2018).
  * Zhou  _et al._ (2023) Y. Zhou, W. Zhang, E. Ma, and V. L. Deringer, Nat. Electr. 6, 746–754 (2023).
  * Carria  _et al._ (2012) E. Carria, A. Mio, S. Gibilisco, M. Miritello, C. Bongiorno, M. Grimaldi, and E. Rimini, J. Electrochem. Soc. 159, H130 (2012).
  * Gourvest  _et al._ (2009) E. Gourvest, S. Lhostis, J. Kreisel, M. Armand, S. Maitrejean, A. Roule, and C. Vallée, Appl. Phys. Lett. 95, 031908 (2009).
  * Raoux  _et al._ (2009a) S. Raoux, B. Muñoz, H.-Y. Cheng, and J. L. Jordan-Sweet, Appl. Phys. Lett. 95, 143118 (2009a).
  * Navarro  _et al._ (2013b) G. Navarro, V. Sousa, A. Persico, N. Pashkov, A. Toffoli, J.-C. Bastien, L. Perniola, S. Maitrejean, A. Roule, P. Zuliani, R. Annunziata, and B. De Salvo, Solid State Electron. 89, 93 (2013b).
  * VandeVondele  _et al._ (2005) J. VandeVondele, M. Krack, F. Mohamed, M. Parrinello, T. Chassaing, and J. Hutter, Comput. Phys. Commun. 167, 103 (2005).
  * Thompson  _et al._ (2022) A. P. Thompson, H. M. Aktulga, R. Berger, D. S. Bolintineanu, W. M. Brown, P. S. Crozier, P. J. in ’t Veld, A. Kohlmeyer, S. G. Moore, T. D. Nguyen, R. Shan, M. J. Stevens, J. Tranchida, C. Trott, and S. J. Plimpton, Comp. Phys. Comm. 271, 108171 (2022).
  * Mickel  _et al._ (2013) W. Mickel, S. C. Kapfer, G. E. Schröder-Turk, and K. Mecke, J. Chem. Phys. 138, 044501 (2013).
  * Bartók  _et al._ (2013) A. P. Bartók, R. Kondor, and G. Csányi, Phys. Rev. B 87, 184115 (2013).
  * De  _et al._ (2016) S. De, A. P. Bartók, G. Csányi, and M. Ceriotti, Phys. Chem. Chem. Phys. 18, 13754 (2016).
  * Musil  _et al._ (2021) F. Musil, A. Grisafi, A. P. Bartók, C. Ortner, G. Csányi, and M. Ceriotti, Chem. Rev. 121, 9759 (2021).
  * Xu  _et al._ (2022) Y. Xu, Y. Zhou, X.-D. Wang, W. Zhang, E. Ma, V. L. Deringer, and R. Mazzarello, Adv. Mat. 34, 2270084 (2022).
  * Laakso  _et al._ (2023) J. Laakso, L. Himanen, H. Homm, E. V. Morooka, M. O. Jäger, M. Todorović, and P. Rinke, J. Chem. Phys. 158 (2023).
  * Himanen  _et al._ (2020) L. Himanen, M. O. J. Jäger, E. V. Morooka, F. Federici Canova, Y. S. Ranawat, D. Z. Gao, P. Rinke, and A. S. Foster, Comp. Phys. Comm. 247, 106949 (2020).
  * Larsen  _et al._ (2017) A. H. Larsen _et al._ , J. Phys.: Condens. Matter 29, 273002 (2017).
  * Urata (2022) S. Urata, J. Phys. Chem. C 126, 21507 (2022).
  * Zhai and Wang (2023) B. Zhai and H. Wang, Comput. Mater. Sci. 216, 111843 (2023).
  * Liu  _et al._ (2023) J. Liu, R. Liu, Y. Cao, and M. Chen, Phys. Chem. Chem. Phys. 25, 983 (2023).
  * Huang  _et al._ (2023) X. Huang, K. Luo, Y. Shen, Y. Yue, and Q. An, Energy and AI 11, 100210 (2023).
  * Laaziri  _et al._ (1995) K. Laaziri, S. Roorda, and J. Baribeau, J. Non-Cryst. Solids 191, 193 (1995).
  * Weber  _et al._ (2017) H. Weber, M. Schumacher, P. Jóvári, Y. Tsuchiya, W. Skrotzki, R. Mazzarello, and I. Kaban, Phys. Rev. B 96, 054204 (2017).
  * Grimme  _et al._ (2010) S. Grimme, J. Antony, S. Ehrlich, and H. Krieg, J. Chem. Phys. 132, 154104 (2010).
  * Sosso  _et al._ (2012b) G. C. Sosso, J. Behler, and M. Bernasconi, Phys. Status Solidi B 249, 1880 (2012b).
  * Klemm and Frischmuth (1934) W. Klemm and G. Frischmuth, Zeitschrift für anorganische und allgemeine Chemie 218, 249 (1934).
  * Schlieper  _et al._ (1999) A. Schlieper, Y. Feutelais, S. Fries, B. Legendre, and R. Blachnik, Calphad 23, 1 (1999).
  * Dünweg and Kremer (1993) B. Dünweg and K. Kremer, J. Chem. Phys. 99, 6983–6997 (1993).
  * Cohen and Grest (1979) M. H. Cohen and G. Grest, Phys. Rev. B 20, 1077 (1979).
  * Angell (1995) C. A. Angell, Science 267, 1924 (1995).
  * Imabayashi and Ishimaru (2017) S. Imabayashi and M. Ishimaru, Materials Transactions 58, 857 (2017).
  * Chen  _et al._ (2017) Y. Chen, G. Wang, L. Song, X. Shen, J. Wang, J. Huo, R. Wang, T. Xu, S. Dai, and Q. Nie, Cryst. Growth Des. 17, 3687 (2017).
  * Deringer  _et al._ (2014) V. L. Deringer, W. Zhang, M. Lumeij, S. Maintz, M. Wuttig, R. Mazzarello, and R. Dronskowski, Angew. Chem. Int. Ed. 53, 10817 (2014).
  * Errington and Debenedetti (2001) J. R. Errington and P. G. Debenedetti, Nature 409, 318 (2001).
  * Mazzarello  _et al._ (2010) R. Mazzarello, S. Caravati, S. Angioletti-Uberti, M. Bernasconi, and M. Parrinello, Phys. Rev. Lett. 104, 085503 (2010).
  * Spreafico  _et al._ (2011) E. Spreafico, S. Caravati, and M. Bernasconi, Phys. Rev. B 83, 144205 (2011).
  * Akola and Jones (2007) J. Akola and R. Jones, Phys. Rev. B 76, 235201 (2007).
  * Goldak  _et al._ (1966) J. Goldak, C. S. Barrett, D. Innes, and W. Youdelis, J. Chem. Phys. 44, 3323 (1966).
  * Chattopadhyay  _et al._ (1987) T. Chattopadhyay, J. X. Boucherle, and H. G. vonSchnering, J. Phys. C: Solid State Phys. 20, 1431 (1987).
  * Fons  _et al._ (2010) P. Fons, A. V. Kolobov, M. Krbal, J. Tominaga, K. S. Andrikopoulos, S. N. Yannopoulos, G. A. Voyiatzis, and T. Uruga, Phys. Rev. B 82, 155209 (2010).
  * Matsunaga  _et al._ (2011) T. Matsunaga, P. Fons, A. V. Kolobov, J. Tominaga, and N. Yamada, Appl. Phys. Lett. 99, 231907 (2011).
  * Xu  _et al._ (2018) M. Xu, Z. Lei, J. Yuan, K. Xue, Y. Guo, S. Wang, X. Miao, and R. Mazzarello, RSC Adv. 8, 17435 (2018).
  * Dangić  _et al._ (2022) D. Dangić, S. Fahy, and I. Savić, Phys. Rev. B 106, 134113 (2022).
  * Sosso  _et al._ (2013) G. C. Sosso, G. Miceli, S. Caravati, F. Giberti, J. Behler, and M. Bernasconi, J. Phys. Chem. Lett. 4, 4241 (2013).
  * Ghezzi  _et al._ (2011) G. E. Ghezzi, J. Y. Raty, S. Maitrejean, A. Roule, E. Elkaim, and F. Hippert, Appl. Phys. Lett. 99, 151906 (2011).
  * Gabardi  _et al._ (2017) S. Gabardi, E. Baldi, E. Bosoni, D. Campi, S. Caravati, G. Sosso, J. Behler, and M. Bernasconi, J. Phys. Chem. C 121 (2017).
  * Gabardi  _et al._ (2019) S. Gabardi, G. C. Sosso, J. Behler, and M. Bernasconi, Faraday Discuss. 213, 287 (2019).
  * Perego  _et al._ (2023) S. Perego, D. Dragoni, S. Gabardi, D. Campi, and M. Bernasconi, Phys. Stat. Solidi RRL 13, 2200433 (2023).
  * Sosso  _et al._ (2015) G. C. Sosso, M. Salvalaglio, J. Behler, M. Bernasconi, and M. Parrinello, J. Phys. Chem. C 119, 6428 (2015).
  * Acharya  _et al._ (2024) D. Acharya, O. Abou El Kheir, D. Campi, and M. Bernasconi, Sci. Rep. 14, 3224 (2024).
  * Wintersteller  _et al._ (2024) S. Wintersteller, O. Yarema, D. Kumaar, F. M. Schenk, O. V. Safonova, P. M. Abdala, V. Wood, and M. Yarema, Nat. Commun. 15, 1011 (2024).
  * Raoux  _et al._ (2009b) S. Raoux, , H.-Y. Cheng, M. Caldwell, and H. Wong, Appl. Phys. Lett. 95, 071910 (2009b).
  * Berthier  _et al._ (2017) R. Berthier, N. Bernier, D. Cooper, C. Sabbione, F. Hippert, and P. Noé, J. Appl. Phys. 122, 115304 (2017).
  * Agati  _et al._ (2020) M. Agati, C. Gay, D. Benoit, and A. Claverie, Appl. Surf. Sci. 518, 146227 (2020).
  * Lee  _et al._ (2020) D. Lee, K. Lee, D. Yoo, W. Jeong, and S. Han, Comp. Mater. Sci. 181, 109725 (2020).

Data availability

The NN potential, the training DFT database, atomic trajectories of the phase
separation and crystallization process, and a video of the relevant processes
will be available in the Materials Cloud repository after acceptance of the
article.

Acknowledgements

The project has received funding from European Union NextGenerationEU through
the Italian Ministry of University and Research under PNRR M4C2I1.4 ICSC
Centro Nazionale di Ricerca in High Performance Computing, Big Data and
Quantum Computing (Grant No. CN00000013). We acknowledge the CINECA award
under the ISCRA initiative, for the availability of high-performance computing
resources and support.

Code availability

LAMMPS, DeePMD and CP2k are free and open source codes available at https://
lammps.sandia.gov, http://www.deepmd.org and http://www.cp2k.org,
respectively.

Competing interests

The authors declare no competing interests.

Authors contributions

D.B. and M.B. designed the research, performed the research and analyzed data.
O.A.E contributed to the generation of the NN potential. D.B. and M.B. wrote
the paper and all the authors edited the manuscript before submission.