Determination of Krypton Diffusion Coefficients in Uranium Dioxide

Dec 16, 2016 - Determination of gaseous fission product behavior near the cerium dioxide Σ 3 (111)/[1 1 ¯ 0] tilt grain boundary via first-principle...
1 downloads 0 Views 3MB Size
Download PDF
Article pubs.acs.org/IC

Determination of Krypton Diffusion Coefficients in Uranium Dioxide Using Atomic Scale Calculations Emerson Vathonne,† David A. Andersson,‡ Michel Freyss,† Romain Perriot,‡ Michael W. D. Cooper,‡ Christopher R. Stanek,‡ and Marjorie Bertolus*,† †

CEA, DEN, DEC, Centre de Cadarache, 13108 Saint-Paul-lez-Durance, France Materials Science and Technology Division, Los Alamos National Laboratory, Los Alamos, New Mexico 87545, United States

‡

ABSTRACT: We present a study of the diffusion of krypton in UO2 using atomic scale calculations combined with diffusion models adapted to the system studied. The migration barriers of the elementary mechanisms for interstitial or vacancy assisted migration are calculated in the DFT+U framework using the nudged elastic band method. The attempt frequencies are obtained from the phonon modes of the defect at the initial and saddle points using empirical potential methods. The diffusion coefficients of Kr in UO2 are then calculated by combining this data with diffusion models accounting for the concentration of vacancies and the interaction of vacancies with Kr atoms. We determined the preferred mechanism for Kr migration and the corresponding diffusion coefficient as a function of the oxygen chemical potential μO or nonstoichiometry. For very hypostoichiometric (or U-rich) conditions, the most favorable mechanism is interstitial migration. For hypostoichiometric UO2, migration is assisted by the bound Schottky defect and the charged uranium vacancy, V4− U . Around stoichiometry, migration assisted by the charged 4− uranium−oxygen divacancy (V2− UO) and VU is the favored mechanism. Finally, for hyperstoichiometric or O-rich conditions, the migration assisted by two V4− U dominates. Kr migration is enhanced at higher μO, and in this regime, the activation energy will be between 4.09 and 0.73 eV depending on nonstoichiometry. The experimental values available are in the latter interval. Since it is very probable that these values were obtained for at least slightly hyperstoichiometric samples, our activation energies are consistent with the experimental data, even if further experiments with precisely controlled stoichiometry are needed to confirm these results. The mechanisms and trends with nonstoichiometry established for Kr are similar to those found in previous studies of Xe.

■

and one more recent by Michel et al.10 Table 1 gives the preexponential factor and the activation energies obtained in these two studies, as well as the results of Michel et al. on xenon for comparison. It can first be seen from the results by Michel et al. that krypton diffuses faster than xenon in uranium dioxide. Second, rather

INTRODUCTION The fission of uranium atoms in uranium dioxide causes the production of large quantities of fission gases such as xenon and krypton, which have a low solubility in the material.1,2 This phenomenon has important consequences for the fuel behavior, since it induces swelling and leads to the decrease of the thermal conductivity and the release of the gases outside the fuel, which are limiting factors for the fuel performance and reactor safety. Further insight is still needed into all stages of the fission gas behavior, starting with the elementary diffusion mechanisms active in segregation and bubble formation. Atomic scale modeling, in particular electronic structure calculations, combined with separate effect experiments, are now invaluable tools to elucidate these mechanisms. Xenon behavior in UO2 has been extensively studied experimentally2−6 and using atomic scale modeling as discussed in the review by Liu et al.7 Experiments and electronic structure calculations were compared in ref 8. Fewer studies have been performed on Kr, which is slightly less abundant in uranium dioxide. Experimentally, only two studies giving Arrhenius diffusion laws for Kr are available, one by Auskern9 from 1985 © XXXX American Chemical Society

Table 1. Experimental Results for the Diffusion of Krypton in UO2, as well as Results on Xe by Michel et al.10 D0 (m2·s−1) Auskern9 Michel et al.10 Michel et al.10 Michel et al.10 Michel et al.10

Krypton Diffusion 4.9 × 10−8 bulk 6.7 × 10−17 surface 2.7 × 10−10 Xenon Diffusion bulk 2.0 × 10−11 surface 3.9 × 10−13

Ea (eV) 3.2 1.4 2.7 3.1 2.1

Received: July 5, 2016

A

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry unexpectedly, the activation energy obtained for krypton in the bulk (1.4 eV), which should correspond to stoichiometric or slightly hyperstoichiometric UO2, is smaller than that for the surface (2.7 eV), which is likely to be more hyperstoichiometric due to the easy oxidation of UO2 even in the presence of very small quantities of oxygen. There are also only a few modeling studies of Kr behavior in UO2, which focus mainly on its incorporation.11−14 The only computed migration data published for Kr are those calculated by Catlow15 at the end of the 1970s using empirical potentials. He obtained activation energies of 5.1, 9.2, and 5.1 eV for migration as interstitial and assisted by neutral and charged trivacancies, respectively. Finally, Kr migration and release was investigated using cluster dynamics, a rate theory method, by Skorek et al.16 using parameters derived among others from the results of atomic scale calculations and separate effect experiments. We present here a systematic study of the diffusion of krypton in UO2 using atomic scale calculations combined with diffusion models adapted to the system studied. The elementary mechanisms for interstitial or vacancy assisted migration are calculated in the DFT+U framework using the nudged elastic band method (NEB) for pathway determination. The migration rates of these mechanisms are calculated from the migration barrier obtained from the NEB method and the attempt frequency obtained from the phonon modes of the defect at the initial and saddle points using empirical potential methods. The diffusion coefficients (effective activation energy and prefactor) of Kr in UO2 are calculated as a function of the oxygen chemical potential μO (or nonstoichiometry) using diffusion models accounting for the concentration of vacancies and the interaction of vacancies with Kr atoms. This paper is organized as follows. We first present the atomic scale calculation methods and the parameters that were used. Second, we describe the diffusion models employed to determine the effective activation energy of Kr considering various elementary mechanisms. Third, we show the results for the elementary mechanisms of vacancy assisted krypton migration in UO2, and the corresponding diffusion coefficient is derived. Fourth, we present the diffusion coefficient for interstitial migration. We then summarize the results for the Kr diffusion coefficient as a function of nonstoichiometry and present our conclusions. Finally, we discuss our results in view of the experimental data and of the electronic structure predictions for Xe by Andersson et al.8,17

as implemented in the Vienna Ab initio Simulation Package (VASP).20−22 A Hubbard-like term (U) is added in order to take into account the strong correlations between the 5f electrons of the uranium atoms according to the Liechtenstein scheme23 of the DFT+U method. The values of the U and J parameters were set to 4.50 and 0.51 eV, respectively, as originally derived by Dudarev24 from XPS results by Kotani and Yamazaki25 and optical spectroscopy measurements by Krupa and Gajek.26 We have used these values extensively since we started to apply the DFT+U method to UO2 in 2009 and have shown that they enable one to obtain a good description of the UO2 properties.27 In order to avoid the convergence to one of the metastable states yielded by the DFT+U method and ensure that the ground state is reached, the occupation matrix control scheme28 was used. We optimized geometries using the generalized gradient approximation (GGA) functional parametrized by Perdew, Burke, and Ernzerhof (PBE)29 in the DFT+U formalism. The energies of the structures obtained were then calculated, still in DFT+U, using the vdW-DF functional, a nonlocal correlation functional developed by Dion et al.30 and implemented in VASP by Klimes et al.31,32 In this functional, the exchange−correlation energy is calculated as follows: Exc = ExDFT + EcDFT + Ecnl

EDFT x

(1)

standard exchange energy, EDFT c energy (LDA functional), and Enlc

