Frédéric Soisson
email: frederic.soisson@cea.fr

Maylise Nastar
email: maylise.nastar@cea.fr

Atomistic simulations of diffusive phase transformations with non-conservative point defects

Keywords: Defects, diffusion, kinetics, phase transformation, simulation

Atomistic simulations of diffusive phase transformations with non-conservative point defects

Introduction

Materials properties and design considerations are essentially related to their composition and microstructure. The action of lattice point defects (PDs) on a microstructure is multiple. The agglomeration of PDs is at the origin of some of the extended defects, their creation/annihilation at pre-existing or nucleated extended defects can induce the displacement, the growth, or the transformation of these [START_REF] Kirkendall | Diffusion of Zinc in Alpha Brass[END_REF][START_REF] Herring | Diffusional viscosity of a polycrystalline solid[END_REF][START_REF] Zinkle | Materials challenges in nuclear energy[END_REF]. Finally, PD concentration greatly determines atomic diffusion in solids, because an atom that shares the same lattice as PDs, diffuses through a series of exchanges with nearest-neighbor PDs [START_REF] Mehrer | Diffusion in Solids[END_REF][START_REF] Van Der Ven | Vacancy mediated substitutional diffusion in binary crystalline solids[END_REF]. The resulting atomic diffusion coefficients and rates of diffusion-controlled phase transformation, are thus proportional to the PD concentration. A change of concentration of one PD population may thus affect the competition between phase transformations whenever they are controlled by different diffusion mechanisms, as for example the vacancy diffusion mechanism and carbon direct-interstitial diffusion mechanism controlling the formation of cementite in iron [START_REF] Schuler | Transport properties of dilute α-Fe(X) solid solutions (X = C, N, O)[END_REF][START_REF] Schuler | Mass-transport properties of ternary Fe(C,O) alloys revealed by multicomponent cluster synergies[END_REF]. Therefore, non-equilibrium PDs may not only change the rates of atomic diffusion and phase transformation, but also the kinetic path, and the resulting steady state microstructure in driven systems [START_REF] Zinkle | Materials challenges in nuclear energy[END_REF][START_REF] Nastar | Radiation-induced segregation[END_REF].

Although lattice PDs are equilibrium defects in crystalline solids, their concentration may be non uniform, as for example in a two-phase equilibrium system. The PD equilibrium concentration corresponds then to a spatial average of local PD equilibrium concentrations within selected unit volumes. Besides, PDs are not always at equilibrium. Departure of PD from its local equilibrium within a given finite volume, happens whenever source-and-sinks (in short SAS) are not active enough and/or PD long-range diffusion is not fast enough, to guarantee instantaneous annihilation and/or creation of PD everywhere. In driven systems such as materials submitted to irradiation, the departure of PDs from their local equilibrium concentration is permanent. Non-equilibrium vacancies are generated upon quenching heat treatment [START_REF] Anthony | Solute segregation in vacancy gradients generated by sintering and temperature changes[END_REF], by mechanical solicitations [START_REF] Militzer | Modelling the effect of deformationinduced vacancies on segregation and precipitation[END_REF], by a Kirkendall's mechanism upon annealing [START_REF] Kirkendall | Diffusion of Zinc in Alpha Brass[END_REF][START_REF] Desgranges | Depletion and Voids Formation in the Substrate During High Temperature Oxidation of Ni-Cr Alloys[END_REF], or consequently to a surface reaction such as oxidation or nitriding [START_REF] Hu | Nonequilibrium point defects and diffusion in silicon[END_REF][START_REF] Bobeth | Modelling of vacancy diffusion and pore formation during parabolic oxide growth[END_REF][START_REF] Nakamura | Formation of oxide nanotubes via oxidation of Fe, Cu and Ni nanowires and their structural stability: Difference in formation and shrinkage behavior of interior pores[END_REF][START_REF] Masse | Point defect creation induced by solid state reaction between nickel and silicon[END_REF][START_REF] Martin | Use of a buried loop layer as a detector of interstitial flux during oxidation of SiGe heterostructures[END_REF]. The response of the system to a non-equilibrium PD concentration is many-fold. It mainly depends on the competition between diffusion, annihilation/creation reactions at SAS, and agglomeration reactions involving PDs [START_REF] Zinkle | Materials challenges in nuclear energy[END_REF] and/or alloying elements [START_REF] Maydet | Precipitate stability under irradiation: Point defect effects[END_REF][START_REF] Russell | The theory of phase stability under irradiation[END_REF][START_REF] Voronkov | Vacancy-type microdefect formation in czochralski silicon[END_REF][START_REF] Nastar | Thermodynamic model for lattice point defect-mediated semi-coherent precipitation in alloys[END_REF]. The agglomeration of non-equilibrium PDs in clusters is recognized to be responsible for damaging electronic and optical properties of semi-conductors [START_REF] Martin | Use of a buried loop layer as a detector of interstitial flux during oxidation of SiGe heterostructures[END_REF][START_REF] Mccluskey | Defects in semiconductors[END_REF], or the dimensional and mechanical properties of nuclear plant components [START_REF] Zinkle | Materials challenges in nuclear energy[END_REF]. On the positive side, the agglomeration of vacancies into voids subsequent to oxidation-induced vacancy injection, is currently used to fabricate hollow nanoparticles for various applications [START_REF] Nakamura | Formation of oxide nanotubes via oxidation of Fe, Cu and Ni nanowires and their structural stability: Difference in formation and shrinkage behavior of interior pores[END_REF], or investigate the diffusion properties of dopants in Si [START_REF] Hu | Nonequilibrium point defects and diffusion in silicon[END_REF].

There have been numerous atomic modeling studies investigating the impact of the variation of PD jump frequencies with the local composition and strain, on diffusion and phase transformations, a lot less on the kinetic impact of a departure from PD local equilibrium concentration. Yet, it may strongly modify the kinetic processes and the resulting microstructure. PDs are non-conservative species, their concentration depends on temperature, applied stress and microstructure. In an evolving microstructure, PDs which diffuse faster than atoms are generally treated as an adiabatic species, that instantaneously reaches local equilibrium with respect to the local atomic composition. When a vacancy is removed (created) a lattice site is removed (created). When a self-interstitials is removed (created) by being transformed into two substitutional atoms on two lattice sites, a lattice site is created (removed). The creation-annihilation reactions of PDs not only occur at free surfaces, but also at extended defects such as climbing dislocations, grain-boundaries, PD clusters, all acting as SAS of PDs. Therefore, PDs not only control diffusion but also the displacement of lattice sites whenever the sources and sinks of PDs are not uniformely distributed.

There are therefore three aspects that we choose to discuss here: the equilibrium PD concentration, its variation in an evolving microstructure, and finally the kinetic processes that produce deviations from equilibrium of PDs and their impact on the microstructure.

In an equilibrium alloy, configurational entropies related to the number of possible permutations of PDs and atoms on their lattice is an essential thermodynamic quantity, that determines the variation of the equilibrium concentration of PDs with temperature. To properly account for the configurational entropy, Van der Ven et al. have developed atomic Monte Carlo methods in the semi-grand canonical ensemble that yield equilibrium PD concentrations of ordered compounds and solid solutions [START_REF] Van Der Ven | Vacancy mediated substitutional diffusion in binary crystalline solids[END_REF][START_REF] Van Der Ven | Vacancies in ordered and disordered binary alloys treated with the cluster expansion[END_REF][START_REF] Belak | Effect of disorder on the dilute equilibrium vacancy concentrations of multicomponent crystalline solids[END_REF], while the present authors have introduced sampling Monte Carlo methods in the canonical ensemble [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF][START_REF] Li | Magnetochemical effects on phase stability and vacancy formation in fcc fe-ni alloys[END_REF].

As long as we may assume PD adiabatically reaches its local equilibrium concentration, we may add to a standard canonical Atomistic Kinetic Monte Carlo (AKMC) method, performed with constant numbers of PDs, an on-thefly correction of the AKMC time to simulates both the physical time and the time evolution of PD concentration in an evolving microstructure [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF]. Note that, although it is a common practice, an homothetic correction of the AKMC time, calibrated on experimentally observed microstructures, may be wrong, because PD concentration is not constant, but a complex, a priory unknown function of time.

