Wang Yuan-Yuan, Ding Jian-Hua, Liu Wen-Bo, Huang Shao-Song, Ke Xiao-Qin, Wang Yun-Zhi, Zhang Chi, Zhao Ji-Jun. Irradiation-induced void evolution in iron: A phase-field approach with atomistic derived parameters. Chinese Physics B, 2017, 26(2): 026102
Permissions
Irradiation-induced void evolution in iron: A phase-field approach with atomistic derived parameters
Wang Yuan-Yuan1, Ding Jian-Hua1, Liu Wen-Bo2, Huang Shao-Song1, Ke Xiao-Qin3, Wang Yun-Zhi3, 4, Zhang Chi5, Zhao Ji-Jun1, †
Key Laboratory of Materials Modification by Laser, Ion and Electron Beams (Ministry of Education), Dalian University of Technology, Dalian 116024, China
Department of Nuclear Science and Technology, Xi’an Jiaotong University, Xi’an 710049, China
Center of Microstructure Science, Multi-Disciplinary Materials Research Center, Frontier Institute of Science and Technology, State Key Laboratory for Mechanical Behavior of Materials, Xi’an Jiaotong University, Xi’an 710049, China
Department of Materials Science and Engineering, Ohio State University, Columbus, Ohio 43210, USA
Key Laboratory of Advanced Materials of Ministry of Education, School of Materials Science and Engineering, Tsinghua University, Beijing 100084, China
A series of material parameters are derived from atomistic simulations and implemented into a phase field (PF) model to simulate void evolution in body-centered cubic (bcc) iron subjected to different irradiation doses at different temperatures. The simulation results show good agreement with experimental observations — the porosity as a function of temperature varies in a bell-shaped manner and the void density monotonically decreases with increasing temperatures; both porosity and void density increase with increasing irradiation dose at the same temperature. Analysis reveals that the evolution of void number and size is determined by the interplay among the production, diffusion and recombination of vacancy and interstitial.
The development of fast reactor and fusion reactor requires advanced structural materials which can tolerate higher temperature and higher irradiation dose under harsh service environment.[1] Owing to their excellent irradiation swelling resistance and good fracture toughness, ferrite/martensite (FM) steel[2,3] has been developed and considered as a candidate structural material for advanced nuclear reactor applications. Body-centered cubic (bcc) iron matrix is the main component (∼90 wt%) for the complex FM steel. Irradiation effects on bcc iron, including swelling, embrittlement and other degradation of the macroscopic mechanical properties,[4–6] have attracted significant attention in the last decades. These effects are microscopically attributed to the increases in vacancy and interstitial content during continuous irradiation, resulting in voids/bubbles, dislocation loops and precipitates in the material.
Numerous efforts have been devoted to investigating the irradiation effects of the structural materials in multiscale.[7–11] As a mesoscale simulation approach, the phase field (PF) method has been widely used to study microstructure evolution in various materials. With respect to irradiation effects, the PF method can well capture the kinetic and morphological evolution of irradiation-induced voids compared with irradiation experiments. The pioneer work on simulating the effect of irradiation on void evolution was reported by Yu and Lu.[12] It was assumed that void formation is caused by spinodal instability and its growth needs continuous supplies of vacancies. In 2009, Rokkam et al. developed a PF model on the basis of Cahn–Hilliard and Allen–Cahn equations to predict the stages of void incubation, nucleation and growth under irradiation.[13] Subsequently, this model was adopted for irradiated metal to simulate void and gas bubble nucleation and annihilation in intergranular and intragranular locations.[14–21] The influence of temperature gradient on void stability and growth kinetics, which is known as the Soret effect, was virtually observed by Li et al., and the void mobility was found to strongly depend on temperature gradient rather than void size.[22] In addition, the effects of grain sizes and intergranular gas bubble sizes on the effective thermal conductivity of the irradiated polycrystalline alloy in a three-dimensional PF model were investigated, and the growth kinetics and stability of interstitial loop were found to be highly dependent on the elastic interaction between the interstitial loops and point defects.[23–28]
So far, most PF simulations of void formation and growth are based on empirical parameters[13–18,22,23] and thus are limited by the lack of real parameters of specific materials. Therefore, real material parameters, such as void/interstitial formation energy, should be taken into account with the inclusion of the temperature effect.[18] This study aims to investigate the irradiation-induced void formation and growth by considering different temperatures and irradiation doses in junction with the atomistic derived material parameters. In Section 2, we describe the PF model for void formation during irradiation. The method of first-principle calculations and the empirical equations for deriving material parameters are presented in Section 3. In Section 4, the results of void evolution influenced by temperature and irradiation dose are presented and discussed. Major findings are summarized in Section 5.
2. Phase field model
As is well known, irradiation can generate vacancies and interstitials at the same time in the solid matrix; however, their kinetic behaviors and production during the microstructure evolution are distinctly different, that is, vacancy clusters form void and interstitial clusters form dislocation loop. Since the diffusion rate of interstitial is much faster than that of vacancy by orders of magnitude,[29] leading to the different kinetic time scales between vacancy and interstitial, the formation of dislocation loop and the effect of lattice strain relating to interstitials are not considered in the present PF model. Irradiation-induced void is a core with rich vacancies and its growth is primarily controlled by absorbing the surrounding irradiation-induced vacancies.[30]
In this work, we consider a system with two phases (i.e., matrix and void) to simulate void evolution during irradiation. Temporal and spatial variations of vacancy and interstitial concentrations in the matrix are characterized by two conserved fields, i.e., an interstitial concentration field ci(x,t) and a vacancy concentration field cv(x,t), where x is the spatial coordinate and t is the time. In addition, a long-range order parameter field, η(x,t), is introduced to distinguish the matrix phase ( and η = 0) from the void phase (cv = 1 and η = 1), with a finite-width diffuse interface in between.
Millett et al.[15] modified the classical Cahn–Hilliard and Allen–Cahn free energy equations for a non-uniform system by considering ci(x,t), cv(x,t), and η(x,t) to distinguish the matrix phase from the void phase. The expression of total free energy function is given by
(1)
The free energy of the matrix phase, fm(cv,ci), includes the enthalpy and entropy terms as a function of the vacancy and interstitial concentrations in the matrix and can be written as
(2)
where and are the formation energies of vacancy and interstitial, respectively; kB is the Boltzmann constant; T is the absolute temperature. The second term of the integrand fv(cv,ci) in Eq. (1) describes the free energy of void and is related to the concentrations of vacancies and interstitials by
(3)
This formula can well describe both the formation of a void due to vacancy cluster and the growth of a void due to absorption vacancies.
The gradient coefficients κv, κi, and κη, which depend on temperature, are related to void/matrix interfacial energy. In addition, the local free energies of the matrix and void phases are connected smoothly in the concentration-order parameter space by the following interpolation functions
(4)
(5)
The time evolutions of the conserved and non-conserved fields are governed by the Cahn–Hilliard and Allen–Cahn equations, respectively, with the consideration of production (source) and recombination (sink) of the irradiation-induced point defects, that is
(6)
(7)
(8)
with
(9)
(10)
where L is the mobility of the matrix and void interface, Mv and Mi are the chemical mobilities of vacancy and interstitial atom, respectively, and Dv and Di are the diffusivities of vacancy and interstitial atom, respectively. The second term ξ (x,t) in Eqs. (6) and (7) and the term ζ (x,t) in Eq. (8) represent the thermal fluctuation of the vacancy and interstitial concentration fields, respectively. In this model, we focus on the growth behavior of void, especially the influences of temperature and irradiation dose on the void morphology evolution. Instead of employing the classical nucleation theory, a noise term is used to trigger the void nucleation for simplicity.
The Frenkel pair production rates of vacancy Pv(x,t) and interstitial Pi(x,t) are used to describe the irradiation damage production. The relationship between Pv(x,t) and Pi(x,t) is given by
(11)
where R1 and R2 are the random numbers between 0 and 1 generated at each time step and each grid point, respectively, and B is the production bias and usually set to be 1.0.[15,18] The parameter Pcasc and VG are the probability of cascades occurring within a unit volume per unit time and the maximum increase in vacancy (or interstitial) concentration due to irradiation, respectively. The critical value of 0.8 assigned for η ensures the formation of void.[14,31] For example, if η < 0.8 and R1 ≤ Pcasc, vacancy only emerges in the matrix phase. Usually, the defects are generated stochastically on both spatial and temporal scales under irradiation.
The units of displacements per atom (dpa) per unit time, K, namely the irradiation dose rate, can be obtained from the following equation:
(12)
The last term Riv(x,t) in Eqs. (6) and (7) is the vacancy–interstitial recombination rate, which is a function of the vacancy and interstitial concentrations
(13)
where Rr depends on both bulk and surface terms with the expression Rr = Rbulk + η2Rsurf.
The above void kinetic equations are solved numerically in a two-dimensional (2D) domain with a uniform square grid and the periodic boundary condition. Three different computational cell sizes are used, i.e., 128 × 128, 256 × 256, and 512 × 512, to examine the influence of spatial resolution. It is found that the simulation results concerning void size, density and porosity are independent of system size. In this work, we present all the simulation results obtained with 256 × 256 grid points. The time step is and the grid spacing is , where and are the reduced units of time and length, respectively. and can be expressed as = l/l* and = t/t* (l* and t* are the real length and time, respectively), respectively. The relationship between l* and t* is described by t* = l*2/Dv. The length scale l* for PF model can be assumed to be about 1.5 times greater than the lattice constant of material,[28] and thus l* is taken to be 0.5 nm here. Hence, the time scale t* can be derived. A similar description of the converting dimensionless parameters into real units can be found in Refs. [15], [30], and [31]. The kinetic equations are solved by using the finite difference method.
3. Calculation of material parameters
Ab initio calculations offer the possibility to qualitatively predict the thermodynamic parameters of the material used in PF simulations. Here ab initio molecular dynamics (AIMD) simulation implemented in the Vienna ab-initio simulation package (VASP),[32,33] which is based on density functional theory (DFT) and plane-wave pseudopotential technique,[34,35] is carried out to determine the formation energy of a vacancy or self-interstitial in iron solid. The interactions between ions and electrons are described by the projector-augmented wave (PAW) potential.[36,37] The electron exchange–correlation interactions are treated by the generalized gradient approximation (GGA) with the Perdew–Burke–Ernzerhof (PBE) functional.[38,39] The plane-wave kinetic energy cutoff is chosen to be 400 eV. A 128-atom bcc supercell formed by 4 × 4 × 4 unit cells is used to model the perfect crystal of iron. The Brillouin zones are sampled by Monkhorst–Pack scheme[40] with a uniform mesh of 2 × 2 × 2 k points. Spin-polarization is considered for all DFT calculations due to the ferromagnetism of bcc iron. The calculated lattice parameter, cohesive energy, magnetic moment and bulk modulus of bcc iron are listed in Table 1, by comparing with experimental data, to validate the present first-principles method.
Table 1.
Table 1.
Table 1.
Values of lattice parameter (in unit Å), cohesive energy (in unit eV), magnetic moment per atom (in unit μB), and bulk modulus (in unit GPa) for bcc iron.
Values of lattice parameter (in unit Å), cohesive energy (in unit eV), magnetic moment per atom (in unit μB), and bulk modulus (in unit GPa) for bcc iron.
.
The MD time step is 1.0 fs and the total simulation time is 15 ps. The AIMD simulations within the canonical ensembles (NVT) use the Nosé–Hoover thermostat to control the temperatures at 293 K, 513 K, 622 K, 722 K, and 822 K, respectively. The temperature-dependent lattice parameters of bcc iron are taken from the experimental data.[45,46] At a given temperature, the formation energy of a vacancy or interstitial is obtained by averaging the free energy of each step of the last 10-ps AIMD simulation after the first 5-ps equilibration. At temperatures of 293 K, 513 K, 623 K, 722 K, and 822 K, the vacancy formation energies are 2.058 eV, 2.162 eV, 2.243 eV, 2.242 eV, and 2.288 eV, respectively, and the interstitial formation energies are 3.384 eV, 3.277 eV, 3.145 eV, 3.035 eV, and 3.002 eV, respectively. According to our DFT calculations, producing an interstitial in the perfect crystal of iron requires more energy than producing a vacancy. More interestingly, the formation energies of vacancy and interstitial show opposite trends of temperature dependence, i.e., increases with T but decreases. The values of and have been determined by previous experiments and ab initio calculations,[47–49] which gave ∼2 eV for and 3 eV .
Diffusivities of vacancy (Dv) and interstitial (Di) in iron host are related to temperature via the Arrhenius diffusion equation
(14)
where D0 is the diffusion constant and Ea is the activation energy. For interstitial diffusion, Ea is equal to the migration energy; while Ea is the sum of formation energy and motion energy for vacancy diffusion.[50] Here we adopt the reported values of and from previous atomistic simulation.[8,49–51] The value of D0 can be estimated from
(15)
where a is the jump length derived by (a0 is the lattice constant of iron), v is the vibration frequency and approximately given by v = (2Ea/ma2)1/2 (m = 9.27 × 10−26 kg is the mass of iron atom). Finally, the diffusivities of vacancy and interstitial in iron from Eq. (14) are obtained as follows:
(16)
(17)
Compared with the numerical simulations on diffusivity coefficients of vacancy and interstitial,[52,53] the calculated Dv and Di are within the same magnitude. For example, at 575 K, the diffusion values of vacancy and interstitial are 7.5 × 10−14 m2⋅s−1 and 7.1 × 10−10 m2⋅s−1 from previous atomistic simulations, respectively, and 2.7 × 10−14 m2⋅s−1 (vacancy) and 2.3 × 10−11 m2⋅s−1 (interstitial) calculated from Eqs. (16) and (17), respectively.
The gradient coefficient κ is used to describe the interface between void and matrix. κv, κi, and κη used in this model are assumed to be the same. In general, κ is related to the interfacial energy γ as
(18)
where Δf is the energy barrier height between two stable phases. Actually, Δf depends on temperature and its value can be derived by minimizing the two free energies of matrix (Eq. (2)) and vacancy (Eq. (3)) under the condition of () and ( (). In a recent paper, Liu et al. pointed out that the interfacial energy anisotropy considered in the PF model plays an important role in determining the void shape.[31] For simplicity, here we assume that the interfacial energy is isotropic. According to the literature, the interface energy for bcc iron is 13.11 eV⋅nm−2.[54]
4. Results and discussion
As demonstrated in previous studies,[14–18] PF model is able to capture void formation and growth in material during irradiation, exhibiting many interesting phenomena on a mesoscale. One significant aspect is the formation of irradiation-induced void, leading to the changes of the material properties. The vacancies and interstitials produced in the form of Frenkel pairs are randomly distributed in bcc iron. The on-going production of vacancies leads to vacancy accumulation, which finally gives rise to void formation.
Figure 1 displays the PF simulation results on the morphology of void, including void number, size and distribution, within a single grain under different irradiation temperatures and doses. The obtained round-shaped voids are consistent with the experimental observations in the neutron irradiated iron.[6] With increasing irradiation temperatures, the number of voids decreases, while the void size increases until 722 K and then decreases at 822 K. At a given temperature, the number of voids and the void size both increase with increasing dose, especially at temperatures higher than 623 K.
Fig. 1. (color online) Evolutions of void morphology in the irradiated iron as a function of temperature and irradiation dose. The gradient color is used to represent the interface between void (red color) and matrix (blue color).
The temperature-dependent mean void diameters of 1 dpa irradiated iron are given in Fig. 2. With increasing temperature, the predicted void size from PF modeling using the temperature-dependent parameters (referred as SL1) increases from 513 K to 722 K and then it shows a slow decrease at 822 K. That is, there is a critical temperature to affect the saturation void size. Interestingly, previous experimental data from Ref. [55] show the same tendency and even the same magnitude of void size. Therefore, the simulation results prove that the atomistic derived material parameters are sufficiently accurate and are promising for further investigating microstructural evolution by irradiation in iron. Hayns and Williams pointed out that the diameter of void is a consequence of the balance between vacancy flux into and emission out of a void.[56] While for the void with more than a critical accumulation number of vacancies, it would collapse into a dislocation loop in order to reach a minimum energy configuration.[57]
Fig. 2. (color online) Comparison between PF and experiment results (red color) at different temperatures after 1-dpa irradiation. The simulation results using the temperature-dependent parameters are marked as SL1 (blue color); those that use a set of temperature-independent parameters are marked as SL2 (green color).
To further elucidate the role of temperature-dependent parameters used in the PF simulations, the values for bcc iron parameters (e.g., Ef, D, and κ) relying on temperature are now assumed to be temperature independent and are the same as those of 293 K (referred to as SL2). As seen in Fig. 2, SL2 yields the same evolution tendency of the average void size as SL1. However, the theoretical void sizes of SL2 are much higher than the experimental data and those of SL1. It clearly indicates that using temperature-independent material parameters is significant to influence the authenticity of the results. This may be partly attributed to the highly temperature-dependent thermal properties of vacancy and interstitial. For example, with temperature increasing, the formation energies of vacancy and interstitial show an increasing and decreasing tendency, respectively, and their diffusivities both increase. The competition among the formation, diffusion and recombination of vacancies and interstitials determines the void morphology.
As for 2D simulation, instead of swelling, porosity defined as the fraction of all voids over the total volume can be used to characterize the irradiation effect.[15,16,18,21]Figure 3 shows that the plots of porosity versus temperature over a wide range of irradiation dose. Overall speaking, higher dose induces larger void fraction, regardless of temperature. At low irradiation dose of 1 dpa, there is no noticeable change in the porosity under different temperatures. At higher doses above 2.5 dpa, void fraction shows obvious increase until 722 K, and then decrease at 822 K. In current PF model, irradiation dose effect is directly related to Frenkel pair production rate P(x,t), see Eq. (11). It can be inferred that the higher irradiation dose, the more vacancies are produced. The obtained tendency of temperature-dependent porosity, which is called bell-shaped curve, is consistent with the experimental observations in irradiated ferrite steels.[57–59] Previous theoretical work explained the influence of temperature on the material porosity as follows:[56,57] at temperatures below peak porosity (T < 722 K), the mobility of vacancy clusters is low, the recombination between vacancies and interstitials has little effect on porosity, resulting in a high vacancy density (i.e., the number of voids over the total simulation size). While when the temperature is higher than 722 K, the vacancy-interstitial recombination becomes a dominant role, leading to a reduction of vacancy density.
Fig. 3. (color online) Plots of porosity versus temperature at different irradiation doses of the saturated voids in the irradiated iron.
The plots of void density versus temperature at different irradiation doses are shown in Fig. 4. The void density monotonically decreases with increasing temperature for a given dose, while it increases with increasing dose as expected. In a dose range from 1 dpa to 10 dpa, there is an obvious increase of the void density in the irradiated microstructure. Meanwhile, the number of irradiation-induced voids shows little increment with irradiation dose above 10 dpa. Since an equal number of vacancies and interstitials are generated under the same dose, the increasing formation energy of vacancy may lead to self-instability and even degradation of the number of void nucleations during increasing temperatures. In short, the temperature effect on an irradiated material can be traced back to the vacancy stability as well as its kinetics, whereas irradiation dose plays a significant role in producing Frenkel pairs.
Fig. 4. (color online) Plots of void density versus temperature at different irradiation doses of the saturated voids in the irradiated iron.
It should be mentioned that the present 2D PF model is insufficient to give a complete description of the initial nucleation of void as well as the continuous growth process that is influenced by dislocation and lattice strain. Nevertheless, it does capture the essential features of void evolution (e.g., void size and density) affected by temperature and irradiation dose. In addition, there is a wide space to extend the current irradiation-induced void PF model for other candidate nuclear materials used in fast reactors and fusion reactors, such as vanadium alloys and SiC, aiming to gain an important understanding of void evolution in the irradiated materials. In future study, the consideration of void nucleation should be improved by introducing the explicit nucleation algorithm[30] and would include the effect of dislocation on the void evolution by using the results from molecular dynamics (MD) or other simulation methods.
5. Conclusions
With the real material parameters derived from atomistic calculations, the PF model is used to study the kinetics of void formation and growth in a single crystal of bcc iron under irradiation. Ab initio molecular dynamics simulations are performed to obtain the temperature-dependent formation energies of vacancy and interstitial. Other key parameters (e.g., vacancy/interstitial diffusivity) used in the PF model are calculated with data from atomistic simulations through using empirical formula.
Void evolution influenced by temperature and irradiation dose is predicted and analyzed. The predicted void size of the irradiated bcc iron increases with increasing temperature and then shows a slow decrease. This evolution tendency is in accordance with the experimentally measured void size. In addition, the void sizes predicted by the temperature-independent parameters in conjunction with the PF model show an obvious deviation from the temperature-dependent cases and the experimental data. Hence, it is necessary to consider the temperature effect on the material properties (Ef, D, and κ) to simulate the irradiation-induced voids in a real material. With increasing temperature, the simulated porosity exhibits a bell-shaped trend and the monotonical decrease of void density, which are observed in the irradiation experiments of ferrite steels. A higher irradiation dose induces more void formation, leading to an increase in porosity and void density. In conclusion, the current PF model with realistic material parameters incorporated is able to describe the void formation and evolution in bcc iron or ferrite steels under irradiation.
Based on the feasibility demonstrated with the present phase-field approach with atomistic derived parameters, our future work will focus on the development of the descriptions of the kinetic behaviors of vacancy and interstitial (e.g., to form dislocation and to induce lattice strain) within the PF framework to more accurately capture the irradiation-induced void evolution in the bcc materials.