where, is the is the standard correlation is the nonlocal energy term based on electron densities interacting via a model response function, whose particular form is still a subject of research.31 We used the optPBE functional for EDFT x . This functional enables one to take approximately into account dispersion interactions, which will ensure a better description of the bonds formed by Kr.33,34 In addition, it has been shown that these functionals yield better energies for ionic solids.35 Table 2 shows the results on bulk UO2 yielded by the Table 2. Results on Bulk UO2 Yielded by the VdW-DF+U (opt-PBE) Functional Compared to the GGA+U (PBE) and Experimental Results

UO2 cell parameters (Å) UO2 bulk modulus (eV) UO2 elastic constants (eV)

■

ATOMIC SCALE CALCULATION METHODS AND PARAMETERS The calculation of diffusion coefficients (see below) implies the calculation of a large number of parameters. On the one hand, defect formation energies, Kr solution energies, and elementary migration energies of defects and Kr must be calculated to obtain the effective activation energies. On the other hand, phonon frequencies must be calculated to evaluate the attempt frequencies and the diffusion coefficient prefactors. The energies needed were calculated using electronic structure methods in the DFT+U framework. Because of the computational cost involved, the calculation of the phonon frequencies were performed using empirical potentials. We describe the details of the methodology used in our work below. Calculation of Defect Formation Energies and Migration Pathways Using DFT. DFT calculations were performed using the projector augmented wave method (PAW)18,19

Formation enthalpy (eV)

a b c11 c12 UO2 U4O9a

GGA +U

VdW-DF +U

experimental

5.57 5.49 194 361 110 −10.54 −42.46

5.58 5.51 195 364 110 −11.23 −45.50

5.47 5.47 209 389 119 −11.23 −46.53

a For the calculations of the U4O9 formation enthalpy, an approximate 58-atom supercell corresponding to a U4O8.889 composition was considered.36

VdW-DF+U functional using the parameters described below, as well as the GGA+U and experimental results for comparison. It can be seen in this table that the structural parameters yielded by the VdW-DF+U are very similar to those obtained in GGA+U. The formation enthalpies of uranium oxides, however, are significantly improved when using the VdW-DF functional: Formation enthalpies of −11.23 eV and −45.50 eV are obtained for UO2 and U4O8.889, respectively, using the vdWDF functional, compared to −10.54 eV and −42.46 eV using PBE and to the −11.23 eV and −46.53 eV experimental values for UO2 and U4O9.37,38 B

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry

atoms in UO2 at 0 K, μe is the electron chemical potential varying from the bottom to the top of the band gap (EVBM), ΔEel is the Madelung term, and ΔV is the potential alignment, which corrects the shift of the band structure between the perfect cell and the defective cell. Using the experimental values for the formation enthalpies of UO2 and U4O9, three stoichiometric regimes can be defined as follows: • hypostoichiometric conditions (or U-rich) with −9.86 < μO < −7.96 eV and −11.12 < μU < −7.33 eV • stoichiometric conditions with μO = −7.96 eV and μU = −11.12 eV • hyperstoichiometric conditions (or O-rich) with −7.96 < μO < −6.06 eV and −14.92 < μU < −11.12 eV The uranium and oxygen chemical potentials in UO2, μU and μO, are linked, and formation energies can be expressed as a function of μO alone (as well as of the electron chemical potential μe for charged defects). From the DFT+U results, we obtain the expressions as follows for the V0U and V4− U vacancies (in eV):

A plane wave cutoff energy of 500 eV and a Gaussian smearing of 0.1 eV were used. All initial point symmetries were removed during calculations. The lattice parameters and atomic positions were relaxed during the geometry optimizations and a convergence of the Hellmann−Feynman forces of less than 0.01 eV/Å was ensured. The total energies of the systems were converged to less than 0.01 meV/atom for all calculations. Uranium dioxide is paramagnetic above 30.8 K and antiferromagnetic below this temperature. Due to the cost and complexity of modeling paramagnetic solids using DFT+U,39 we considered a 1k antiferromagnetic ordering, which is an approximation of the noncollinear 3k order observed experimentally at low temperature.40 This is justified by the fact that the two AFM orderings exhibit only small differences in energy41 and that the 1k antiferromagnetic order is very good approximation of the paramagnetic state for the energetic properties of UO2.42 Spin−orbit coupling (SOC) was neglected because of the computational costs and since previous studies accounting for spin−orbit coupling showed only small differences in the bulk properties of UO2, including energetics, compared to those obtained without spin−orbit coupling.41,43,44 In addition, a preliminary study of the influence of SOC on the UO2 defect energies calculated in DFT+U was performed using 12-atom supercells.45 The results showed only small differences between the results with and without SOC: 0.1 and 0.02 eV for the neutral oxygen vacancy and interstitial, respectively. This is negligible compared to the energy scale that governs the cation and Kr diffusion properties. Migration pathways were determined using the standard nudged elastic band method (NEB).46 Several migration pathways were calculated in VdW-DF+U, and the results were compared to the calculation of the energy using VdW-DF of the image configurations yielded by GGA+U. Given the very slight differences in the geometries and energies between the two calculation methods, most elementary migration pathways were determined using the latter method. Particular care was taken to avoid the convergence of the image configurations to metastable states. Two types of supercells were used. Defect energies for the monovacancies, U−O divacancy, and bound Schottky defects were calculated in supercells containing 96 atomic sites (2 × 2 × 2 repetitions of the fluorite unit cell) using a 2 × 2 × 2 Monkhorst−Pack k-point mesh.47 For all NEB calculations, we used supercells with 144 atomic sites (2 × 2 × 3 repetitions of the unit cell) and a 2 × 2 × 1 Monkhorst−Pack k-point mesh. Defect Formation Energies and Stoichiometry Regimes. Formation energies of defects in solids, and consequently activation energies for migration, vary with the material composition and the electron chemical potential. This variation is all the more significant in UO2 since the range of stoichiometries accessible is quite large. We calculated the formation energies of the various vacancies involved in the Kr migration mechanisms using the formalism presented in ref 48 derived from the work by Zhang and Northrup:49

(3)

Ef (V U4−) = 21.96 − 2μO − 4(E VBM + μe )

(4)

At the middle of the band gap, for μe = 1.2 eV Ef (V U4−) = −13.24 − 2μO

(5)

Defect formation energies in UO2 therefore depend on the oxygen chemical potential, μO, that is, on nonstoichiometry. For comparison with experiments, formation energies should be calculated using the experimental oxygen chemical potential. Similarly to defect formation energies, the Kr solution energy, that is, the energy necessary to incorporate Kr in a perfect UO2 crystal, is calculated as Esol(Kr, VX , q) = Etot(VX , Kr, q) − nOμO − nUμU − nKrμKr + ΔEel(q) + q(E VBM + μe + ΔV )

(6)

where Etot(VX, Kr, q) is the total energy of the supercell with the charge q containing the defect X and Kr, nKr is the number of Kr atoms in the cell containing the defect, and μKr is the chemical potential of atomic Kr at 0 K. Calculation of Attempt Frequencies for Elementary Migration Mechanisms. The attempt frequencies, νi, of the elementary mechanisms ωi involved in the diffusion models were obtained using the Vineyard method:50 νi =

3N − 3

3N − 4

∏

νndefect / ∏ νnsaddle

n=1

n=1

(7)

where ν and ν are the normal modes of vibration for the initial configuration and the migration saddle. The three zerofrequency acoustic modes are excluded, as well as one extra negative mode corresponding to the saddle mode. The modes are obtained by diagonalization of the dynamical force matrix of the system (Hessian). These costly calculations were performed using a rigid-ion empirical potential of the Buckingham form:51 defect

Ef (VX , q) = Etot(VX , q) − nOμO − nUμU + ΔEel(q) + q(E VBM + μe + ΔV )

Ef (V 0U) = −8.76 − 2μO

(2)

where Etot(VX,q) is the total energy of the supercell with the charge q containing the defect X, nU and nO are the numbers of uranium and oxygen atoms in the cell containing the defect, μU and μO are the chemical potentials of uranium and oxygen

saddle

⎛ r⎞ C V (r ) = A exp⎜ − ⎟ + 6 ⎝ ρ⎠ r C

