OKMC simulation of vacancy-enhanced Cu solute segregation affected by temperature/irradiation in the Fe-Cu system

OKMC simulation of vacancy-enhanced Cu solute segregation affected by temperature/irradiation in the Fe-Cu system 1 2 first-author 3 4 corresp Sha-Sha Lv lvss@bnu.edu.cn Sha-Sha Lv lvss@bnu.edu.cn lvss@bnu.edu.cn https://orcid.org/0000-0002-5219-7349 5 6 2 corresp De-Sheng Ai ai_desh@tsinghua.edu.cn De-Sheng Ai ai_desh@tsinghua.edu.cn ai_desh@tsinghua.edu.cn 1 corresp Zheng-Cao Li zcli@tsinghua.edu.cn Zheng-Cao Li zcli@tsinghua.edu.cn zcli@tsinghua.edu.cn Key Lab of Advanced Materials (MOE), School of Materials Science and Engineering, Tsinghua University, Beijing 100084, China Key Lab of Advanced Materials (MOE), School of Materials Science and Engineering, Tsinghua University, Beijing 100084, China Nuclear and New Energy Institute, Tsinghua University, Beijing 100084, China Nuclear and New Energy Institute, Tsinghua University, Beijing 100084, China Department of Engineering Physics, Tsinghua University, Beijing 100084, China Department of Engineering Physics, Tsinghua University, Beijing 100084, China Key Laboratory of Beam Technology, Ministry of Education, College of Nuclear Science and Technology, Beijing Normal University, Beijing 100875, China Key Laboratory of Beam Technology, Ministry of Education, College of Nuclear Science and Technology, Beijing Normal University, Beijing 100875, China School of Materials Science and Engineering, Shanghai Jiao Tong University, Shanghai 200240, China School of Materials Science and Engineering, Shanghai Jiao Tong University, Shanghai 200240, China Department of Nuclear Engineering and Management, The University of Tokyo, Tokyo 113-8656, Japan Department of Nuclear Engineering and Management, The University of Tokyo, Tokyo 113-8656, Japan The effects of annealing and irradiation on the evolution of Cu clusters in α-Fe are investigated using object kinetic Monte Carlo simulations. In our model, vacancies act as carriers for chemical species via thermally activated diffusion jumps, thus playing an important role in solute diffusion. At the end of the Cu cluster evolution, the simulations of the average radius and number density of the clusters are consistent with the experimental data, which indicates that the proposed simulation model is applicable and effective. For the simulation of the annealing process, it is found that the evolution of the cluster size roughly follows the 1/2 time power law with the increase in radius during the growth phase and the 1/3 time power law during the coarsening phase. In addition, the main difference between neutron and ion irradiation is the growth and evolution process of the copper-vacancy clusters. The aggregation of vacancy clusters under ion irradiation suppresses the migration and coarsening of the clusters, which ultimately leads to a smaller average radius of the copper clusters. Our proposed simulation model can supplement experimental analyses and provide a detailed evolution mechanism of vacancy-enhanced precipitation, thereby providing a foundation for other elemental precipitation research. Object kinetic Monte Carlo Irradiation effect Solute segregation Reactor pressure vessel Introduction 1 The formation of solute-decorated clusters is one of the three main damage mechanisms in reactor pressure vessel (RPV) steels under irradiation, and this process is sensitive to the structure and composition of RPV steels. Generally, the alloying elements in RPV steels of A508-III mainly comprise C (0.19 wt.%), Mn (1.20–1.43 wt.%), Ni (0.73–0.79 wt.%), Cr (0.06–0.12 wt.%), Mo (0.48–0.51 wt.%), and Cu (0.034–0.070 wt.%). In the past, the formation of Cu-rich clusters in iron has been widely studied experimentally [ 1 1 - 2 2 1-2 ]. Many studies have found that the formation of these clusters under service temperatures and neutron irradiation is a major cause of mechanical property degradation in RPV steels, such as hardening and embrittlement, which limits the lifespan of nuclear power plants. In earlier studies, it was observed that the precipitation of Cu solute in RPV steel could induce enrichment of the Mn/Ni/Si content; thus, the Cu content in RPV steel is strictly controlled. However, the influence of Cu on the evolution of matrix defects is of great significance for studies of the radiation-induced segregation process under irradiation; more importantly, it provides a reference for the subsequent study of the solute cluster precipitation of other elements. Therefore, vacancy-enhanced Cu solute segregation affected by temperature/irradiation in Fe-Cu systems has often been studied. Using the positron annihilation technique, the enhanced diffusion of Cu was evidenced; meanwhile, vacancy dissociation from the Cu-vacancy cluster and the combination of clusters was observed [ 3 3 - 4 4 3-4 ]. Furthermore, Jin et al. investigated the effects of the irradiation dose and irradiation depth on the Cu cluster formation process using positron annihilation spectroscopy [ 5 5 ]. However, considering the difficulty of neutron irradiation experiments, computer simulations have been used as an alternative to neutron irradiation experimental studies [ 6 6 - 7 7 6-7 ], and these simulations can demonstrate the evolution of nanostructures under irradiation conditions and evaluate the mechanical properties of materials. Among common simulation methods, the Monte Carlo method is widely used for the analysis of reactor physics owing to its advantages in solving complex problems. Moreover, with the advancement of computer technology, its application for simulations of various physical reactor processes is increasing [ 8 8 - 10 10 8-10 ]. In fact, atomic kinetic Monte Carlo (AKMC) and object kinetic Monte Carlo (OKMC) methods allow large time-scale simulations over large spatial regions, up to the life span of a reactor vessel (30–40 y), and the microstructural evolution under irradiation in ferritic steels can be analyzed. AKMC is a type of rigid lattice model that can describe atomic motion, but crystallographic changes cannot be simulated. However, for OKMC, the atomic structures of the lattice and clusters are simplified and replaced by abstract objects, which are used to describe the properties of different types of defect structures through corresponding parameters. The AKMC model is based on the first neighbor jumps in α-Fe, thereby reproducing in detail the diffusion of defects or solutes, which leads to the formation of voids and precipitates. An introduction to the main features and different applications of the AKMC model is available in reference [ 11 11 ]. Correspondingly, OKMC simulations treat all defects as objects, including point defects (i.e., vacancies and self-interstitial atoms (SIAs)) and solute-decorated matrix defects (i.e., Cu-vacancy clusters). The characteristics of these objects, such as their shape, type, size, and position, are fully abstracted and parameterized. During the simulation, these objects follow several possible events, mainly migrating and dissociating or reacting with other objects, i.e., merging into larger clusters or annihilating. Soisson and Fu obtained AKMC simulation results that were consistent with experiments [ 12 12 ]. Considering the strong interaction between vacancies and copper clusters, their model simulated a series of vacancy transitions on the surface of copper clusters, which resulted in movement of the copper clusters. It was confirmed that the movement of clusters was an important mechanism for the precipitation of copper clusters in α-Fe. Furthermore, Castin and Pascute combined AKMC and OKMC simulations and also obtained good results consistent with those obtained experimentally: when the size of one copper cluster was larger than a certain threshold, the atomic structure was ignored, and the cluster simulation performed as an object [ 13 13 ]. This led them to extend the simulation scale to the coarsening state, allowing the simulation results to be compared with broader experimental results. Based on these previous studies, an OKMC simulation is proposed to describe the evolution of copper cluster segregation under irradiation in Fe-Cu systems. In this model, we further parameterize the properties of small clusters and apply a new method to define the factor for rescaling the simulation annealing time when dealing with the annealing process. This method enables us to simulate Cu segregation in α-Fe under irradiation. The simulation model in this study is built based on the framework of an open-source code (MMonCa) [ 14 14 ], and the fundamentals of our algorithm are described in the first part of this paper. Then, the evolution of Cu clusters in the Fe-Cu system during the annealing process is successfully simulated using the proposed model. Finally, we conduct a comparative study of the precipitation of Cu-vacancy clusters under neutron irradiation and ion irradiation, and the main evolutionary characteristics of these clusters are analyzed under different conditions. Object kinetic Monte Carlo simulation method 1 OKMC model 2 The kinetic Monte Carlo algorithm assumes that the system transitions randomly between different states. The transition is a Markovian process in which the transition probability, r ij , depends only on the properties of the initial state, i, and the final state, j. The transition probability, r ij , is the basic input parameter of the Monte Carlo algorithm. In this study, the transition probability, r ij , is calculated using the Arrhenius formula [ 15 15 ]: <math display="block" id="M1"><mrow><msub><mi>r</mi><mrow><mi>i</mi><mi>j</mi></mrow></msub><mo>=</mo><msub><mi>P</mi><mrow><mi>i</mi><mi>j</mi></mrow></msub><mo>×</mo><mi>exp</mi><mo stretchy="false">(</mo><mo>−</mo><msub><mi>E</mi><mrow><mi>i</mi><mi>j</mi></mrow></msub><mo>/</mo><msub><mi>k</mi><mtext>B</mtext></msub><mi>T</mi><mo stretchy="false">)</mo><mo>,</mo></mrow></math> where the pre-factor, P ij , is the frequency of attempts, and E ij is the activation energy when the system transitions from state i to state j. After the transition probabilities, r ij , of all possible states in the system are defined, we can perform an OKMC simulation using the residence time algorithm, which is the most commonly used algorithm in OKMC simulations [ 16 16 , 17 17 ]. Basic assumptions 2 Generally, we assume that copper clusters migrate in the iron matrix through a vacancy mechanism, which implies a series of vacancy jumps at the surface of the copper clusters. However, because of the foundational assumptions of the OKMC simulation, the atomic structure of the cluster is simplified and the clusters are treated as independent objects with different properties, such as their size, trapping radius, attempt frequency, activation energy, and binding energy functions. Thus, a series of macro-equivalent descriptions of the properties of the cluster are needed, among which the two most important parameters are the diffusion coefficient and the migration coefficient; both of these parameters can be calculated from the AKMC simulation results of previous studies. For example, for the migration capacity of copper clusters, a vacancy and a pure copper cluster were placed in the simulation at a certain temperature, and the migration process of the vacancy in the lattice site was simulated; the calculation was stopped when the cluster produced an overall movement through the vacancy transition after a long period of operation. Using this process, a sufficient amount of data was obtained through repeated independent simulations, and the macroscopic equivalent diffusion coefficient of copper clusters could be obtained through statistical calculations [ 13 13 ]. The other basic parameters required in this study were also derived from the results of first-principles calculations and AKMC simulations in previous studies [ 13 13 , 18 18 , 19 19 ]. In addition, although the internal structures of defects and clusters are ignored, defect clusters in body-centered cubic (BCC) iron have two possible shapes in the OKMC simulation: spheres (3D) or dislocation loops (2D), where vacancy-type clusters usually exist as 3D spheres and SIA-type clusters usually exist as 2D dislocation loops. All objects move discretely in a given space according to their physical laws. For example, the mobile direction of the dislocation loops is along the direction of the Burgers vector. Moreover, defect clusters in the system interact during the evolution process to produce a defect reaction network with a corresponding reaction distance. When the distance between the two instances is less than the sum of their respective reaction radii, a reaction can occur. The reaction radius for each particle object and cluster object is defined as 0.55 nm in this model. In particular, although there is a strong interaction between vacancies and Cu clusters, a stable Cu-SIA coupling is difficult to form in the Fe matrix; thus, copper atoms are considered to only interact with single vacancies or vacancy clusters but not with SIA clusters (dislocation loops) in this study. However, SIA clusters can still react with the vacancies in Cu-vacancy composite clusters, which would lead to elimination of the same number of vacancies in the Cu-vacancy composite cluster. The interactions of the defect clusters in this simulation are shown in Fig. 1 Fig. 1 . Although the copper lattice has a face-centered cubic (FCC) structure, experiments have confirmed that when the diameter of the clusters is less than 5 nm, the copper clusters will maintain the same BCC structure as the iron matrix during their growth. When the diameter of the clusters is greater than 12 nm, the crystal structure of the copper clusters is transformed into the FCC structure [ 20 20 , 21 21 ]. This conclusion has also been verified by molecular dynamics simulations, which further revealed that the existence of vacancy clusters in the Cu-vacancy composite cluster enhanced the structural stability of the cluster and inhibited the occurrence of this transition [ 21 21 ]. Overall, although the OKMC simplification ignores the description of FCC-structured copper clusters formed in BCC-structured iron, it is still acceptable to simulate the precipitation of BCC-structured copper clusters in iron, especially when the copper clusters are less than 5 nm in diameter. Parameterization 1 In this study, the main parameters required for a cluster object include its mobility and stability, which are functions of the number of defects contained in the cluster. The abbreviations of these parameters are listed in Table 1 Table 1 , referring to vacancy-type clusters and SIA-type clusters, and the selection of their numerical values is described below. 1 1 <math id="M2"><mrow><msub><mi>N</mi><mi>δ</mi></msub></mrow></math> 1 1 Number of defects contained in the cluster 1 1 <math id="M3"><mrow><msubsup><mi>v</mi><mrow><mtext>mig</mtext></mrow><mi>δ</mi></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>δ</mi></msub><mo> </mo></mrow><mo>)</mo></mrow></mrow></math> 1 1 Attempt frequency for migration of the cluster 1 1 <math id="M4"><mrow><msub><mi>E</mi><mtext>p</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>δ</mi></msub></mrow><mo>)</mo></mrow></mrow></math> 1 1 Potential energy for defects to form in the cluster 1 1 <math id="M5"><mrow><msub><mi>E</mi><mtext>f</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>δ</mi></msub></mrow><mo>)</mo></mrow></mrow></math> 1 1 Formation energy of the cluster 1 1 <math id="M6"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mi>δ</mi></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>δ</mi></msub></mrow><mo>)</mo></mrow></mrow></math> 1 1 Binding energy of the defects to the mother cluster 1 1 <math id="M7"><mrow><msubsup><mi>E</mi><mrow><mtext>dis</mtext></mrow><mi>δ</mi></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>δ</mi></msub><mo> </mo></mrow><mo>)</mo></mrow></mrow></math> 1 1 Dissociation energy of a defect Vacancy-type clusters 2 For a single vacancy, the migration energy has been calculated to be 0.64 eV [ 22 22 ] and 0.67 eV [ 23 23 ] by density functional theory (DFT); it has been experimentally determined to be 0.55–0.57 eV [ 24 24 ], while the calculated value from the classical Mendelev potential is 0.63 eV [ 25 25 ]. In our study, the migration energy of the vacancy is set to 0.63 eV, between the above theoretical calculation and experimental values, to ensure good reliability. At the same time, we use the Debye frequency as the attempt frequency of the migration, <math id="M8"><mrow><msub><mi>v</mi><mn>1</mn></msub></mrow></math> = 6×10 12 s -1 . For vacancy clusters containing <math id="M9"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> vacancies, when <math id="M10"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 2–10, the AKMC simulation results [ 26 26 ] are used as the migration energies and corresponding transition attempt frequencies. For vacancy clusters with <math id="M11"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> > 10, the calculated migration energy fluctuates between 1.15 and 1.50 eV, and thus an intermediate value of 1.30 eV is selected in this study. Similarly, it is assumed that the diffusion mechanism is surface diffusion [ 27 27 , 28 28 ]; that is, the movement of these vacancy clusters is the result of a series of surface vacancy transition events. Thus, the corresponding migration attempt frequency can be calculated as follows: <math display="block" id="M12"><mrow><msubsup><mi>v</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><msub><mi>v</mi><mn>1</mn></msub></mrow><mrow><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow></mrow><mrow><mn>4</mn><mo>/</mo><mn>3</mn></mrow></msup></mrow></mfrac><mo>.</mo></mrow></math> Before parameterizing the vacancy cluster stability, we must define the formation energy of vacancy clusters. In MMonCa, the cluster formation energy is the sum of the formation energies of the constituent particles minus the potential energy E p ( N v )oftheclusters, E n ( N v ) <math display="block" id="M13"><mrow><msub><mi>E</mi><mtext>f</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow><mo>=</mo><msub><mi>N</mi><mtext>v</mtext></msub><mo>×</mo><msubsup><mi>E</mi><mtext>f</mtext><mtext>v</mtext></msubsup><mo>−</mo><msub><mi>E</mi><mtext>p</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow><mo>,</mo></mrow></math> where <math id="M14"><mrow><msubsup><mi>E</mi><mtext>f</mtext><mtext>v</mtext></msubsup></mrow></math> is the vacancy formation energy in the iron matrix, which is set as 2.02 eV in this study according to the value calculated in the AKMC simulation. Then, <math id="M15"><mrow><msub><mi>E</mi><mtext>p</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow><mo> </mo></mrow></math> is the accumulation of binding energies for each unit vacancy bound to this cluster: <math display="block" id="M16"><mrow><msub><mi>E</mi><mtext>p</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow><mo>=</mo><msub><mi>E</mi><mtext>p</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><mo>−</mo><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><mo>,</mo></mrow></math> where <math id="M17"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow></mrow></math> is the binding energy of a single vacancy to the cluster with <math id="M18"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> -1. In this study, the corresponding AKMC results are used for the binding energies of vacancy clusters with <math id="M19"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 2–4. For vacancy clusters with <math id="M20"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> > 4, the fitting function in the literature [ 19 19 ] is used: <math display="block" id="M21"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow><mo>=</mo><mn>1.71</mn><mo>+</mo><mn>3.39716</mn><mrow><mo>[</mo><mrow><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow></mrow><mrow><mfrac><mn>2</mn><mn>3</mn></mfrac></mrow></msup><mo>−</mo><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo>+</mo><mn>1</mn></mrow><mo>)</mo></mrow></mrow><mrow><mfrac><mn>2</mn><mn>3</mn></mfrac></mrow></msup></mrow><mo>]</mo></mrow><mo>.</mo></mrow></math> The stability of vacancy clusters can then be defined to determine the rate of decomposition of vacancy clusters. Because the event in which a vacancy cluster emits another vacancy cluster has a very low probability [ 22 22 ], the main means of decomposing vacancy clusters is the continuous emission of single vacancies into the matrix. The dissociation energy, <math id="M22"><mrow><msubsup><mi>E</mi><mrow><mtext>dis</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow></mrow></math> , of vacancy clusters with <math id="M23"><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow></math> is the sum of the binding energy of a single vacancy to the mother cluster, <math id="M24"><mrow><mo> </mo><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow><mo>,</mo></mrow></math> and the single vacancy migration energy, <math id="M25"><mrow><msubsup><mi>E</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> : <math display="block" id="M26"><mrow><msubsup><mi>E</mi><mrow><mtext>dis</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow><mo>=</mo><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow><mo>+</mo><msubsup><mi>E</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup><mo>.</mo></mrow></math> It should be noted that the vacancy migration energy during dissociation from the mother cluster is equal to its migration energy when transitioning in the matrix, 0.63 eV, which is independent of the <math id="M27"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> of the vacancy cluster. In addition, the dissociation attempt frequency of the cluster should be proportional to the size, <math id="M28"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> , which is also the condition required for the cluster absorption and emission point defects to reach a balance. Then, the attempt frequency is defined in this study as follows: <math display="block" id="M29"><mrow><msubsup><mi>v</mi><mrow><mtext>dis</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow><mo>=</mo><msub><mi>v</mi><mn>1</mn></msub><mo>×</mo><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow></mrow><mrow><mfrac><mn>1</mn><mn>3</mn></mfrac></mrow></msup><mo>.</mo></mrow></math> As mentioned previously, the movable copper clusters will contain at least one vacancy, whereas pure copper clusters without vacancies are regarded as immovable objects. Then, the Cu-vacancy clusters with <math id="M30"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 1 are defined first. For Cu-vacancy clusters with <math id="M31"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> = 2–5 and <math id="M32"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 1, the frequency of transition attempts is based on the AKMC results [ 18 18 ]. For clusters with <math id="M33"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> > 5 and <math id="M34"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 1, the mobility decreases with increasing <math id="M35"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> following the surface diffusion mechanism [ 29 29 ]: <math display="block" id="M36"><mrow><msubsup><mi>v</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>+</mo><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow><mo>=</mo><mo> </mo><mfrac><mrow><msub><mi>v</mi><mn>1</mn></msub></mrow><mrow><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>+</mo><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow></mrow><mrow><mn>4</mn><mo>/</mo><mn>3</mn></mrow></msup></mrow></mfrac><mo>.</mo></mrow></math> At the same time, the migration energy of the Cu-vacancy pair is close to that of a single vacancy [ 19 19 ] according to DFT calculations, and thus it is also assigned a value of 0.63 eV. For clusters with <math id="M37"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> = 2–4 and <math id="M38"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 1, we use the AKMC simulation results mentioned previously [ 18 18 ]. For clusters with <math id="M39"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> = 15–6000 and <math id="M40"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 1, the migration energy is assigned a value of 1.2 eV, as calculated in reference [ 13 13 ]. For clusters with <math id="M41"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> = 5–14 and <math id="M42"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 1, the corresponding value is taken as 1 eV. Figure 2 Figure 2 shows the decrease in <math id="M43"><mrow><msub><mi>D</mi><mi>N</mi></msub></mrow></math> with increasing <math id="M44"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> , which is similar to the AKMC results obtained in a previous study. For the stability of Cu-vacancy clusters with <math id="M45"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> = 1, only two types of dissociation events are considered in the simulation: the emission of single vacancies or Cu-vacancy pairs from their mother clusters. The corresponding dissociation energy is equal to the migration energy of the emitted defect object plus its binding energy to the mother cluster. For example, the dissociation energy, <math id="M46"><mrow><msubsup><mi>E</mi><mrow><mtext>dis</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> , for a vacancy emission can be calculated as follows: <math display="block" id="M47"><mrow><msubsup><mi>E</mi><mrow><mtext>dis</mtext></mrow><mtext>v</mtext></msubsup><mo>=</mo><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>,</mo><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow><mo>+</mo><msubsup><mi>E</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup><mo>,</mo></mrow></math> where <math id="M48"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>,</mo><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow></mrow></math> is the binding energy of one vacancy to a cluster with <math id="M49"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> and <math id="M50"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> ; <math id="M51"><mrow><msubsup><mi>E</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> is the migration energy of the vacancy, and the attempt frequency for emission is defined as 6×10 12 . The migration energy of the emitted object is also considered to be independent of the size of its mother cluster, and <math id="M52"><mrow><msubsup><mi>E</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> = <math id="M53"><mrow><msubsup><mi>E</mi><mrow><mtext>mig</mtext></mrow><mrow><mtext>Cu</mtext><mo>−</mo><mtext>vacancy</mtext></mrow></msubsup></mrow></math> = 0.63 eV for all <math id="M54"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> . We then calculate <math id="M55"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>,</mo><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow><mo>=</mo><msubsup><mi>E</mi><mrow><mtext>Dis</mtext></mrow><mtext>v</mtext></msubsup><mo>−</mo><msubsup><mi>E</mi><mrow><mtext>mig</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> from the dissociation energies, <math id="M56"><mrow><msubsup><mi>E</mi><mrow><mtext>Dis</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> . For clusters with <math id="M57"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> = 2–5, the dissociation energies are obtained from reference [ 18 18 ]. For <math id="M58"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> ≥ 6, the values are extrapolated by fitting a function to the AKMC results using the following equation: <math display="block" id="M59"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>,</mo><mn>1</mn><mo> </mo></mrow><mo>)</mo></mrow><mo>=</mo><mn>0.9</mn><mo>+</mo><mn>1.6</mn><mo>×</mo><mo stretchy="false">(</mo><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow><mo>)</mo></mrow></mrow><mrow><mfrac><mn>2</mn><mn>3</mn></mfrac></mrow></msup><mo>−</mo><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>+</mo><mn>1</mn></mrow><mo>)</mo></mrow></mrow><mrow><mfrac><mn>2</mn><mn>3</mn></mfrac></mrow></msup><mo stretchy="false">)</mo><mo>.</mo></mrow></math> This function is assumed to have the same functional form as Eq. (7), and the corresponding coefficients are then trained using the AKMC values. For Cu-vacancy clusters with <math id="M60"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> ≥ 2, the cluster mobility also decreases as the number of vacancies increases, following the same surface diffusion mechanism described by Eq. (10). To determine the stability of these clusters, we use the results of molecular dynamics simulations that provide the binding energy of copper–vacancy composite clusters through computational simulations of metropolis Monte Carlo (MMC) and molecular dynamics binding [ 19 19 , 29 29 ]; the stability is a binary function of <math id="M61"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> and <math id="M62"><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub></mrow></math> . In this study, the partial differential form is used as the input parameter, <math id="M63"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>,</mo><msub><mi>N</mi><mtext>v</mtext></msub><mo> </mo></mrow><mo>)</mo></mrow></mrow></math> . SIA-type clusters 2 The method for defining the mobility of SIA-type clusters is similar to that for vacancy-type clusters, and the diffusion pre-factor is determined using the same equations. For SIA clusters of size <math id="M64"><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow></math> < 7, the migration energy and migration attempt frequency used in this study are obtained from molecular dynamics studies [ 19 19 , 29 29 ] and experimental studies [ 25 25 ] based on Mendeleev’s potential. It should also be noted that when the SIA cluster grows beyond a certain threshold, it will transform into 1/2 < 111 > loops, and <math id="M65"><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow></math> = 7 is selected as the threshold for this transition. At the same time, we do not distinguish between < 111 > dislocation loops and < 100 > dislocation loops when simulating their migration, and thus we use an effective migration energy and set a dislocation loop distinction threshold, <math id="M66"><mrow><msub><mi>N</mi><mrow><mtext>th</mtext></mrow></msub></mrow></math> . For dislocation loops with <math id="M67"><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow></math> = 7 - <math id="M68"><mrow><msub><mi>N</mi><mrow><mtext>th</mtext></mrow></msub></mrow></math> , it is assumed that only part are < 100 > dislocation loops. Therefore, we use the effective migration energy for the dislocation loops whose size is in this interval. It is assumed that the < 111 > dislocation loops are transformed into < 100 > dislocation loops with a higher migration energy when the number of atoms in the dislocation loop is greater than <math id="M69"><mrow><msub><mi>N</mi><mrow><mtext>th</mtext></mrow></msub></mrow></math> , and the corresponding value is 0.90 eV [ 29 29 ]. Jansson et al. used <math id="M70"><mrow><msub><mi>N</mi><mrow><mtext>th</mtext></mrow></msub></mrow></math> = 90 as the transition threshold, and the corresponding migration energy was 0.2 eV [ 29 29 ]. We selected the same values for these two effective parameters in this study. For interstitial clusters with <math id="M71"><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow></math> ≥ 7, the attempt frequency is calculated as follows: <math display="block" id="M72"><mrow><msubsup><mi>v</mi><mrow><mtext>mig</mtext></mrow><mi>i</mi></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><msub><mi>v</mi><mn>1</mn></msub></mrow><mrow><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow></mrow><mi>e</mi></msup></mrow></mfrac><mo>,</mo></mrow></math> where the exponential part, e , of <math id="M73"><mrow><msub><mi>N</mi><mtext>i</mtext></msub></mrow></math> varies theoretically from 0.50 to 1.00, corresponding to the SIA cluster diffusion model [ 30 30 ] and kink movement model of dislocation loop edges [ 31 31 ], respectively. In this study, the value is taken as 0.8, as in reference [ 29 29 ]. Similar to the definition of stability for vacancy-type clusters mentioned previously, the dissociation energy of an SIA cluster is the sum of the binding energy of an SIA to its mother cluster plus and the migration energy. For the binding energy of the SIA to its mother cluster, we use the fitting function given by Domain et al. [ 28 28 ]: <math display="block" id="M74"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mi>i</mi></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><mo>=</mo><mn>4</mn><mo>+</mo><mn>5.107249</mn><mrow><mo>[</mo><mrow><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow></mrow><mrow><mfrac><mn>2</mn><mn>3</mn></mfrac></mrow></msup><mo>−</mo><msup><mrow><mrow><mo>(</mo><mrow><msub><mi>N</mi><mi>i</mi></msub><mo>+</mo><mn>1</mn></mrow><mo>)</mo></mrow></mrow><mrow><mfrac><mn>2</mn><mn>3</mn></mfrac></mrow></msup></mrow><mo>]</mo></mrow><mo>.</mo></mrow></math> The final assumption of this model is the use of a simple stochastic description to simulate SIA defect traps. The simulation does not explicitly introduce trap objects but defines the probability of an SIA-type defect transforming into a pinned defect. The transition probability is calculated using the relevant sink strength formula, as in reference [ 32 32 ]. The main trap type considered in the iron matrix is a C-vacancy trap; thus, its concentration is set to 50 appm. Results and discussion 1 Evolution of Cu clusters during thermal aging 2 Based on this model, we simulated the thermal aging process of copper clusters under the same conditions as the experiments performed by Chen et al. [ 33 33 ] and Vincent et al. [ 34 34 ] in a simulation box with dimensions of 72 × 72 × 72 nm 3 . The experimental alloys with Cu contents of 0.5 at.% and 0.6 at.% were annealed at 450 ℃ and 500 ℃, respectively, and the precipitation of Cu clusters was characterized by atom probe tomography (APT). Because the concentration of vacancies in the simulation box is much higher than that in the real physical environment, the simulation time must be rescaled before comparison with the experimental data on the same time scale. Generally, the ratio of the Monte Carlo simulation time, <math id="M75"><mrow><msub><mi>t</mi><mrow><mtext>MC</mtext></mrow></msub></mrow></math> , to real physical time, <math id="M76"><mrow><msub><mi>t</mi><mrow><mtext>real</mtext></mrow></msub></mrow></math> , is defined as follows [ 13 13 ]: <math display="block" id="M77"><mrow><msub><mi>t</mi><mrow><mtext>real</mtext></mrow></msub><mo>=</mo><msub><mi>t</mi><mrow><mtext>MC</mtext></mrow></msub><mfrac><mrow><msubsup><mi>C</mi><mrow><mtext>MC</mtext></mrow><mtext>v</mtext></msubsup></mrow><mrow><msubsup><mi>C</mi><mrow><mtext>real</mtext></mrow><mtext>v</mtext></msubsup></mrow></mfrac><mo>,</mo></mrow></math> where <math id="M78"><mrow><msubsup><mi>C</mi><mrow><mtext>real</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> is the experimental vacancy concentration, and <math id="M79"><mrow><msubsup><mi>C</mi><mrow><mtext>MC</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> is the simulation vacancy concentration. Considering that the interaction between solute atoms and vacancies will change the local vacancy formation energy during the precipitation of solute clusters, the effective concentration of vacancies in the simulation needs to be corrected [ 35 35 ]. In this study, we approximate <math id="M80"><mrow><msubsup><mi>C</mi><mrow><mtext>MC</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> using the mobile concentration [ 14 14 ], <math id="M81"><mrow><msubsup><mi>C</mi><mrow><mtext>mobile</mtext></mrow><mtext>v</mtext></msubsup></mrow></math> , of vacancies as follows: <math display="block" id="M82"><mrow><msubsup><mi>C</mi><mrow><mtext>mobile</mtext></mrow><mtext>v</mtext></msubsup><mo>=</mo><mfrac><mrow><msub><mi>N</mi><mrow><mtext>jump</mtext></mrow></msub></mrow><mrow><mtext>Δ</mtext><mi>t</mi><mtext>Δ</mtext><mi>V</mi><mi>ν</mi><mrow><mo>(</mo><mrow><mi>p</mi><mi>t</mi></mrow><mo>)</mo></mrow></mrow></mfrac><mo>,</mo></mrow></math> where <math id="M83"><mrow><mtext>Δ</mtext><mi>V</mi><mi>ν</mi><mrow><mo>(</mo><mrow><mi>p</mi><mi>t</mi></mrow><mo>)</mo></mrow></mrow></math> is the sum of the capture volumes of movable particle objects at each moment (i.e., the vacancies in this case), and <math id="M84"><mrow><msub><mi>N</mi><mrow><mtext>jump</mtext></mrow></msub></mrow></math> is the number of vacancy jumps in the Fe matrix in time interval <math id="M85"><mrow><mi>�</mi><mi>t</mi></mrow></math> ; this value can be output directly from the simulation in MMonCa. It is worth mentioning that its ratio to the absolute vacancy concentration, <math id="M86"><mrow><msub><mi>f</mi><mi>V</mi></msub></mrow></math> , is positively related to the sink strength of Cu clusters to vacancies. For the real experimental vacancy concentration, the equilibrium estimation in a pure Fe environment can be obtained as follows: <math display="block" id="M87"><mrow><msubsup><mi>C</mi><mrow><mrow><mo>(</mo><mrow><mtext>real</mtext></mrow><mo>)</mo></mrow></mrow><mtext>v</mtext></msubsup><mo>=</mo><mi>e</mi><mi>x</mi><mi>p</mi><mrow><mo>(</mo><mrow><mfrac><mrow><mo>−</mo><msubsup><mi>h</mi><mtext>f</mtext><mtext>v</mtext></msubsup></mrow><mrow><msub><mi>k</mi><mi>B</mi></msub><mi>T</mi></mrow></mfrac></mrow><mo>)</mo></mrow><mo>,</mo></mrow></math> where <math id="M88"><mrow><msubsup><mi>h</mi><mtext>f</mtext><mtext>v</mtext></msubsup></mrow></math> is the enthalpy of formation of the vacancy in a pure Fe environment, which is taken to be 1.7 eV. Figure 3 Figure 3 shows a comparison of the thermal aging of the Fe-Cu alloy under different conditions, in which the OKMC simulation results are plotted as lines and correspond to the experimental results plotted as triangles. It is clear that our model is consistent with the experimental results and successfully shows the evolution of the density and size of the clusters. At 450 ℃, as shown in Fig. 3(a) and (b) Fig. 3(a) and (b) , the average number density and radius of clusters versus time in the Fe-0.5 at.% Cu alloy increase gradually and then reach a steady state after 200 h; this behavior is closely reproduced at the temperature investigated. The APT results show that the micro-hardness reaches maximum after 200 h of annealing, which may indicate the largest value of the cluster number density multiplied by the average radius [ 33 33 ]. At 500 ℃, the number density increases rapidly at the initial stage of thermal aging and then gradually decreases, indicating obvious coalescence between the clusters, as shown in Fig. 3(c) Fig. 3(c) . The average radius of copper clusters was also measured during thermal aging at 500 ℃ for up to 20 h, from which three different cluster growth processes were observed. In the first nucleation stage, the continuous precipitation of small-sized clusters is observed. Then, in the growth stage of the clusters, the average growth rate of the clusters is close to the classical time law t 1/2 . Finally, the growth of the average cluster size is further reduced and approaches the classical time law t 1/3 , corresponding to the coarsening stage [ 36 36 ]. Simultaneously, the number density of clusters decreases in the later stage of cluster growth. For this result, the prediction of the OKMC model is in good agreement with the experimental data. Through comparison with the experimental results, we verified the validity of the model in this study. Moreover, to explain the evolution mechanism of copper clusters under different conditions, the parameters of the simulation process are further analyzed. As mentioned previously, in the simulated annealing process, we defined a recalculation for the simulation time; then, there is a fraction, <math id="M89"><mrow><msub><mi>f</mi><mi>V</mi></msub></mrow></math> , of the vacancy concentration between the mobile concentration and absolute concentration, which should be a variable of temperature and time. The changes in <math id="M90"><mrow><msub><mi>f</mi><mi>V</mi></msub></mrow></math> over the Monte Carlo time are shown in Fig. 4 Fig. 4 . It is clear that the fraction decreases with increasing simulation time, which indicates that the growing Cu clusters become a stronger sink for vacancies, and vacancies will then spend more time inside the clusters. It is also noticeable that there is a minimum value over the curing period at 500 ℃, which may be related to the decrease in the average cluster number density according to the evolution of the copper clusters shown above. Furthermore, the cluster density during thermal aging decreases because of the merging of large clusters and cluster dissociation through the emission of Cu-vacancy pairs, which can be observed by monitoring the simulation process. The one-by-one emission of Cu-vacancy pairs from mother Cu-vacancy clusters is important for the formation of large Cu clusters; when the Cu clusters maintain proper size and mobility, the merging of large clusters is possible. Considering the temperature effect, the Cu cluster grows much faster at 500 ℃ than at 450 ℃. Because the cluster mobility is highly related to temperature, clusters have a better chance of merging together at higher temperatures. At the same time, it is assumed that the emission rate of Cu-vacancy couples reaches a balance with the absorption rate, after which the entire microstructure evolution process appears very slow; thus, there some of Cu atoms remain dispersed in the matrix. Evolution of Cu clusters under irradiation 2 In this section, OKMC simulations were conducted for neutron irradiation under the representative conditions reported by Liu et al. [ 37 37 ] and ion irradiation experiments by Shu et al. [ 38 38 ]; the simulations were performed in a simulation box with dimensions of 72 × 72 × 72 nm 3 . Chen33] obtained the size distribution of Cu clusters using histogram statistics, and the results showed a unimodal distribution with a peak at approximately 1.6 nm. In the neutron irradiation experiment performed by Liu[ 37 37 ], the peak size was approximately 1.45 nm. We carried out a neutron irradiation simulation in Fe-0.5 at.% Cu and compared the results with those from the annealing simulation, as shown in Fig. 5 Fig. 5 . It should be noted that the kernel density estimation (KDE) is used to obtain the cluster size distribution instead of a histogram because of the small sample size in the simulation box; correspondingly, its ordinate values are the probability density. This approach provides an estimation with excellent reliability, which can quickly stabilize with increasing samples by repeating the simulation. The results show that the size distribution obtained by the simulation of the thermal aging process has a peak value of 1.6 nm, whereas for the neutron irradiation simulation, the peak value is approximately 1.4 nm. This result is consistent with the experimental results, indicating that the cluster segregation produced by annealing at 450 ℃ for 200 h is comparable to that under neutron irradiation at 290 ℃ with a flux of 5.4×10 18 n/cm 2 . Furthermore, we compared the neutron irradiation and ion irradiation in Fe-0.8Cu using a damage rate of 3×10 -5 dpa/s to simulate ion irradiation and a damage rate of 3×10 -7 dpa/s to simulate neutron irradiation; the simulation results are shown in Fig. 6 Fig. 6 . The simulation results for neutron irradiation at lower dpa levels are shown as black lines on the left, corresponding to the experiments plotted as blue dots. The simulation results for ion irradiation at higher dpa levels are shown as red lines and correspond to the ion irradiation experiments plotted as purple triangles. Under irradiation, the number density of Cu-vacancy clusters formed in the alloy initially increases and then decreases, accompanied by an increase in the mean radius of the clusters, which is similar to the annealing simulation and corresponds to nucleation, growth, and coarsening stages. Meanwhile, the mean radius and number density ratios of precipitates in Fe-Cu alloys under ion irradiation are much smaller and larger, respectively, compared with those under neutron irradiation. For the irradiation simulations, it is found that under the same irradiation damage dose, the mean radius and number density of clusters under neutron and ion irradiation are significantly different; thus, we further analyze the vacancy clusters in the solute clusters. Figure 7 Figure 7 shows the 2D distribution of copper and vacancies in solute clusters; it should be noted that this distribution is also obtained by the 2D KDE. Under neutron irradiation, from 0.0001dpa to 0.001dpa, the increase in the number of vacancies is small, as shown in Fig. 7(a) and (b) Fig. 7(a) and (b) . Thus, we can infer that the cluster maintains a high mobility, and the cluster merger mechanism occupies the later stage of evolution. The average radius of the clusters still increases rapidly, as shown in Fig. 6 Fig. 6 . Under ion irradiation, the content of vacancies in the cluster is high even at 0.03 dpa ( Fig. 7(c) Fig. 7(c) ), which also leads to weakening of the cluster mobility. It can be inferred that the absorption of Cu-vacancy pairs should be dominant in the later stage of copper cluster precipitation such that the cluster size reaches equilibrium, and most of the clusters contain less than 1000 Cu atoms, as shown in Fig. 7(d) Fig. 7(d) . Moreover, we compared the distributions of the binding energy, <math id="M91"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>,</mo><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow></mrow></math> , as shown in Fig. 8 Fig. 8 . It is clear that with an increase in <math id="M92"><mrow><msub><mi>N</mi><mtext>v</mtext></msub></mrow></math> in the Cu-vacancy cluster, the binding energy, <math id="M93"><mrow><msubsup><mi>E</mi><mrow><mtext>bind</mtext></mrow><mtext>v</mtext></msubsup><mrow><mo>(</mo><mrow><msub><mi>N</mi><mrow><mtext>Cu</mtext></mrow></msub><mo>,</mo><msub><mi>N</mi><mtext>v</mtext></msub></mrow><mo>)</mo></mrow></mrow></math> , gradually increases. Under neutron irradiation, because of the slow damage rate, the vacancy accumulation rate is close to the vacancy dissociation rate in the cluster; therefore, the distribution in Fig. 7(b) Fig. 7(b) is very similar to this binding energy distribution map. However, in the lower area of the map, a mutation range is observed. Therefore, when the irradiation rate is fast, i.e., corresponding to ion irradiation, the rapid accumulation of vacancies will lead to a rapid increase in the stability of the cluster. A similar tendency can be observed in the distribution shown in Fig. 7(c) Fig. 7(c) . Conclusion 1 In this study, we establish an OKMC model under the MMonca framework to simulate the solute precipitation process in Fe-Cu systems. The proposed model can be used to simulate the formation process of Cu-vacancy composite clusters under neutron/ion irradiation conditions, in addition to enabling prediction of the annealing process. Based on a modified description of the reaction network in the iron matrix, this model enables the simulation of sufficiently large volumes with long-term irradiation. The corresponding energies used for the parameterization are calculated from published DFT and AKMC results and provide a complementary simulation result to APT experiments. Our model includes all of the important mechanisms driving Cu precipitation in Fe. It is demonstrated that the merging of mobile Cu clusters containing thousands of Cu atoms is the main mechanism responsible for the growth and coarsening of Cu clusters during thermal aging. A more complex reaction network will occur in the Fe-Cu system under irradiation than under thermal aging, but there is still a mapping relationship in the cluster size distribution between both processes. In addition, the evolution behavior of the copper clusters differs under neutron and ion irradiation. The thermal aging results confirm that the merging of clusters is an important mechanism for the growth of copper clusters. However, because of the strong reaction between vacancy clusters and copper clusters, both clusters easily become nucleation sites for each other under ion irradiation, thus inhibiting the migration and coarsening of the copper clusters, ultimately leading to the smaller average size of copper clusters observed by APT experiments. References R.E. Smallman , K.H. Westmacott , Structure of quenched and irradiated metals . J. Appl. Phys. 30 ( 5 ), 603 - 616 ( 1959 ). doi: 10.1063/1.1735203 http://doi.org/10.1063/1.1735203 Structure of quenched and irradiated metals P.D. Styman , J.M. Hyde , K. Wilford et al. , Characterisation of interfacial segregation to Cu-enriched precipitates in two thermally aged reactor pressure vessel steel welds . Ultramicroscopy 159 , 292 - 298 ( 2015 ). doi: 10.1016/j.ultramic.2015.05.013 http://doi.org/10.1016/j.ultramic.2015.05.013 Characterisation of interfacial segregation to Cu-enriched precipitates in two thermally aged reactor pressure vessel steel welds T. Onitsuka , M. Takenaka , E. Kuramoto et al. , Deformation-enhanced Cu precipitation in Fe-Cu alloy studied by positron annihilation spectroscopy . Phys. Rev. B 65 ( 1 ), 012204 ( 2001 ). doi: 10.1103/physrevb.65.012204 http://doi.org/10.1103/physrevb.65.012204 Deformation-enhanced Cu precipitation in Fe-Cu alloy studied by positron annihilation spectroscopy S. Jin , X. Cao , G. Cheng et al. , Thermally promoted evolution of open-volume defects and Cu precipitates in the deformed FeCu alloys . J. Nucl. Mater. 2018, 501 , 293 - 301 ( 2018 ). doi: 10.1016/j.jnucmat.2018.01.061 http://doi.org/10.1016/j.jnucmat.2018.01.061 Thermally promoted evolution of open-volume defects and Cu precipitates in the deformed FeCu alloys S. Jin , P. Zhang , E. Lu , L. Guo et al. , Correlation between Cu precipitates and irradiation defects in Fe–Cu model alloys investigated by positron annihilation spectroscopy . Acta Mater. 103 , 658 - 64 ( 2016 ). doi: 10.1016/j.actamat.2015.10.051 http://doi.org/10.1016/j.actamat.2015.10.051 Correlation between Cu precipitates and irradiation defects in Fe–Cu model alloys investigated by positron annihilation spectroscopy V. Slugeň , A. Kryukov , Microstructural study of WWER reactor pressure vessel steels . Nucl. Eng. Des. 263 , 308 - 312 ( 2013 ). doi: 10.1016/j.nucengdes.2013.06.013 http://doi.org/10.1016/j.nucengdes.2013.06.013 Microstructural study of WWER reactor pressure vessel steels T. Muroga , N. Sekimura , New insights into the temperature effects on neutron irradiation of structural materials1 . Fusion Eng. Des. 41 ( 1-4 ), 39 - 46 ( 1998 ). doi: 10.1016/s0920-3796(97)00141-5 http://doi.org/10.1016/s0920-3796(97)00141-5 New insights into the temperature effects on neutron irradiation of structural materials1 X. Wang , J.L. Li , Z. Wu et al. , CMGC - A CAD to Monte Carlo Geometry Conversion Code . Nucl. Sci. Tech. 31 ( 8 ), 82 ( 2020 ). doi: 10.1007/s41365-020-00793-8 http://doi.org/10.1007/s41365-020-00793-8 CMGC - A CAD to Monte Carlo Geometry Conversion Code T.Y. Huang , Z.G. Li , K. Wang et al. , Hybrid windowed networks for on-the-fly Doppler broadening in RMC code . Nucl. Sci. Tech. 32 ( 6 ), 62 ( 2021 ). doi: 10.1007/s41365-021-00901-2 http://doi.org/10.1007/s41365-021-00901-2 Hybrid windowed networks for on-the-fly Doppler broadening in RMC code Z.G. Li , K. Wang , Y.C. Guo et al. , Forced propagation method for Monte Carlo fission source convergence acceleration in the RMC . Nucl. Sci. Tech. 32 ( 3 ), 27 ( 2021 ). doi: 10.1007/s41365-021-00868-0 http://doi.org/10.1007/s41365-021-00868-0 Forced propagation method for Monte Carlo fission source convergence acceleration in the RMC C. S. Becquart , C. Domain . Malerba Introducing chemistry in atomistic kinetic Monte Carlo simulations of Fe alloys under irradiation . Phys. Status Solidi B 247 ( 1 ), 9–22 ( 2010 ). doi: 10.1002/pssb.200945251 http://doi.org/10.1002/pssb.200945251 Malerba Introducing chemistry in atomistic kinetic Monte Carlo simulations of Fe alloys under irradiation F. Soisson , C.C. Fu , Cu-precipitation kinetics in Fe from atomistic simulations: Vacancy-trapping effects and Cu-cluster mobility . Phys. Rev. B 76 , 214102 ( 2007 ). doi: 10.1103/physrevb.76.214102 http://doi.org/10.1103/physrevb.76.214102 Cu-precipitation kinetics in Fe from atomistic simulations: Vacancy-trapping effects and Cu-cluster mobility N. Castin , M. I. Pascuet , L. Malerba , Modeling the first stages of Cu precipitation in α-Fe using a hybrid atomistic kinetic Monte Carlo approach . J. Chem. Phys. 135 , 064502 ( 2011 ). doi: 10.1063/1.3622045 http://doi.org/10.1063/1.3622045 Modeling the first stages of Cu precipitation in α-Fe using a hybrid atomistic kinetic Monte Carlo approach I. Martin-Bragado , A. Rivera , G. Valles et al. , MMonCa: An Object Kinetic Monte Carlo simulator for damage irradiation evolution and defect diffusion . Comput. Phys. Comm. 184 ( 12 ), 2703 - 2710 ( 2013 ). doi: 10.1016/j.cpc.2013.07.011 http://doi.org/10.1016/j.cpc.2013.07.011 MMonCa: An Object Kinetic Monte Carlo simulator for damage irradiation evolution and defect diffusion P. J. Othen , M. L. Jenkins & G. D. W. Smith , High-resolution electron microscopy studies of the structure of Cu precipitates in α-Fe . Phil. Mag. A 70 ( 1 ), 1 - 24 ( 1994 ). doi: 10.1080/01418619408242533 http://doi.org/10.1080/01418619408242533 High-resolution electron microscopy studies of the structure of Cu precipitates in α-Fe W.M. Young & E.W. Elcock , Monte Carlo studies of vacancy migration in binary ordered alloys: I . Proc. Phys. Soc. 89 ( 3 ), 735 - 746 ( 1966 ). doi: 10.1088/0370-1328/89/3/329 http://doi.org/10.1088/0370-1328/89/3/329 Monte Carlo studies of vacancy migration in binary ordered alloys: I A. Kristen , Fichthorn and W. H. Weinberg . Theoretical foundations of dynamical Monte Carlo simulations . J. Chem. Phys. 95 , 1090 ( 1991 ). doi: 10.1063/1.461138 http://doi.org/10.1063/1.461138 Theoretical foundations of dynamical Monte Carlo simulations M.I. Pascuet , N. Castin , C.S. Becquart et al. , Stability and mobility of Cu–vacancy clusters in Fe–Cu alloys: A computational study based on the use of artificial neural networks for energy barrier calculations . J. Nucl. Mater. 412 ( 1 ), 106–115 ( 2011 ). doi: 10.1016/j.jnucmat.2011.02.038 http://doi.org/10.1016/j.jnucmat.2011.02.038 Stability and mobility of Cu–vacancy clusters in Fe–Cu alloys: A computational study based on the use of artificial neural networks for energy barrier calculations A.T. Al-Motasem , M. Posselt , F. Bergner , Nanoclusters in bcc-Fe containing vacancies, copper and nickel: Structure and energetics . J. Nucl. Mater. 418 ( 1-3 ), 215–222 ( 2011 ). doi: 10.1016/j.jnucmat.2011.07.002 http://doi.org/10.1016/j.jnucmat.2011.07.002 Nanoclusters in bcc-Fe containing vacancies, copper and nickel: Structure and energetics S. Pizzini , K. J. Roberts , W. J. Phythian et al. , A fluorescence EXAFS study of the structure of copper-rich precipitates in Fe–Cu and Fe–Cu–Ni alloys . Phil. Mag. Lett. 61 ( 4 ), 223 - 229 ( 1990 ). doi: 10.1080/09500839008202362 http://doi.org/10.1080/09500839008202362 A fluorescence EXAFS study of the structure of copper-rich precipitates in Fe–Cu and Fe–Cu–Ni alloys J. J. Blackstock & G. J. Ackland , Phase transitions of copper precipitates in Fe–Cu alloys . Phil. Mag. A 81 ( 9 ), 2127 - 2148 ( 2001 ). doi: 10.1080/01418610108217139 http://doi.org/10.1080/01418610108217139 Phase transitions of copper precipitates in Fe–Cu alloys C.S. Becquart , C. Domain , Ab initio contribution to the study of complexes formed during dilute FeCu alloys radiation . Nucl. Instr. and Meth. B 202 , 44 - 50 ( 2003 ). doi: 10.1016/s0168-583x(02)01828-1 http://doi.org/10.1016/s0168-583x(02)01828-1 Ab initio contribution to the study of complexes formed during dilute FeCu alloys radiation C.C. Fu , J. Torre , F Willaime , et al. , Multiscale modelling of defect kinetics in irradiated iron . Nature Mater. 4 , 68 - 74 ( 2005 ). doi: 10.1038/nmat1286 http://doi.org/10.1038/nmat1286 Multiscale modelling of defect kinetics in irradiated iron L. Malerba , M.C. Marinica , N. Anento et al. , Comparison of empirical interatomic potentials for iron applied to radiation damage studies . J. Nucl. Mater. 406 ( 1 ), 19 - 38 ( 2010 ). doi: 10.1016/j.jnucmat.2010.05.017 http://doi.org/10.1016/j.jnucmat.2010.05.017 Comparison of empirical interatomic potentials for iron applied to radiation damage studies M. I. Mendelev , S. Han , D. J. Srolovitz , Development of new interatomic potentials appropriate for crystalline and liquid iron . Phil. Mag. Vol. 83 ( 35 ), 3977 – 3994 ( 2003 ). doi: 10.1080/14786430310001613264 http://doi.org/10.1080/14786430310001613264 Development of new interatomic potentials appropriate for crystalline and liquid iron M.I. Pascuet , N. Castin , C.S. Becquart et al. , Stability and mobility of Cu-vacancy clusters in Fe-Cu alloys: A computational study based on the use of artificial neural networks for energy barrier calculations . J. Nucl. Mater. 42 ( 1 ), 106 - 115 ( 2011 ). doi: 10.1016/j.jnucmat.2011.02.038 http://doi.org/10.1016/j.jnucmat.2011.02.038 Stability and mobility of Cu-vacancy clusters in Fe-Cu alloys: A computational study based on the use of artificial neural networks for energy barrier calculations F.A. Nichols , Kinetics of diffusional motion of pores in solids: A review . J. Nucl. Mater. 30 ( 1-2 ), 143 - 165 ( 1969 ). doi: 10.1016/0022-3115(69)90176-7 http://doi.org/10.1016/0022-3115(69)90176-7 Kinetics of diffusional motion of pores in solids: A review C. Domain , C.S. Becquart , L. Malerba , Simulation of radiation damage in Fe alloys: an object kinetic Monte Carlo approach . J. Nucl. Mater. 335 ( 1 ), 121 - 145 ( 2004 ). doi: 10.1016/j.jnucmat.2004.07.037 http://doi.org/10.1016/j.jnucmat.2004.07.037 Simulation of radiation damage in Fe alloys: an object kinetic Monte Carlo approach V. Jansson , L. Malerba , Simulation of the nanostructure evolution under irradiation in Fe–C alloys . J. Nucl. Mater. 443 ( 1-3 ), 274–285 ( 2013 ). doi: 10.1016/j.jnucmat.2013.07.046 http://doi.org/10.1016/j.jnucmat.2013.07.046 Simulation of the nanostructure evolution under irradiation in Fe–C alloys A.V. Barashev , S.I. Golubov , H. Trinkaus , Reaction kinetics of glissile interstitial clusters in a crystal containing voids and dislocations . Phil. Mag. A 81 ( 10 ), 2515 - 2532 ( 2001 ). doi: 10.1080/01418610108217161 http://doi.org/10.1080/01418610108217161 Reaction kinetics of glissile interstitial clusters in a crystal containing voids and dislocations B.D. Wirth , G.R. Odette , D. Maroudas et al. , Dislocation loop structure, energy and mobility of self-interstitial atom clusters in bcc iron . J. Nucl. Mater. 276 ( 1-3 ), 33 - 40 ( 2000 ). doi: 10.1016/s0022-3115(99)00166-x http://doi.org/10.1016/s0022-3115(99)00166-x Dislocation loop structure, energy and mobility of self-interstitial atom clusters in bcc iron N. Castin , G. Bonny , A. Bakaev et al. , Object kinetic Monte Carlo model for neutron and ion irradiation in tungsten: Impact of transmutation and carbon impurities . J. Nucl. Mater. 500 , 15 - 25 ( 2018 ). doi: 10.1016/j.jnucmat.2017.12.014 http://doi.org/10.1016/j.jnucmat.2017.12.014 Object kinetic Monte Carlo model for neutron and ion irradiation in tungsten: Impact of transmutation and carbon impurities L. Chen , K. Nishida , Kenta Murakami et al. , Effects of solute elements on microstructural evolution in Fe-based alloys during neutron irradiation following thermal ageing . J Nucl. Mater. 498 , 259 - 268 ( 2018 ). doi: 10.1016/j.jnucmat.2017.10.026 http://doi.org/10.1016/j.jnucmat.2017.10.026 Effects of solute elements on microstructural evolution in Fe-based alloys during neutron irradiation following thermal ageing E. Vincent , C. S. Becquart , C. Pareige et al. , Precipitation of the FeCu system: a critical review of atomic kinetic Monte Carlo simulations . J. Nucl. Mater. 373 ( 1-3 ), 387 - 401 ( 2008 ). doi: 10.1016/j.jnucmat.2007.06.016 http://doi.org/10.1016/j.jnucmat.2007.06.016 Precipitation of the FeCu system: a critical review of atomic kinetic Monte Carlo simulations Y. L. Bouar , F. Soisson , Kinetic pathways from embedded-atom-method potentials: Influence of the activation barriers . Phys. Rev. B 65 ( 9 ), 094103 ( 2002 ). doi: 10.1103/physrevb.65.094103 http://doi.org/10.1103/physrevb.65.094103 Kinetic pathways from embedded-atom-method potentials: Influence of the activation barriers I.M. Lifshitz , V.V. Slyozov , The kinetics of precipitation from supersaturated solid solutions , J. Phys. Chem. Solids 19 ( 1–2 ), 35 - 50 ( 1961 ). doi: 10.1016/0022-3697(61)90054-3 http://doi.org/10.1016/0022-3697(61)90054-3 The kinetics of precipitation from supersaturated solid solutions L. Liu , K. Nishida , K. Dohi et al. , Effects of solute elements on hardening and microstructural evolution in neutron-irradiated and thermally-aged reactor pressure vessel model alloys . J. Nucl. Sci. Tech. 53 ( 10 ): 1546 - 53 ( 2016 ). doi: 10.1080/00223131.2015.1136902 http://doi.org/10.1080/00223131.2015.1136902 Effects of solute elements on hardening and microstructural evolution in neutron-irradiated and thermally-aged reactor pressure vessel model alloys S. Shu , N. Almirall , P.B. Wells et al. , Precipitation in Fe-Cu and Fe-Cu-Mn model alloys under irradiation: Dose rate effects . Acta Mater. 157 , 72 - 82 ( 2018 ). doi: 10.1016/j.actamat.2018.07.017 http://doi.org/10.1016/j.actamat.2018.07.017 Precipitation in Fe-Cu and Fe-Cu-Mn model alloys under irradiation: Dose rate effects The authors of Zi-Qin Shen, Sha-Sha Lv, De-Sheng Ai and Zheng-Cao Li contributed to the OKMC simulation processing and design of the study. The detailed OKMC simulation was performed and analyzed by Zi-Qin Shen and Jie Gao. The discussion and analyses were performed by Sha-Sha Lv, Liang Chen and Dong-Yue Chen. The original draft of the manuscript was written by Zi-Qin Shen, and all the authors have reviewed and commented on the manuscript. 10.1007/s41365-022-01122-x This work was supported by the National Natural Science Foundation of China (Nos. 11975135 and 12005017) and China Postdoctoral Science Foundation (No. 2021M701829). 2022-05-13 2022-09-08 2022-09-11 2022-11-19 2022-12-16 © China Science Publishing & Media Ltd. (Science Press), Shanghai Institute of Applied Physics, the Chinese Academy of Sciences, ChineseNuclear Society and Springer Nature Singapore Pte Ltd. 2022 This is a post-peer-review, pre-copyedit version of an article published in Nuclear Science and Techniques. The final authenticated version is available online at: http://dx.doi.org/10.1007/s41365-022-01122-x http://dx.doi.org/10.1007/s41365-022-01122-x . Springer Nature or its licensor (e.g. a society or other partner) holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law. 2022 Nuclear Science and Techniques 11 33