In dilute alloys, simplified SAS models have been explicitly introduced in AKMC simulation box [26,[START_REF] Hin | Kinetics of heterogeneous dislocation precipitation of NbC in alpha-iron[END_REF]. To keep constant the number of lattice sites, a removed PD is systematically replaced by a solvent atom and vice versa. The PD creation/removal reaction rates was assumed to be constant, independant of the local composition. Apart from that case, models accounting for the creation mechanisms of non-equilibrium point defects, their impact on diffusion and lattice motion, and the resulting microstructure, are mostly mesoscopic, phenomenological rate theory models relying on a coupled set of reactiondiffusion equations, such as the Kirkendall model of interdiffusion [START_REF] Kirkendall | Diffusion of Zinc in Alpha Brass[END_REF], or the Nabarro-Herring model of diffusional creep [START_REF] Herring | Diffusional viscosity of a polycrystalline solid[END_REF]. Cluster dynamics methods have been developed to treat the competition between PD and atom longrange diffusion, PD clustering, and vacancy-solute clustering reactions under thermal [START_REF] Jourdan | Influence of cluster mobility on Cu precipitation in α-Fe: A cluster dynamics modeling[END_REF], and irradiation conditions [START_REF] Maydet | Precipitate stability under irradiation: Point defect effects[END_REF][START_REF] Russell | The theory of phase stability under irradiation[END_REF][START_REF] Mathon | Experimental study and modelling of copper precipitation under electron irradiation in dilute FeCu binary alloys[END_REF]. However, they are limited to dilute alloys and a large part of flux couplings driven by concentration gradients of non-equilibrium PDs is ignored.

Despite the intense modeling activity in this field, we are not aware of atomic scale kinetic methods tackling the competition between diffusion and the annihilation/creation mechanisms of PDS at SAS in concentrated alloys. We therefore propose an AKMC simulation method including these mechanisms at point-like SAS. We present a simple application of the method, by simulating the vacancy-mediated coherent precipitation of a B-rich phase in an unmixing A-B solid solution, under isothermal and non-isothermal annealing treatments. Then, we discuss the capabilities of the existing mesoscopic models and the present AKMC method regarding various kinetic phenomena that are driven by non-equilibrium PD concentrations.

2 Atomic kinetic Monte Carlo simulation including PD sources and sinks

Equilibrium vacancy concentration

We first present two atomic Monte Carlo methods that have been developed in the canonical ensemble, to evaluate the equilibrium vacancy concentration of an equilibrium or a slowly evolving non-equilibrium microstructure: an extension of the Widom insertion method (2.1.1) and one based on the measurements of local PD concentrations and the correction of the simulated AKMC time (2.2).

Widom-based method

The standard Widom insertion method can be used to determine the chemical potential of an alloy component, and its equilibrium mixing properties (see, e.g. [START_REF] Sindzingre | Partial enthalpies and related quantities in mixtures from computer simulation[END_REF]). An extension of this method has been deployed [START_REF] Li | Magnetochemical effects on phase stability and vacancy formation in fcc fe-ni alloys[END_REF][START_REF] Sindzingre | Partial enthalpies and related quantities in mixtures from computer simulation[END_REF], to yield the equilibrium atomic fraction of vacancy as an ensemble average in the canonical ensemble,

X eq V = n P n X A exp - E f,A→V (n) k B T , (1)
where P n is the probability of occurrence of the atomic configuration n, which is naturally produced by a Monte Carlo simulation. X A accounts for the number of possible insertions of a single vacancy in place of an atom A in configuration n. The energy of vacancy formation from one atom of species A in configuration n, is given by

E f,A→V (n) = ∆E A→V n + µ A , (2)
with the substitution energy, ∆E A→V n , equal to the difference of enthalpies of the atomic configuration n with vacancy V in place of atom A on a given site, and the atomic configuration with atom A on the same site. X A is the atomic fraction of species A, equal to (X A = N A /N ), with N and N A respectively the number of atoms and atoms A. This formulae accounts for the fact that to form a vacancy in place of atom A, we have to put atom A in a reservoir. The energy of this atom in the reservoir, which is equal to its chemical potential, µ A , contributes to the vacancy formation energy. This way, the vacancy formation energy in a given atomic configuration, n, accounts for nominal composition and temperature of the alloy. The probability of forming this vacancy in configuration n is (X A exp -E f,A→V (n)/k B T ). The equilibrium vacancy concentration results then from an ensemble average of such probability over the atomic configurations. Note that, this formulae implies PDs are created as non-interacting isolated PDs.

To determine the equilibrium vacancy concentration of an equilibrium binary alloy AB, we first perform equilibrium Monte Carlo simulations in the semi-grand canonical ensemble (in short, GCMC). We impose the alloy chemical potential ∆µ = µ B -µ A , then by using a Metropolis algorithm to obtain an equilibrium sampling of the low-energy atomic configurations, we deduce the equilibrium alloy composition from an ensemble average of the atomic configuration compositions, then extract the atom chemical potentials from the calculation of the Gibbs free energy, equal to the integral of ∆µ with respect to X B . Finally, we compute the equilibrium atomic fraction (also called concentration) of vacancy by means of the Widom insertion method. As shown by Eq. 1, the equilibrium vacancy fraction is an ensemble average of exponential of vacancy formation energy. Note that for a given atomic configuration, instead of randomly choose one atom A, we may run the sum over all atoms A. In that case, X A in RHS of Eq. 1, should be removed when computing the vacancy equilibrium concentration. Note that, the same equilibrium vacancy concentration is obtained if instead of atoms A, atoms B are virtually replaced by a vacancy. This method will be referred as GCMC-Widom. This method is similar to the one first introduced by Belak et al [START_REF] Belak | Effect of disorder on the dilute equilibrium vacancy concentrations of multicomponent crystalline solids[END_REF] for their study of the effect of disorder on the equilibrium vacancy concentration in Ti-Al alloys, and by Barnard et al. in [31] for their AKMC simulations of Ni 2 Cr ordering.

Time average as an ensemble average

In the canonical ensemble, the usual AKMC simulation box contains a fixed number, N V , of vacancies (in most cases, N V = 1). At equilibrium and if the AKMC run is long enough to guarantee an ergodic exploration of all accessible local microstates, the residence time of a vacancy in any local environment i, should be proportional to exp[-G for V (i)/k B T ], with G for V (i) the Gibbs free energy of vacancy formation of an infinite system at composition equal to the one of environment i [START_REF] Nastar | Radiation-induced segregation[END_REF][START_REF] Bouar | Kinetic pathways from embedded-atom-method potentials: Influence of the activation barriers[END_REF][START_REF] Soisson | Cu-precipitation kinetics in α-Fe from atomistic simulations: Vacancy-trapping effects and cu-cluster mobility[END_REF]. Although, one vacancy in a simulation box generally yields a vacancy concentration higher than the equilibrium one, the relative residence times of the vacancy between the various local environments are the thermodynamic ones. Therefore, it is only needed to correct the vacancy concentration in a chosen environment i -for which we know its equilibrium concentration. From the time average concentration of vacancy in environment i, X MC V (i), we then easily deduce the equilibrium fraction of vacancy,

X eq V = X MC V X MC V (i) X eq V (i) . (3)
Note that, to account for vacancy cluster contributions, we may need to introduce several vacancies in the simulation box.

On-the-fly rescaling time method

The equilibrium methods presented above, can also be employed to simulate kinetics of slowly evolving microstructures. Although, most of the time we cannot work with the physical PD concentration, we can correct the simulated time to reproduce the correct diffusion coefficients. Indeed, in cases where it is possible to ignore contributions of vacancy clusters, the diffusion coefficients of substitutional atomic species are proportional to the vacancy concentration.

Hence, as previously proposed (see e.g. [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF][START_REF] Soisson | Cu-precipitation kinetics in α-Fe from atomistic simulations: Vacancy-trapping effects and cu-cluster mobility[END_REF]), we rescale the time of the Monte Carlo simulations (t MC ), to get a physical time scale t:

t = t MC X MC V X V , (4)
where X MC V = N V /N is the vacancy atomic fraction in the simulation box, and X V the real time-dependent fraction. Using the fact that vacancies diffuse much faster than atoms (in a ratio that is approximately the vacancy atomic fraction), we compute X V from the hypothesis, that vacancy adiabatically remains at equilibrium during the transformation, i.e. for each local environment i, the vacancy fraction is the equilibrium one, X eq V (i). Note that, the local equilibrium hypothesis not only applies for PDs but also for the atoms, because the local PD concentration is the equilibrium one only if the atomic shortrange order is the thermodynamic one. Even though, the system itself did not reach a macroscopic equilibrium, one can define a microstructure-constrained equilibrium vacancy concentration, X eq V , and set X V = X eq V in eq. 4.

To rescale the time, a first method -referred hereafter as the AKMC-Widom method -consists in performing an AKMC simulation for a given N V and X MC V , and to evaluate X eq V from time to time on chosen atomic configurations by means of the Widom method described in the previous section (using the equilibrium chemical potentials previously given by semi-grand canonical Monte Carlo simulations). A second method -referred hereafter as the "local environment time rescaling" method (AKMC-LETR) -has been proposed previously (see e.g. [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF][START_REF] Soisson | Cu-precipitation kinetics in α-Fe from atomistic simulations: Vacancy-trapping effects and cu-cluster mobility[END_REF]). To emphasize the similarities and differences between both methods, we give a brief presentation of this second method just below.

In the AKMC-LETR method, we rescale the time by measuring the average vacancy fraction of a selected local environment i:

t = t MC X MC V (i) X eq V (i) , (5)
where X MC V (i) is the vacancy concentration for the environment i in the AKMC simulation (which can be easily measured). In principle, any choice of the reference environment (i), must give the same result, provided the exploration of the energetic landscape by the vacancy is ergodic. For a phase separation in a demixing binary alloy AB, convenient choices for i is pure A or pure B environment volumes. Since the equilibrium fraction, X eq V , depends on the atomic configuration (in the proportions of the various i environments), it changes during the phase transformation and its variation must be evaluated "on-thefly". The time of exploration must be sufficiently long to get good statistics and sufficiently short to avoid a significant evolution of the microstructure.

Note that the comparison between eq. 4 and eq. 5 leads to eq.3. A time rescaling method can also be employed to estimate the time evolution of X eq V and G for V = -k B T ln X eq V during the transformation, using either the AKMC-LETR or the AKMC-Widom method. Details are given in ref. [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF].

Vacancy sources and sinks

The above methods are based on an adiabatic approximation and a local equilibrium hypothesis: the time to reach the equilibrium vacancy concentration, in a given microstructure, is negligible compared to the evolution time of the microstructure. As mentioned in the introduction, while this assumption is often reasonable, there are clearly situations where it is not necessarily true. For instance, in the beginning of an isothermal annealing after a rapid quench from high temperature, or in alloys under irradiation (where high supersaturations of PDs are permanently sustained). For such situations, we propose here a new Monte Carlo simulation method with non-conservative vacancies (the NCV-AKMC method), which we describe in detail in this section. It includes an explicit treatment of the thermally activated creation and elimination mechanisms of vacancies in concentrated alloys, and allows to perform simulations with real time-dependent vacancy concentrations (a similar method was first proposed for dilute alloys [26,[START_REF] Hin | Kinetics of heterogeneous dislocation precipitation of NbC in alpha-iron[END_REF]). Note that in the above canonical methods, vacancy was implicitly assumed to be an infinite dilute species, while the NCV-AKMC method allows, in principle, a proper modeling of PD clusters and agglomeration processes.

We apply the NCV-AKMC method (as well as the AKMC-Widom and AKMC-LETR ones) to a model alloy AB with a phase separation tendency, using the same diffusion model and parameters as in ref. [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF] (given in Appendix A). We first recall the main features of this diffusion model, before to present the details of the NCV-AKMC method. The alloy has a bodycentered cubic (bcc) structure and its energy is written as a sum of constant pairwise interactions between first nearest neighbors ε AA , ε BB and ε AB on a rigid lattice. Its phase diagram displays a symmetrical miscibilty gap below T c ≃ 928 K, with a phase separation between a A-rich phase α and a B rich phase β (with respective B atomic fractions X α B and X β B ). Interactions between atoms and vacancies (ε AV and ε BV ) control the formation energies. For example, G for V (A) = -(z/2)ε AA + zε AV in pure A (where z = 8 is the coordination number). With the chosen parameters, G for V (A) = 1.4 eV in pure A and is larger than G for V (B) = 1.0 eV in pure B. Diffusion of atoms A an B, occurs by vacancy jumps: one vacancy can exchange with one of its first nearest neighbors A or B. The jump frequency for an A-V exchange writes

Γ mig AV = ν A exp - E mig AV k B T , (6)
where the migration barrier of an A-V exchange, E mig AV , depends on the local atomic environment through the pair interactions, which are broken around the A-V pair at the initial position and created at the saddle-point [START_REF] Bouar | Kinetic pathways from embedded-atom-method potentials: Influence of the activation barriers[END_REF][START_REF] Soisson | Cu-precipitation kinetics in α-Fe from atomistic simulations: Vacancy-trapping effects and cu-cluster mobility[END_REF][START_REF] Girifalco | Vacancy concentration and diffusion in order-disorder alloys[END_REF]. With the AKMC parameters of Appendix A, the migration barriers are of the order of 1 eV (for example: 1.15 and 1.025 eV respectively for the A-V and B-V exchanges in pure A, 1.025 and 1.25 eV for the A-V and B-V exchanges in pure B).

AKMC simulations are performed with a residence time algorithm [START_REF] Bortz | A new algorithm for monte carlo simulation of ising spin systems[END_REF]. In simulations performed with a constant number of vacancies (section 2.2), all the possible jump frequencies are computed at each Monte Carlo step. One jump is chosen, using a random number, with a probability proportional to its jump frequency. The residence time associated with this Monte Carlo step is t MC = 1/ i Γ iV , where the sum runs over all the neighbors A or B of the vacancies. In the AKMC-LETR and AKMC-Widom methods, this residence time is then rescaled to take into account the difference between the vacancy atomic fraction of the simulation box and the real one, using eq. 4 or 5, as explained in 2.2.

In the NCV-AKMC method, to avoid the adiabatic assumption and the vacancy local equilibrium constraint (corresponding to setting X V (t) = X eq V (t)), and to explicitly simulate the dynamics of vacancies during the transformation, we include the vacancy creation and annihilation mechanisms in the set of kinetic events. We introduce N s sites in the lattice, which act as SAS for vacancies. In the following, we introduce a random spatial distribution of the SAS sites, although they could also be arranged in lines or planes, to mimic the geometry of dislocations or grain boundaries.

Q sp = (1 -X s B )Q sp XA + X s B Q sp
XB determines the height of the formation and elimination barriers, but not the formation energy. The migration barrier for the exchange between a vacancy on site i and a atom Y on site i is also represented. The height of the barriers depends on the nature of the atom to be replaced (Y = A, B) and on the local composition.

At each MC step, one vacancy can be created on a SAS lattice site initially occupied by an atom A, with a frequency equal to

Γ + AV = ν s A exp - E + AV k B T , (7)
where

E + AV = (1 -X s B )Q sp AA + X s B Q sp AB -n s A ε AA -n s B ε AB . (8)
X s B is the B atomic fraction around the SAS site, n s A and n s B the numbers of A and B atoms among the SAS first nearest-neighbors (see Fig. 1). Q sp AA , Q sp AB are additional saddle-point parameters that control the rate of vacancy creation. A similar equation applies if the SAS site is occupied by an atom B. When the vacancy is created in replacement of an atom A or B, this atom is placed in a reservoir that contains n res A and n res B atoms A and B. Inversely, when a vacancy V is on a SAS lattice site, it can be removed to be replaced by an atom A with a frequency

Γ - VA = 2(1 -X res B )ν s A exp - E - VA k B T , (9)
where X res B = n res B /(n res A + n res B ) and

E - AV = (1 -X s B )Q sp AA + X s B Q sp AB -n s A ε AV -n s B ε BV -µ A (X s B ), (10)
where µ A (X s B ) is the chemical potential of atoms A, at the local composition X s B (defined as the B atomic fraction among the first and second nearestneighbors of the SAS site). The vacancy can also be replaced by an atom B, with a frequency similar to that given by eq. 9.

In a stable solid solution, above the miscibility gap, the equilibrium values of µ A and µ B can be estimated for any composition using an equilibrium Monte Carlo simulation performed in the semi grand canonical ensemble.

In case of a phase separation, the local configuration within a finite volume surrounding a SAS site can change. The most probable local environment in the initial solid solution is a disordered configuration with the nominal composition X B . When the phase separation proceeds, the most probable environments will correspond to A-rich or B-rich solid solutions (with local composition evolving towards the solubility limits, X α B and X β B ), or to A-rich/B-rich interface (the latter less and less frequently as the A-rich and B-rich domains grow). The actual chemical potential should therefore be somewhere between the values of a mean-field Bragg-Williams approximation (associated with a fully homogeneous and disordered solid solution), and the equilibrium values

µ A (X α B ) = µ A (X β B ) and µ B (X α B ) = µ B (X β B )
of the final two-phase equilibrium state, as given by the GCMC method (see Appendix B). In all the results presented here, we assume that the chemical potentials have always these equilibrium values. One can consider that this corresponds to a local equilibrium approximation, justified by the fact that the atomic mobility near the SAS sites is fast enough to reach a local equilibrium state. The main advantage of this choice is that it ensures that the simulated vacancy concentration tends toward the proper equilibrium vacancy concentration of the final equilibrium state. Note that, although not shown here, the same AKMC simulations performed with µ A (X s B ) and µ B (X s B ) set to their corresponding Bragg-Williams approximated values, lead to very small differences in the evolution of X eq V (t). From the equation (7-9), we easily demonstrate that the Gibbs free energy of vacancy formation only depends on the difference, E + AV -E - AV . The same difference appears in eq. 2 of the Widom insertion method, that we may use to compute the vacancy formation free energy. In AKMC simulations of dilute A-B alloys, one may assume that a SAS always operates in pure A, and this formation energy is a constant [26,[START_REF] Hin | Kinetics of heterogeneous dislocation precipitation of NbC in alpha-iron[END_REF]. The saddle-point parameters, Q sp AA and Q sp AB , are not involved in the average differences. Hence, they only affect the kinetics of vacancy concentrations, but not their equilibrium values. The kinetics of vacancy concentration also depends on the exponential attempt frequencies, ν s A and ν s B , and on the number N s of SAS sites. In principle, these parameters can be chosen to simulate a given SAS microstructure; that of a given dislocation network, for example. In the present study, we have simply chosen parameters (Appendix A) that ensure a rapid evolution of the vacancy concentration and a good statistics. If Q sp AA and Q sp AB are too small, the energies E - AV are significantly lower that the migration barriers and the vacancies are eliminated immediately after their formation: the simulation is then stuck in a rapid succession of formation and annihilation events and is enable to generate the time evolution of the microstructure in a reasonable CPU time. If Q sp AA and Q sp AB are too high, the formation and annihilation events are sparse, and X V cannot follow X eq V (t). Note also that the frequencies Γ - VA and Γ - VB are proportional to the atomic fractions of atoms A and B in the reservoir. This ensures that the composition of the system remains almost constant. Without this constraint, a significant shift of the composition may happen -especially at high temperatures, where the vacancy concentration is high. The factor 2 in eq. 9 comes from the fact that the replacement of the vacancy by atoms A or B, is mutually exclusive. Without this factor, the vacancy concentration evolves towards a steady-state value two times larger than the equilibrium one.

These additional events are handled by the residence time algorithm, in addition to the vacancy jumps and contribute to the residence time, t MC , according to their respective frequency. No time rescaling is required, since these events drive the vacancy concentration towards its physical value.

To assess our NCV-AKMC method and its SAS model, we first check that it yields the proper vacancy equilibrium concentrations, in homogeneous solid solutions as well as in two phase states.

Let us consider the steady-state values of vacancy concentrations in alloys AB, with various compositions, at 600 and 1000 K, i.e. respectively well below and just above T c . To determine the equilibrium vacancy atomic fraction, we first perform GCMC-Widom equilibrium Monte Carlo simulations (section 2.1.1). The steady-state values of the NCV-AKMC simulations are measured after a run sufficiently long to reach steady-state values of the short range order in the solid solution at 1000 K, or steady-state compositions and volume fractions of the phases α and β, at 600 K.

At 1000 K, the steady-state values of the NCV-AKMC simulations are in very good agreement with the equilibrium values predicted by the GCMC-Widom method (Fig. 2). Both give vacancy concentrations larger than those predicted by the following Bragg-Williams approximation [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF],

G for V (X B ) = (1 -X B )G for V (A) + X B G for V (B) + (zv/2)X B (1 -X B ), (11)
which neglects the short range order effects on the occupation probabilities of the lattice sites. The vacancy atomic fraction derived from the Lomer approximation ( [START_REF] Mehrer | Diffusion in Solids[END_REF], p. 82), which is valid in dilute solid solutions, is also shown for comparison (in the A-rich side, for example:

X eq V (X B ) = X eq V (0)[(1 -zX B ) + zX B exp(E b BV /k B T )]
, where E b BV is the B-V binding energy in pure A). At 600 K, the solid solution is stable between X α B = 0.03 and X β B = 0.97. The GCMC-Widom simulations only give the equilibrium values of X eq V in the stable or metastable solid solutions (between X B ≃ 0.05 and X B ≃ 0.95). For Fig. 2 Equilibrium vacancy concentration in A-B solid solutions at 600 and 1000 K: comparison between the GCMC-Widom and the NCV-AKMC methods. At 600 K, the red line is deduced from the equilibrium vacancy concentrations in the two coexisting phases and the lever rule (eq. 12)

these compositions, the steady state values extracted from the NCV-AKMC simulations are again in good agreement with the GCMC-Widom ones. Within the miscibility gap, the steady-state values of the NCV-AKMC simulations are the ones of the two phase-state. There are found to be in good agreement with the ones obtained from the lever rule applied to the equilibrium GCMC-Widom atomic fractions (the full red line in Fig. 2), as written below

X eq V (X B ) = X β B -X B X β B -X α B X eq V (X α B ) + X B -X α B X β B -X α B X eq V (X β B ). ( 12
)
At the temperature of the simulation, we observe that using the Bragg-Williams values of µ A (X s B ) and µ B (X s B ) in the AKMC simulations with sources, instead of the equilibrium values in the coexisting α and β phases, as discussed above, leads to a small underestimation of X eq V (by approximately 15%).

Simulation of precipitation upon isothermal annealing

Let us consider now the evolution of vacancy concentration during the precipitation of the B-rich phase β in a A-10%B alloy at 600 K. We perform two AKMC simulations, both in a system of N = 2 × 64 3 bcc sites, with periodic boundary conditions.

The first simulation is performed with N V = 1 vacancy, and the time is rescaled assuming that the vacancy concentration remains at equilibrium, using either the Widom or the local environment time rescaling methods, with the pure A or pure B environment as a reference. The respective results are referred as AKMC-Widom, AKMC-LETR(A), AKMC-LETR(B). Note that they all correspond to the same AKMC run and only differ by the time correction.

The second simulation is performed with the NCV-AKMC method, with a random distribution of N s = 64 SAS sites (this corresponds approximately to a dislocation network with a density ρ d ≃ 10 15 m -2 , if one considers that all the sites along the dislocation lines act as SAS).

Fig. 3 shows the evolution of the number per unit volume d p and the average radius R of the β precipitates (defined as clusters with more than 30 atoms B). The AKMC simulations performed with N V = 1 and with NCV-AKMC give very similar results, in which one can roughly distinguish the usual successive nucleation-growth and coarsening regimes. The minor differences are related to the estimation of the vacancy atomic fractions. Note the very long time scale (solid state kinetics are often very slow): one simulation corresponds to series of approximately 5 × 10 11 vacancy jumps, lenghty times are due to the fact the both the vacancy jump frequencies and the equilibrium vacancy concentrations are small at this temperature. This is of course very system dependant. To give a simple example, a decrease of all the migration barriers by 0.3 eV (giving an average value of approx. 0.7 eV instead of 1 eV) would accelerate the kinetics by a factor of 331 at 600 K. Irradiation can also accelerate the precipitation by orders of magnitudes, due to point-defects supersaturations (see e.g. the case of Cr precipitation in α-iron, under neutron irradiation [START_REF] Soisson | Radiation-accelerated precipitation in Fe-Cr alloys[END_REF]).

Figure 4 shows the evolution of X V in the NCV-AKMC simulation, and the one deduced from the time rescaling methods, AKMC-LETR(A) and AKMC-LETR(B). The results are similar (Fig. 2): X V starts from the equilibrium value in the initial random solid solution, X eq V ≃ 10 -11 , then increases when β precipitates nucleate and grow (due to the lower vacancy formation energy in the β phase), and stabilizes at X V ≃ 2.4 × 10 -10 after t > 10 8 h, in the coarsening stage, where the α and β phases reach their equilibrium composition and volume fraction. During the precipitation, X V increases by more than one order of magnitude. A realistic time scale is obtained, only if such increase of X V is not ignored. The simulated kinetics obtained from the NCV-AKMC method, shows that with the SAS density and parameters used here, X V is indeed close to equilibrium from the beginning of the precipitation (t > 2 × 10 4 h), as it is assumed in the conserved-vacancy simulation methods such as the time rescaling ones, although the values are significantly scattered up to t ≃ 10 3 h. The values are very close to those obtained with the AKMC-LETR(B) method and slightly above those obtained with the AKMC-Widom and AKMC-LETR(A) methods (by a factor ≃ 2 at the beginning of the precipitation).

The AKMC-LETR(B) method gives some results for t > 3 × 10 4 h only, because it requires a sufficient proportion of pure B local environments, which is significant only when the precipitation starts. It also gives more scattered results, than AKMC-LETR(A) or AKMC-Widom, because this proportion remains small (much smaller than that of pure A environments) in the early stages of precipitation. The small difference between AKMC-TR(A) and AKMC-TR(B) is due to a Gibbs-Thomson effect, described in details in ref. [START_REF] Nastar | Atomistic modeling of phase transformations: Point-defect concentrations and the time-scale problem[END_REF]: the curvature of the precipitate matrix interface affects the vacancy formation energies. It disappears at longer times, when the size of precipitates is sufficiently large (approx. R > 2 nm).

Fig. 4 AKMC simulations of precipitation in A-10%at.B at 600 K: evolution of the vacancy concentration.

Simulation of precipitation after a rapid quenching

A rapid quenching from high temperature is an usual method to produce nonequilibrium vacancies. We now simulate, using the NCV-AKMC method, the precipitation in the same alloy as before, again at T 1 = 600 K, but following an infinitely rapid quenching from a higher temperature T 0 . The simulation box initially contains 100 or 1000 vacancies (X V = 1.9 × 10 -4 or X V = 1.9 × 10 -3 ), corresponding to the equilibrium values at respectively T 0 = 1373 and T 0 = 1707 K.

Fig. 5 NCV-AKMC simulation of precipitation in A-10%at.B 600 K, after an infinitely rapid quench from 1707 K (with N V = 1000 vacancies) or 1373 K (with N V = 100): evolution of the density and radius of B-rich precipitates.

Fig. 6 NCV-AKMC simulation of precipitation in A-10%at.B at 600 K: evolution of the vacancy atomic fraction after a quench from 1707 K to 600 K Figures 5 and6 show two very distinct regimes. At the beginning of the simulation, the excess vacancies rapidly annihilate on the SAS sites: most of the decrease occurs for t < 10 -1 h. During this early stage, the vacancy concentration remains well above its 600 K equilibrium value, and diffusion is sufficiently rapid to observe the beginning of the the β-precipitate nucleation (after a quenching from 1373 K) or even the complete nucleation (after a quenching from 1707 K). The nucleation rate is thus increased by approximately 6 or 7 orders of magnitude with respect to an isothermal annealing with no prior rapid quenching. After 10 -1 h the vacancy concentration has almost reached its equilibrium value. Then nothing happens, because of the now much slower diffusion of atoms A and B, up to 10 4 -10 5 h, at which time the precipitation restarts, with growth and coarsening regimes very similar to the ones observed upon isothermal annealing. The coarsening regime in particular, is practically unaffected by the quenching thermal treatment, because it occurs over much larger times. The quenched-in excess vacancies have somehow dissociated the nucleation stage from the growth-coarsening stage.

Point defect induced kinetic phenomena

Although, PDs diffuse much faster than atoms, we have shown above that there are situations where kinetics of PDs toward equilibrium may interplay with the evolution of the microstructure. It happens whenever the characteristic time of PDs toward local equilibrium is close to the atom ones. PDs could then be locally or globally in supersaturation, leading to various nonequilibrium mechanisms, first one being the enhancement of atomic diffusion as investigated in section 2.5, second one being the activation of dislocation climbing and subsequent moving of the atomic planes, last ones involving PD supersaturation as an additional diffusion or precipitation driving force of PD sinks in the microstructure.

We present below, three kinetic phenomena involving the non-conservative nature of PDs: (a) interdiffusion with non-conservative PDs, (b) radiationinduced segregation at PD sinks, and (c) non-coherent precipitation of second phases driven by an excess of PDs. We will provide a short description of these phenomena, then try to answer the question about the capabilities of existing models, how they break and what an atomic approach could bring models in such circumstances.

Diffusion with non-conservative PDs

Kirkendall's matter transport in interdiffusion

In the famous Kirkendall's interdiffusion experiment, a Cu-Zn solid solution was sandwiched between two layers of pure copper. The sample was then placed in a furnace to be annealed at T = 785 • C for 45 days [START_REF] Kirkendall | Diffusion of Zinc in Alpha Brass[END_REF]. After thermal annealing, the initially pure layers became a solid solution of Cu-Zn; and the initial concentration of Cu in the solid solution of the central zone increased. So all the layers tend towards the same composition of a uniform solid solution.

From the decreasing rate of the distance between the Cu/Cu-Zn interfaces, Kirkendall could evaluate the speed at which the atomic planes are moving inward. He has then demonstrated that only a diffusion mechanism involving non-conservative vacancies, could explain the displacement of these planes. The net flow of vacancies caused by the difference of the diffusion velocities of the Zn and Cu atoms, generates local PD supersaturations or undersaturations, that are almost instantaneously suppressed by creation/annihilation reactions at uniformly distributed PD SAS. From this simplifying hypothesis of an infinite sink strength that guarantees the vacancy local equilibrium everywhere, the Kirkendall speed could be written as a difference of Cu and Zn intrinsic diffusion coefficients, and the resulting interdiffusion coefficient as a sum of these intrinsic diffusion coefficients weighted by the atomic fractions. Creation-annihilation reactions of point defects are most often on dislocations and/or grain-boundaries [START_REF] Lothe | Dislocation Climb Forces[END_REF][START_REF] Weertman | The Peach-Koehler equation for the force on a dislocation modified for hydrostatic pressure[END_REF], but when these extended defects are scarce or not uniformly distributed, or the PD/defect reactions are prevented or too slow, it happens that the relaxation rate of PDs toward local equilibrium becomes comparable to that of the microstructure. Therefore, when sink strength cannot be considered to be infinite, PDs may take other reactive paths, that may lead to different equilibrium states. Backward diffusion of PDs is an alternative path. By slowing down the interdiffusion process, it gives time to PDs to reach a local equilibrium. When this reactive path operates alone, PDs behave like a conservative species with its own diffusion driving force. Indeed, the PD chemical potential gradient generates a second diffusion driving force acting against the interdiffusion mixing driving force, leading to an effective interdiffusion mobility controlled by the slowest species [START_REF] Qin | Computer simulation of chemical diffusion in a binary alloy with an equilibrium concentration of vacancies[END_REF][START_REF] Gusak | Nonequilibrium vacancies in nanosystems[END_REF][START_REF] Nastar | Atomic diffusion theory challenging the Cahn-Hilliard method[END_REF]. In nanoscale interdiffusion processes such as the ones operating in the early stage of a spinodal decomposition, we expect the kinetics to be controlled by the conservative PD mobility [START_REF] Gusak | Nonequilibrium vacancies in nanosystems[END_REF]. Another reaction path is the clustering of PDs into voids or dislocation loops, who themselves become SAS. This reaction can be classified as a mechanism of precipitation of PD SAS. In many interdiffusion experiments as the one shown in Fig. 7, voids are formed. This is a direct proof that vacancies do not fully behave as a nonconservative species. A fraction of these vacancies, instead of being eliminated, agglomerate to form voids.

There have been rate theory models accounting for partial elimination/creation of PDs at sinks [START_REF] Martin | Diffusion in crystals with nonconservative defects[END_REF][START_REF] Svoboda | Diffusion and creep in multicomponent alloys with non-ideal sources and sinks for vacancies[END_REF]. The solid is modeled as a 1D-stack of slabs along the diffusion direction, whose thicknesses are governed by a continuum variation of the slab's lattice site concentrations, resulting from local PD elimination/reactions. The diffusion/reaction equations modeling the evolution of PDs and atomic species concentrations together with the slab thickness, include a PD reaction term equal to the product of a local density of sinks, PD diffusion coefficient, and an elimination/creation driving force written as a difference between the PD concentration and its equilibrium value. Such a rate theory model when employed to investigate interdiffusion in alloys, leads to interdiffusion profiles depending on the density of PD sinks [START_REF] Martin | Diffusion in crystals with nonconservative defects[END_REF], and moving alloy/oxide interfaces [START_REF] Desgranges | Depletion and Voids Formation in the Substrate During High Temperature Oxidation of Ni-Cr Alloys[END_REF]. In these studies, interdiffusion is treated at the atomic scale but the internal stress field generated by non uniform elimination/creation reactions of PDs is ignored. PDs with their relaxation volumes, local changes of composition, together with the accumulation or depletion of matter resulting from elimination/creation of lattice sites, may lead to an internal strain field affecting diffusion and vice versa [START_REF] Svoboda | Diffusion and creep in multicomponent alloys with non-ideal sources and sinks for vacancies[END_REF][START_REF] Garikipati | A lattice-based micromechanical continuum formulation for stress-driven mass transport in polycrystalline solids[END_REF].

High temperature diffusional creep

The mechanism of diffusional creep occurring at high temperature, is the transport of matter resulting from the coupled fluxes of lattice sites and PDs induced by an applied stress. In case of vacancy-mediated diffusion, the flow of lattice sites (which are infinitely small matter pieces) is opposite to the flux of vacancies. Creation/annihilation reactions of PDs accommodate the applied stress by inducing a net flow of matter mediated by PD diffusion, that leads to macroscopic deformation [START_REF] Herring | Diffusional viscosity of a polycrystalline solid[END_REF]. Whenever the applied stress generates a gradient of internal stress, in particular between different SAS, the PD diffusion driving force is non zero, potentially leading to a net flux of atoms and PDs between SAS. This mechanism operates at sufficiently high temperatures for PD long-range diffusion to be significant. The modeling of PD diffusion under applied stress, requires an explicit treatment of the diffusion/mechanical coupling and associated motion of PD sinks following the matter flow, such as surfaces and grain-boundaries [START_REF] Garikipati | A lattice-based micromechanical continuum formulation for stress-driven mass transport in polycrystalline solids[END_REF][START_REF] Kovacevic | Diffusion-induced stress concentrations in diffusional creep[END_REF], or dislocations [START_REF] Mesarovic | Dislocation creep: Climb and glide in the lattice continuum[END_REF][START_REF] Wang | Vacancy diffusion coupled discrete dislocation dynamic modeling of compression creep of micro-pillars at elevated temperature[END_REF]. First model was the analytic model of Nabarro-Herring: the sources and sinks of vacancies were grain boundaries of a rectangular grain [START_REF] Herring | Diffusional viscosity of a polycrystalline solid[END_REF]. Finite element methods solving the coupled diffusion/mechanical equations with moving boundaries, were able to simulate the formation of nanoscale heterogeneities of the stress-strain field and PD concentration, and investigate the effect of stress, temperature, grain size and the density and mobility of dislocations, on creep.

Constitutive laws and diffusion equations reside in a local equilibrium hypothesis. The reaction/diffusion driving forces are the thermodynamic ones with strain and composition set to their local values in every unit volume. Continuous diffusion equations implicitly impose a constraint of homogeneity of composition and stress in unit volumes, and ignore the thermal fluctuations of unit volume composition, which prevents simulating the nucleation of PD clusters such as voids. Hybrid approach consisting of simulating the evolution of PD clusters within unit volumes by means of a cluster dynamic method, and the composition field by means of a phase field method, could be the approach to follow in the near future [START_REF] Badillo | A phase field model for segregation and precipitation induced by irradiation in alloys[END_REF]. Besides, elasto-diffusion contributions to PD diffusion should not be neglected, while not only they affect the diffusion rates but also may generate anisotropies of PD fluxes, leading for instance to preferential directions of faceted-void growth [START_REF] Carpentier | Effect of saddle point anisotropy of point defects on their absorption by dislocations and cavities[END_REF].

As long as we may ignore the heterogeneities of the internal stress/strain field, atomic Monte Carlo simulation methods including an explicit treatment of non-conservative PDs, could be developed to simulate these diffusioncontrolled phenomena. However, a big challenge remaining in this field, lies in the difficulty of dealing with diffusion mechanisms and reactions of PDs in strongly distorted lattices nearby climbing dislocations or wide angle grainboundaries. Besides, atomistic simulations can hardly take into account the flow of matter resulting from the creation/annihilation of lattice sites, unless SAS is a free surface that is explicitly treated in the simulation. In pure solid, off-lattice molecular dynamics method is the appropriate method whenever PD diffusion is fast enough to be simulated within pico-seconds [START_REF] Kabir | Predicting dislocation climb and creep from explicit atomistic details[END_REF]. Analytical models can be used to extrapolate the simulation results to realistic dislocation densities [START_REF] Clouet | Predicting dislocation climb: Classical modeling versus atomistic simulations[END_REF]. Note that molecular dynamics methods cannot tackle alloying effects mediated by a PD diffusion mechanism, because the atomic diffusion timescale (inversely proportional to the concentration of PDs) is generally not accessible.

Diffusion and precipitation induced by non-equilibrium PDs

Atomic diffusion mediated by PDs is enhanced whenever PDs are in excess, because atomic diffusion coefficients are proportional to the PD concentration. In bulk materials supersaturated in PDs, the chemical potential of PD is non-zero and positive. Whereas, at PD sinks, due to the creation-annihilation reactions ensuring local equilibrium, the local chemical potential of PD is null.

Non-uniform PD chemical potential generates a diffusion driving force equal to the gradient of PD chemical potential divided by temperature [START_REF] Mehrer | Diffusion in Solids[END_REF]53]. Below, we give an example of local chemical redistribution driven by PD chemical potential gradients, and an example of phase transformation triggered by a non-zero PD chemical potential.

Radiation-induced segregation

Irradiation creates excess PDs as Frenkel pairs (vacancies and self-interstitial atoms), which are eliminated by mutual recombination, clustering or annihilation/creation reactions at pre-existing defects such as surfaces, grain boundaries, and dislocation lines, or nucleated defects such as dislocation loops or voids. As a result, permanent irradiation generates gradients of PD chemical potential, leading to net fluxes of PDs toward PD SAS and atom fluxes by an inverse Kirkendall mechanism. In case of any preferential transport of one of the alloy components proceeding by exchanges with nearby PDs, fluxes of PDs lead to a local chemical redistribution at SAS. Under irradiation, this is the so-called radiation-induced segregation (RIS) phenomenon. It almost systematically occurs at SAS, and has important technological implications [START_REF] Nastar | Radiation-induced segregation[END_REF][START_REF] Russell | Phase stability under irradiation[END_REF][START_REF] Ardell | Radiation-induced solute segregation in metallic alloys[END_REF]. RIS was first predicted by Anthony from the analysis of non-equilibrium segregation of various solute elements on pyramidal cavities formed in aluminum after quenching from high temperature [START_REF] Anthony | Solute segregation in vacancy gradients generated by sintering and temperature changes[END_REF]. Anthony suggested that similar flux coupling phenomena should produce non-equilibrium segregation in alloys under irradiation. He predicted that the segregation should be much stronger than after quenching, because under irradiation, not only the excess of PDs is sustained for very long times but also, in addition to vacancies, self-interstitials are in excess. Fig. 8 shows how an atom probe tomography analysis is able to highlight a self-interstitial clustering mechanism induced by irradiation leading to the formation of dislocation loops and a positive radiation-induced segregation of Ni around the loops [START_REF] Huang | Impact of the local microstructure fluctuations on radiation-induced segregation in dilute feni and ni-ti model alloys: A combined modeling and experimental analysis[END_REF].

Finite-difference rate theory models [START_REF] Huang | Impact of the local microstructure fluctuations on radiation-induced segregation in dilute feni and ni-ti model alloys: A combined modeling and experimental analysis[END_REF][START_REF] Huang | Multiscale modeling of the effects of temperature, radiation flux, and sink strength on point-defect and solute redistribution in dilute Fe-based alloys[END_REF] or phase field methods [START_REF] Piochaud | Atomicbased phase-field method for the modeling of radiation induced segregation in Fe-Cr[END_REF][START_REF] Thuinet | Multiscale modeling of radiation induced segregation in iron based alloys[END_REF], and atomic kinetic Monte Carlo (AKMC) methods have become efficient tools to simulate RIS [START_REF] Senninger | Modeling radiation induced segregation in iron-chromium alloys[END_REF][START_REF] Martínez | Role of sink density in nonequilibrium chemical redistribution in alloys[END_REF]. Thanks to high throughput ab initio investigations of PD formation energies and jump frequencies with respect to the local chemical environment, these simulations could provide a fine description of the PD thermodynamic and diffusion properties, as well as flux couplings of a large number of dilute Fe-, Zr-and Ni-based binary alloys [START_REF] Messina | Systematic electronicstructure investigation of substitutional impurity diffusion and flux coupling in bcc iron[END_REF][START_REF] Jain | First-principles calculations of solute transport in zirconium: Vacancy-mediated diffusion with metastable states and interstitial diffusion[END_REF][START_REF] Toijer | Solute-point defect interactions, coupled diffusion, and radiation-induced segregation in fcc nickel[END_REF]. Radiation induced precipitation (RIP) of particles subsequent to RIS, even though the parent phase is not supersaturated, currently occurs. RIP is triggered by the formation of a local supersaturation in solute at SAS, which brings the system into a two-phase domain of the phase diagram. Stochastic AKMC method accounting for the fluctuations of lattice occupancy, is the appropriate method to simulate radiation-enhanced precipitation, as for example the coherent precipitation of Cr-rich α ′ particles in Fe-Cr [START_REF] Soisson | Radiation-accelerated precipitation in Fe-Cr alloys[END_REF][START_REF] Soisson | Atomistic modeling of α precipitation in Fe-Cr alloys under charged particles and neutron irradiations: Effects of ballistic mixing and sink densities[END_REF], or radiationinduced precipitation subsequent to RIS at grain-boundary [START_REF] Senninger | Modeling radiation induced segregation in iron-chromium alloys[END_REF]. Note that annihilation of PDs at SAS was modeled as an irreversible reaction. Such treatment is valid as long as the thermal equilibrium PD concentration is negligible compared to the bulk PD concentration, which is generally the case at temperatures below 0.6 times the melting point [START_REF] Russell | Phase stability under irradiation[END_REF]. Above 0.6 times the melting point, RIS and other radiation-induced phenomena are significantly reduced, because the equilibrium concentration of vacancy is too high to generate a significant vacancy supersaturation leading to PD fluxes. Methods accounting for the competition between annihilation and creation reactions at sinks and subsequent suppression of RIS at high temperature, are mostly mesoscopic kinetic methods [START_REF] Huang | Impact of the local microstructure fluctuations on radiation-induced segregation in dilute feni and ni-ti model alloys: A combined modeling and experimental analysis[END_REF][START_REF] Huang | Multiscale modeling of the effects of temperature, radiation flux, and sink strength on point-defect and solute redistribution in dilute Fe-based alloys[END_REF][START_REF] Piochaud | Atomicbased phase-field method for the modeling of radiation induced segregation in Fe-Cr[END_REF][START_REF] Thuinet | Multiscale modeling of radiation induced segregation in iron based alloys[END_REF].

The present modeling of PD SAS should allow the study of RIS or any diffusion-controlled phenomena at temperatures where thermal creation reactions of PDs at SAS cannot be ignored. Besides, the modeling of competitive diffusion mechanisms such as vacancy, self-interstitial and direct interstitial diffusion mechanisms, requires an explicit modeling of the thermal creation/annihilation reactions of every PD at SAS, since we cannot proceed to two different scaling of time in the same simulation.

PD-induced non-coherent precipitation

In alloys super-saturated in PDs, concurrent nucleation paths compete with formation of PD clusters. An excess of PDs by generating supplementary precipitation driving forces, may trigger phase transformations in undersaturated parent phase, or change the selection of phases, or even lead to the formation of new phases. It happens for instance in systems submitted to irradiation [START_REF] Russell | The theory of phase stability under irradiation[END_REF][START_REF] Nastar | Thermodynamic model for lattice point defect-mediated semi-coherent precipitation in alloys[END_REF][START_REF] Garner | Spinodal-like decomposition of fe-35ni and fe-cr-35ni alloys during irradiation or thermal aging[END_REF], mechanically solicited [START_REF] Militzer | Modelling the effect of deformationinduced vacancies on segregation and precipitation[END_REF], or thermally quenched-in [START_REF] Khachaturyan | Theory of Structural Transformations in Solids[END_REF][START_REF] Bourgeois | Transforming solid-state precipitates via excess vacancies[END_REF], but also as a consequence of surface reactions. Experimental studies have clearly established that the complex behavior of oxide precipitation in Si was involving both the injection of self-interstitials and subsequent oxide precipitation at the substrate/oxide interface triggered by the excess of self-interstitials [START_REF] Hu | Nonequilibrium point defects and diffusion in silicon[END_REF].

In case of a large attractive PD-solute binding energy and solute concentration in the order of magnitude of the PD concentration, irradiation-induced PDs may act as a chemical component [START_REF] Schuler | Vacancy-induced dissolution of precipitates in out-of-equilibrium systems: A test case of FeX (X = C, N, O) alloys[END_REF]. In iron, due to the large attraction between vacancy and oxygen, we expect an excess of vacancy stabilizes more oxygen atoms in the parent phase, thereby a vacancy-induced dissolution mechanism of oxide particles at intermediate temperature and high radiation flux [START_REF] Schuler | Vacancy-induced dissolution of precipitates in out-of-equilibrium systems: A test case of FeX (X = C, N, O) alloys[END_REF]. Growing semi-coherent or incoherent particles may act as efficient PD sinks. To accommodate the precipitate eigenstrain, radiation-induced selfinterstitials may be transformed into substitutional atoms of the precipitate when negative eigenstrain (undersized particle) [START_REF] Nastar | Thermodynamic model for lattice point defect-mediated semi-coherent precipitation in alloys[END_REF], or vacancy are removed when positive eigenstrain (oversized particle) [START_REF] Maydet | Precipitate stability under irradiation: Point defect effects[END_REF][START_REF] Russell | The theory of phase stability under irradiation[END_REF]. For instance, a joint agregation of self-interstitials and Ni atoms in bcc Fe-Ni, yields the precipitation of the unexpected face-centered-cubic phase γ [START_REF] Nastar | Thermodynamic model for lattice point defect-mediated semi-coherent precipitation in alloys[END_REF][START_REF] Belkacemi | Radiation-induced bcc-fcc phase transformation in a Fe -3Ni alloy[END_REF]. Fig. 9 is an illustration of the combined modeling-experimental investigation that could highlight the role of radiation-induced self-interstitials on the precipitation of the undersized austenite particles γ.

The transformation of PDs into substitutional atoms is both a powerful mechanism of PD annihilation and an accommodation mechanism of the precipitate-matrix volume-mismatch. PD contribution to the precipitation driving force is directly related to the excess energy produced by nonequilibrium PDs, e.g. the PD chemical potential. The energy released upon the growth of a non-coherent precipitate, is equal to the fraction of PDs required to compensate the volume mismatch times the PD chemical potential. We may derive non-equilibrium Gibbs free energy of the precipitating phase by just adding to the standard precipitation driving force, this PD released energy, then compute constrained phase diagrams including additional coordinate axes for the PD supersaturation [START_REF] Nastar | Thermodynamic model for lattice point defect-mediated semi-coherent precipitation in alloys[END_REF].

Such a formulation can be easily combined with thermodynamic CAL-PHAD databases commonly used to compute phase diagrams [START_REF] Nastar | Thermodynamic model for lattice point defect-mediated semi-coherent precipitation in alloys[END_REF], or else to cluster dynamics methods [START_REF] Maydet | Precipitate stability under irradiation: Point defect effects[END_REF][START_REF] Russell | The theory of phase stability under irradiation[END_REF]. The modeling of the competition between void and oxide precipitation in PD supersaturated metallic or semiconductor alloys could be performed by relying on the classical nucleation theory [START_REF] Desgranges | Depletion and Voids Formation in the Substrate During High Temperature Oxidation of Ni-Cr Alloys[END_REF][START_REF] Bobeth | Modelling of vacancy diffusion and pore formation during parabolic oxide growth[END_REF][START_REF] Voronkov | Vacancy-type microdefect formation in czochralski silicon[END_REF]. However, these phenomenological methods do not inform on the nucleation mechanism of semi-coherent or incoherent precipitates in interaction with an excess of PDs. As for the simulation of climbing dislocations, molecular dynamics methods should provide some insight on the role of excess PDs in non-coherent off-lattice precipitation.

Conclusion

We have presented a short review of technologically important non-equilibrium phenomena that may be affected by the fact that PDs are non-conserved species, which are not always at equilibrium. Even in the simple situation where one can assume that they remain at local equilibrium during a phase transformation, it is essential to take into account the evolution of their concentration to predict physical time scales. Taking into account the evolution of PD concentrations, their formation and elimination mechanisms, becomes even more important in systems that are temporarily or permanently far from equilibrium. Only an explicit treatment of the PDs concentrations allows to investigate complex kinetic competitions between concurrent clustering reactions in highly PD supersaturated materials, including multiple diffusion mechanisms.

We proposed a first attempt to explicitly include the modeling of creation and annihilation reactions of point defects in atomistic kinetic simulation methods. We have shown simple applications to dilute and concentrated alloys, from the simulation of the equilibrium vacancy concentration in equilibrium and evolving non-equilibrium microstructures and the kinetics of a coherent phase separation in a model A-B alloy upon isothermal annealing after a rapid quench. Like the previous methods, that were based on the assumption that X V = X eq V , this NCV-AKMC method is able to predict a strong increase of X V during the precipitation when the vacancy formation of vacancies is significantly smaller in B-rich than in A-rich metals. But this method can also deal with the evolution of far-from-equilibrium vacancy concentrations arising after a rapid quench. The model still relies on strong approximations concerning the sink/source spatial distribution, the modeling of the formation and elimination barriers, and the definition of the local unit volume and associated chemical potentials. Even though, the strain effects resulting from the non-conservation of lattice sites are not easy to simulate at the atomic scale, the present NCV-AKMC method should allow investigating the competition of the many reaction paths of the point defects in alloys, including diffusion mechanisms and PD clustering events. Such simulation methods should be even more useful when multiple diffusion mechanisms control the microstructure. Therefore, the present Monte Carlo method could be directly applied to simulate PD induced kinetic phenomena such as radiation-induced segregation and coherent precipitation.

Apart from these direct applications, developments of the method to introduce a simplified modelling of the effect of PD creation/annihilation reactions on the variation of the lattice site number, the motion and transformation of SAS, should be possible, at least when SAS is a free surface or a planar grain-boundary delimiting the simulation box. These developments would allow first atomic-scale studies of the contribution of non-equilibrium PDs to interdiffusion and diffusional creep.

Fig. 1

 1 Fig. 1 chematic energetic landscape of the vacancy near a SAS site s (see Eq.7 to 10). The various energies contributing to the formation of a vacancy in replacing atom X (X = A, B) are represented: E + XV is the formation barrier, E - XV the elimination barrier, E for V (X) the local formation energy. The saddle-point contributionQ sp = (1 -X s B )Q sp XA + X s B Q spXB determines the height of the formation and elimination barriers, but not the formation energy. The migration barrier for the exchange between a vacancy on site i and a atom Y on site i is also represented. The height of the barriers depends on the nature of the atom to be replaced (Y = A, B) and on the local composition.

Fig. 3

 3 Fig. 3 AKMC simulations of precipitation in A-10%at.B at 600 K: evolution of the density and radius of B-rich precipitates (clusters containing more than 30 B atoms).

Fig. 7

 7 Fig. 7 Backscatter electron micrograph of Kirkendall voids formation resulting from the clustering of non-equilibrium vacancies, in the Cu vs. Ni interdiffusion couple annealed at 1000 °C for 48 h [37]

Fig. 8

 8 Fig. 8 Atom Probe Tomography analysis of nickel (in green) distribution in ion-irradiated Fe-Ni at T = 670 K (in the middle) and the uni-dimensional depletion profiles crossing segments of the potential dislocation loops formed upon the clustering of non-equilibrium self-interstitials [56]. The square symbols of the uni-dimensional profiles represent the measured local Ni composition, and the shaded areas indicate the uncertainty domain.

Fig. 9

 9 Fig.9High-Resolution Transmission Electron Microscopy image of a face-centered cubic precipitate γ in bcc Fe-Ni ion-irradiated at the JaNNuS -CEA Saclay facility, and simulated atomic planes of the Kurdjumov-Sachs oriented precipitate. Red lines are predicted misfit dislocations. Self-interstitials (colored disks) are transformed into substitutional atoms of the precipitate to accommodate the variation of volume between the initial bcc and the secondary fcc phases[START_REF] Nastar | Thermodynamic model for lattice point defect-mediated semi-coherent precipitation in alloys[END_REF].

Acknowledgements

We thank Chu-Chun Fu, Kangming Li and Thomas Schuler for fruitful discussions.

Data availability

The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.

Conflict of interest statement

On behalf of all authors, the corresponding author states that there is no conflict of interest.

Appendix A: AKMC parameters

The interatomic pair interactions of the A-B system are between first-nearest neighbor lattice sites only. They are equal to ε AA = -1.05, ε BB = -1.05 and ε AB = -1.025 eV, which gives an ordering energy equal to v = ε AA + ε BB -2ε AB = -0.05 eV. The atom-vacancy pair interactions are set to ε AV = -0.40 and ε BV = -0.35 eV.

When an atom is at the saddle-point, during its exchange with a vacancy, it interacts with with the 6 nearest neighbors. The corresponding pair interactions are equal to ε sp AA = -1.500 and ε sp AB = -1.55 eV for a A atom, and ε sp BA = -1.500 and ε sp BB = -1.55 eV for a B atom. These interaction parameters control the migration barriers. For instance, they give 1.15 eV for an A-V exchange in pure A, 1.025 eV for a B-V exchange in pure A, 1.25 eV for a B-V exchange in pure B, and 1.025 eV for a A-V exchange in pure B.

The saddle-point parameters that control the kinetics of vacancy creation/annihilation reactions at source/sink sites are set to Q sp AA = Q sp BA = -6.6, Q sp AB = Q sp AB = 6.7 eV for the simulations at 1000 K. The same parameters are all set to -6.35 eV for the simulations at 600 K. The chemical potentials of atoms A and B in binary alloys AB of different compositions at T = 600 K are plotted in Fig. 10. The circles are the values extracted from Monte Carlo simulations in the semi-grand canonical ensemble (GCMC). In GCMC simulations, the difference ∆µ = µ B -µ A is imposed, the values of µ A and µ B are then deduced from the Helmholtz free energy: F = (1 -X B )µ A + X B µ B , obtained by integration of µ B -µ A with respect to the atomic fraction of B. The values predicted by the Bragg-Williams approximation (regular solution model, see e.g. [START_REF] Porter | Phase Transformations in Metals and Alloys[END_REF]) are shown for comparison.

Appendix B: Chemical potentials

The GCMC simulations only give the chemical potentials of stable or metastable solid solutions (X B ≲ 0.05 and X B ≳ 0.95 in Fig. 10). The values of µ A and µ B used in the NCV-AKMC simulations are set equal to the GCMC measured values as long as the composition of the local finite volume centered on the SAS site, X s B (fraction of B atoms among the first and second nearest neighbors of the SAS site) is the one of a stable solid solution (

, the chemical potentials are set to their equilibrium two-phase values, i.e.:

) (corresponding to the dashed-dotted line in Fig. 10). For the equilibrium compositions X α B and X β B , the differences between the Monte Carlo and the Bragg-Williams chemical potentials are relatively small (∼ 0.008 eV), but leads to differences of approx. 15% on the equilibrium vacancy concentration, X eq V , at 600 K.