(8) DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry Each mechanism is thermally activated so that

Table 3. Parameterization of the Buckingham Potential Used in This Study interaction

A (eV)

ρ (Å)

C (eV/Å6)

ref

O2−−O2− U4+−O2− Kr0−O2− Kr0−U4+

9547.96 1761.78 800.38 5912.78

0.2192 0.3582 0.3888 0.3191

32.0 12.3 52.48 71.84

54 55 56 56

ωi = νi e−Ei /(kT )

(9)

where Ei is the activation energy for the jump ωi and νi is its attempt frequency. It is important to note that in the absence of interaction between atoms further than first nearest neighbors, the association and dissociation rates are related to the interaction between the impurity and the vacancy, so the impurity−vacancy binding energy can be calculated as EB = E4 − E3. Using this definition, EB is negative when the vacancy and impurity are bound. In the general case, the effective mechanism for vacancyassisted migration of the impurity is the impurity-trap exchange ω2.57,58 The diffusion coefficient can then be written as

The parametrizations for A, ρ, and C for the O−O, U−O, Kr−O, and Kr−U are shown in Table 3. In addition, the formal +4 and −2 charges were used for U cations and O anions. A short-range cutoff of 8 Å and a 4 × 4 × 4 supercell (768 atoms for stoichiometric UO2) were used. The lateral dimensions (∼20 Å, more than twice rc) were chosen big enough to avoid the influence of the image cells in the periodic boundary condition formalism. It is to be noted that, since the charges are fixed with the Buckingham potential, the empirical potential calculations consider exclusively charged defects. The saddle configurations were obtained by performing climbing image nudged elastic band calculations52 using the molecular dynamics code CLSMAN.53 Eleven or nine images were used to determine the saddle point for all mechanisms except for the VU migration in bulk UO2, where more images were needed to resolve the saddle point, which is located near a very shallow minimum. The NEB calculations and optimizations were performed at constant volume, fixed at the equilibrium one for stoichiometric UO2.

D=

1 2 d 2 pf2 ω2 6

(10)

where d2 is the jump distance of the impurity from U to U site, equal to a/√2, where a is the UO2 cell parameter; p is the probability of the presence of the assisting defect in a neighboring site of the impurity. This probability is expressed as a function of the number of possible adjacent sites n2 (equal to 12 in a fcc structure) and of the formation energy of a vacancy at that site, that is, the formation energy of the assisting vacancy in the bulk, Ef , plus the binding energy between the impurity and the vacancy, EB. This is due to the fact that a second assisting defect must be present next to the krypton sitting in a vacancy to enable migration. If the vacancy and impurity attract (EB < 0), the probability of finding a vacancy in a first nearest neighbor site is enhanced with respect to the bulk concentration. This expression of p assumes that the defect cluster involved in migration is also the lowest energy trap. In the cases when the cluster must change from the lowest energy trap to the migrating defect cluster, an additional reconfiguration energy must be included, ER, which is the difference in Kr solution energy between the defect cluster and the lowest energy defect. Finally, p can be expressed as

■

DIFFUSION MODELS Five-Frequency Model. Vacancy-assisted migration of a solute involves several subsequent and correlated steps. In uranium dioxide, the uranium ions sit on a fcc lattice in the fluorite structure. In a fcc lattice, the vacancy-assisted migration is commonly described by the diffusion model developed by Lidiard57,58 called the five-frequency model. This model assumes that there is no interaction between atoms further than first neighbors. The five distinct elementary mechanisms of this model, shown in Figure 1, are as follows: • ω0, uranium vacancy migration between two first nearest neighbor sites of the lattice without the presence of impurity. • ω1, uranium vacancy migration between two first nearest neighbor sites around the impurity. • ω2, krypton migration toward a neighboring uranium vacancy (krypton−uranium vacancy exchange). • ω3, uranium vacancy migration from a first nearest neighbor site of the impurity to a second nearest neighbor site (dissociation). • ω4, reverse transition of ω3 (association).

p = n2 e−(Ef + EB + ER )/(kT )

(11)

f 2 is the correlation factor between two ω2 jumps, which reflects the fact that the various atomic hops are not independent and that the relative directions of two subsequent steps determine if there will be net migration. It has been calculated by Manning,60 and for a fcc structure, it is equal to f2 =

ω1 + 7F3ω3/2 ω1 + ω2 + 7F3ω3/2

(12)

Figure 1. Five frequency model (adapted from Mehrer59). D

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry where F3 is the vacancy escape probability, that is, the likelihood that if the vacancy dissociates from the impurity, it will not ω return. It can be calculated as a function of α = ω4 :

Application of the Five-Frequency Model to UO2. Uranium dioxide is a binary compound, so the application of the five-frequency model to this material calls for approximations. Given the size of the impurity considered (Kr) and the large incorporation energy in the O vacancy (6.37 eV using vdW-DF+U), we apply the five-frequency model to the U sublattice only, which is also a fcc lattice. The reorganization of the O sublattice is taken into account in the DFT+U calculations of the elementary pathways. This implies, however, that mechanisms involving only the O sublattice are not taken into account. In this approximation, U atoms are considered to be first nearest neighbors, whereas they are formally second nearest neighbors in the UO2 lattice. The distance between second nearest neighbors in the U lattice is large enough to make the assumption that defect and impurity are noninteracting as second nearest neighbors, which is the main assumption of the five-frequency model. In the following, we call trap sites the vacancies that can accommodate the Kr atoms and mobile vacancy the second vacancy assisting the migration. We then need to determine the most probable trap sites for Kr migration. The formation energies of the trap sites involving one U vacancy are reported in Table 4 for various stoichiometric regimes, and the

0

10α 4 + 180.5α 3 + 927α 2 + 1341α 7F3(α) = 7 − 2α 4 + 40.2α 3 + 254α 2 + 597α + 436 (13)

F3 takes different values depending of the ratio between ω0 and ω4. F3 is often described as the escape probability when the vacancy is in a second nearest site. Hence it depends on the rate of the association versus bulk diffusion. Three limiting cases can be noted: • For ω4 ≫ ω0, the vacancy moves faster from second to first neighbor position of the impurity than in the bulk, that is, the latter is the limiting step; then α → ∞ and F3 = 2/7. • In contrast, for ω4 ≪ ω0, the vacancy moves faster far from the impurity than from second to first neighbor of the impurity; then α = 0 and F3 = 1. • For ω4 = ω0, then α = 1 and F3 = 0.7375. Combining eqs 9, 10, 11, and 12, we obtain D = a 2ν2

ω1 + 7F3ω3/2 e−(E2 + Ef + EB + ER )/(kT ) ω1 + ω2 + 7F3ω3/2

(14)

Table 4. Formation Energies of the Uranium Vacancy, Uranium−Oxygen Divacancy, and Bound Schottky Defect Yielded by vdW-DF+U Calculations at the Middle of the UO2 Band Gap

The diffusion factor does not in general follow an Arrhenius law since the ωi rates, and therefore the prefactor, depend on the temperature. The non-Arrhenius behavior, however, would only occur in rather narrow temperature regions where transitions between mechanisms take place. We list below the cases for which D follows an Arrhenius behavior and the resulting expressions of D. • For ω3 ≪ ω1, ω2, dissociation jumps are very unlikely. Then ω1 f2 = ω1 + ω2 (15)

trap site

UO2−x

UO2

UO2+x

V0U V4− U V0UO V2− UO V0UO2

10.96 6.48 7.09 4.44 3.76

7.17 2.68 5.19 2.54 3.76

3.37 −1.12 3.30 0.65 3.76

• If ω1 ≪ ω2, the vacancy−impurity exchange is much faster than the vacancy jump to the vicinity of the impurity. The impurity often moves back and forth between two adjacent sites. The correlation factor is then ω ν f2 = 1 = 1 e−(E1− E2)/(kT ) ω2 ν2 (16) In this case the diffusion coefficient can be expressed as D = a 2ν1 e−(E1+ Ef + EB + ER )/(kT )

(17)

• If ω1 ≫ ω2, the vacancy−impurity exchange is much slower than the vacancy jump to the vicinity of the impurity. Then f 2 = 1 and

Figure 2. Solution energies of Kr in the uranium vacancy, uranium− oxygen divacancy, and bound Schottky defect, as well as interstitial, yielded by vdW-DF+U calculations at the middle of the UO2 band gap.

D = a 2ν2 e−(E2 + Ef + EB + ER )/(kT ) (18) • For ω3 ≫ ω1, ω2, then dissociation of the impurity− vacancy pair is the fastest mechanism, f 2 = 1, and D = a 2ν2 e−(E2 + Ef + EB + ER )/(kT )

corresponding Kr solution energies are represented as a function of the oxygen chemical potential μO in Figure 2. For the bound Schottky defect, three configurations are possible with oxygen vacancies in the (100), (110), and (111) direction (BSD1, BSD2, and BSD3, respectively). We consider the BSD1 here, even though BSD2 is slightly more stable in the calculations. It is the most compact one and the most favorable for Kr incorporation (less than 0.01 eV incorporation energy

(19)

1

• For ω1 = ω2 ≪ ω3, then f2 = 2 1 D = a 2ν1 = 2 e−(E1=2 + Ef + EB + ER )/(kT ) 2

(20) E

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry in vdW-DF+U), as well as the Kr trap observed experimentally using X-ray absorption spectra (XAS) in Kr implanted UO2.61 In addition, this Schottky configuration leads to the most simple migration mechanism for Kr. It can be seen from Table 4 and Figure 2 that the most stable trap vacancies and the most favorable ones for Kr incorporation are • the BSD1 Schottky defect for UO2−x • the V2− UO divacancy in stoichiometric UO2 • the V4− U monovacancy in UO2+x We will therefore study the migration of Kr in these three traps, shown in Figure 3. For the VU monovacancy and the VUO

divacancy, we will also study the neutral defects. In all cases, we will consider the U monovacancy (neutral or charged) as second assisting vacancy (mobile vacancy). The Catlow Model. The five-frequency model must be rearranged for the special case when a stable bound complex is formed between the impurity and the two assisting defects, as observed for Xe with the various clusters containing two U vacancies in UO2.15,17 In this case, there is no vacancy−solute exchange mechanism, and migration is due to the relocation of the impurity in the middle of the vacancy cluster after the ω1 jump of the mobile U vacancy around the Kr. This mechanism is shown in Figure 4.

Figure 4. Mechanism for impurity migration in the presence of a bound complex between the impurity and the two assisting defects: ω1 jump of the second U vacancy and relocation of the impurity in the middle of the vacancy cluster.

The effective mechanism for migration is now ω1, and the diffusion coefficient can be written as 1 D = d12p′f1 ω1 (21) 6 where d1 is the jump distance of the impurity from its position in the middle of two adjacent U sites to the nearest equivalent position and is therefore equal to a/2√2; p′ is the probability of the presence of the assisting defect in a neighboring site of the impurity. p′ = n1′ e−(Ef + EB + ER )/(kT )

(22)

n′1 is the number of nearest neighbors of Kr in its trap (U, UO, or UO2 vacancy) and is equal to 12 in the fcc structure. f1 is the correlation factor between two ω1 steps. Following the reasoning presented by Mehrer,59 we determined f1, which can be expressed as a function of the jump probability averaged value of the cosine of the angle between two consecutive jumps, θ: Figure 3. Possible assisting vacancies for the migration of krypton in UO2: (a) uranium vacancy, (b) U−O divacancy, and (c) neutral Schottky defect (BSD1).

f1 = F

1 + ⟨cos θ ⟩ 1 − ⟨cos θ ⟩

(23) DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry Z

⟨cos θ ⟩ =

∑ Pl cos θl l=1

and independent of the prior jump. This will be true especially in hypostoichiometric and stoichiometric regimes. ωint is the rate of the interstitial jump:

(24)

where θl is the angle between two possible subsequent jumps with destination the site l and Pl is the corresponding probability. A computation of Pl should take into account the infinity of possible trajectories leading from the initial site to the site l. If we suppose that a first jump (in our case ω1) has taken place, the simplest estimate can be obtained by considering only the trajectory leading to a jump back to the initial impurity site. In this case, ⟨cos θ ⟩ = Pback cos(180°)

ωint = νint e−Eint /(kT )

Similarly to the vacancy-assisted case, an additional reconfiguration energy, ER, which is the difference in solution energy of Kr in interstitial position and in the lowest energy defect, must be included. For low interstitial concentrations, the interstitial diffusion coefficient can then be expressed as 1 2 D= a νint e−(Eint + ER )/(kT ) (33) 12

(25)

■

This probability can also be expressed in terms of the jump rates involved in the model. ω1 Pback ≈ n1ω1 + n3ω3 (26)

ELEMENTARY MECHANISMS OF KR VACANCY-ASSISTED MIGRATION Elementary Activation Energies. We report in Table 5 the activation energies obtained in vdW-DF+U using the

where n1 is the number of possible ω1 jumps and n3 the number of possible ω3 jumps for the impurity. In the particular case when the impurity is in the middle of the two uranium vacancies in UO2, there are only 8 possible ω1 jumps, so that n1 = 8. If we insert this in eq 25 and then in the expression of f1 (eq 23): f1 ≈

(n1 − 1)ω1 + n3ω3 (n1 + 1)ω1 + n3ω3

Table 5. Activation Energies Yielded by vdW-DF+U (in eV) for the Elementary Mechanisms Involved in the VacancyAssisted Migration of Kr, as well as Resulting Binding Energy between Impurity and Assisting Vacancy for the Various Assisting Defects

(27)

And the diffusion coefficient can then be calculated as D≈

1 2 (n1 − 1)ω1 + n3ω3 −(E1+ Ef + EB + ER )/(kT ) a ν1 e 4 (n1 + 1)ω1 + n3ω3

(28)

Two special cases are worth noting: • For ω3 ≪ ω1, f1 = 7/9 7 2 −(E1+ Ef + EB + ER )/(kT ) D= a ν1 e 36

(29)

• For ω3 ≫ ω1, f1 = 1 D=

1 2 −(E1+ Ef + EB + ER )/(kT ) a ν1 e 4

1 2 d int p″fint ωint 6

Kr in V0U/V0U

Kr in V0UO/V0U

Kr in V0UO2/V0U

E0 E1 E2 E3 E4 EB

4.39 2.98 1.06 5.74 4.90 −0.84 4− Kr in V4− U /VU

4.39 4.59 N/A 4.86 2.98 −1.88 4− Kr in V2− UO/VU

4.39 5.16 N/A 5.20 3.04 −2.16 Kr in V0UO2/V4− U

E0 E1 E2 E3 E4 EB

3.91 2.39 1.57 3.22 2.68 −0.54

3.91 4.08 N/A 5.26 3.77 −1.49

3.91 4.30 N/A 4.62 3.38 −1.24

NEB method for the elementary mechanisms ωi involved in the diffusion models for the various traps and mobile vacancies: • Kr in neutral VU, VUO, and Schottky assisted by an additional neutral VU 2− • Kr in a charged V4− U , VUO, and neutral Schottky combined 4− with an additional VU . It is worth noting that a single NEB calculation enables one to determine ω3 and ω4, as well as the krypton−uranium vacancy binding energy (EB = E4 − E3), as shown in Figure 5 in 4− the case of Kr in a V4− U assisted by a second VU . It can be seen in Table 5 that the ω2 frequency only exists for the migration of Kr assisted by two uranium vacancies (neutral or charged). For the VUO and VUO2 (neutral or charged) combined to the mobile vacancy, the krypton atom forms a stable complex with the vacancies and sits in the middle of the two uranium vacancies, similarly to what happens for Xe,8,17 which means that there is no ω2 mechanism and justifies the use of the Catlow model. If we compare the various elementary mechanisms, we can see that the E0 and E1 migration energies are lower for charged defects while E2, E3, and E4 are lower for neutral ones.

(30)

Diffusion through Interstitial Mechanism. In the case of diffusion through interstitial mechanism, D can be expressed as

D=

(32)

(31)

where dint is the jump distance of the impurity from an interstitial position to the next one in the (110) direction and is therefore equal to a/√2. p″ is the probability of an empty interstitial being present in the vicinity of the impurity. If the concentration of interstitials is low enough, that is, the number of interstitial atoms is much less than the number of the available sites, p″ = 1.59 Given the very high formation energy of U interstitials,48 we are mostly concerned with O interstitials. This approximation might therefore not be valid for highly hyperstoichiometric conditions. f int is the correlation factor between 2 subsequent interstitial jumps and will be equal to 1 if again the interstitial concentration is low enough. In this case, each interstitial atom has a high probability of being surrounded by empty interstitial sites. All directions for jumps of an interstitial solute to neighboring empty sites will be equally probable G

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry

Table 7. Attempt Frequencies (in Hz) Associated with the Five-Frequency and Catlow Models for the Three Trap Sites Considered for the Kr Impurity (VU, VUO, and VUO2) Yielded by the Empirical Potential Calculations 4− Kr in V4− U /VU

ν0 ν1 ν2 ν3 ν4

Figure 5. Elementary ω3/ω4 pathways for Kr in a V0U assisted by a second V0U.

8.99 1.00 1.10 2.93 8.72

× × × × ×

12

10 1012 1013 1011 1011

4− Kr in V2− UO/VU

8.99 × 4.61 × N/A 2.09 × 1.72 ×

12

10 1011 1011 1012

Kr in V0UO2/V4− U 8.99 × 2.09 × N/A 1.97 × 2.42 ×

1012 1012 1012 1012

frequency for the Kr migration assisted by 2V4− U . The attempt frequency of this mechanism is significantly higher (1.10 × 1013 Hz) and is of the same order of magnitude as the frequency for the V4− U vacancy jump in the bulk. Rates of Elementary Mechanisms. The ωi rates of the elementary migration mechanisms involving charged defects calculated from the Ei and the νi presented in the sections above are shown in Figure 6.

We can also observe that the binding energy EB is always negative, which means that it is always favorable for a uranium vacancy to be close to the Kr atom in its trap site, which favors vacancy assisted migration. EB is larger for neutral vacancies than for charged ones. This can be due to the fact that the interaction of the lattice with Kr, which has a quite far-reaching electron cloud, is slightly more repulsive in the case of negative vacancies than with neutral ones. Finally, here are the relative orders obtained for the migration energies of the various elementary mechanisms: • For Kr in V0U assisted by V0U, E2 < E1 < E0 < E4 < E3 4− • For Kr in V4− U assisted by VU , E2 < E1 < E4 < E3 < E0 • For the four other mechanisms involving Kr trapped in neutral or charged VUO and VUO2, E4 < E0 < E1 < E3. Attempt Frequencies. To evaluate the accuracy of the empirical potential used to calculate the attempt frequencies, we first show in Table 6 the migration energies yielded by the Table 6. Activation Energies Yielded by the Empirical Potential (in eV) for the Three Trap Sites Considered for the Kr Impurity (VU, VUO, and VUO2) and Comparison with DFT+U Results 4− Kr in V4− U /VU

E0 E1 E2 E3 E4 EB

4− Kr in V2− UO/VU

Kr in V0UO2/V4− U

EEP

EDFT+U

EEP

EDFT+U

EEP

EDFT+U

6.91 3.90 0.13 4.12 6.92 2.80

3.91 2.39 1.57 3.22 2.68 −0.54

6.91 6.10 N/A 4.72 4.38 −0.34

3.91 4.08 N/A 5.26 3.77 −1.49

6.91 7.95 N/A 6.94 3.93 −3.01

3.91 4.30 N/A 4.62 3.38 −1.24

empirical potential calculations compared to the DFT+U results. It can be observed that the empirical energies are often higher than the DFT ones and that the biggest disagreement occurs for the highest charged defect (migration of V4− U and Kr migration assisted by 2V4− ). The binding energies are not well reproduced, U which comes as no surprise considering the rather simple nature of the interatomic potential. Similar discrepancies were found for other solid and gaseous fission products in UO2.8,62 The empirical potential, however, is able to capture the main qualitative physical trends and most of the relative orders between Ei, which justifies its use here. The attempt frequencies νi obtained from the empirical calculations, which are necessary to compute the ωi migration rates, are summarized in Table 7. It can be seen that they are between 2.09 × 1011 and 2.42 × 1012 Hz, except the ν2

2− Figure 6. Elementary rates for the migration of Kr in V4− U , VUO, and 4− neutral Schottky defects assisted by VU .

H

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry The relative orders between the rates of the various elementary mechanisms, which determine the equation to use to calculate the diffusion coefficient, are as follows: 4− • For Kr in V4− U assisted by VU , ω0 ≪ ω3 ≪ ω4 ≪ ω1 ≪ ω2, 2− 0 • For Kr in VUO and VUO2 assisted by V4− U , ω3 ≪ ω1 ≪ ω0 ≪ ω4. These are the same as those obtained from the activation energies alone, so we can also conclude on the relative order of the elementary rates and on the equations to use for the neutral defects, even without attempt frequencies. Finally, we see in Table 6 for neutral assisting vacancies and Figure 6 for charged ones that ω0 ≪ ω4 for all assisting vacancies, except for the migration of Kr assisted by two neutral uranium vacancies for which ω0 ≫ ω4. This means that the vacancy moves faster far from the impurity than from second to first neighbor of the impurity.

Table 8. For the Various Assisting Vacancies Studied, Synthesis of the Information and Parameters Needed to Calculate the Activation Energiesa and Effective Activation Energies Obtained for Various Stoichiometric Conditionsb model rate inequalities equation E1 EB Ef O-rich Ef stoich. Ef U-rich ER O-rich ER stoich. ER U-rich Ea O-rich Ea stoich. Ea U-rich

■

RESULTS ON KR VACANCY-ASSISTED DIFFUSION Effective Activation Energies. For the various assisting vacancies (traps and mobile vacancies) studied, we synthesize in Table 8 the data needed to calculate the effective activation energy: the diffusion model used, the inequalities between the rates of the elementary mechanisms, the equation to use to calculate D and Ea from the relative orders between the elementary rates, the formation energy of the neutral or charged VU, the binding energy between the impurity and the assisting vacancy, the reconfiguration energy, as well as the effective activation energies predicted by the diffusion models for various nonstoichiometry conditions. In addition, from the expressions of the formation energies of the assisting uranium vacancies (eq 5), Ea can be expressed as a function of the oxygen potential. Figure 7 shows the evolution of Ea as a function of μO. It can be seen from the results shown in Table 8 and Figure 7 that due to the high formation energy of V0U compared to V4− U on the one hand and to the high solution energy of Kr in neutral VU and VUO on the other hand, the mechanisms involving only neutral defects are not favorable. Kr migration will mainly be assisted by charged vacancies. Then, because of the decreasing formation energy of all considered vacancies except the Schottky defect when going from U-rich to O-rich conditions, the activation energies of all migration mechanisms decrease from hypo- to hyperstoichiometric conditions. Kr migration will be favored in hyperstoichiometric UO2 with low activation energies (0.73 eV for the mechanism assisted by two V4− U ). In contrast, in U-rich UO2, there are very few neutral or charged VU or VUO vacancies, and the activation energies of the vacancy-assisted migration of Kr are quite high, between 9.54 and 12.90 eV depending on the Kr trap. Finally, we observe in Figure 7 that the vacancy-assisted mechanism with lowest activation energy depends on UO2 stoichiometry. For U-rich conditions (−9.86 < μO < −8.99 eV) the mechanism with the lowest Ea is the migration assisted by the Schottky defect and the V4− U vacancy. For μO between −8.99 and −7.37 eV, that is, for slightly hypo- to slightly hyperstoichiometric UO2, the migration assisted by V2− UO and V4− U exhibits the lowest activation energy. Finally for O-rich conditions (−7.37 < μO < −6.06 eV), the lowest activation energy is obtained for the migration assisted by two V4− U . Diffusion Coefficient Prefactors. The D0 prefactors of the diffusion coefficients for the Kr migration assisted by the

model rate inequalities equation E1 EB Ef O-rich Ef stoich. Ef U-rich ER O-rich ER stoich. ER U-rich Ea O-rich Ea stoich. Ea U-rich

Kr in V0U/V0U

Kr in V0UO/V0U

Kr in V0UO2/V0U

5-frequency ω0 = ω4 ω3 ≪ ω1 ≪ω2 17 2.98 −0.84 3.36 7.16 10.96 4.28 5.60 8.85 9.78 14.91 21.95

Catlow ω0 ≪ ω4 ω1 ≪ ω3 30 4.59 −1.88 3.36 7.16 10.96 2.97 2.40 3.73 9.04 12.27 17.41

Catlow ω0 ≪ ω4 ω1 ≪ ω3 30 5.16 −2.16 3.36 7.16 10.96 3.03 0.56 0.00 9.39 10.72 13.96

4− Kr in V4− U /VU

4− Kr in V2− UO/VU

Kr in V0UO2/V4− U

5-frequency ω0 ≪ ω4 ω3 ≪ ω1 ≪ω2 17 2.39 −0.54 −1.12 2.68 6.48 0.00 1.33 4.57 0.73 5.86 12.90

Catlow ω0 < ω4 ω3 ≪ ω1 29 4.08 −1.49 −1.12 2.68 6.48 0.57 0.00 1.34 2.04 5.27 10.41

Catlow ω0 ≪ ω4 ω1 ≪ ω3 30 4.30 −1.24 −1.12 2.68 6.48 3.03 0.56 0.00 4.97 6.30 9.54

a

Diffusion model used; inequalities between the rates of the elementary mechanisms; equation used to calculate D and Ea; formation energy of the VU; binding energy between the impurity and the vacancy; reconfiguration energy. All energies are given in eV. bO-rich, μO = −6.06 eV; stoichiometric UO2, μO = −7.96 eV; U-rich, μO = −9.86 eV

Figure 7. Evolution of Ea as a function of μO for the vacancy assisted and interstitial mechanisms. For the charged defects, μe is taken equal to 1.2 eV (middle of the band gap).

various charged defects calculated using the equations indicated in Table 8 are listed in Table 9. For the sake of simplicity I

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry Table 9. Prefactors (in m2·s−1) of the Kr Vacancy-Assisted Migration for the Various Charged Assisting Defects 4− V4− U /VU

4− V2− UO/VU

V0UO2/V4− U

3.07 × 10−7

7.08 × 10−9

5.35 × 10−8

and given the small variation of a with nonstoichiometry, we always use the cell parameter of the perfect crystal. In this approximation, the D0 parameters are independent from μO.

■

INTERSTITIAL DIFFUSION Figure 8 shows the energy along the migration pathway in the ⟨111⟩ direction of Kr as an interstitial yielded by the NEB method.

Figure 8. Migration pathway of Kr as an interstitial in the ⟨111⟩ direction obtained using the vdW-DF functional.

The activation energy for the interstitial mechanism calculated according to eq 33 from the migration energy shown in Figure 8 and the solution energy shown in Figure 2 (6.05 eV) is given in Table 10 as a function of μO. Table 10. Activation Energy (in eV) for the Interstitial Mechanism As a Function of μO μO Ea (eV)

−9.86 8.01

−8.52 8.01

−7.96 8.57

−6.63 9.90

−6.06 11.04

It is high, in particular in hyperstoichiometric conditions. From Figure 7, this mechanism will probably only come into play in hypostoichiometric conditions when there are no or very few neutral or charged VU or VUO vacancies to assist Kr migration. Under these conditions, the approximation made on the low interstitial concentration is valid. The attempt frequency for the interstitial mechanism was calculated to be 2.50 × 1012 Hz, yielding a diffusion prefactor of 6.40 × 10−8 m2·s−1. We have not considered in this work the indirect mechanism assisted by O interstitials that was investigated by Liu et al.63 using empirical potentials and standard GGA because of the significant instability of Xe as interstitial or on the O sublattice in DFT+U, contrary to what is observed in GGA. Despite this instability, this mechanism could be significant in hypostoichiometric UO2, and it should be investigated in the future.

Figure 9. Diffusion coefficients for the interstitial and vacancy-assisted migration of Kr for 4 different oxygen chemical potentials corresponding to hypostoichiometric (2 values of μO, −9.86 and −8.99 eV), stoichiometric (μO = −7.96 eV), and hyperstoichiometric UO2 (μO = −6.63 eV).

■

KRYPTON DIFFUSION COEFFICIENT AS A FUNCTION OF STOICHIOMETRY The diffusion coefficients for the interstitial and charged vacancyassisted migration of Kr are shown in Figure 9 for 4 different oxygen chemical potentials corresponding to hypostoichiometric (2 different μO), stoichiometric, and hyperstoichiometric UO2.

We also present in Figure 10 the Kr diffusion coefficient for the most favorable mechanism (interstitial or vacancy assisted one) for various stoichiometric conditions. J

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry

with fission gas atoms also seem to follow similar trends. Direct comparison between the Xe and Kr diffusion rates would require application of the same nonstoichiometry models and DFT methodology, which is left as future work.

■

CONCLUSIONS We performed a study of the Kr diffusion in UO2 using a combination of diffusion models, DFT+U (VdW-DF functional) and empirical potentials. We determined the diffusion coefficient (activation energy and prefactor) for the migration of Kr either by interstitial mechanism or assisted by uranium vacancies as a function of the oxygen chemical potential μO, that is, for various nonstoichiometries. We found that the vacancy− impurity exchange exists only for the migration assisted by two monovacancies. For VUO and VUO2 traps associated with the mobile VU vacancy, the Kr is positioned in the middle of the two uranium vacancies, so there is no vacancy−impurity exchange. The Catlow diffusion model was used in this case. It was first observed that the binding energy EB between Kr in its traps and the assisting vacancy is negative for all considered traps, which means that it is always favorable for a vacancy to be close to the Kr impurity, which results in vacancyassisted migration. In addition, we show that the mechanisms involving neutral vacancies are not favorable because of the high formation energies of the uranium neutral vacancy compared to the charged one and of the high solution energy of Kr in neutral VU and VUO. Finally, we determined the most favorable mechanism for Kr migration and the corresponding diffusion coefficient as a function of nonstoichiometry or μO. For very U-rich conditions (−9.86 < μO < −9.06 eV), the most favorable mechanism is the interstitial migration with an activation energy equal to 8.01 eV. For μO between −9.06 and −8.99 eV, that is, for hypostoichiometric UO2, migration assisted by the Schottky defect and the V4− U vacancy is the most favorable, and Ea is between 7.80 and 8.01 eV. Around stoichiometry (−8.99 < μO < 4− −7.96 eV), migration assisted by V2− UO and VU dominates and 7.80 < Ea < 4.09 eV. For O-rich conditions (−7.96 < μO < −6.06 eV), migration assisted by two V4− U is the most favorable mechanism and in this case Ea will be between 4.09 to 0.73 eV depending on nonstoichiometry. Kr diffusion is thus enhanced at higher μO. The experimental values presented in Table 1 are in the latter interval. Since it is very probable that these values were obtained for at least slightly hyperstoichiometric samples given the easy oxidation of UO2, our Ea results are consistent with the experimental data and the variations observed between the activation energies measured are probably due to differences in the stoichiometry of the samples. Precise measurement of the Kr diffusion coefficient in UO2 as a function of nonstoichiometry is necessary to confirm these results. Finally, the mechanisms and trends with nonstoichiometry established for Kr are similar to previous studies of Xe, although direct comparison of activation energies and pre-exponentials is difficult due the application of different energy models for the nonstoichiometry and different exchange−correlation functionals.

Figure 10. Kr diffusion coefficient for the most favorable mechanism (interstitial or vacancy assisted one) for various oxygen potentials.

It can be seen that the most favorable migration mechanism for Kr also depends on the nonstoichiometry when considering the diffusion coefficient (and not only the activation energy) and all possible mechanisms (interstitial as well as vacancy assisted) and that Kr diffusion is enhanced when μO increases. • For very U-rich conditions (−9.86 < μO < −9.06 eV) the most favorable mechanism is the interstitial migration. The diffusion coefficient is then (in m2·s−1 with the activation energy in eV): D = 6.40 × 10−8 e−8.01/(kT )

(34)

• For μO between −9.06 and −8.99 eV, that is, for slightly hypostoichiometric UO2, the migration assisted by the Schottky defect and the V4− U vacancy is the most favorable and (35) D = 5.35 × 10−8 e−Ea(hypo)/(kT ) with 8.01 < Ea(hypo) < 7.80 eV. • Around stoichiometry (μO between −8.99 and −7.37 eV), 4− the migration assisted by V2− UO and VU is the most favorable mechanism and (36) D = 7.08 × 10−9 e−Ea(stoich)/(kT ) with 7.80 < Ea(stoich) < 4.09 eV. • Finally, for O-rich conditions (−7.37 < μO < −6.06 eV) (37) D = 3.07 × 10−7 e−Ea(hyper)/(kT ) with 4.09 < Ea(hyper) < 0.73 eV. The experimental values for Ea presented in Table 1 are in this latter interval. Since it is very probable that these values were obtained for at least slightly hyperstoichiometric samples because of the very easy oxidation of UO2, our Ea results are consistent with the experimental data, and the variations observed between the activation energies measured are probably due to differences in the stoichiometry of the samples. Further experiments with precisely controlled sample stoichiometry are needed to confirm these results. Finally, the mechanisms and trends with nonstoichiometry established for Kr are similar to those found in previous studies of Xe.8,17 Direct comparison of activation energies and preexponentials, however, is difficult due the application of slightly different energy models for the nonstoichiometry and the resulting uranium vacancy concentration, as well as the use of different exchange−correlation functionals. The binding energies and migration barriers for the different vacancy complexes

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Marjorie Bertolus: 0000-0001-6598-5312 K

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry Notes

(20) Kresse, G.; Furthmüller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 11169−11186. (21) Kresse, G.; Hafner, J. Ab initio molecular dynamics for liquid metals. Phys. Rev. B: Condens. Matter Mater. Phys. 1993, 47, 558−561. (22) Kresse, G.; Furthmüller, J. Efficiency of Ab initio total energy calculations for metals and semiconductors using a plane-wave basis set. Comput. Mater. Sci. 1996, 6, 15−50. (23) Liechtenstein, A. I.; Anisimov, V. I.; Zaanen, J. Densityfunctional theory and strong interactions: Orbital ordering in MottHubbard insulators. Phys. Rev. B: Condens. Matter Mater. Phys. 1995, 52, R5467−R5470. (24) Dudarev, S.; Manh, D.; Sutton, A. Effect of Mott-Hubbard correlations on the electronic structure and structural stability of uranium dioxide. Philos. Mag. B 1997, 75, 613−628. (25) Kotani, A.; Yamazaki, T. Systematic Analysis of Core Photoemission Spectra for Actinide Di-Oxides and Rare-Earth Sesqui-oxides. Prog. Theor. Phys. Suppl. 1992, 108, 117−131. (26) Krupa, J.; Gajek, Z. Crystal field levels of tetravalent actinide ions in actinide dioxides UO2, NpO2 and PuO2. Eur. J. Solid State Inorg. Chem. 1991, 28, 143−146. (27) Dorado, B.; Freyss, M.; Amadon, B.; Bertolus, M.; Jomard, G.; Garcia, P. Advances in first-principles modelling of point defects in UO2: f electron correlations and the issue of local energy minima. J. Phys.: Condens. Matter 2013, 25, 333201. (28) Dorado, B.; Amadon, B.; Freyss, M.; Bertolus, M. DFT+U calculations of the ground state and metastable states of uranium dioxide. Phys. Rev. B: Condens. Matter Mater. Phys. 2009, 79, 235125. (29) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (30) Dion, M.; Rydberg, H.; Schröder, E.; Langreth, D.; Lundqvist, B. Van der Waals Density Functional for General Geometries. Phys. Rev. Lett. 2004, 92, 246401. (31) Klimes, J.; Bowler, D. R.; Michaelides, A. Chemical accuracy for the van der Waals density functional. J. Phys.: Condens. Matter 2010, 22, 022201. (32) Klimes, J.; Bowler, D. R.; Michaelides, A. Van der Waals density functionals applied to solids. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 83, 195131. (33) Teobaldi, G.; Ohnishi, H.; Tanimura, K.; Shluger, A. L. The effect of van der Waals interactions on the properties of intrinsic defects in graphite. Carbon 2010, 48, 4145−4161. (34) Silvestrelli, P. L.; Ambrosetti, A.; Grubisiĉ, S.; Ancilotto, F. Adsorption of rare-gas atoms on Cu(111) and Pb(111) surfaces by van der Waals corrected density functional theory. Phys. Rev. B: Condens. Matter Mater. Phys. 2012, 85, 165405. (35) Zhang, G.-X.; Tkatchenko, A.; Paier, J.; Appel, H.; Scheffler, M. Van der Waals interactions in ionic and semiconductor solids. Phys. Rev. Lett. 2011, 107, 245501. (36) Andersson, D.; Baldinozzi, G.; Desgranges, L.; Conradson, D.; Conradson, S. Density Functional Theory Calculations of UO2 Oxidation: Evolution of UO2+x, U4O9−y, U3O7, and U3O8. Inorg. Chem. 2013, 52, 2769−2778. (37) Guéneau, C.; Baichi, M.; Labroche, D.; Chatillon, C.; Sundman, B. Thermodynamic assessment of the uranium-oxygen system. J. Nucl. Mater. 2002, 304, 161−175. (38) Guéneau, C.; Dupin, N.; Sundman, B.; Martial, C.; Dumas, J.-C.; Gosse, S.; Chatain, S.; Bruycker, F. D.; Manara, R. J.; Konings, D. Thermodynamic modelling of advanced oxide and carbide nuclear fuels: Description of the U-Pu-O-C systems. J. Nucl. Mater. 2011, 419, 145−167. (39) Abrikosov, I.; Ponomareva, A.; Steneteg, P.; Barannikova, S.; Alling, B. Recent progress in simulations of the paramagnetic state of magnetic materials. Curr. Opin. Solid State Mater. Sci. 2016, 20, 85− 106. (40) Ikushima, K.; Tsutsui, S.; Haga, Y.; Yasuoka, H.; Walstedt, R. E.; Masaki, N. M.; Nakamura, A.; Nasu, S.; O̅ nuki, Y. First-order phase transition in UO2: 235U and 17O NMR study. Phys. Rev. B: Condens. Matter Mater. Phys. 2001, 63, 104404.

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS Work at Los Alamos National Laboratory was sponsored by the U.S. Department of Energy, Office of Nuclear Energy, Nuclear Energy Advanced Modeling and Simulation (NEAMS) program. This work was partly performed using HPC resources from GENCI-CCRT (Grants x2014086922 and x2015086922). The authors express their gratitude to G. Carlot and C. Sabathier for fruitful discussions. This research contributes to the joint programme on nuclear materials (JPNM) of the European energy research alliance (EERA).

■

REFERENCES

(1) Hayns, M.; Wood, M. On the rate theory model for fission-gas behaviour in nuclear fuels. J. Nucl. Mater. 1976, 59, 293−302. (2) Matzke, H. J. Gas release mechanisms in UO2 - a critical-review. Radiat. Eff. 1980, 53, 219−242. (3) Matzke, H. J. Atomic transport properties in UO2 and mixed oxides (U, Pu)O2. J. Chem. Soc., Faraday Trans. 2 1987, 83, 1121− 1142. (4) Matzke, H. J. In Diffusion Processes in Nuclear Materials; Agarwala, R. P., Ed.; North-Holland: Amsterdam, 1992. (5) Miekeley, W.; Felix, F. W. Effect of stoichiometry on diffusion of xenon in UO2. J. Nucl. Mater. 1972, 42, 297−306. (6) Lindner, R.; Matzke, H. J. Diffusion von Xe-133 in Uranoxyd verschiedenen Sauerstoffgehaltes. Z. Naturforsch., A: Phys. Sci. 1959, 14, 582−584. (7) Liu, X.-Y.; Andersson, D. A.; Uberuaga, B. P. First-principles DFT modeling of nuclear fuel materials. J. Mater. Sci. 2012, 47, 7367−7384. (8) Andersson, D.; Garcia, P.; Liu, X.-Y.; Pastore, G.; Tonks, M.; Millett, P.; Dorado, B.; Gaston, D.; Andrs, D.; Williamson, R.; et al. Atomistic modeling of intrinsic and radiation-enhanced fission gas (Xe) diffusion in UO2±x: Implications for nuclear fuel performance modeling. J. Nucl. Mater. 2014, 451, 225−242. (9) Auskern, A. The diffusion of krypton-85 from uranium dioxide powder. US Report WAPDTM-185, 1960. (10) Michel, A. Etude du comportement des gaz de fission dans le dioxyde d’uranium: mécanismes de diffusion, nucléation et grossissement de bulles. Ph.D. thesis, Université de Caen, 2011. (11) Petit, T.; Jomard, G.; Lemaignan, C.; Bigot, B.; Pasturel, A. Location of krypton atoms in uranium dioxide. J. Nucl. Mater. 1999, 275, 119−123. (12) Crocombette, J.-P. Ab initio energetics of some fission products (Kr, I, Cs, Sr and He) in uranium dioxide. J. Nucl. Mater. 2002, 305, 29−36. (13) Thompson, A.; Wolverton, C. First-principles study of noble gas impurities and defects in UO2. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 84, 134111. (14) Tian, X.; Gao, T.; Jiang, G.; He, D.; Xiao, H. The incorporation and solution of krypton in uranium dioxide: Density functional theory calculations. Comput. Mater. Sci. 2012, 54, 188−194. (15) Catlow, C. Fission gas diffusion in uranium dioxide. Proc. R. Soc. London, Ser. A 1978, 364, 473−497. (16) Skorek, R.; Maillard, S.; Michel, A.; Carlot, G.; Gilabert, E.; Jourdan, T. Modelling Fission Gas Bubble Distribution in UO2. Defect Diffus. Forum 2012, 323−325, 209−214. (17) Andersson, D.; Uberuaga, B.; Nerikar, P.; Unal, C.; Stanek, C. U and Xe transport in UO2±x: density functional theory calculations. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 84, 054105. (18) Kresse, G.; Joubert, D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 1758−1775. (19) Blöchl, P. E. Projector augmented-wave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 50, 17953−17979. L

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX

Article

Inorganic Chemistry (41) Laskowski, R.; Madsen, G. K. H.; Blaha, P.; Schwarz, K. Magnetic structure and electric-field gradients of uranium dioxide: An ab initio study. Phys. Rev. B: Condens. Matter Mater. Phys. 2004, 69, 140408. (42) Dorado, B.; Garcia, P. First-principles DFT+U modeling of actinide-based alloys: Application to paramagnetic phases of UO2 and (U,Pu) mixed oxides. Phys. Rev. B: Condens. Matter Mater. Phys. 2013, 87, 195139. (43) Gryaznov, D.; Heifets, E.; Sedmidubsky, D. Density functional theory calculations on magnetic properties of actinide compounds. Phys. Chem. Chem. Phys. 2010, 12, 12273−12278. (44) Boettger, J. Predicted spin-orbit coupling effect on the magnetic ordering of crystalline uranium dioxide. Eur. Phys. J. B 2003, 36, 15− 20. (45) Lei, S. Theoretical study using electronic structure calculations of uranium and cerium dioxides containing defects and impurities. Ph.D. thesis, Aix-Marseille Université, 2016. (46) Henkelman, G.; Jónsson, H. Improved tangent estimate in the nudged elastic band method for finding minimum energy paths and saddle points. J. Chem. Phys. 2000, 113, 9978−9985. (47) Monkhorst, H.; Pack, J. Special points for Brillouin-zone integrations. Phys. Rev. B 1976, 13, 5188−5192. (48) Vathonne, E.; Wiktor, J.; Freyss, M.; Jomard, G.; Bertolus, M. DFT+U investigation of charged point defects and clusters in UO2. J. Phys.: Condens. Matter 2014, 26, 325501. (49) Zhang, S. B.; Northrup, J. E. Chemical potential dependence of defect formation energies in GaAs: Application to Ga self-diffusion. Phys. Rev. Lett. 1991, 67, 2339−2342. (50) Vineyard, G. H. Frequency factors and isotope effects in solid state rate processes. J. Phys. Chem. Solids 1957, 3, 121−127. (51) Buckingham, R. A. The Classical Equation of State of Gaseous Helium, Neon and Argon. Proc. R. Soc. London, Ser. A 1938, 168, 264− 283. (52) Henkelman, G.; Uberuaga, B. P.; Jonsson, H. A climbing image nudged elastic band method for finding saddle points and minimum energy paths. J. Chem. Phys. 2000, 113, 9901−9904. (53) Voter, A. F. CLSMAN program; Los Alamos National Laboratory: Los Alamos, NM, 1991. (54) Grimes, R. W. J. Am. Ceram. Soc. 1994, 77, 378. (55) Stanek, C. R.; Bradford, M. R.; Grimes, R. W. Segregation of Ba2+, Sr2+, Ce4+ and Zr4+ to UO2 surfaces. J. Phys.: Condens. Matter 2004, 16, S2699−S2714. (56) Grimes, R.; Catlow, C. The stability of fission products in uranium dioxide. Philos. Trans. R. Soc., A 1991, 335, 609−634. (57) Lidiard, A. B. Impurity Diffusion in Crystals (Mainly Ionic Crystals with the Sodium Chloride Structure). Philos. Mag. 1955, 46, 1218−1235. (58) Lidiard, A. B. The influence of solutes on self-diffusion in metals. Philos. Mag. 1960, 5, 1171−1180. (59) Mehrer, H. Diffusion in solids: fundamentals, methods, materials, diffusion-controlled processes; Springer series in solid-state sciences; Springer: Berlin, 2007; Vol. 155. (60) Manning, J. Correlation Factors for Impurity Diffusion. bcc, Diamond, and fcc Structures. Phys. Rev. 1964, 136, A1758−A1766. (61) Martin, P.; Vathonne, E.; Carlot, G.; Delorme, R.; Sabathier, C.; Freyss, M.; Garcia, P.; Bertolus, M.; Glatzel, P.; Proux, O. Behavior of fission gases in nuclear fuel: XAS characterization of Kr in UO2. J. Nucl. Mater. 2015, 466, 379−392. (62) Perriot, R.; Liu, X.-Y.; Stanek, C.; Andersson, D. Diffusion of Zr, Ru, Ce, Y, La, Sr and Ba fission products in UO2. J. Nucl. Mater. 2015, 459, 90−96. (63) Liu, X.-Y.; Uberuaga, B. P.; Andersson, D. A.; Stanek, C. R.; Sickafus, K. E. Mechanism for transient migration of xenon in UO2. Appl. Phys. Lett. 2011, 98, 151902.

M

DOI: 10.1021/acs.inorgchem.6b01560 Inorg. Chem. XXXX, XXX, XXX−XXX