1 DEFECT ENERGETICS IN URANIUM DIOXIDE BY ELECTRONIC STRUCTURE AND ATOMISTIC LEVEL SIMULATIONS By PANKAJ VILAS NERIKAR A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2009
2 2009 Pankaj Vilas Nerikar
3 To my family
4 ACKNOWLEDGMENTS There are a lot of people whom I wish to acknowledge in thi s work. First and foremost, my supervisory committee chair, Prof. Susan B. Sinnott. She has always been very encouraging, supportive and really cares about her students. She also has a very enthusiastic approach towards research and science. I respect her as a scientist and an individual. I am also grateful to Prof. Phillpot, Prof. Baney, Prof. Nino and Prof. Tulenko for their guidance and willingness to help. I would like to thank Dr. Blas P. Uberuaga, with whom I had the fortune of working as a student intern at the Los Alamos National Laboratory. I have learned a lot from his deep insight into the subject and he has helped me carry out some of the calculations in this work. The people at MST 8, especially Dr. Chris Stanek have made my time productive and enjoyable at Los Alamos and have impressed me as friends and scientists. On the personal side, I will be ever grateful to my family and I dedicate this work to them. They have always supported and encouraged me at every stage of life and have given me t his opportunity. Finally, I would like to thank the past and present members of the Computational Materials Science Group and other friends at the Department of Materials Science and Engineering. They have been very supportive and helpful and this has had a positive impact on my life both inside and outside the school. Partial funding from DOE NERI, DOE CMSN and DOE -AFCI is gratefully acknowledged.
5 TABLE OF CONTENTS page ACKNOWLEDGMENTS .................................................................................................................... 4 LIST OF TABLES ................................................................................................................................ 7 LIST OF FIGURES .............................................................................................................................. 8 ABSTRACT ........................................................................................................................................ 12 CHAPTER 1 GENERAL INTRODUCTION .................................................................................................. 14 Uranium Oxide ............................................................................................................................ 14 Point Defects ................................................................................................................................ 15 Fission Products .......................................................................................................................... 20 Grain Boundaries ......................................................................................................................... 22 Objective and Scope of this Dissertation ................................................................................... 27 2 GENERAL INTRODUCTION TO SIMULATION METHODS ........................................... 36 Density Functional Theory ......................................................................................................... 36 Local Density and Ge neralized Gradient Approximations ............................................... 38 Periodic Supercells .............................................................................................................. 40 Electron Ion Interactions ..................................................................................................... 42 Strong Correlation ............................................................................................................... 43 Accuracy of Calculations ............................................................................................................ 46 Lattice Statics .............................................................................................................................. 46 Defective Lattice .................................................................................................................. 48 Computational Software ............................................................................................................. 50 3 ENERGETICS OF INTRINSIC POINT DEFECTS IN URANIUM DIOXIDE ................... 54 Introduction ................................................................................................................................. 54 Review of Literature ................................................................................................................... 55 Computational Methodolog y ...................................................................................................... 57 Electronic Structure Calculations ....................................................................................... 57 Defect Formation Energy .................................................................................................... 58 Thermodynamic Component ............................................................................................... 59 Results and Discussion ............................................................................................................... 60 Equilibrium Properties of Bulk UO2 and Uranium Metal ................................................. 60 Formation Energy of Neutral Defects ................................................................................ 62 Effect of Charged States ...................................................................................................... 66 Effect of Temperatu re and Oxygen Partial Pressure ......................................................... 69
6 4 THERMODYNAMICS OF FISSION PRODUCTS IN UO2 X ............................................... 87 Introduction ................................................................................................................................. 87 Computational Methodology ...................................................................................................... 90 Electronic Structure Calculations ....................................................................................... 90 Incorporation Energy ........................................................................................................... 90 Solution Energy ................................................................................................................... 91 Oxide Solution Energy ........................................................................................................ 93 Results and Discussion ............................................................................................................... 94 Validation of Approach ....................................................................................................... 94 Incorporation Energy ........................................................................................................... 97 Solution Energy ................................................................................................................... 97 Oxide Solution Energy ........................................................................................................ 99 Comparison with other Theoretical Approaches ..................................................................... 101 Non -Spin Polarized GGA .................................................................................................. 101 Empirical Potentials ........................................................................................................... 103 5 SEGREGATION OF FISSION PRODUCTS TO GRAIN BOUNDARIES ........................ 116 Introduction ............................................................................................................................... 116 Review of Literature ................................................................................................................. 117 Computational Methodology .................................................................................................... 121 Results ........................................................................................................................................ 123 Structural Analysis ............................................................................................................. 123 Xe Segregation to the 5 Tilt Grain Boundary ............................................................... 124 Xe Segregation to the 5 Twist Grain Boundary ........................................................... 127 Xe Segregation to the Amorphous Grain Boundary ........................................................ 128 Discussion .................................................................................................................................. 131 Comparison of Potentials .................................................................................................. 131 Comparison of Different Boun daries ............................................................................... 133 6 CONCLUSIONS ....................................................................................................................... 160 LIST OF REFERENCES ................................................................................................................. 164 BIOGRAPHICAL SKETC H ........................................................................................................... 170
7 LIST OF TABLES Table page 3 1 The calculated properties of bulk uranium oxide are shown using different approximations. These values are compare d to published experimental and theoretical data. ...................................................................................................................... 72 3 2 Relative stabilities of the three different phases of bulk uranium metal, the orthorhombic -phase, the body centered tetragonal phase, and the body centered cubic -phase. ......................................................................................................................... 73 3 3 Calculated formation energies of neutral point defects compared to reported work using various approxi mations. ............................................................................................... 74 3 4 The values for Frenkel energies are given for different compounds having the fluorite structure. .................................................................................................................... 75 4 1 Comparison of equilibrium and defect properties of UO2 from GGA, SP GGA+U and experiment. .................................................................................................................... 1 04 4 2 Incorporation energies of xenon, cesium and strontium calculated using SP -GGA+U and GGA and compared to pr evious empirical potential results. For each of the fission product, the lowest incorporation energy is shown in bold. .................................. 105 4 3 Solution energies of xenon, cesium and strontium calculated using S P GGA+U and GGA and compared to previous empirical potential results. For each of the fission product and stoichiometry, the lowest solution energy is shown in bold. ........................ 106 4 4 Solution energies of oxides of cesium and strontium calculated using SP GGA+U and GGA and compared to previous empirical potential results. A positive energy implies insolubility in the fuel matrix phase. ..................................................................... 107 4 5 Binding energies of strontium -strontium interaction calculated using GGA. .................. 108 5 1 Average segregation energy for the 5 tilt boundary as a function of the supercell size. ........................................................................................................................................ 135 5 2 Average segregation energy and the grain boundary energy for the three different types of boundaries considered. .......................................................................................... 136
8 LIST OF FIGURES Figure page 1 1 The fluorite structure of UO2; uranium atoms are illustrated in blue and oxygen atoms in red. ........................................................................................................................... 27 1 2 The U O binary phase diagram, subscripts S and L indicate the solid and liquid phases. The numbers correspond to the partial pressure of oxygen present in that particular region of the phase diagram. ................................................................................. 28 1 3 Point defects and defect complexes ..................................................................................... 29 1 4 The free energy change as a function of the number of defects created. The top part represents the chan ge in enthalpy to create a defect and the bottom part is the change in the configurational entropy. .............................................................................................. 30 1 5 The yield of fission products for a 3.2% enriched pressurized water reactor pin rated after 2.9% burn up. ................................................................................................................ 31 1 6 A schematic representation of a low angle tilt grain boundary. It is made of a series of noninteracting edge dislocations with burgers vector b. .............................................. 32 1 7 Coincident site lattice model for a 5 (100) grain boundary in a cubic crystal. The left side of the figure shows a tilt boundary and the right side illustrates a twist boundary. ................................................................................................................................ 33 1 8 Relative energy for a Cu (110) tilt grain boundary as function of the misorientation. .... 34 1 9 Segregation of a free atom to different surfaces with changes in the Gibbs free energy. ..................................................................................................................................... 35 2 1 Illustration of the DFT approach. .......................................................................................... 51 2 2 A schematic of the supercell method for a point defect in a solid, the supercell being the area enclosed by the dashed lines. .................................................................................. 52 2 3 Illustration of the Mott Littleton method. ............................................................................. 53 3 1 Summary of the oxygen diffusion and metal diffusion data in fluorite type oxides, showing the temperature scale normalized to the melting point of the solid. .................... 76 3 2 Illustration of the 96 atoms supercell used in the present work. ......................................... 77 3 3 Total DOS projected to the U 5f and 6d orbitals as the conduction band and to the U 5f and O 2p as the valence band for antiferromagnetic UO2 ............................................... 78
9 3 4 Normalized formati on energy of neutral defect complexes calculated as function of system size by atomic level simulations using an empirical potential within GULP. ....... 79 3 5 Stability transitions for the neutral, + 1 and +2 charged oxygen vacancies as a function of Fermi level over the entire band gap. Only the most dominant defect with reference to a particular Fermi level is shown for clarity. ................................................... 80 3 6 Ca lculated defect formation energies for variously charged oxygen and uranium point defects as function of reference state. The empty symbols denote uranium while the filled symbols are oxygen. .................................................................................... 81 3 7 Calculated formation energy of the anti Frenkel defect as a function of the Fermi level; the formation energies of the neutral complexes are compared to charged individual defects. .................................................................................................................. 82 3 8 Calculated formation energy of the Schottky defect as a function of the Fermi level; the formation energies of the neutral complexes are compared to charged individual defects. .................................................................................................................................... 83 3 9 Predict ed effect of temperature on the formation of neutral point defects at standard atmospheric pressure for the temperature range of 300 to 1400 K. The vertical dashed lines indicate a phase change of bulk uranium metal. ............................................. 84 3 10 Predicted effect of oxygen partial pressure on the formation of neutral point defects in UO2 at 800 K; the partial pressure of oxygen is varied from pO2 = 1 x 1010 atm. to 1 atm. This corresponds to highly reducing and oxidizing conditions. .............................. 85 3 11 Illustration of the Brouwer diagram for UO2+x, the diagram assumes that the intrinsic electronic reaction dominates and neglects defects on the uranium sublattice. ................ 86 4 1 Lattice parameter of UxM1 xO2 pseudo-binary solid solutions. ......................................... 109 4 2 Illustrations of the various trap sites for fission produ cts considered in this work. Uranium atoms are shown in blue and oxygen atoms in red. ............................................ 110 4 3 Normalized incorporation energy of cesium in pre -existing trap sites calculated as function of syste m size by atomic level simulations using an empirical potential within GULP. ....................................................................................................................... 111 4 4 Formation energy of a uranium vacancy calculated using SP -GGA+U method as a function of temperature and stoic hiometry using the point defect model. The temperatures taken into account are 300 K, 1000 K and 2000 K. ..................................... 112 4 5 Solution energies of Xe, Cs and Sr in UO2 relative to the lowest solution energy site for the hypostoichiometric (UO2 x) case. ............................................................................ 113 4 6 Solution energies of Xe, Cs and Sr in UO2 relative to the lowest solution energy site for the stoichiometric (UO) case. ........................................................................................ 114
10 4 7 Solution energies of Xe, Cs and Sr in UO2 relative to the lowest solution energy site for the hyperstoichiometric (UO2+x) case. .......................................................................... 115 5 1 A schema tic of the grain boundary energy versus tilt angle for low energy orientations. .......................................................................................................................... 137 5 2 The three types of grain boundaries used in this work a) 5 (310)/ (001) symmetric tilt boundary, b) 5 twist boundary and c) amorphous grain boundary. .......................... 138 5 3 Radial distribution functions for the uranium uranium interaction .................................. 139 5 4 Radial distribution functions for the uranium -oxygen interaction .................................... 140 5 5 Segregation of Xe ................................................................................................................. 141 5 6 Effect of oxygen vacancies on segregation ........................................................................ 142 5 7 Effect of reconstruction. ...................................................................................................... 143 5 8 Illustration of the optimized struc tures after performing rigid body translations ............ 144 5 9 Grain boundary energy as a function of UO2 formula units in the supercell. .................. 145 5 10 Relative segregation energy as a function of distance perpendicular to the grain boundary plane using the Basak potential for different number of Xe atoms as illustrated in the caption of each figure. ............................................................................. 146 5 11 Average first nearest neighbors for the uranium atoms as a function of the relative energy of each uranium atom for the 5 tilt grain boundary. .......................................... 148 5 12 Relative energy o f xenon as a function of distance for the 5 twist grain boundary, relative energies are illustrated in blue while average energies are shown in red. .......... 149 5 13 Average first nearest neighbors for the uranium atoms as a function of the relative energy of each atom for the 5 twist grain boundary. ...................................................... 150 5 14 Plot of segregation energy as a function of distance for the amorphous grain boundary calculated using the Basak potential. ................................................................. 151 5 15 Average first nearest neighbors for the uranium atoms as a function of the relative energy of each atom for the amorphous grain boundary. .................................................. 152 5 16 Relative segregation energy of Xe atoms as a function of the relative energy of each atom for the amorphous grain boundary. ............................................................................ 153 5 17 Virial stress on each atom as a function of the distance along the direction perpendicular to the grain boundary plane for the amorphous grain boundary. .............. 154
1 1 5 18 Virial stress on each atom as a function relative energy of each uranium atom for the amorphous grain boundary. ................................................................................................. 155 5 19 Virial stress on each atom as a function of the relative Xe segregation energy for the amorphous grain boundary. ................................................................................................. 156 5 20 Relative segregation energy as a function of distance perpendicular to the grain boundary plane using the Busker potential. ........................................................................ 157 5 21 Schematic of the variation in segregation energy with change in the position of solute according to the local arrangement of atoms around the solute. ....................................... 158 5 22 Effect of oxygen vacancies on segregation: plot of normalized energy versus position. ................................................................................................................................. 159
12 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy DEFECT ENERGETICS IN URANIUM DIOXIDE BY ELECTRONIC STRUCTURE AND ATOMISTIC LEVEL SIMULATIONS By Pankaj Vilas Nerikar May 2009 Chair: Susan B. Sin nott Major: Materials Science and Engineering Defect energetics in uranium dioxide is investigated using electronic structure and atomistic level simulations. The stability range of intrinsic point defects in uranium dioxide is determined as a function of temperature, oxygen partial pressure, and non-stoichiometry. In particular, the density functional theory (DFT) calculations are performed at the level of the spin polarized, generalized gradient approximation and includes the Hubbard U term; as a result they predict the correct anti -ferromagnetic insulating ground state of uranium oxide. The predicted equilibrium properties and defect formation energies for neutral defect complexes match trends in the experimental literature quite well. Next, the stabi lities of selected fission products Xe, Cs, and Sr are investigated as a function of non-stoichiometry x in UO2 using DFT and empirical potentials. In general, higher charge defects are more soluble in the fuel matrix and the solubility of fission prod ucts increases as the hyperstoichiometry increases. The solubility of fission product oxides is also explored. These observations mirror experimentally observed phenomena. Finally, the segregation of Xe to different grain boundaries is studied using empir ical potentials. The segregation behavior is observed to be similar for the different boundaries
13 qualitatively but the ease with which segregation occurs is found to depend on the type of boundary.
14 CHAPTER 1 GENERAL INTRODUCTION The outlook towards nuclea r energy has generally brightened over the years with progressive improvement in the operating performance of existing reactors. This ensures the economic competitiveness of electricity produced through nuclear methods and the trend towards greater use of nuclear energy is projected to grow in the future. At the end of 200 4 approximately 441 nuclear power plants, with a total installed capacity of 358 gigawatts of electricity GW (e), were in operation worldwide, generating rough ly 17% of global electricity .1 Uranium Oxide Uranium dioxide (UO2) is an oxide of uranium that naturally occurs in the mineral uraninite.2 It is also known as urania, and is a black, radioactive, crystalline powder. Uranium dioxide (UO2) is the standard fuel utilized in light -water nuclear reactors. Some of the reasons it is used are: (i) it exists as a single phase at all temperatures up to its high melting point of 3125 K, (ii) it can accommodate substantial quantities of fission products and point defects without significant changes in the lattice, (iii) it also maintains its useful pr operties in a variety of environments and is also chemically stable. Uranium oxide has a variety of other applications. UO2 is mainly as a fuel for nuclear power plants. At the same time, depleted uranium compounds (DUO2) are well known to be good candidat es as potential catalysts for various chemical processes because of the ability of their delocalized f and d -electrons to interact with numerous organic and inorganic compounds.3, 4 UO2 manifests itself in the fluorite structure (space group Fm 3m ) with a l attice parameter of 5.47 see Figure 1 1 for an illustration .5 In the fluorite structure, the anions (oxygen) fill all the tetrahedral interstices (Wyckoff position 8c) of the close -packed cation sub lattice (uranium, Wyckoff position 4a). The resulting compound is MX2. The oxides of large quadrivalent cations
15 (Pu, Zr, Hf, Th, Ce) and the fluorides of large divale nt cations (Ca, Sr, Ba) both crystallize in this structure. Another way to view th e structure is to focus on the oxygen atoms, which are in a simple cubic arrangement with the uranium atom s l ying in between the cubes The ionic radii of uranium (U4+) and oxygen (O2-) are 1.00 and 1.38 respectively.6 Based of radius radio (rcation/ranion) arguments, the bonding in fluorite is predicted to be cubic in nature and hence this structure obeys Paulings rules of chemical bonding.7 The uranium -oxygen system exhibits other crystalline phases that are preferred at temperatures and partial pressure s of oxygen other than ambient condit ions These phases can have important impacts on the performance of the UO2 within a nuclear reactor because the fuels are subjected to a variety of temperatures (8001600 K) and partial pressures and hence t he knowledge of these phases becomes important. The calculated phase diagram is shown in Figure 1 2.8 As can be observed, there is a wide non-stoichiometry range owing to the creation of point d efects in the lattice. The fluorite phase is stable up to an oxygen to uranium ratio of 2.25 after which the U4O9 phase forms.810 Further increase in the oxygen partial pressure leads to the formation of U3O8 and UO3 phases. These are not important in fuel chemistry but are encountered during reprocessing of the fuel. The hypostoichiometric UO2x phase exists only at elevated temperatures. Point Defects A point defect is defined as any point in the lattice that is not occupied (va cancy) or occupied (interstitial) by the proper ion or atom needed to preserve the long-range periodicity of the crystal. It is so called, as it is not extended in space in any dimension. There is no strict limit as to how small a "point" defect should be. However, typically the term is used to define a few extra or missing atoms without an ordered structure of the defective positions. Other defects such as dislocations (one dimensional) and grain boundaries (two dimensional) may be viewed as
16 extensions of point defects. Point defects can be classified into different categories see Figure 1 3 for an illustration : (i) Vacancy: site where an atom or ion is missing. This can occur on both cation and anion sublattices. (ii) Interstitial: An atom occupying a site that is generally empty. This depends on the availability of empty sites in the lattice. (iii) Anti -site: An atom occupying a site on a second atoms sub-lattice and vice a versa. This can happen in covalent ceramics where the atoms are not charged and in ceramics which have one or more cations (iv ) El ectronic defects: extra electrons which ar e present in the system The location of these electrons depends on the type on the band gap of the system. After this basic introduction to the types of point de fects, it is important to understand the thermodynamics of point defect formation. Consider a perfect crystal consisting of N atoms on N lattice sites. T he Gibbs free energy Gperfect, of this crystal is given as: Gperfect Hperfect T Sperfect (1 1 ) where Hperfect is the change in enthalpy with reference to the initial state Sperfect is the change in entropy of the system and T is the temperature. The entropy contribution Sperfect is a combination of the arrangement of atoms (configurational) and the thermal vibrations of these atoms (vibrational). Following the methodology given by Barsoum,11 t he configurational part, in the absence of defects, does not contribute, as there is only one way in which N atoms can be arranged on N lattice sites. The vibrational part is given by: Sv Nk ln kT h 1 (1 2)
17 where k is the Boltzmanns constant, h is the Planks constant, is the characteristic vibrational frequency of the bond. Now consider a crystal in which a vacancy has been created. At equilibrium, there will be nV concentration of vacancies and if hd is the cost of creating one vacancy, the change in enthalpy for a defective crystal is given by: Hdefect Hperfect hdnV (1 3) The entropy Sdefect will now consist of a configurational part and a vibrational part. First considering the configurational part, it can be shown that the number of ways of arrang ing n vacant sites and N atoms on n+N lattice sites is given by: n N n N (1 4) Using this equation and applying Stirlings approximation, the configurational part is given by: Sconfiguation k N ln N N nv nvln nvnv N (1 5) The above equation assume s no interaction between vacancies and can be seen as an ideal solution model of two component mixing. When considering the vibrational part, it is important to remember that the atoms near the vacancy will vibrate at a different frequency ( ) than the atoms far away ( ). Further assuming that this frequency change is limited to first nearest neighbors, the total number of atoms affected is nV where is the coordination number of the vacancy. The vibrational part is given as: Svibration k ( N nv) ln kT h 1 knvln kT h' 1 (1 6)
18 Combining E quations 1 3 and 1 6 and subtracting them from Equation 11 yields the change in Gibbs energy upon the introduction of n vacancies: G = Gdefect Gperfect nvhd kTnvln' kT N ln N N nv nvln nvnv N (1 7) If G is plotted against nV at constant T, w e find that the equation goes through a minimum as seen in Fig ure 1 4.11 This is the reason why the formation of vacancies is thermodynamically stable. We have to take the minimization of the Gibbs free energy into account rather than just the minimization of the enthalpy. The change in entropy becomes impor tant when there are many defects in the system. In this work, only one defect is taken into account each time and hence the Gibbs free energy can be accounted for using the ideal solution model. Future work would include multiple defects and the above form alism would be adopted. However, the situation in ionic ceramics is more complicated. The defect formation, whether interstitial or vacancy, has to proceed in such a manner as to preserve the electro neutrality of the system. Such formation reactions can b e broadly classified into stoichiometric and non -stoichiometric reactions. Stoichiometric reactions are those in which the composition of the crystal does not change. There are two major types of stoichiometric reactions in ceramics: 1) Frenkel Disorder: T his disorder occurs when an atom gets displaced from its regular lattice position, thus creating a vacancy and interstitial in the process see Figure 1 -3 for an illustration This can occur on either sub lattice. In Kr ger -Vink notation, for the cation su b -lattice: MM x VM Mi (1 8) For the anion su b -lattice, it is called the anion Frenkel defect. i O x OO V O (1 9)
19 The anti -Frenkel defect is characteristic of the fluorite structure with UO2, ThO2, CaF2 all demonstrating this disorder.12 The reason is that a large amount of energy is required to displace a heavy cation atom such as uranium. 2) Schottky defect: This disorder consists of charge compensating vacancies. The thermodynamics of this disorder is slightly complicated as both sublattices are involved. For a binary oxide MO2, null VM ' 2 VO (1 10) The equilibrium concentrations of these sto ichiometric defects can be calculated in a similar manner as illustrated for a single vacancy. The other broad classification of defect reactions is the non -stoichiometric reaction where the composition of the crystal changes. These reactions depend more on the partial pressure of oxygen and the band gap rather than thermal activation.13 For example, at high oxygen partial pressures, extra oxygen is incorporated into the crystal interstitially through the reaction: 1 2 O2( g ) Oi x (1 11) Since the ext ra oxygen cannot exist as a neutral species, it has to obtain two electrons. This can be achieved by creating holes in the valence band. However, a more probable reaction for transition elements and heavy metals is to donate two electrons to the oxygen int erstitial since the energy cost associated with changing the valence state is not too large relative to the creation of holes in the valence band. Point defects are especially important for UO2 since they are created during fuel operation, a significant nu mber of which persist in the fuel although many recombine or otherwise annihilate each other. Over time, these defects are detrimental to the structure and performance of the fuel, since they form defect clusters that cause fuel swelling and ultimately cha nge the
20 crystal structure of UO2; this phase change is accompanied by a significant change in volume, which has highly negative effects on the microstructure and mechanical stability of the fuel pellet. These intrinsic point defects also affect the transpo rt properties of the fission gases formed during operation, allowing them to escape from the fuel matrix and damage the fuel cladding. 14 The point defec ts produced include isolated vacancies and interstitials and defect complexes. These individual point defects can form these complexes in a variety of possible charged states that can influence the relative stability of a cluster. The number of possible combinations of point defects and charge states is quite large, especially since UO2 is known to be non-stoichiometric over a wide range of compositions.15 Any hyperstoichiometry (UO2+x) or hypostoichiomtery (UO2 x) will stabilize some point defects in preference to others. Thus, an improved knowledge of the stability of these point defects is required for improved understanding of fuel performance. Fission Products Fission products (FP) can simply be defined as the atomic species obtained after a fission reaction. A typical fission reaction is as follows: n Xe Sr U U n1 0 136 54 96 38 236 92 235 92 1 03 (1 12) Here a heavy atom like uranium is hit with a neutron; the atom becomes unstable after absorbing the neutron and decomposes The end products of this reaction are lighter elements, extra neutrons to carry forwar d the reaction and about 200 MeV of energy. The yield or the types of fission products produced in % per parent atom depend on the reactor and the type of atom subjected to fission. The yield for a U 235 fission atom in a pressurized water reactor (PWR) i s shown in Figure 1 5.14 While fission products include elements from zinc through the lanthanides, the majority of the fission products occur in two
21 peaks. One peak occurs at about strontium to ruthenium while the other peak occurs at about tellurium to n eodymium. FP are produced during neutron irradiation of the nuclear fuel.16 Although the amount of a particular fission product produced depends on the type of reactor and the composition of the fuel, the qualitative trend of FP produced is the same. Since each product has very different physical and chemical interac tions with the fuel matrix, it is not advisable to assign relative importance based on yield alone.17 Some fission products do not react, some are as found as metallic precipitates, and others react with the fuel matrix to form separate phase s. Although there are many ways to classify the behavior of fission products, the physical state of the product is the most important for fuel applications The first class is volatile products. Fission gases xenon (Xe) and krypton (Kr) are essentially in soluble in the UO2 matrix and are present in the UO2 matrix only by virtue of the fission process. Nevertheless, they can migrate to grain boundaries where they form bubbles. This bubble formation can lead to considerable swelling of the fuel and severely degrades its mechanical properties.18 These gases can also escape to the plenum region between the fuel rod and the cladding, which can contribute to internal stresses on the cladding,19 which can shorten the lifetime of the cladding. Halogens such as bromine (Br) and Iodine (I) are also volatile fission products that are chemically reactive and may react with the zircal loy cladding to cause stress corrosion cracking.20 In addition, 129I has a large half life and this is one of the major radioactive elements that enter the sea from reprocessing plants. The second category includes elements like cesium (Cs), rubidium (Rb) and tellurium (Te) and is termed oxi de precipitates A solid -state fission product, such as Cs, can assume various
22 roles depending on conditions prevailing in the matrix. It is found as a minor component of the grey phase, which is a separate phase consisting of a number of fission products.21, 22 Cs is found as a metallic inclusion found in the fuel, though implantation studies have found it to be present on lattice sites in hyperstoichiometric fuels.23 It can also react with iodine (I) and cause corrosion of cladding materials such as zircalloy and stainless steel.24 Although the yield of Rb is high, it is only found in minor quantities due to its volatile and reactive nature. The behavior of Te is complex as its yield is only 16% of Cs. Te with platinum metals such as rubidium (Rb) and palladium (Pd) forms hexagonal phases that appear as metallic inclusions in irradiated fuels that appear both in the fuel and the fuel -cladding interface.25 In addition, it can precipitate out of solution in the form of complex oxide phases. The final classification according to the physical state is when the fission product forms a separate oxide phase collec tively known as the grey phase. This phase consists of the alkaline earth metals barium (Ba) and strontium (Sr), some amount of zirconium (Zr) and occasionally cerium (Ce). The oxides BaO and SrO are generally regarded as s oluble in the fuel matrix. Howeve r, i n the presence of another fission product zirconium (Zr), they precipitate out in small amounts and tend to form a grey phase perovskite BaZrO3 and SrZrO3. 26 Thi s phase accumulates at the boundary between the columnar and equiaxed grains of the irradiated fuel.27 The fission products affect the fuel matrix. Therefore, in order to be able to predict the system behavior at a macroscopic level, it becomes very import ant to understand the fission product behavior at the microscopic level. Grain Boundaries A grain boundary (GB) is defined as the interface between two grains in a polycrystalline material. A simple way of distinguishing between different types of bounda ries is through the misorientation angle between two grains. The misorientation angle is given by:
23 d b sin (1 13) Here b is the Burgers vector and d is the spacing between dislocations. Low angle grain boundaries have misorientation angles less than 11 and their properties and structure depend on the angle of misorientation. The two major classes of low angle grain boundaries are tilt and twist boundaries. A tilt GB is formed by gradually bending a single grain so that the rotation axis is parallel to the boundary plane see Fig ure 1 6 for an example .11 Another way to understand this boundary is to visualize it as consisting of noninteracting edge dislocations, which are extra half plane s of atoms in a perfect crystal.2 8 If the atomic planes parallel to the GB on either side are the same, then the GB is symmetric. If the atomic planes are not the same, it is an asymmetric GB. A twist boundary is formed when one half of the crystal is rotated with respect to the other half perpendicular to the boundary plane. It consists of a series of screw dislocations. T hese concepts of tilt and twist boundarie s represent somewhat idealized cases.28 The majority of general boundaries are of a mixed type, containing dislocations of different types with different Burgers vector s in order to create the best fit between the neighboring grains. High angle grain boundaries have misorient ation angles greater than 11 The properties of thes e boundaries are considered independent of the misorientation. When the dislocations in the crystal increase due to continued deformation, they begin overlapping and the ordered structure of the GB breaks down. Indeed, these boundaries are more open and disordered. Previously they were visualized as liquid or amorphous layer s between two crystalline grains.29 However, this could not explain the strength of the grain boundaries With the invention of electron microscopy, it is possible to directly observe t he boundaries and it is now accepted that a grain boundary is
24 the region between two grains whose properties depend on the relative orientation of the two grains. T hese boundaries are made up of special structural units. An important phenomenon when discu ssing grain boundaries is the grain boundary energy. It can be defined as the energy required to create a boundary: EGB() E () E A (1 14) Here E ( ) is the energy of structure with the grain boundary, E is the energy of the perfect crystal a nd A is the area of the grain boundary. The energy of the GB depends on orientation relationship of the two grains and also on the Miller indices of GB plane.30 These phenomena can be related to the geometrical criterion of close packing of atoms in a crys tal. Coincident site lattice theory (CSL)31 is one of the most important models proposed to describe the structure of GB, especially high angle grain boundaries. The CSL is a geometrical construction the two grains is described by the reciprocal of the ratio of coincidence sites to the total number of sites. area of coincident lattice cell area of original lattice cell (1 15) Lattices cannot actually overlap. Since a fraction of lattice sites are coincident, boundary structure is more regular than a general boundary. = 1 corresponds to complete coincidence of the two lattices The rotation of the second lattice is limited to those values that bring a (lattice) point into coincidence with a different point in the first lattice see Fig ure 1 -7 for an illustration .31 Atomic positions are not accounted for in this theory. The value is always odd for cubic lattices.
25 The grain boundary energy depends on the values of Since the GB energy is proportio 32 it has been found that the relation between energy and is not straightforward. The plot of GB energies for different values of is shown in Figure 1 8.33 As can be seen, there are larger cusps at lower values of as compared to higher values. Even for the same values of different energies exist due the different orientation planes of the two crystals. However, there exist s a maximum value of beyond that the GB energy value is independent. Deviations from the ideal CSL orientation may be accommodated by local atomic relaxation or the inclusion of dislocations into the boundary. Grain Boundaries play an important role in a variety of material phenomena. Grain Boundary Engineering (GBE) is the practice of obtai ning microstructures with a high fraction of boundaries with desirable properties. In general, desirable properties are associated with boundaries that have simple, low energy structures. Such low energy structures are, in turn, associated with CSL boundar ies. Increasing the fraction of grain boundaries in a material is one of the common methods of strengthening as they impede the motion of dislocations. However, since ultimately they are defects, they tend to decrease the electrical and thermal conductivit y of the material. They are also areas of high interfacial energy and relatively weak bonding that often makes them preferred sites for the onset of corrosion. They are also important for many of the creep mechanisms.29 However, the most relevant applica tion for nuclear application is atomic segregation. Atomic s egregation can be described as selective enrichment of a material at a grain boundary, dislocation or an interface between two materials and even a free surface This produces a phase that has dif ferent properties as compared to the parent matrix and often has deleterious effects. It
26 can lead to grain boundary fracture as a result of stress corrosion cracking, hydrogen embrittlement, environmentally assisted fatigue, and grain boundary corrosion. I n addition, it can vastly influence optical, mechanical and dielectric properties.34 There are various theo ries proposed to describe equilibrium and non-equilibrium segregation. The most famous are the Langmuir Mclean, BET adsorption isotherm and the Fowler isotherm. A general equation for segregation is given by: CGBCBulk exp GsegregationkT (1 16) Here the l eft hand side is the ratio of the concentration of solute at the grain boundary to that in the bulk. G is the free energy change due to segregation. This equation is only valid at low concentration of solute at the grain boundary and assumes no interaction between solute particles. Segregation to different surfaces is different as illustrated in Fig ure 1 9.34 One of the contributions to the decrease in free energy due to segregation comes from the decrease in strain energy.11 It has been proven that the decreas e in strain energy is approximately proportional to [(r2r1)/r1] 2 where r2 and r1 are the atomic radii of the solute and matrix atoms. Hence, the larger the difference between the radii, the greater is the segregation. The absolute size of the solute atom does not matter, just the difference. The grain boundary is a region of misfit in the lattice; hence both smaller and larger solute atoms can be accommodated. Segregation is especially important in the context of nuclear materials. As previously mentione d fission gases xenon and krypton are insoluble in the UO2 matrix and migrate to grain boundaries. This delays the fission gas release to the surface.35 These gases can also migrate via grain boundary diffu sion. It is generally accepted that this mechanism occurs at a higher rate than lattice diffusion.35 All these phenomena can lead to bubble format ion inside the fuel, which is
27 deleterious to the normal operation of the reactor. Hence, it becomes critical to study the segregation phenomena of fission gases at grain boundaries in UO2. Objective a nd Scope of t his Dissertation The objective of the work contained in this dissertation is to elucidate the mechanisms responsible for the production of point defects and fission products in UO2. The influence of microstructural features is understood on the segregation behavior of fission products. Electronic l evel calculations, in particular density functional theory, and atomic level simulations in combination with thermodynamics are used. These methods are discussed in detail in the next chapter. Figure 1 1 The fluorite structure of UO2; uranium atoms are illustrated in blue and oxygen atoms in red.
28 Figure 1 2 The U O binary phase diagram, subscripts S and L indicate the solid and liquid phases. The numbers correspond to the partial pressure of oxygen present in that particular region of the phase diagram. [Ada pted from reference 8 (Figure 4).]
29 Figure 1 3 Poin t defects and defect complexes a ) Vacancy, b) Interstitial, c) Frenkel pair, and d ) Schottky pair. (a) (b) (c) (d)
30 Figure 1 4 The free energy change as a func tion of the number of defects created. The top part represents the change in enthalpy to create a defect and the bottom part is the change in the configurational entropy. [Adapted from reference 1 1 (Figure 6.2).]
31 Figure 1 5 The yield of fission product s for a 3.2% enriched pressurized water reactor pin rated after 2.9% burn up [Adapted from reference 1 4 (Figure 1).]
32 Figure 1 6 A schematic representation of a low angle tilt grain boundary. It is made of a series of non interacting edge dislocations with burgers vector b. [Adapted from reference 1 1 with minor modifications (Figure 6.11).]
33 Figure 1 7 Coincident site lattice model for a 5 (100) grain boundary in a cubic crystal. The left side of the figure shows a tilt boundary and the right si de illustrates a twist boundary. [Adapted from reference 3 2 ]
34 Figure 1 8 Relative energy for a Cu (110) tilt grain boundary as function of the misorientation. [Adapted from reference 3 3 (Figure 4).]
35 Figure 1 9 Segregati on of a free atom to different surfaces with changes in the Gibbs free en ergy. [Adapted from reference 34 (Figure 5. 1).]
36 CHAPTER 2 GENERAL INTRODUCTION TO SIMULATION METHODS Density Functional Theory T he description of the interactions between e lectrons and the nuclei of atoms requires quantum mechanics. These interactions may include interactions between electrons, between electrons and nuclei the kinetic energy of electrons and nuclear nuclear interactions. A difficult problem in the quantum m echanical description of a system consisting of electrons and nuclei is to take into account the electron -electron interaction.36 Two important concepts while discussing electron -electron interactions are exchange and correlation. Electrons have both wave -like and particle like nature. The wave function (a mathematical tool used in quantum mechanics to describe the wave particle duality) of a many electron system is anti -symmetric si nce the electrons are fermions. This anti -symmetry produces a spatial separation between two elec trons with the same spin and reduces the Coulomb energy of the system. This reduction in the energy is called the exchange energy. This energy is taken into account in a Hartree -Fock type calculation.37 The other contribution in electron -electron interactions is the correlation energy. The Coul omb energy of a system can further be reduced if the electrons having opposite spins are also spatially separated. The difference in the total energy of the many electron system and the energy calculated using the Hartree -Fock approximation is the correlat ion energy. It is extremely difficult to calculate the correlation energy of a many body system. Density Functional Theory (DFT) is one of the most powerful approaches to deal with the problem of exchange and correlation. It was proposed in 1964 by Hohenberg and Kohn38 and further developed by Kohn and Sham39 in 1965. It is based on two theorems:
37 Theorem 1: The total energy of a system of interacting electrons (including exchange and correlation) in the presence of a static external potential is a unique functional of the electron densit y. Theorem 2: It is possible to formally replace a system of many electrons by an exactly equivalent set of self -consistent single electron equations. The first theorem reduces the problem of solving a system of 3N spatial coordinates with into a much simpler solution of electron density, see Figure 21 for an illustration.40 The minimum value of the KohnSham functional is the ground state energy of the system and the density that yields this energy is the exact electron density of single particles.36 The Schrdinger equation can be written in terms of the set of KohnSham equations as: ) r ( ) r ( ) r ( V ) r ( V ) r ( V m hi i i XC H ion 2 22 (2 1) Here Vion (r) is the static electron ion potential, i is the Kohn-Sham eigenvalue and i is the wave function of state i. The Hartree potential is given by: r d r r ) r ( n e ) r ( VH 3 2 (2 2) where n(r) is the electron density. The exchange -correlation potential is given as: ) r ( n )] r ( n [ E ) r ( VXC XC (2 3) Here Exc is the exchange -correlation energy. The Kohn Sham equations are a set of eigen equations that must be solved self -consistently. The terms within the bracket of Equation 2 1 can be seen as a Hamiltonian. The highest occupied eigenvalu e is approximately equal to the unrelaxed ionization energy of that system. The limitation of this approach is that it is nearly
38 impossible to calculate the exchange correlation energy of a many electron system. If Exc functional were exactly known, then t he use of Equation 2 3 would result in a potential that would take into account the effects of exchange and correlation exactly. Local Density and Generalized G radient Approximations There have been many approximations proposed to deal with the problem of inexact description of the exchange -correlation functional. The most common and almost universal is the local density approximation (LDA). Here the exchange -correlation energy per electron in the inhomogeneous electron system [ xc(r)] is assumed to the same as the exchange -correlation energy per electron for a homogeneous gas having the same density at r. Thus the expression becomes: r d ] ) r ( n [ ) r ( n ] ) r ( n [ EXC LDA XC (2 4) ) r ( n ] ) r ( ) r ( n [ ] ) r ( n [ VXC XC (2 5) ) r ( n n o hom XC XC| ) n ( ) r ( (2 6) LDA ignores the corrections to the exchange -correlation energy due to nearby inhomogeneities and is purely local in nature. Despite this, calculations performed using LDA are fairly successful in reproducing key properties. The exchange part is the exact expression derived for a uniform gas. The available versions of LDA differ only in their representation of correlation.40 One of the reasons LDA is successful is because it satisfies the sum rule of the exchange -correlation hole.37 It gives a well-defined energy minimum for a system where the spin of the electron is not taken into account. But for systems where the electron density is changing
39 rapidly or for a magnet ic system, which has several local minimas, LDA fails to give accurate results and one has to look beyond this simple expression. The other functional in increasing hierarchy takes into account the gradient of the electron density as an independent variab le and called the generalized gradient approximation (GGA). The gradient introduces non locality into the expression for exchange and correlation. The expression now becomes: r d ] n n n n [ f ] n n [ EGGA XC (2 7) To facilitate practical calculations, ] n n n n [ f has to be determined. This search has evolved in two separate directions. One method is where the expression is determined from known expansion coefficients and other theoretical conditions. This is done by taking the second order density gradient expansion of the exchange correlation hole of the electron and cutting off the longrange components to satisfy the properties of the exchange correlation hole. This method is called the parameter free approach and is mostly used in physics. In the other approach, parameters are determined by fitting to experimental data or accurately calculating atomic and molecular properties. Most GGAs used in chemistry applications are empirical. In comparison to LDA, GGA improves total energies, barriers and structural properties. However, it is not advisable to assume that the more complicated functional would always give better properties. It is well known that LDA works better for certain classes of properties and even systems, especially for surface pr operty calculations.40 Selection of a proper functional is particularly challenging for problems bridging physics and chemistry such as reactions on surfaces. Preliminary tests have to be performed using different functionals to determine the best approximations for that system and particular set of properties one is interested in.
40 Periodic Supercells Alt hough many properties of a many -electron system can be calculated using a series of single particle equations, there still remains the problem of calculating the wave functions for these particles. A wave function ( ) is a mathematical tool used in quantum mechanics to describe the wave particle duality of the electron. The square of the wave function | |2 gives the probability of finding an electron at a certain time and position and is a complex conjugate The second problem is that since the wave func tion for each electron extends over the entire solid, so the basis set becomes infinite. Both these tasks can be overcome by using the Blochs theorem41 which states that in a periodic solid, the wave function of an electron can be written as a product a wave like part and a cell periodic par t: ) r k ( i exp( ) r ( u ) r (i i (2 8) The cell periodic part )r ( ui can be expanded using a set of discrete plane waves whose wave vectors are just reciprocals of the lattice vectors of the solid: ) r G i exp( C ) r ( uG i G i (2 9) where the reciprocal lattice vectors G are represented as rG =2 m, l being the lattice vectors of the solid and m being an integer. Therefore, the electronic wave function becomes: ] r ) k G ( i exp[ C ) r (G k i G i (2 10) The Bloch theore m reduces the problem of calculating infinite electronic wave functions to a problem of calculating the wave functions of finite particles at special points in the reciprocal space called k points.42 The density of these points is proportional to the volume of the solid. The
41 occupied electronic states at these k points contribute to the total electronic potential and the wave functions at closely spaced points are almost identical. Hence, only a few k points are needed to cover the entire reciprocal space. A higher accuracy of calculations can be reached by using a tighter k point mesh. The total energy calculated will converge with increasing k points and the magnitude of error approaches zero. In general, the total energy and electronic potential of semi -conductors and insulators can be determined by using a smaller number of k points. In con trast, to define the Fermi surface of a metallic system, a denser k point mesh has to be used.42 The other component of the Bloch theorem is the wave like part. To expand the electronic wave functions in terms of a wave -like part, one would need an infinite set of plane waves. However, it is known that the plane waves having smaller coeffici ents (CiG), corresponding to lower kinetic energy, have a more significant contribution to the calculation of wave functions. Hence, it is reasonable to truncate the plane waves based on a kinetic energy cut -off. The relative energies are more important in many cases as compared to the total energies and the total energy difference converges rapidly with cut -off. There is obviously some error associated with this truncation. The error can be reduced by systematically increasing the energy cut -off. However, one significant error associated with plane waves can occur when the volume of the cell is changing and the energy does not converge with respect to that volume. As discussed previously, the k point mesh is proportional to the volume of the solid. If the v olume is changing, discontinuities will occur for different k point meshes at different cut -offs. This problem can be resolved by performing constant volume calculations and then testing for the energy convergence with respect to k points. However, most materials are non -periodic, containing defects and surfaces. In its original form, the Bloch theorem cannot be applied to such systems. Hence, we have to use the supercell
42 approach while performing defect calculations. A schematic of the supercell is shown in Figure 2 2.37 Here, the cell with th e defect is repeated infinitely in space by applying periodic boundary conditions (PBCs). Thus, the energy per unit cell of a periodic solid with an array of defects is calculated as opposed to the energy of cell with a single defect. Care has to taken to include enough bulk solid around the defect so that the neighboring defect images do not interact.43, 44 This can be done by employing larger and larger supercell and test ing for energy convergence. For a surface, a large enough vacuum region has to be included in a direction perpendicular to the surface. Electron -Ion Interactions The plane wave basis set is impractical to expand the electronic wave functions despite the a pplication of the Bloch theorem and kinetic energy cut -off. The core orbitals of an atom are tightly bound to its nucleus and the valence electrons oscillate rapidly in the core region. These two phenomena would require a large basis set to perform an all electron calculation. This limitation can be overcome by the application of the pseudopotential approximation. This theory assumes the frozen core approximation which in other words, takes advantage of the fact that most of the properties of materials are determined by their valence electrons. The strong ionic potential of the core electrons is replaced by a weaker pseudopotential which acts on pseudo electrons rather than the true valence electrons. Thus, the oscillation of the valence electrons disappears and the pseudo wave functions have no radial nodes in the core region. The most general form of the pseudopotential is: | lm V lm | Vlm l NL (2 11) where | lm > are the spherical harmonics and Vl is the pseudopotential for angular momentum l.
43 The pseudopotential is constructed in such a manner that the scattering properties of the pseudo electrons are similar to those of the true valence electrons. The two potentials are the same outside the core region. Since the scattering properties are diff erent for each ion, the pseudopotential is angular momentum dependent. Even though these norm -conserving pseudopotentials45 have the advantage of formal simplicity, the computational cost of treating transition metals and actinid e elements is quite high. The two most popular pseudopotentials for treating these elements are the ultrasoft pseudopotentials (US -PP) by Vanderbilt46 and the projector a ugmented wave (PAW) method by Blchl.47 Generally, PAWs are expected to be more accurate than US -PPs. There are two reasons for it: the radial cut -offs or core radii are smaller for P AWs and they exactly reconstruct the exact valence function with all nodes in the core region. Since the core radii are smaller, the energy cut offs are larger for PAWs. But in practice, the computational cost is slightly high. PAWs have been used in this study. Strong Correlation Electronic correlation in a quantum system refers to the interaction between electrons. The term correlation comes from statistics, which means that two distribution functions, m and n, are not independent of each other. UO2 is a very special class of material. It is an insulator at 0 K48 but changes into a semi cond uctor as it is heated up. It is also anti -ferromagnetic below its Neel temperature at 30 K.49 For these reasons, it is s ometimes called a Mott Hubbard insulator. Here, electrons are viewed as occupying their normal orbitals and then hopping between orbitals during conduction. This is represented as a hopping integral mathematically. However, the general crystalline band theories do not include the interaction between electrons. By setting up conduction in terms of hopping integrals, this model is able to include the on-site repulsion between electrons which is
44 Coulombic in nature. The hopping integral is a function of distance between atoms and the Coulombic repulsion is not. The Hubbard model can thus explain the shift from conductor to insulator for certain transition metal oxides as the temperature increases because this increases the nearest neighbor distance and the on -site repulsion becomes dominant. Mott proposed a similar model for nickel oxide (NiO). The band gap of the material is given as: EGap = U 2 zt. (2 12) where, U is the Coulombic repulsion, t is the hopping integral and z is the atomic number. This model can also explain the transition from conductor to insulator for rare -earth elements as the atomic number increases. The strong correlation is the reason why conventional DFT methods such as LDA and GGA, which are one -electron methods with an orbital dependent potential, fail to capture the ground state behavior of Mott insulators. This was discovered by Kelly and Brooks in 1987,50 where they performed spin dependent LDA (LSDA) and found the resulting structure to be unstable with a lattice parameter of 90 % the experimentally determined one. The reason suggested was partial fillin g of the 5f orbitals of uranium (U) and the resulting strong correlation. The problem of metallic ground states predicted with conventional DFT has attracted much attention over the years and there have been many attempts to correct this. However, the mos t practical (in terms of the computational cost associated and the size of the supercell) is the LSDA+U approach proposed by Anisimov et al.5153 It is based on the Mott model and includes the on -site repulsion. In general, the second quantized Hamiltonian for the f electrons can be written as: jl ^ jl ^ l l j jl ^ jl ^ l l l l j il ^ il ^ l l j i l l j i ^n n U n n ) J U ( c c t H 2 2 1 (2 13)
45 Here th e first term represents the hopping of electron between lattice sites i and j. l and l are indices of the f orbitals localized on one particular atom and jn^ = ^ ^c c The second and the third terms of the equation descr ibe the on-site repulsion between electrons with U being the onsite Coulomb interaction parameter and J being the exchange interaction parameter. In conventional LSDA, which is a mean field approach, the orbital occupation numbers njl, are independent of the magnetic quantum number, l. The energy EMF, which corresponds to a DFT -LSDA treatment of the potential energy terms of equation 2 13, is represented as: j j j j j j j j MFn n U )] n ( n ) n ( n [ ) J U ( E 1 1 2 1 (2 14) where l j jn n and l j jn n Taking into account the terms not described by equation 2 14 and no longer assuming jn to be independent of l, one gets: ) n n ( ) J U ( E, l j l j l j 22 1 (2 15) Differentiating E with the occupation number of a particul ar orbital gives the correction to the potential acting on that orbital: ) n )( J U ( V, l j l j 2 1 (2 16) The inclusion of this term into the conventional DFT leads to the so -called LSDA+U or SP GGA+U approach:54 ) n n ( ) J U ( E E, m m GGA SP U GGA SP 22 (2 17) Here nm, is the occupation number of the mth f state and the second term on the right hand side of the equation accounts for double counting.
46 Accuracy of Calculations It is important to know the error bar on the calculated properties. In general, DFT is regarded to be very accurate in predicting various properties, espe cially total energies, as it does not depend on any fitted parameters. By this argument, all properties that can be related to the total energies can also be determined accurately.40 For example, the lattice parameter is determined by performing a series o f total energy calculations as a function of different lattice parameters and finding the minimum energy. The most stable defect is one which has the least total energy of the system. For defect formation energy calculations, it is differences in the total energy that are most important and hence the systematic errors cancel out to some degree. However, there are possible sources of error which should be taken into account. One component of the error is purely numerical in nature, associated with the size of the supercell which is important while performing defect calculations (electrostatic interactions) the plane wave cut -off, the k point mesh, the force and the energy convergence criteria Well converged calculations are typically regarded as being very precise, with an error of about 0.1 eV or so Another source of error is DFT itself. This may include the choice of pseudopotential, which approximates an all electron calculation and the exchange -correlation functional In some cases, DFT underestimates the band gap which becomes important while considering charged defects. In all cases, a conservative error bar of 0.1 eV must be assumed.37 Regardless of any error, the physi cal trends are well reproduced. Use of different approximations would lead to dif ferent values of absolute energies. However, for solid state systems, DFT reproduces the physically observed phenomena. Lattice Statics The second approach adopted in this work is the energy minimization through the use of empirical potentials. The advanta ges of this approach are that larger supercell sizes can be
47 treated at a small fraction of the computational cost of DFT calculations. The calculated values are less accurate than those obtained with DFT, but the effects of supercell size are not expected to change significantly with method. These energy minimization calculations are referred to as static because they do not explicitly take into account the lattice vibrations or the configurational entropy. The internal energy of a solid is a many body fu nction that depends on the position and momenta of all electrons and nuclei present. This is difficult to solve by any level of theory and hence the assumption of ionic solid is made where the crystal consists of regularly spaced positively and negatively charged ions. The total energy due to N atoms is divided into multiple components using Madelung expansion: U Ui i 1 N 1 2 i 1 NUij j 1 N 1 6 i 1 NUijk k 1 Nj 1 N (2 18) where the first term is the self -energy of the atom; the second term is the energy due to pair -wise interaction and so on. The contribution of the higher terms becomes progressively smaller and beyond a certain point, these terms are ignored. For ionic, cubic solids, generally the expansion to the second term is sufficient. The energy minimization in these calcul ations is carried out by optimizing the geometry of the structure. The internal energy U may be written as a function of the position of the ion (x) and its displacement ( x) as: .... dx u d dx du ) x ( U ) dx x ( Ux x 0 2 2 0 0 02 1 (2 19) Here, the first derivatives are collectively written as the gradient vector, g, and the second derivatives are the Hessian matrix H. The common methods of structural optimization are the methods of steepest desc ents, conjugate gradients and Newton Raphson. All three methods have been used in this study.
48 Empirical potentials are functions that describe the interaction of atoms at close range. The general form of the potential consists of a van der Waals attractiv e term and a short range repulsive term to prevent like charged ions coming in close proximity of each other. For solids, one type of potential used is the Buckingham potential which is given by, 6r C Ae Ur ij (2 20) where A is related to the hardness of the ions, is the size of the ions while C represents the van der Waals attraction. These parameters are fit to ab initio data or experiments to give the best properties. Other forms of potentials, in which the ions are partially charge d and have a Morse term, to represent covalence have become very popular for UO2 as they give better match with experiments.55 One of the limitations of empirical potentials is that they are fit to certain properties of the material and hence the transferability of the potential to simulate other properties is always a challenge. Defective Lattice The calculation of defect energies is always more difficult and approximate that the calculation of bulk properties. The defect, in theory, can have long range i nteraction with the lattice. This is especially true if the defect is charged, as is the case of point defects in ionic solids. In addition to the supercell approach described above for defect calculations in DFT, the Mott Littleton56, 57 approach is another method of choice for modeling of defects. The basis of this approach is to divide the space around the defect in three separate regions, known as I, IIa and IIb, see Figure 2 3 for an illustration.56 In region I, the atoms around the defect are allowed to relax explicitly and all interactions are treated at an atomistic level. However, it is computationally inefficient to keep treating all interactions explicitly to achieve the desired
49 degree of convergence. Therefore, some allowance is made for these interactions in region IIa, in a more economical manner. In this re gion, the ions are assumed to be in a harmonic well and they respond to the perturbations caused by the defect. It is important to optimize the positions of the ions before a defect calculation since this assumption is only valid for small disturbances. In region IIb, only the implicit polarization of the sub lattices is considered, not the actual ions. This region is assumed to be sufficiently away from the defect and does not contribute significantly towards the defect energy. The total energy of the sy stem can be written in terms of the position of the ions (x) in region I and the displacement of the ions ( ) in region IIa as: ) ( E ) x ( E ) x ( E E 2 12 1 (2 21) Here E1 and E2 are the energies of the regions I and II respectively and E12 is the intera ction energy between the two regions. Another important point in the discussion of defect calculations is how the energy of the defective lattice is minimized. There are two approaches of doing this: constant pressure and constant volume. The constant vol ume approach minimizes the strain of individual ions while the constant pressure approach determines the ground state energy by accounting for the strain on the ions as well as minimizing the forces on the unit cell itself by adjusting the shape and volume of the unit cell. There is a lot of debate as to which method is to be adopted. As pointed out by Catlow,57 these calculations must be done at constant volume with a fixed lattice parameter. At high temperatures, there may be significant differences between the two approaches. Allowing the lattice parameter of the cell to change (constant pressure) corresponds to a hig h concentration
50 of defects.58 The modeling of high concentration was not the objective of this study and hence the constant volume approach has been adopted. Computational Software VASP (Vienna Ab Initio Simulation Package)59 60 is a program for performing ab initio quantum -mechanical molecular dynamics (MD) using pseudopotentials and a plane wave basis set. The approach implemented in VASP is a self -consistency cycle t o evaluate the solution of the Kohn-Sham functional. The interaction between ions and electrons is described using ultrasoft Vanderbilt pseudopotentials (US -PP) or the projector augmented wave method (PAW). Forces and stress can be easily calculated with V ASP and used to relax atoms into their instantaneous ground state. The strong correlation of actinide elements is included within this scheme. GULP (General Utility Lattice Program)6 1, 6 2 is a program for performing a variety of types of simulation on materials using periodic boundary conditions. These include 0 D (molecules and clusters), 1 D (polymers), 2 D ( surfaces, slabs and grain boundaries), or 3 D (periodic solids). The capabilities of this program include energy minimization, calculation of transition states, crystal properties and fitting of empirical potentials. A variety of force fields can be used w ithin GULP that includes the shell model and rigid ion model for ionic materials, molecular mechanics for organic systems, and the embedded atom method for metals.
51 Figure 2 1 Illustration of the DFT app roach. [Adapted from reference 4 0 (Figure 1).]
52 Figure 2 2 A schematic of the supercell method for a point defect in a solid, the supercell being the area enclosed by the dashed lines. [Adapted from reference 3 7 (Figure 3).]
53 Figure 2 3 Illustration of the Mott Littleton me thod. [Adapted from reference 5 6 (Figure 1).]
54 CHAPTER 3 ENERGETICS OF INTRIN SIC POINT DEFECTS IN URANIUM DIOXIDE1 Introduction UO2 contains a variety of point defects, which are created during its use as nuclear fuel in a pressuriz ed water reactor (PWR). These may include intrinsic point defects such as vacancies and interstitials and extrinsic defect like fission products that are produced during irradiation. The focus of this chapter is on intrinsic defects and fission products ar e considered in Chapter 4. In stoichiometric UO2, the possible point defects are cation and anion vacancies and interstitials. These point defects can also combine to form neutral defect complexes such as anion Frenkel, cation Frenkel and Schottky defects. A majority of these defects persist in the fuel although many recombine or otherwise annihilate each other. Over time and higher burnups, these point defects form clusters that cause fuel swelling and ultimately change the crystal structure of UO2; this p hase change is accompanied by a significant change in volume, which has highly detrimental effects on the microstructure and mechanical stability of the fuel pellet. These intrinsic point defects also affect the diffusion of the fission gases formed during operation by forming traps and allowing them to escape from the fuel matrix and damage the fuel cladding. 14 They also have negative effects on the thermal transport properties of the fuel. Thus, an improved knowledge of the sta bility of these point defects is required for improved understanding of fuel performance. The point defects produced include isolated vacancies and interstitials and defect complexes. These individual point defects can also occur in a variety of charge sta tes, especially uranium which can easily change its valence state from U4+ to U5+ and U6+ depending on the 1 In this chapter, I describe the electronic structure and atomic level simulations of point defects in UO2. This work was published in Ref. 63 P. Nerikar, T. Watanabe, J. Tule nko, S. Phillpot, and S. Sinnott, Energetics of intrinsic point defects in uranium dioxide from electronic -structure calculations, Journal of Nuclear Materials 384, 61 (2009)
55 partial pressure of oxygen present. The effect of charge can influence the relative stability of the defect cluster. Since UO2 is known to be non-sto ichiometric over a wide range of temperatures and partial pressures, the possible combinations of point defects and charge states is quite large. The relative stability of defect is affected depending on if the fuel is hypostoichiometric or hyperstoichiome tric. Review of Literature There have been numerous experimental studies on defect properties of UO2. Matzke et al.64 have presented a comprehensiv e review of the available experimental data and critically analyzed it. The defect properties are usually measured using diffusion data, as shown in Figure 3 1. Generally oxygen diffusion is much faster than metal diffusion. The diffusion coefficient is de rived from measurements of the rate of thermodynamic equilibrium achieved by the migration of point defects when the temperature and partial pressure are changed. Unfortunately, this equilibrium is not always achieved. In the case of oxygen diffusion, the situation is complicated by numerous factors. The first limitation is that conventional diffusion techniques cannot be used as oxygen does not have a naturally occurring isotope. The mechanism of oxygen migration in UO2 also changes with stoichiometry. In hypostoichiometric and stoichiometric UO2, the migration mechanism is through oxygen vacancies while in UO2+x, it is through oxygen interstitials. Thirdly, in UO2+x, oxygen forms cluster, as shown in neutron diffraction studies by Willis et al.65 These difficulties in measurement of diffusion properties lead to a wide range of formation energies of defect complexes, as will be discussed later in this work. Theoretical calculations are complementary to exper imental methods and are well positioned to provide fundamental insight into defect chemistry. 66, 6 7 This is due to the fact that systems examined computationally have well -defined compositions and the influence of each change in condition can be analyzed separately. Catlow12, 15, 57 undertook the first theoretical
56 survey of the defect properties in UO2 and other fluorite materials. He used an empirical potential that was fitted to a Hartree -Fock type calculation. His work confirmed that the experimental observation that the anti -Frenkel pair, consisting of an oxygen vacancy and interstitial, was the dominant defect in UO2. His work also showed that the oxygen vacancy was the dominant migratin g defect in stoichiometric UO2 instead of the oxygen interstitial, as was previously believed. With increases in computational resources and both more powerful and more efficient algorit hms, DFT can be applied to larger systems and has been successfully us ed to examine defect formation in numerous materials, 40, 6871 including UO2, with greater accuracy. For example, Petit and co -workers 717 3 used the LDA and GGA to predict the energetics associated with point defect formation. Their resul ts agree well with experimentally determined values. 6 4 However, LDA and GGA predict bulk UO2 to be a metal while it is experimentally known to be a n insulator. This restricts the applicability of their analysis to neutral systems 73 and prevents them from exploring the effect of charge on defect stability. In addition to being limited to neutral defects, earlier studies were restricted to small supercell sizes. Petit et al. 72 and Freyss et al. 73 both used a 24atom supercell owing to the high computational expense of modeling UO2. Despite these limitations, the importance of their work cannot be underestimated as these were the first calculations on UO2 and showed that DFT could predict experime ntal trends successfully. With recent increases in computational power, larger system sizes have be considered. Iwasawa 74 and Gupta et al. 75 each utilized the spin polarized GGA with Hubbard U correction term (SP -GGA+U) method for a 96 atom supercell and predicted the correct insulating structure of bulk UO2. Moreover, they correctly predicted its anti -ferromagnetic ground state and found that including the correct magnetic state changes the absolute values of the formation energies of neutral point defects.
57 In this study, the earlier work on neutral defects has been extended to include charged po int defects in UO2. A combination of experimental thermodynamic data and DFT calculations has been used. The DFT calculations themselves should give the most accurate results to date, since they use the currently most sophisticated approximation within GGA the PAW method, and 96 atom supercells, as in the previous work of Gupta et al. Most significantly, a new bridge is made between the DFT calculations and temperature, oxygen partial pressure and the charge associated with a defect, which are key paramete rs for controlling the type and concentration of defects. When the current and previous results on defects are considered as a whole, a coherent and consistent picture of the energetics of neutral intrinsic defects emerges The introduction of charge leads to significant changes in the formation energies that vary with Fermi energy. In addition, the predicted defect formation energies of Schottky and anti -Frenkel defect complexes that consist of charged components are significantly lower than the formation energies of complexes that consist of neutral components. Computational Methodology Electronic Structure C alculations The DFT calculations were performed with spin polarization by the projector augmentedwave (PAW) method 47, 76 utilizing the GGA exchange correlation functional as implemented in the Vienna Ab -Initio Simulation Package (VASP) 59 6 0 The U 6s2 6p6 7s2 5f2 6d2 and O 2s2 2p4 electrons are treated as valence electrons. The Dudarevs simplified scheme (the so called SP GGA+U method) 54 was used to include the effect of the strong correlation of 5f electrons in uranium. In addition, t he values of U=4.5 eV and J=0.54 eV give the best match with experiment for equilibrium properties (see Table 3 1). Here the 5f Coulomb correlation energy (U) has been determined with a combination of X ray photoelec tron spectroscopy and Bremsstahlung isochromat spectroscopy and the onsite exchange energy (J) was determined using X ray
58 photoemission spectroscopy .48 This leads to a Ueff value of 3.96 eV and hence this value was chosen These values are similar to those previously used by others 75, 76, 77 and also to experimental measurements. 48 UO2 is experimentally known to be anti -ferromagnetic below its Neel temperature of 30 K, 77 and this state has been taken for these calculations. With regard to other calculation parameters, a 1x1x1 unit cell was used for the initial validation The Brillouin zone sa mpling used a 4x4x4 Monkhorst -Pack k -point mesh, 42 and the cut -off energy for the plane waves was 400 eV. The combination of the above parameters resulted in a lattice parameter of 0.549 nm and a band gap of 1.92 eV, both of which are in good agreement with experiment (see Table 3 1). A 2x2x2 unit cell consisting of 96 atoms was then used for the pristine and defect containing structures for the defect energetics calculations. In the case of the pristine supercell, cell shape and atomic positions were relaxed to their equilibrium positions. However, in the case of the supercells that contained point defects, on ly the atomic positions were relaxed while the cell shape was kept fixed. This is critical for defect calculations in the limit of dilute defect concentrations, 58 since allowing the cell volume to relax corresponds to calculating the lattice constant of a system with a high concentrat ion of defects (the actual defect concentration within the supercell), which was not the objective of this work. It should however be noted that if the system size if large enough, both constant pressure and constant volume calculations should converge to the same result. Defect Formation Energy The formation energy of a defect as a function of temperature, partial pressure, defect species and charge state q is given as: 58 ] V E E [ q ) P T ( n ) perfect ( E ) q ( E ) P T q ( GV F i i total total f (3 1)
59 Here Etotal ( q) is the optimized energy of the supercell containing the defect of charged state q; Etotal (perfect) is the optimized energy of the perfect supercell. Both these values are obtained directly from the DFT calculations. In Equation 3 1, ni represents the number of atoms of species i added to (ni> 0) or subtracted from (ni < 0) the system; i is the chemical potential of species i, which can be defined as the change in free energy when a species is ad ded or taken out of a system. The last term is seen as an electronic chemical potential that represents the change in energy associated with charged defects. In this term, EF is the Fermi level in the bulk with reference to the valence band maximum and EV is the energy of the bulk val ence band maximum. Since the bulk value of the valence band maximum to estimate the energy of a defect supercell was used the electrostatic potential of the defective structure has to be aligned with the bulk value; the term V represents this alignment. The Fermi energy is treated as a variable in this approach and is dependent upon the charge associated with the cumulative effect of defects and dopants in the system. A physically meaningful range up to the band gap of the m aterial, around the Fermi energy of pristine UO2 is considered. The vibrational entropy is neglected; it has been previously demonstrated in the context of TiO2 that such vibrational contributions are small. 78 The configurational entropy is treated within the ideal solution model. 7 9 Thermodynamic Component The chemical potential i (T, P) is also treated as a variable in these calculations. Following the approach of He et al., 78, 79 the oxygen chemical potential in UO2 is expressed as: 0 0 0 2 0 0 22 1 2 1 P P ln T k ) T ( ] G [ ) P T (B O UO f U UO O (3 2) where 0 UO2 and 0 U are the chemical potentials of bulk UO2 and uranium metal as calculated with DFT. Uranium metal forms three phases with increasing temperature (see Figure
60 1 2).8 The parameters U and J were kept constant while calculating the chemical potentials of these phases, in order to maintain the same reference state while calculating the defect properties. G0 f,UO2 is the Gibbs energy per mole of UO2 under standard conditions obtained from thermodynam ic data, 80 and 0 O (T) is the change in oxygen chemical potential with change in temperature obtained from thermodynamic data, 8 0 Thus, Equation 3 2 can be used to calculate a chemical potential of oxygen that varies in a physically meaningful way with temperature and oxygen partial pressure; this would not have been possible if this c hemical potential had been calculated directly from a DFT calculation of an oxygen molecule that is, by definition, carried out at zero Kelvin in perfect vacuum. Thus, with the combination of Equations 3 1 and 3 2, the effect of system conditions on the de fect formation energies obtained through DFT calculations can be considered. Results and Discussion Equilibrium Properties of Bulk UO2 and Uranium Metal The first step to validate the approach adopted in this work was to calculate the lattice parameter an d band gap of bulk UO2 using different pseudopotentials. Table 3 1 compares the lattice parameter, the band gap and the magnetic moment of the uranium ion from published computational and experimental work with values calculated using various approximation s within DFT. The agreement with the experimental lattice parameter improved as these approximations become more sophisticated. For example, pure LDA and GGA give the poorest agreement with experiment, but the agreement significantly improved with the addi tion of spin polarization. Of all the methods considered, the PAW -SPGGA+U method gave the best match with experimental data, even though it tends to overestimate the lattice parameter slightly. These results are similar to the PAW SPGGA+U results of Gupt a et al.; 75 the small discrepancy
61 between the two sets of results most probably arises from the slight differences in the value of Ueff (U J) used in the calculations (in this work a value of 3.96 eV for Ueff has been used and Gupta et al. used a value of 4.00; both values fall within the experimentally determined range. 48 Uranium dioxide is an electrical insulator, 48 however, conventional DFT predicts it to be a metal unless the 5f electron on site repulsion is included using in the +U term (refer to Table 3 1). Figure 3 3 shows the density of states (DOS) plots calculated from PAW -SP GGA and PAW SP GGA+U calculations. While the former predicts UO2 to be metallic, the latter gives a band gap of 1.92 eV, which is only 0.08 eV smaller than the experimental value. More detailed analysis showed that the valence band consists mainly of U 5f and O 2p orb itals while the conduction band as shown in Figure 33 b) was mainly derived from 6d and 5f orbitals of uranium. In the PAW -SP GGA method, all of these orbitals are clustered together leading to the metallic state shown in Figure 3 3 a). The spin polariza tion of the U 5f orbital controls the anti -ferromagnetic ground state. The calculations take this effect into account. The magnetic moment of each uranium ion is predicted to be 1.93 Bohr Magneton ( B) a value that is about 10% larger than the experimenta l value. Of course, the net magnetization is 0 as should be the case for an anti -ferromagnetic material. This correct description of the electronic state is an indication that this method will be able to correctly describe charged defects. The different ch arges associated with a defect affect the Fermi level (EF) and therefore influence the relative stability of the charged defect relative to the neutral defect; this is discussed in more detail later in the work. One of the goals of this work was to predic t the role of temperature on point defect formation. Bulk uranium metal exists in three different phases and depending on temperature. The -phase is the most stable phase under ambient conditions and has an
62 orthorhombic crystal structure with four atoms per unit cell. The phase, which is stable at intermediate temperatures, has a body -centered tetragonal structure; the highest temperature phase has a body -centered cubic structure. The PAW -SP GGA+U calculations predict this phase order correctly as indicated in Table 3 2. This provides further validation for this approach and is also critical for predicting the effect of temperature on defect formation energies from Equation 3 1. However, it should be pointed out that the corresponding calculati ons for the other DFT approximations were not performed. These calculations have been carried out previously by Freyss et al. 72 who also found the correct phase order. These initial calculations of lattice parameter, band gap, magnet ic moment and phase order agree well with experiment and provide an indication that the PAW -SPGGA+U approach will be appropriate for defect energy calculations. Formation Energy of Neutral Defects Although DFT has been found to accurately describe the energetics of neutral defects, it is still limited to a relatively small number of atoms and electrons in the system. UO2 is especially challenging for DFT because it has a large number of electrons, strong electron -electron correlations, and is anti -ferroma gnetic, all of which must be included to properly model the material. The calculations performed, therefore use a 2x2x2 supercell consisting of 96 atoms, see Figure 3 2 for an illustration. This is effectively the current system -size limit for such calcula tions. This limitation on system size is important as there could be interactions between periodic images of the defect in neighboring supercells. To estimate the effects of the limited supercell size used, defect calculations were also carried out using a n empirical potential within the General Utility Lattice Program (GULP). 61 6 2 The advantages of such an approach were that larger supercell sizes could be treated at a small fraction of the computational cost of DFT
63 calculations. The calculated values are known to be less accurate than those obtained with DFT, but the effects of supercell size are not expected to change significantly with method. The Yamada potential 8 1 was chosen for these static energy calculations. This potential includes parti al charges and a Morse term to model covalent character. It yields more accurate values for energetics of point defects than most other Buckingham type potentials. 55 The normalized defect formation energies of the anti -Frenkel, Frenkel and the Schottky defects are plotted i n Figure 3 4. The defect energy converged the system size was increased. For the anti -Frenkel and the Schottky defects, the energy is almost converged as the size of the supercell is increased to 2x2x2. For the Frenkel defect, the change in the defect ene rgy is 12% between the 2x2x2 and 3x3x3 supercells and changes very little ( 5%) beyond this as the system size increases further. It was also noted that the empirical potential predicted the experimentally observed trend for the defect energies. After the validation of system size using empirical potentials, DFT was used for the defect properties. All the point defects whose formation energies have been reported so far in the literature from DFT calculations have been neutral. For the sake of completene ss, and to unify the sometimes dissimilar literature values, the formation energies associated with neutral uranium and oxygen vacancies and interstitials were calculated. These findings are compared to published results in Table 3 3. The qualitative trend predicted by all the calculations is the same. In particular, they all predict that the oxygen interstitial has negative formation energy. This suggests that UO2 is unstable in the presence of oxygen vapor, which is consistent with the fact that UO2 is we ll known to be hyperstoichiometric. 65 The large positive formation energy for O
64 vacancies is also consistent with the very narrow range of stabilit y of hypostoichiometric UO2, experimentally observed at only high temperatures. Although the various methods agree on the order of stability of the defects, they do not agree on the absolute values of defect formation energies. There could be a number of possible reasons for this. Some calculations 72 used pure GGA as their exchange correlation functional while others used SP GGA+U. 75 In the SP GGA+U methods, a spin polarized, anti ferromagnetic state was consi dered while in the GGA calculations a non -spin polarized and non magnetic state was considered. Even when the same exchange correlation functional was used, differences in the results are most likely due to the fact that different types of pseudopotentials were used. For example, Iwasawa 74 used the PBE -SP GGA+U pseudopotential while Gupta et al. 75 utilized the PAW SPGGA+U pseudopotential; the latter is expected to be the more accurate. Importantly, Iwasawa allowed the defect -containing supercell to relax under constant pressure (approximating the high defect concentration limit), while Gupta et al performed a constant volume calculation that approximates the low concentration limit. In this work, constant volume calculations have been performed. Considering the similarities between the approach adopted in this work and the approach of Gupta et al. it is not surprising that the defect formation energies are in reasonable agreement. Table 3 3 further indicates that the choice of the reference state for the defect calculations is important. Freyss et al. 72 used a reference stat e of uranium for uranium defects and an oxygen molecule for oxygen defects, where in each case the reference state chemical potential was calculated directly with DFT. This approach is valid when one is not considering distributed Schottky defects or the effect of temperature and pressure on defect formation energies is not being taken into account. In this work, however, the goal was to investigate the effect of
65 temperature and partial pressure on defect energies. Therefore a slightly different approach was adopted embodied in Equations 3 1 and 32. Thus, when uranium is the reference state, its chemical potential was calculated with DFT. The chemical potential of oxygen ( O) was then determined using Equation 3 2, through which the effect of temperature and oxygen partial pressure on defect formation energies was considered. The justification for this choice is that temperature and atmospheric pressure have a greater influence on the chemical potential of gaseous oxygen than on solid uranium. However, when an oxygen molecule was used as the reference state, its chemical potential calculated with DFT, and obtain the chemical potential of uranium ( U) from Equation 3 2, somewhat different values for defect formation energies, were obtained as indicated in Table 3 3. Although this choice of reference state chang es the values of the formation energies, it does not affect the qualitative picture that the formation energy for neutral oxygen interstitials is negative, the formation energy for oxygen vacancies is large, and those for uranium defects are even larger. The individual point defect formation energies in Table 3 3 can be combined to explore the relative stability of neutral defect complexes in UO2, which can be compared more directly with experiment. Such a calculation corresponds to all of the point defe cts being dissociated from one another and at a sufficient distance apart that they do not interact; this thus ignores association energies. The experimental data is also based on the assumption that point defects do not interact. The experimental values i n Table 3 4 are reported as a range, indicative of the difficulty associated with measuring absolute values. The calculated values are in good agreement with experiment for the anti -Frenkel pair and are quite reasonable for the Schottky defect. However, th e calculated values for the formation energy of the Frenkel defect are significantly larger than the experimental value. This might be attributable to the fact that the Frenkel formation energy
66 was determined indirectly in experiment from oxygen diffusion data; indeed it has be en suggested that the experimentally measured value is an underestimate 55 A similar calculated and experimental trend is observed in other fluorites such as PuO2, ThO2 and CaF2 where the anti Frenkel defect is found to dominate (see Table 3 4). The c alculations were carried out using empirical potentials and oxygen diffusion experiments were carried out for the anti -Frenkel formation energies. The data cited in the table is limited to the dominant defects in each system. The actual values for the anti Frenkel defect are also rather similar to that for UO2. Effect of Charged States Point defects in metal oxides are usually thought of as being charged rather than neutral. For example, when an oxygen vacancy is created it is considered to carry a charge of +2 in traditional Krger -Vink notation, which represents the loss of two electrons. However, most calculations of the formation energies of point defects, such as those discussed previously, treat them as neutral. The next step was to explore the effect s of charge on the formation energy of the defect. This was done by explicitly applying charge to the supercell, which is mostly localized when using the SP GGA+U method. In conventional DFT, the charge is spread around in a manner representative of deloca lization. The influence of charge on the defect formation energy was taken into account in two places in Equation 3 1. The first is through the DFT calculation of the charged, defect containing supercell where explicitly electrons are added or taken out of the system and the other is through the electron chemical potential term, which effectively shifts the Fermi energy. For example, for an n type defect, which donates electrons to the system, the Fermi level will be near the conduction band, while for a p type defect it will be near the valence band. The most stable defect is that for which the defect energy is lowest for a specific Fermi level.
67 This is illustrated in Figure 3 5, which compares the energies associated with forming V O x, V O and V O oxygen vacancies. A transition was observed as the Fermi level is increased from the valence band (0.0 eV) to the conduction band (1.92 eV), spanning the predicted band gap of the material Close to the v alence band, t he +2 charged oxygen vacancy is favored. This i n turn means that oxygen vacancies have a tendency to donate electrons or behave as an ntype defect. With an increase in Fermi level up to 0.6 eV, the +1 charged oxygen vacancy becomes energetically favorable. As the Fermi level moves towards the conduction band the tendency of the oxygen vacancy to donate electrons dimini shes and the neutral oxygen vacancy becomes dominant as th e Fermi level approaches the con duction band. The energies of all the individual point defects considered are shown in Figure 3 6 as a function of the position of the Fermi energy. When considering only the vacancy defects, the +2 charged oxygen vacancy is predicted to be the most stable defect near the valence band. Thi s in turn predicts that even in the presence of a uranium vacancy, which is a p -type defect (accepting electrons from the system), the system will still donate electrons. However, a s the Fermi level approach es the conduction band, the 4 charged uranium va cancy becomes increasingly favorable and remains the dominant defect over the range from 0.3 to 2.0 eV. When considering only differently charged uranium vacancies, the 4 charged vacancy is the most favorable The formation energy of this charged vacancy is significantly lower than that of the neutral vacancy. The effect of charge on the oxygen interstitial was also examined (see Figure 3 6). The effect of charge on uranium interstitials was not taken into account because the formation energy of the Frenke l pair is so high. Finally, it was observed that the 2 charged oxygen interstitial (Oi ) dominated over the entire range of the Fermi energies considered, but increased in stability as the Fermi energy approaches the conduction band. Using the same for malism for chemical
68 potentials the effect of stoichiometry can also be simulated by choosing different reference states. Choosing uranium as the reference state (Fig ure 3 6 a)) would correspond to a hypostoichiometric (UO2x) case while choosing an oxygen molecule as a reference (Fig ure 3 6 b )) corresponds to a hyperstoichiometric (UO2+x) case. The former reference state would shift the Fermi level towards the conduction band due to the presence of oxygen vacancies and the latter would shif t it towards t he valence band. It was observed that t he changes in environmental conditions do not affect the general conclusions. However, there are changes in quantative trends when comparing both cases. The formation energy of an oxygen vacancy is significantly lower in the hypostoichiometric case as it is easier to form an oxygen vacancy in an oxygen deficient scenario. On the other hand, it becomes easier to form an oxygen interstitial and a uranium vacancy when an oxygen molecule is considered to be the reference, as both these defects contribute to the hyperstoichiometry in the system. This analysis of the relative stabilities and behavior of individual point defects is informative. However, point defects do not occur in isolation but rather in combinations, such a s Frenkel, anti -Frenkel and Schottky complexes. 55 Therefore the energies associated with complexes that consist of either a combination of neutral or charged defects were next considered. B ecause the complexes are charge neutral, the position of the Fermi level does not enter the calculation s, even though the energy of the individual charged defects varies with the Fermi level. The formation en ergy of the anti -Frenkel complex is shown in Figure 3 7 It is clear that the charges on the individual defects that make up the complex influence its overall formation energy, an effect previously predicted for TiO2. 78 Specifically, the combination of an oxygen anti -Frenkel pair of charged defects ( OV + O i) has a lower formation energy than the corresponding combination of neutral point defects ( X OV + Ox). However, a s previ ously
69 discussed the calculated formation energy for Frenkel complex made from neutral defects agrees reasonably well with experimental measurements. Thus t he combination of charged components yields a formation energy that is considerably lower than the e xperimental value. Probable reasons are the neglect of charge transfer in experiment and the distribution of charge in DFT A similar trend is seen for the Schottky defects, which are illustrated in Figure 3 8. The formation energy of the Schottky defect a s a combination of charged vacancies ( U V+ 2 OV ) is significantly lower ( 50 %) than the combination of neutral vacancies. As for the anti -Frenkel, the calculated energy of the Schottky defect of neutrals components a lso agreed reasonably well with experiment. The difference in energy between the anti -Frenkel and Schottky defects is predicted to vary considerably with charge. This is because the uranium vacancy is affected to a greater degree than the oxygen point def ects (see Figure 3 6 ) when charges associated with the constituent point defects are considered. This phenomenon can be seen as a consequence of the strong electron correlation of the uranium atom, in addition to the factors listed above for the charged an ti Frenkel defect complex. Effect of Temperature and Oxygen Partial Pressure The typical operating conditions of a nuclear reactor span a variety of temperatures and oxygen partial pressures. The influence of temperature on the defect formation energies is discussed first in Figure 3 9. The temperature range considered is from 300 to 1400 K at standard atmospheric pressure of 1 atm. This range was specifically chosen after examining the phase diagram (see Figure 1 2) in order to observe the effect of temper ature on the stoichiometry of the fuel. 8 In general, UO2 is hyperstoichiometric but this tendency decreases as the
70 temperatur e increases; hypostoichiomtery exists in measurable quantities only at high temperatures (~1600 K). The oxygen interstitial is still predicted to be the dominant defect present. However, the formation energy of this defect increases as the temperature increases. This observation in turn means that the formation energy of an oxygen vacancy should decrease with temperature, which is consistent with the narrowing of the hyperstoichiometric range. This trend is indeed seen in the calculations, where the oxygen point defects are more energetically favorable than their uranium counterparts. The calculations also indicate that the uranium defects are affected to a greater degree by the phase transitions of bulk uranium with temperature. In particular, the formatio n energy of the uranium vacancy decreases while that of the uranium interstitial increases as temperature increases. The other significant issue is the equilibrium partial pressure of oxygen within the fuel which, to a substantial degree, determines wheth er the fuel will oxidize the metallic cladding surrounding it 82 and contributes to the non -stoichiometry observed in the fuel. This leads to an undesirable thinning of the cladding. In Figure 3 10, the effect of oxygen partial pressure on the defect f ormation energies at 800 K is shown. This temperature was specifically chosen because it is the rim temperature of a typical UO2 pellet in an LWR. Additionally, many experiments have been carried out in the partial pressure range used here, which allows th ese predictions to be readily correlated with experimental results.8 This thermodynamic analysis predicts the oxygen interstitial to be the most stable point defect. However, it becomes increasingly difficult to form this defect as the partial pressure is reduced. This is expected as at highly reduced conditions there is little free oxygen in the system and the tendency to lose oxygen to the atmosphere is significant. On the other hand, it becomes
71 energetically feasible to form oxygen vacancies in reducing conditions. This contrasting nature of oxygen can be understood on the basis of the Brouwers diagram for UO2 in Figure 3 11 where the concentration of defects is plo tted as a function of oxygen partial pressure. 83 As indicated in the figure, the concentration of oxygen vacancies is proportional to pO2 1/2 which accounts for their being energetically unfavorable at high oxygen partial pressures. The decrease in oxygen vacancy concentration increases the uranium vacancy concentration. Furthermore, the formation energy of the cation vacancy drops dramatically as the oxygen partial pressure increases. This is due to the fact that Schottky equilibrium 11 or charge neutrality through vacancies must be maintained in the system
72 Table 3 1 T he calculated properties of bulk uranium oxide are shown using different approximations. These values are compared to published experimental and theoretical data. Method Lattice parameter (nm) Band Gap (eV) Magnetic moment ( B) Experim ent 5, 48 0.547 2.00 1.74 Freyss 7 2 GGA 0.540 0.00 0.00 D udarev 54 LMTO LSDA+U 0.537 2.10 1.70 Gupta 7 5 PAW SP GGA+U 0.552 1.80 1.94 This Work LDA GGA LSDA SP GGA LSDA+U PAW SP GGA+U 0.526 0.533 0.530 0.542 0.543 0.549 0.00 0.00 0.00 0.00 1.68 1.92 0.00 0.00 1.93
73 Table 3 2 Relative stabilities of the three different phases of bulk uranium metal, the orthorhombic -phase, the body centered tetragonal phase, and the body centered cubic -phase. Phase Temperature range (K) U (eV/atom) -U up to 941 8.08 -U 941 1049 7.93 U 1049 1408 7.57
74 Table 3 3 Calculated formation energies of neutral poin t defects compared to reported work using various approximations Method Oxygen Vacancy (eV) Oxygen Interstitial (eV) Uranium Vacancy (eV) Uranium Interstitial (eV) System size Freyss 7 2 GGA 6.1 2.5 4.8 7.0 2x1x1 Crocombette 71 LDA 6.7 2.9 3.3 7.3 2x1x1 Gupta 7 5 PAW SP GGA+U 5.6 1.6 6.0 8.2 2x2x2 Iwasawa 7 4 PBE SP GGA+U 4.46 0.44 8.45 4.70 2x2x2 This Work PAW SP GGA+U Ref. state: -U Ref. state: O 2 4.32 5.29 0.37 1.34 8.96 7.04 6.12 8.04 2x2x2
7 5 Table 3 4 The values for Frenkel energies are given for different compounds having the fluorite structure. The comparison indicates that the anti -Frenkel is the dominant defect for a ll compounds with the fluorite structure and is not limited to actinide oxides. Method Anti Frenkel (eV) Frenkel (eV) Schottky (eV) Melting Temperature (K) Experiment 6 4 3.0 4.6 8.5 9.6 6.0 7.0 Freyss 7 2 GGA 3.6 11.8 5.6 3125 Crocombette 7 1 LDA 3.8 10.6 5.8 Gupta 7 5 PAW SP GGA+U 4.0 14.2 7.2 This Work PAW SP GGA+U 3.95 15.08 7.6 CaF 2 37 Oxygen diffusion 2.7 1691 PuO 2 84 Oxygen diffusion 2.7 2.9 2740 ThO 2 84 Oxygen diffusion 2.3 4.7 3663
76 Figure 3 1 Summ ary of the experimental oxygen diffusion and metal diffusion data in fluorite type oxides, showing the temperature scale normalized to the melting point of the solid. [Ada pted from reference 6 4 (Figure 1).]
77 Figure 3 2 Illustration of the 96 atom su percell used in the present work.
78 Figure 3 3 Total DOS projected to the U 5f and 6d orbitals as the conduction band and to the U 5f and O 2p as the valence band for antiferromagnetic UO2 a) SP GGA, and b) SP GGA+U method. (a) (b) 0 10 20 30 40 50 60 70 80 -7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 8 Energy (eV) UO2: SP-GGA Fermi energy = 0 eV Band Gap = 0 eV 0 5 10 15 20 25 30 35 40 45 50 55 60 -7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 Energy (eV) UO2:SPGGA+U Band Gap = 1.92 eV Experimental = 2eV
79 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 1 2 3 4 5 System Size Normalized Formation energy (eV) Frenkel Anti-Frenkel Schottky 1x1x1 2x2x2 3x3x3 5x5x5 4x4x4 Figure 3 4 Normalized formation energy of neutral defect complexes calculated as function of system size by atomic level simulations using an empirical potential within GULP. The energies were normalized with respect to the 5x5x5 supercell.
80 3 3.5 4 4.5 5 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Fermi Level (eV) Defect Formation energy (eV) Figure 3 5 Stability transitions for the neutral, +1 and +2 charged oxygen vacancies as a function of Fermi level over the entire band gap. Only the most dominant defect with reference to a particular Fermi level is shown for clarity. O V OV OV
81 Figu re 3 6 Calculated defect formation energies for variously charged oxygen and uranium point defects as function of reference state a) uranium and b) oxygen molecule. The empty symbols denote uranium while the filled symbols are oxygen. -6 -4 -2 0 2 4 6 8 10 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Fermi Level(eV) Defect Formation Energy (eV) UV OV OV OV '' UV Oi X '''' UV '' iO (a) -8 -6 -4 -2 0 2 4 6 8 10 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Fermi Level(eV) Defect Formation Energy (eV) UV OV OV OV '' UV Oi X '''' UV '' iO (b)
82 -6 -4 -2 0 2 4 6 8 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Fermi Level (eV) Defect Formation energy (eV) Fig 3 7 Calcu lated formation energy of the anti -Frenkel defect as a function of the Fermi level; the formation energies of the neutral complexes are compared to charged individual defects. Anti -Frenkel ( OV + '' iO ) Anti -Frenkel ( X OV + X iO ) O V '' iO
83 0 2 4 6 8 10 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Fermi Level (eV) Defect Formation energy (eV) Fig 3 8 Calculated formation energy of the Schottky defect as a f unction of the Fermi level; the formation energies of the neutral complexes are compared to charged individual defects. Schottky ( VU' + 2 OV ) Schottky ( X UV +2 X OV ) '' UV OV
84 -2 0 2 4 6 8 10 300 400 500 600 700 800 900 1000 1100 1200 1300 1400 Temperature (K) Formation Energy (eV) O Vacancy O Interstitial U Vacancy U Interstitial pO2=1 atm Figure 3 9 Predicted effect of temperature on the formation of neutral point defects at standard atmospheric pressur e for the temperature range of 300 to 1400 K. The vertical dashed lines indicate a phase change of bulk uranium metal
85 -2 0 2 4 6 8 10 12 -23 -21 -19 -17 -15 -13 -11 -9 -7 -5 -3 -1 1 Oxygen Partial Pressure, ln (P/Po) Formation Energy (eV) O Vacancy O Interstitial U vacancy U interstitial 800 K Fig. 3 10 Predicted effect of oxygen partial pressure on the formation of neutral point defects in UO2 at 800 K; the partial pressure of oxygen is varied from pO2 = 1 x 1010 atm. to 1 atm. This corresponds to highly reducing and oxidizing conditions.
86 Figure 3 11 Illustration of the Brouwer diagram for UO2+x, the diagram assumes that the intrinsic electronic reaction dominates and ne glects defects on the uranium sublattice. [Ada pted from reference 83 (Figure 1).]
87 CHAPTER 4 THERMODYNAMICS OF FI SSION PRODUCTS IN UO2 X I ntroduction Fission products (FP) are produced during the burning of UO2.16 Although the specific amount of a particular fission product generated depends on the compositi on of the fuel and the type of reactor, the general trends are similar, see Figure 1 5 for an illustration. Since each FP has very different physical and chemical interactions with the fuel matrix, an assessment of its impact on the structural evolution of the fuel cannot be based on yield alone.17 As discussed in Chapter 1, some fission products do not react with the fuel matrix, some are as found as metallic precipitates, and others react with the fuel matrix to form separate phases. Therefore, in order to be able to pr edict the system behavior at a macroscopic level, it becomes very important to understand the behavior of fission products at the microscopic level. Fission gases xenon (Xe) and krypton (Kr) are essentially insoluble in the UO2 matrix and are present onl y by virtue of the fission process. Nevertheless, they can migrate to grain boundaries where they form bubbles. The migration of fission gas atoms has been studied theoretically by Catlow84 and later by Grimes.85 They have proposed different mechanisms for migration which depend on the stoichiometry of the fuel owing to the different trap sites at different stoichiometries. This bubble for mation can lead to considerable swelling of the fuel and severely degrades mechanical properties.18 These gases can also escape to the plenum region between fuel rod and the cladding, which can contribute to internal stresses on the cladding.19 Solid -state fission products, such as cesium (Cs), can assume various roles depending on conditions prevailing in the matrix. They are found as a minor component of the grey phase, which is a separate phase consisting of a number of fission products.21, 22 Cs in particular is found as a metallic inclusion in the fuel, though implantation studies have observed it to be
88 present on lattice sites in hyperstoichiom etric fuels.23 It can also react with iodine (I) and cause corrosion of cladding materials such as zircalloy and stainless steel.24 The fission product strontium (Sr) is generally considered to be soluble in UO2, 21, 87 although the extent of solubility critical ly depends on the oxygen to metal ratio. In the presence of another fission product zirconium (Zr), it tends to form a grey phase perovskite SrZrO3.26 Kleykamp et al .86 have studied the solubilities of selected fission products in UO2. In particular, they have considered the solubilities of Cs, Sr, Ba, Y, La, Ce, Pr, Nd, Zr, Nb, Mo and Te in UO2. Based on the solubility data, they have measured the change in lattice parameter of UO2 with change in the solid solution composition, as seen in Figure 4 1. It was found that for all the fission product oxides except LaO2, the lattice parameter decreased as the so lubility of the fission product increased. It is to be noted that solubility of Cs2O was negligible in stoichiometric while that of SrO was 12 mol % at 1773 K. These results were used as validation for the calculations performed in this work. There is a s ignificant need for an in -depth analysis of the interaction of fission products with the various defects present in UO2 from reliable calculations and simulations. Grimes and Catlow14 have performed a comprehensive study of the variety of fission products f ormed and developed a theoretical model that provides atomic level explanations to experimentally observed phenomena. In particular, they used the Mott -Littleton56, 57 approach and empirical potentials to calculate all the relevant energies. However, their work used f ixed charge empirical potentials to calculate defect energies, which neglect charge transfer effects that may be important in some situations. Petit and co -workers87 90 performed similar calcul ations using DFT. They used LDA and GGA for the calculations However, their calculations were limited to 24 atom supercells. Moreover, it has be en shown63 7 4, 7 5 that the localization of the 5f orbital of
89 uranium plays a dominant role in determin ing the correct electronic structure of UO2 and consequently on the defect formation energies, which conventional LDA and GGA fail to capture; this effect can be captured using the Hubbard +U on -site repulsion. Recently, Brillant et al.91 used the spin polarized SP GGA+U method and a 96 atom supercell to calculate the incorporation and solution energies of barium (Ba) and zirconium (Zr) and their respective oxides in UO2 x. It was found that the most stable inc orporation site for the fission product depended on the type of exchange correlation functional used and the effects of functional form were more pronounced for higher charged defects. A very recent paper by Gupta et al.92 discusses solution and migration of cesium in UO2. They observed the solution of cesium to be very low and the migration to be highly anisotropic. In this work, the SP GGA+U and non-spin polarized GGA method has been used to analyze the incorporation and solubility of Xe, Cs, and Sr in the UO2 matrix. These fission products have been selected due to their very different characteristics and interactions with the fuel matrix. In contrast to GGA, the SP GGA+U method correctly predicts the experimentally observed insulating behavior of UO2. 48 The DFT calculations were supported by empirical potential calculations to test the interaction of periodic images. This becomes important when considering large -sized fission products and extended defects such as the bound Schottky defect. Finally, the SP GGA+U results were also used in a point defect model that can predict the change in formation energies of point defects with variations in temperature and stoichiometry.
90 C omputation al M ethodology Electronic Structure C alculations The DFT calculations were performed with the projector augmented -wave (PAW)47 method. The SP G GA+U54 was utilized to include the effect of the strong correlation of 5f electrons in uranium and non spin polarized GGA as implemented in the Vienna Ab Initio Simulation Package (VASP).59, 60 A value of 3.96 eV Ueff (U J) was used in this work This value is similar to those previously used by others54, 7 4, 7 5 and also to experimental measurements.48 With regard to other calculation parameters, a 2x2x2 unit cell consisting of 96 atoms was utilized for the structural optimizations where the cell volume was kept constant and the atomic positions relaxed. This is important for simulating low concentration of defects as discussed in Chapter 3. The Brillouin zone sampling used a 2x2x2 Monkhorst -Pack k -point mesh;42 the cut -off energy for the plane waves was 400 eV. A strict force and energy convergence criteria was used for these calculations (see Nerikar et al.63 for details) Th e combination of the above parameters resulted in good agreement with experiments for bulk and defect properties (see Table I). All the fission products considered were charge neutral. Incorporation E nergy The incorporation energy is defined as the energy required taking a fission product from infinity and placing it at a pre -existing trap site. The incorporation energy of a fission product is defined by Grimes et al.14 as: i total total incE defect E E E ) ( ) (. (4 1) Here E total ( ) is the total energy of the cell with the fission product at a particular defect site, Etotal (defect) is the total energy of the cell with a particular defect, and i is the energy of a single isolated fission product. This energy does not account for t he formation of the trap site and
91 assumes that there is always an excess of available sites. This assumption is only valid at low burn ups, i.e., when the concentration of fission products is low This definition also does not account for changes in fuel s toichiometry. Nevertheless, it is a useful quantity that measures the energetics of fission products at different trap sites. A positive value of the incorporation energy means th at energy is required to place a fission product at a particular trap site. I n this work, the trap sites considered we re the empty octahedral interstitial site, oxygen vacancy, uranium vacancy, a neutral divacancy (consisting of a uranium and oxygen vacancy in close proximity) and finally a neutral trivacancy (bound Schottky defect ). See Figure 4 2 for illustrations of these trap sites. It is to be noted that there may be other trap sites depending on the stoichiometry which may be more energetically favorable than the ones considered. However, in equilibrium, these are the most lik ely trap sites without taking the effects of clustering of defects into account. Solution E nergy To compensate for limitations of incorporation energy, the solution energy of fission products has to be considered. T he solution energy is defined as:14 trap inc solutionE E E (4 2) where trap is the trap site formation energy, which is defined as the energy required to form a particular trap site for the incorporation of a fission product. Thus, the solution energy accounts for both the formation o f the trap site as well as the incorporation of the fission product into that trap site. This energy will be a function of stoichiometry and hence burnup. For example, it will be harder to form an oxygen vacancy in a hyperstoichiometric case (UO2+x) comp ared to a hypostoichiometric case (UO2x). One advantage of the solution energy compared to the incorporation energy is that it can be correlated to experimental trends. As with the
92 incorporation energy, this definition assumes no interaction, such as clus tering, between fission products. Experimentally, the concentration of fission products is determined by the fission reactions, whereas computationally the fission products are introduced manually. To evaluate the trap site formation energy as a function of stoichiometry, a point defect model (PDM) was implemented. This model helps relate the experimental and computational findings to one another. Using such a model, the apparent trap site formation energy can be calculated, effectively taking into account several defect formation reactions at once. As the relative importance of different reactions will change with temperature and stoichiometry, the effect of temperature and stoichiometry on the formation energy can also be determined. The apparent trap si te formation energy is obtained through Etrap kT ln([ X ]) (4 3) where [X] is the concentration of the type of defect X considered. As discussed by Brillant et al.,92 the evaluation of the appare nt defect formation energies can be expressed in the framework of a PDM, originally introduced by Matzke64 and Lidiard.94 In the PDM, the concentrations of point defects in UO2 x are governed by of the kinetic equations: ) exp( ] ][ [ kT E I VF FP O OO (4 4) ) exp( ] ][ [ kT E I VF FP U Uu (4 5) [ VO]2[ VU] exp( ES FkT ) (4 6) 2 [ VU] [ IO] 2 [ IU] 2 [ VO] x (4 7)
93 Here, [VO], [IO] [VU] and [IU] are the concentrations of the oxygen vacancy, oxygen interstitial, uranium vacancy and uranium interstitial respectively; EF FPO, EF FPU and EF S are the formation energies of the anti -Frenkel, Frenkel and Schottky complexes. Equation 4 7 is the electroneutrality equation that also accounts for non -stoichiometry. In addition, the concentrations of divacancy and bound Schottky defects are given by [ DV ] [ VO][ VU] exp( BDVkT ) (4 8) [ Sch ] [ VO]2[ VU] exp( BSchkT ) (4 9) Here, BDV and BSch are the binding energies of divacancy and Schottky defects. To solve Equations 4 4 to 47, an assumption has to be made that the uranium interstitial concentration is negligible compared to other defect concentrations, an assumption which is experimentally valid and is also verified in these calculations. With this assumption, Equations 4 4 to 4 7 lead to a cubic equation for [VO] or [IO]. By solving the cubic equations using standard techniques with complex solutions, the concentration of oxygen vacancies or oxygen interstitials and, in turn, their formation energies as a function of stoichiometry are obtained. The concentrations of other defects ar e then determined with Equations 4 4 to 4 -9. Oxi de Solution E nergy Solid fission products such as cesium and strontium can react with oxygen in the fuel, and form secondary oxide phases.16 Any oxide (for example, strontium oxide) before solution into the UO2 matrix will decompose into the fission product and oxygen by the following reaction: ) ( ) ( g O Sr s SrO (4 10) Thus, although the calculation of the solution energy is useful in predicting the most stable trap site as a function of stoichiometry, its definition is limited to isolated fission products. F or
94 example, it does not provide information about whether the fission product would remain in the UO2 matrix or whether it would react with oxygen to form a stable second phase. Therefore, a third definition of energy is required; the fission product oxide solution energy is given as: ESrO solution ESr solution EO solution ESrO formation (4 11) The first term on the right hand side of the equation is the solution energy of the particular fission product at the most stable trap site based on stoichiometry; the second term is the soluti on energy of the oxygen from the fission product oxide into any vacant oxygen site in the UO2 matrix. This vacant site also depends on stoichiometry. In UO2x, the oxygen will be soluble in an oxygen vacancy site, while in UO2+x it will be present as an ox ygen interstitial site, and in UO2 the most stable solution site will be a mixture of vacancy and interstitial sites. The third term is the formation energy of the oxide. If the fission product oxide solution energy is negative, it means that it is energet ically favorable for the oxide to be soluble in the fuel. This in turn means that the fission product will not form a stable second phase. This definition neglects contributions from interfacial energy that would be required to create a stable second phase Results and D iscussion Validation of Approach The first step in validating this approach was to calculate equilibrium and defect properties for UO2. Table 4 1 summarizes the earlier63 comparison of the lattice parameter and formation energies of defect complexes in UO2 calculated using SP GGA+U and GGA with published experimental results. It was observed that the agreement with experiment improves if SP GGA+U method is used. As expected, the GGA method underestimates the lattice parameter of UO2 slightly. The calculated defect energies using SP GGA+U are in good agreement with experiment for the case of the anti -Frenkel defect and agree reasonably well for the unbound
95 Schottky defect as compared to experimental data. Howev er, the SP GGA+U method predicts a higher Frenkel formation energy than experiment as discussed in Chapter 3. This value has been suggested to be an underestimate in experimental results.55 Regardless of whether it is an underestimate or not, the Frenkel defect has very high formation energy and is ignored in the further calculations. The focus of the DFT results was on the SP GGA+U results, a further comparison to GGA and published empirical potential data was made to understand the effect of different levels of theory in predicting the experimentally observed trends. T he DFT calculations we re all performed in a 2x2x2 supercell size with 96 atoms The consideration of limitations of system size is important, as there could be interactions between periodic images of the defect in neighboring supercells especially for the case where a fission product is placed in a vacant site or the cell is not charge neutral. To estimate the effects of the limited supercell size used, fission product incorporation calculations were carried out using an empirical potential within the Gene ral Utility Lattice Program (GULP).6 1, 62 The advantage of this approach is that larger supercell sizes can be treated with negligible computational cost. Because of the intrinsi c limitations of empirical potential approaches, t he se values can be expected to be less accurate than those obtained with DFT However, the effects of supercell size which arise primarily from electrostatic and elastic interactions, are not expected to c hange significantly with method.63 T he Grimes potential14 was chosen for these GULP calculations which has the advantage of being parameterized for a variety of fission product interactions with the UO2 host. T he effect of system size was analyzed for Cs incorporation (see Figure 4 3 ) i nto the empty octahedral interstitial site, oxygen vacancy site, and uranium vacancy site. As expected, it was found that th e incorporation energy converged as the system size increased A s imilar trend was observed for defect complexes as discussed in Chapter 3. For the
96 octahedral and the oxygen vacancy sites, the energy is almost converged at a system size of 2x2x2. For the uranium vacancy site, the energy changes by about 16 % as the syste m size is increased from 2x2x2 to 3x3x3. This kind of analysis provides an indication that there are no strong interactions between neighboring periodic images and the supercell considered in the DFT calculations is sufficient ly large The third step in t he proposed approach was to validate the point defect model that was implemented Grimes and Catlow14 develop ed a series of equations (hereafter referred to as the Grimes model) to calculate the trap site formation energy ( trap) based on the anti -Frenkel and Schottky defect formation energies This model provides a convenient way of estimating the solution energy of fission products as a function of stoichiometry. However, since it considers only the most dominant reaction for the formation of a given trap site, it does not account for temperature variations. The approach adopted in this work using Equations 44 to 4 7 compares the trap site formation energy (Etrap) calculated using both the Grimes model and the PDM and anal yze s the scenarios by which these simplified expressions are valid to determine when a more complex assessment is required. The formation energy of a uranium vacancy is plotted as a function of stoichiometry in Figure 4 4 at 300 K, 1000 K and 2000 K. The stoichiometry is varied from UO1.90 to UO2.1. The temperature effect is included through the PDM described by Equations 4 4 to 4 9, while the results using the Grimes model correspond to 0 K. The inputs (defect formation energies) for both models are from SP -GGA+U calculations. The results of the DFT calculations confirm the conclusion of the Grimes model on the low temperature defect formation reactions, where there is good agreement. The agreement gets worse as the temperature is increased. The reason for this trend is that at lower temperature s only the anti -Frenkel mechanism is dominant which is similar
97 to the prediction of the Grimes model. However, as the temperature increases, other mechanisms such as the Schottky defect formation become important. T hus, at lower temperatures, it is valid to use the Grimes model to calculate the trap site formation energy but at higher temperatures such as those found in a light water reactor more complicated expressions have to be considered This effect is even m ore pronounced when the deviations from stoichiometry are small. Incorporation E nergy The incorporation energies of Xe, Cs and Sr calculated using DFT are reported in Table 4 2 The lowest value for each fission product is shown in bold. First considering the inert fission gas (Xe) the bound Schottky defect is predicted to be the most stable incorporation site. Moreover, t he trend can be related to the relative size of the incorporation sites. The bound Schottky defect is the largest of the defects consi dered here, even though it has a high formation energy. However, the incorporation ene rgy, by definition does not account for the formation of the defect. The SP GGA+U method predicts the uranium vacancy to be the most stable incorporation site for cesiu m. This can be understood by considering Coulombic effects. The Cs+ ion would prefer to reside on a uranium vacancy (U4+) since the Coulombic forces are stronger. The trend is in agreement with the calculations of Gupta et al.92 However, the incorporation energy values are higher even though the calculation method is similar. This can be attributed to the differences in o rdering of magnetic moments of the uranium ions. Differences in ordering can lead to different magnetic ground states and thus different energies. For similar reasons, t he uranium vacancy is predicted to be the most stable incorporation site for a Sr2+ ion Solution E nergy The calculated solution energies are presented in Table 4 3; for each stoichiometry, the lowest solution energy is shown in bold. The results in Table 4 3 are also illustrated in Figures 4 -
98 5 to 4 7. The apparent trap site formation energ y values used are those found with the PDM at 300 K. The results can be understood based on the relationship between the incorporation energy and the apparent formation energy of a particular defect as a function of stoichiometry. The first point to notice is that the solution energy of xenon is independent of stoichiometry for the interstitial and the bound Schottky site. The interstitial site is empty and therefore no defects have to be created to incorporate xenon into that site. The apparent formation of the Schottky site is defined as the energy difference between the actual Schottky formation energy (when the constituent defects are far apart) and the binding energy of the neutral trivacancy (when the defects are in close proximity). For UO2x, the bou nd Schottky site is predicted to be the most stable solution site. This is true even though it is easier to form an oxygen vacancy than a bound Schottky defect, since the incorporation of xenon in the latter is energetically favorable compared to the former. This is due to the larger size of the bound Schottky defect. In UO2, the divacancy formation energy is less than that of bound Schottky defect since it is more difficult to form oxygen vacancies in UO2 compared to UO2x. Hence there is a competition be tween the most stable solution sites. It becomes even more difficult to form oxygen vacancies and easier to form uranium vacancies in UO2+x. Hence, xenon is predicted to occupy the single cation site. As previously mentioned, solution energy results can be correlated directly with experimental trends. The calculated solution energy of xenon is positive and high for all stoichiometries. Therefore, it is insoluble in the UO2 matrix in accordance with experimental findings18 and is only found as a consequence of th e fission process. The results for fission gases can be understood as a balance between the size of the trap and the energy to form the trap. However, for solid fission products such as cesium and strontium, the effect of charge must also be taken into ac count. In UO2+x, the uranium vacancy
99 has the lowest formation energy and hence solution of cesium is favored at this site. In UO2, the divacancy has a larger size and slightly overcomes the more highly charged uranium vacancy as the most stable solution site even though both have similar formation energies. In hypostoichiometric conditions, due to the availability of oxygen vacancies the most stable site is the bound Schottky defect. In agreement with an implantation study by Matzke,24 the calculations predict that cesium will occupy a lattice site (uranium vacancy in this case) in hyperstoichiometric fuels. The divacancy trap site is predicted to be more stable for the solut ion of strontium in sub stoichiometric (UO2 x) environments relative to the bound Schottky defect. However, as the matrix becomes hyperstoichiometric, the formation energy of a uranium vacancy drops significantly and it becomes the dominant solution site f or UO2+x. A similar trend has been observed for barium by Brilliant et al.91 This is expected as both elements belong to the same group in the periodic table and the calculation techniques are quite similar. S trontium solubility increases in the hyperstoichiometric regime as observed experimentally.86 It must be noted, however, that the solution energies have been calculated at 300 K and the most st able site can change with temperature. This is observed for the case of Cs solution in the hypostoichiometric case where the solution site changes from the Schottky site to the divacancy site with increase in temperature. This can have important ramificati ons on the prediction of fuel behavior as different parts of the fuel are at different temperatures. Oxide Solution E nergy The oxygen solution energies are reported in Table 44. A positive value corresponds to the oxide being insoluble in the fuel and f orming a second stable phase. It was found that oxides with higher oxygen-to -metal ratio were more soluble in the fuel. This can be understood on the basis of their ability to donate oxygen to the lattice and the energy cost associated with charge
100 compensa tion when a uranium atom is replaced with a fission product. For similar charged ions, there will be no energy cost for charge compensation and hence these will be more stable compared to lesser charged ions. The solubility also increases as the hyperstoic hiometry of the UO2 increases. Considering individual oxides, Cs2O is predicted to be stable at all stoichiometries, though in UO2+x, it has borderline stability. With regard to UO2+x, it is observed experimentally21 that the solubility of Cs increases as we go from hypostoichiometric to hyperstoichiome tric fuels, in agreement with the calculations performed in this work. This can again be related to the argument that the more oxygen an oxide donates to the fuel matrix, the more soluble it is. 57 The maximum solubility of SrO in UO2 was estimated to be 12 mol% at 1773 K by Kleykamp et al.86 and the solubility has to decrease with decreasing temperature as some degree of that solubility has to be driven by entropy. However, these calculations pr edict SrO to be soluble in UO2 for all the three stoichiometries at the dilute limit. If no Sr -Sr interaction is assumed, this implies complete solubility of Sr, which is clearly at odds with the experimental result. To resolve this contradiction, calculations were performed where two Sr atoms are placed 3.88 apart (the closest cation-cation distance in a fluorite structure) in the same supercell. The calculated interaction energy is 1.35 eV, which suggests strong Sr -Sr repulsion and is of the same order of magnitude as the oxide solution energy. This result strongly suggests that there is a limit to the solubility of Sr in UO2. The results in Table 4, valid only for a dilute limit of strontium, predict Sr is soluble at the 3 mol% level. However, as more S r is added, the Sr ions would begin to repel, and the solubility would peak. Assuming that this repulsion is only for nearest -neighbor Sr ions, this would predict a maximum solubility of 25%, which clearly indicates that the repulsion is even longer ranged The calculations neglect entropic contributions
101 but take charge transfer into account (four U5+ species to compensate for the two strontium atoms; these species are separated by at least 6 to minimize repulsion). In addition, the solubility of Sr was observed to be much higher in UO2 than in UO2x by Kleykamp et al., which mirrors the calculated trend. This can be explained by the fact that the solution of strontium becomes increasingly favorable as one approaches hyperstoichiometric conditions owing to the strong Coulombic interactions. The higher charged uranium vacancy (4+) is the most stable solution site in UO2 as compared to the divacancy (2+) site in UO2 x and hence has a strong preference for charged cations. The calculations agree well with Ba O solubility calculations of Brillant et al.91 and experimental observations of Thomas et al.,98 in which they found that barium does not precipitate in low burn up hypostoichiometric fuels and mid burn up hyperstoichiometric fuels. This is again expected since barium and strontium are quite similar in their chemical and physical make up. Comparison with o ther Theoretical Approaches Theoretical progress has led, and continues to lead, to DFT formalisms of ever increasing materials fidelity. It is thus valuable to assess the applicability of various formalisms to this problem. Thus, in this section the results from GGA and SP GGA+U ar e compared. The results of the DFT calculations are also compared with the results obtained with empirical potentials. Non -Spin P olarized GGA The incorporation energies of Xe, Cs and Sr calculated us ing GGA are reported in Table 4 2 and compared with the S P GGA+U results. When considering the incorporation of Xe, there is qualitative and quantitative agreement between both approaches. Both methods predict the bound Schottky site to be the most stable. This is not unexpected, as the stability of this particu lar fission product is a direct result of the relative sizes of the incorporation sites.
102 When considering the incorporation of cesium, there is a discrepancy between SP GGA+U and GGA. GGA predicts the most stable site to be the bound Schottky defect while SP GGA+U predicts it to be the uranium vacancy site. Since GGA does not predict an ionic or insulating ground state for UO2, it only partially captures the relevant Coulombic effects. These effects become even stronger for Sr and hence both methods predic t the uranium vacancy to the most stable site, although the stability of this site with GGA is only marginal. The solution energies of these defects are reported in Table 4 3 as a function of stoichiometry and the most stable site is again shown in bold. F or Xe, GGA and SP GGA+U predict the bound Schottky to be the most stable solution site for UO2x and UO2. However, for UO2+x, GGA predicts the divacancy site to be energetically favorable in contrast to SP GGA+U, which predicts the uranium vacancy site to be the most stable. As discussed previously, the solution energy is a balance of the incorporation and apparent formation energies. There is a significant difference in the incorporation energy of Xe at the uranium vacancy site between the two methods and this becomes dominant in the hyperstoichiometric case. For Cs, GGA does not predict any solubility for any stoichiometry owing to the incorrect treatment of the Coulombic effects, as opposed to SP GGA+U that finds Cs to be soluble in UO2 and UO2+x. For Sr there is qualitative agreement between the two methods. The solution of Sr becomes more energetically favorable as the matrix becomes hyperstoichiometric. Thus, in general, it was found that the non -spin polarized GGA results agreed with the SP GGA+U re sults when the Coulombic interactions are absent, as in the case of Xe, or when they are present and strong (Sr). However, the strong repulsion for Sr -Sr interaction is not observed for GGA, unlike in the case for SP GGA+U. This can be seen in Table 4 5. H ere the binding energy is shown as a function of the Sr -Sr nearest neighbor distance. As can be observed, the
103 maximum repulsion predicted is 0.23 eV, which is almost an order of magnitude less than that predicted by SP GGA+U. Even though GGA predicts the c orrect trends, a higher approximation is required if experimental trends are to be observed. Empirical Potentials The incorporation energies for Xe, Cs and Sr from reported empirical potential calculations14 are also compared with the calculated SPGGA+U r esults in Table 4 2 There is an overall qualitative agreement in trends between the two methods. However, the absolute values of incorporation energies predicted by empirical potentials are significantly higher, especially for higher charged ions. A proba ble cause for this could be the use of formal charges in the empirical potentials, which tends to lead to overestimations of binding energies.55 With regard to the solution energies ( see Table 4 3 ), the trends are similar for both methods for all stoichiometries. Both methods predict that the solution of all the three fi ssion products considered becomes energetically favorable as the matrix becomes hyperstoichiometric. In particular, they agree on the most stable solution sites except for Cs in the stoichiometric case.
104 Table 4 1 Comparison of equilibrium and defect pr operties of UO2 from GGA, SP GGA+U and experiment. Method Lattice Parameter () Anti Frenkel (eV) Frenkel (eV) Schottky (eV) Experiment 5.47 5 3.0 4.6 6 4 8.5 9.6 6 4 6.0 7.0 6 4 GGA 5.35 3.77 9.09 5.02 SP GGA+U 5.49 6 3 3.95 6 3 15.08 6 3 7.6 6 3
105 Table 4 2 Incorporation energies of xenon, cesium and strontium calculated using SP GGA+U and GGA and compared to previous empirical potential results. For each of the fission product, the lowe st incorporation energy is shown in bold. SP GGA+U GGA Grimes 14 Xe Interstitial 11.11 12.8 17.23 Oxygen vacancy 9.5 9.71 13.34 Uranium vacancy 2.5 6.04 4.99 Divacancy 2.45 3.29 2.84 Schottky 1.38 2.12 1.16 Cs Interstitial 10 10.1 9. 93 Oxygen vacancy 8.4 8.1 9.1 Uranium vacancy -3.4 0.75 -6.08 Divacancy 1.99 0.23 5.63 Schottky 0.8 -0.38 5.47 Sr Interstitial 4.68 4.3 11.04 Oxygen vacancy 7.18 5.3 8.87 Uranium vacancy 9.66 5.4 27.09 Divacancy 7.53 4.97 25.31 Schottky 4.55 4.74 23.36
106 Table 4 3 Solution energies of xenon, cesium and strontium calculated using SP -GGA+U and GGA and compared to previous empirical potential results. For each of the fission product and stoichiometry, the lowest sol ution energy is shown in bold. SP GGA+U GGA Grimes 14 Trap Site UO 2 x (UO 1.9 ) UO 2 UO 2+x (UO 2.1 ) UO 2 x (UO 1.9 ) UO 2 UO 2+x (UO 2.1 ) UO 2 x UO 2 UO 2+x Xe Interstitial 11.10 11.10 11.10 12.75 12.75 12.75 17.23 17.23 17.23 Oxyge n vacancy 9.57 11.48 13.26 9.79 11.69 13.47 13.34 16.75 20.16 Uranium vacancy 9.94 6.13 2.57 13.48 9.67 6.12 18.32 11.50 4.68 Divacancy 6.30 4.40 2.62 7.14 5.24 3.46 12.93 9.52 6.11 Schottky 3.88 3.88 3.88 4.62 4.62 4.62 9.57 9.57 9.57 Cs Interstitial 10.00 10.00 10.00 10.10 10.10 10.10 9.93 9.93 9.93 Oxygen vacancy 8.48 10.38 12.16 8.18 10.08 11.86 9.10 12.5 15.92 Uranium vacancy 4.04 0.23 3.32 8.19 4.38 0.83 7.26 0.43 6.39 Divacancy 1.86 0.04 1.82 4.08 2.18 0.40 4.47 1.06 2. 35 Schottky 1.70 1.70 1.70 2.12 2.12 2.12 2.94 2.94 2.94 Sr Interstitial 4.68 4.68 4.68 4.30 4.30 4.30 11.04 11.04 11.04 Oxygen vacancy 7.26 9.16 10.94 5.38 7.28 9.06 8.87 5.46 2.05 Uranium vacancy 2.22 6.03 9.58 2.04 1.77 5.323 13.76 20.58 27.41 Divacancy 3.68 5.58 7.36 1.12 3.02 4.8 15.22 18.63 22.04 Schottky 2.05 2.05 2.05 2.24 2.24 2.24 14.95 14.95 14.95
107 Table 4 4 Solution energies of oxides of cesium and strontium calculated using SP GGA +U and GGA and compared to previous empirical potential results. A positive energy implies insolubility in the fuel matrix phase. SP GGA+U Grimes 14 Oxide UO 2 x (UO 1.9 ) UO 2 UO 2+x (UO 2.1 ) UO 2 x UO 2 UO 2+x SrO 0.85 1.3 3.09 2.43 0.48 2.93 Cs2O 2.07 2.23 0.73 10.58 8.98 1.25
108 Table 4 5 Binding energies of strontium -strontium interaction calculated using GGA. Sr Sr Distance () Binding energy (eV) 3.87 0.13 5.49 0.23 6.72 0.02 7.76 0.06
109 Figure 4 1 Lattice parameter of UxM1xO2 pseudo-binary solid solut ions. [Adapted from reference 86 (Figure 1.)]
110 Figure 4 2 Illustrations of the various trap sites for fissio n products considered in this work. Uranium atoms are shown in blue and oxygen atoms in red. A sample fission product is shown in green. Oxygen vacancy Uranium vacancy Divac ancy Bound Schottky defect
111 -0.3 -0.1 0.1 0.3 0.5 0.7 0.9 1.1 1.3 1.5 1.7 1.9 1 2 3 4 5 System Size Cs Incorporation Energy (eV) Interstitial Oxygen _Vacancy Uranium _Vacancy 1x1x1 2x2x2 3x3x3 4x4x4 5x5x5 Figure 4 3 Normalized incorporation energy of cesium in pre -existing trap sites calculated as function of system si ze by atomic level simulations using an empirical potential within GULP.
112 -1.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 -0.10 -0.08 -0.06 -0.04 -0.02 0.00 0.02 0.04 0.06 0.08 0.10 Stoichiometry Defect formation energy (eV) 300 K 2000 K 1000 K Grimes -x +x Figure 4 4 Formation energy of a uranium vacancy calculated using SP -GGA+U method as a function of temperature and stoichiometry using the point defect model. The temperatures taken into account are 300 K, 1000 K and 2000 K.
113 UO2-x 0 2 4 6 8 10 12 Xe_int. Xe_O Xe_U Xe_divacancy Xe_schottky Cs_int. Cs_O Cs_U Cs_divacancy Cs_schottky Sr_int. Sr_O Sr_U Sr_divacancy Sr_schottky Relative Solution Energy (eV) SP-GGA+ U GGA Grimes F igure 4 5 Solution energies of Xe, Cs and Sr in UO2 relative to the lowest solution energy site for the hypostoichiometric ( UO2 x) case.
114 UO2 0 2 4 6 8 10 12 14 16 Xe_int. Xe_O Xe_U Xe_divacancy Xe_schottky Cs_int. Cs_O Cs_U Cs_divacancy Cs_schottky Sr_int. Sr_O Sr_U Sr_divacancy Sr_schottky Relative Solution Energy (eV) SP-GGA+ U GGA Grimes F igure 4 6 Solution energies of Xe Cs and Sr in UO2 relative to the lowest solution energy site for the stoichiometric ( UO ) case.
115 UO2+x 0 2 4 6 8 10 12 14 16 18 20 22 24 26 28 Xe_int. Xe_O Xe_U Xe_divacancy Xe_schottky Cs_int. Cs_O Cs_U Cs_divacancy Cs_schottky Sr_int. Sr_O Sr_U Sr_divacancy Sr_schottky Relative Solution Energy (eV) SP-GGA+ U GGA Grimes F igure 4 7 Solution energies of Xe, Cs and Sr in UO2 relative to the lowest solution energy site for the hyperstoichiometric ( UO2 + x) case.
116 CHAPTER 5 SEGREGATION OF FISSI ON PRODUCTS TO GRAIN BOUNDARIES Introduction Segregation can be understood as the interaction between extrinsic impurity and structural defects that includes point defects, line defects such as dislocations and surface defects such as grain boundaries and free surfaces. The driving force for this process is the free energy difference between the impurity at a bulk site and at a defect site. If the energy at the defect site is less that in the bulk, the impurity segregate s to the defect site. An accurate description of segregation is the basis for understanding a variety of surface and interfacial phenomena. Surface segregation influences such properties as oxidation, wetting, and catalysis, while segregation at grain boundaries affects diffusion, creep, corrosion and electrical properties.34 Segregation is especially importa nt in the case of nuclear fuels due to the complex behavior of fission gases xenon (Xe) and krypton (Kr).16 Both these gases are virtually insoluble in the fuel matrix and in their pure states and under ambient conditions they exist as gases not solids. The consequences of this are that either the gases are released fro m the fuel or they contribute to the internal atmosphere of the fuel matrix. Both these routes are detrimental to fuel performance. If the gases are released from the fuel, they contribute to the gaseous atmosphere within the fuel pin and the pressure insi de correspondingly increases leading to stresses which can ultimately lead to the failure of the fuel pin. If these gases are retained inside the fuel, they migrate to grain boundaries where they form bubbles.19 This leads to swelling of the fuel matrix. Swelling is detrimental to fuel performance as it causes fuel -cladding contact and the resulting stresses can shor ten the lifetime of the pin. Swelling and release are complimentary phenomena and need to be studied in depth to gain a fundamental understanding of the complex behavior of
117 fission gases. This work has concentrated on the swelling effects of fission gases, in particular Xe, on UO2. Review of Literature A very large body of segregation work has been devoted to the study of high angle grain boundaries in metals and alloys. Most of these studies96 have been performed on special grain boundaries, in particular symmetrical twist and tilt boundaries. However, the results of these studies are applica ble to more general and random grain boundaries since these boundaries are composed, in part, by special grain boundaries. Such special boundaries are usually described by their low values, where is the reciprocal density of coincident sites in the CSL model, which is described in Chapter 2.31 Segregation to grain boundaries is similar in concept to the segregation that occurs in surfaces. However, due to the wide variety of boundaries observed, the segregation phenomena is much more complex and is found to be dependent on the type of boundary. In sphere sintering experiments by Gleiter and co -workers,97 they observed the occurrence of pointed minima, so called cusps, in the grain boundary energy versus misorientation curves, see Figure 5 1 for an illustration. As can be observed in the figure, the addition of a solu te changes the grain boundary energy. The driving force for segregation is the decrease in interfacial energy. Thus, boundaries with orientations between cusps are more prone to segregation compared to shallow cusps (B) and deeper cusps are observed to ha ve the have the least influence of segregation (A, A). However, it must be pointed out that the energy versus orientation curve has only been established for low values of and other parameters, beyond the misorientation angle, will influence segregation Experimental evidence of the variation of the amount of segregation with grain boundary orientation has been found in many studies of metals and alloys. Watanabe et al.98 showed that
118 above a tilt angle of 15, the segregation of tin in oriented bicrystals of iron is high with a 20% variation between dif ferent orientations. It was also found in segregation studies of bismuth in copper that lower energy boundaries were less prone to segregation than higher angle grain boundaries.99 The ionic character of the bonding of ceramics leads to additional complications and makes segregation studies difficult in ionic interfaces as compared to metallic interfaces. The formation energies of cations and anions within ceramic are, in general, different from one another and are lower at the boundary relative to the bulk. Therefore an electrical double layer or potential will form between the bo undary and bulk.34 This affects both the formation energies of intrinsic defects as well as impurity segreg ation to boundaries. This potential must be taken into account, along with bond-breaking and strain -energy relief, in the determination of segregation behavior. Noboru and Takagi100 have examined fission product and plutonium distribution in irradiated UO2 samples. In their work, UO2 pellets were ground s o as to fracture them along grain boundaries and the distribution of cesium (Cs), ruthenium (Ru) and cerium (Ce) was determined in the grains and grain boundaries for three different irradiated samples. It was observed that while Ru and Ce were uniformly d istributed in the grains, Cs accumulated along the grain boundaries. However, at higher temperatures, Cs also exhibited inhomogeneous distribution. A scanning electron microscopy (SEM) and transmission electron microscopy (TEM) study has been carried out b y Sonoda et al.101 to observe th e microstructural evolution of irradiated UO2 samples as a function of burnup and temperature. It was observed through TEM that the original grains restructured into sub -grains and these sub grain boundaries were low angle and heavily decorated with fissi on gas bubbles. Kubo et al.102 doped UO2 with gadolinium (Gd) to study the elec trical properties of UO2 grain boundaries. They observed that the boundary
119 conductivity decreased with increasing Gd content, owing to the segregation of Gd to the boundaries and forming a potential barrier to the migration of electron holes. There also ha s been a large amount of theoretical work devoted to the segregation of dopants and impurities to grain boundaries. A theoretical understanding of the relationship between orientation and segregation must be based on the atomic structure of the grain bound ary. There are primarily two theoretical approaches employed in the study of segregation. One is based on interatomic potentials, which are used both by molecular statics and molecular dynamics to find the minimum energy structures of the grain boundary in the presence of an impurity. This approach deals with the strain energy contribution to the segregation enthalpy. The other method concentrates on the change in bonding characteristics at the boundary using molecular orbital approaches. Sutton and Vi tek103 examined the segregation behavior of bismuth (Bi) and silver in a copper (Cu) matrix using molecular statics. The have used 5 (210) and 17 (530) symmetric tilt grain boundaries in their study. They observed no significant difference in the segregation behavior of the two boundaries. However, the segregation energy for the seven differe nt sites considered varied from 1.1 eV to 2.3 eV. As stated above, segregation only occurs at sites associated with negative segregation energy. The suggested reason for this trend was the atomic size difference between Bi and Cu. This atomistic model pre dicted a significant change in the grain boundary structure upon segregation. The atomistic structure of ordered metal alloys, particularly transition metal aluminides, has received special attention due their inherent brittleness at the grain boundaries. Farkas et al.96 simulated grain boundaries in NiAl and their interaction with point defe cts using the embedded atom (EAM) potential. The 5 (310)/ grain boundary has also been examined in significant
120 detail in fluorite structures. Mao et al.104 studied the segregation of yttrium (Y) to symmetric tilt grain boundaries in cubic zirconia (ZrO2) using DFT. Five different models of grain boundaries were constructed based on Z -contrast scanning transmission electron microscopy (STEM) and electron energy loss spectra (EELS) measurements.105 The calculated lowest energy structure was observed to be in agreement with experiments. The Y3+ segregation to Zr4+ sites in the grain boundary was observed to be energetically favorable compared to the bulk. There have been many simulation studies on grain boundaries in UO2 but none examining the segregation of fission products to the grain boundaries, thus providing the motivation for this work. Van Brutzel et al.106 carried out collision damage simulations on different models of symmetric tilt boundaries. It was found that the cascade morphologies were dependent on the misorientation angle. Watanabe et al.107 1 08 performed thermal conductivity studies on pol ycrystalline UO2 using molecular dynamics simulations with Busker and Yamada potentials.81 As the thermal conductivity is influenced by the presence of fission products, the study by Watanabe et al. becomes the basis of future work. Stanek et al.109 examined the segregation of barium (Ba2+), strontium (Sr2+), cerium (Ce4+) and zirconium (Zr4+) to the stable (100), (110) and (111) low index surfaces of UO2 using empirical potentials. In the case of Ba2+ and Sr2+, charge compensatio n was required in the form of oxygen vacancies. It was observed that the segregation energies strongly depended on the type of surface and the orientation of the defect cluster. Owing to the similarity between surface and grain boundary segregation, a simi lar trend is expected for grain boundary segregation of fission products in UO2. Recently, Geng et al.110 have carried out molecular dynamics simulations of Xe (111) planar clustering in UO2. They used the Basak potential1111 and predicted the planar configuration to be stabl e.
121 Computational Methodology The segregation of Xe to grain boundaries is examined in this work. The goals were to understand how Xe segregates to grain boundaries in UO2 and to determine how sensitive this segregation behavior is to the details of the gr ain boundary structure. T hree different grain boundaries were considered based on the relative orientation of the two grains. The first is the 5 (310)/  tilt grain boundary which is a low energy boundary commonly found in fluorite materials. The second one is a 5 twist boundary to analyze the influence of orientation of grains on the segregation energies. Both these boundaries were constructed using the software GBstudio.124 Different sizes of supercells were explored to understand the effect of system size and to remove any artificial interactions betwee n images of atoms by the use of periodic boundary con ditions. The third structure considered wa s an amorphous grain boundary. Interfaces in real materials are more general and random compared to the symmetric crystalline boundaries used in simulation studi es. With advances in computational power, it is possible to compare these more realistic interfaces for simulation studies. This grain boundary was constructed using two grains of equal size and which were annealed using molecular dynamics (MD) simulations with the Basak potential. Care was taken to ensure that a w ell -equilibrated grain boundary with no external strain was obtained. It has a 45 tilt and the grain boundary plane is normal to the z direction. The three structures are shown in Figure 5 2. Af ter the satisfactory structures were achieved for the boundaries, Xe segregation studies were performed. In particular, every uranium atom was replaced by a Xe atom in the bulk of the structure as well as in the grain boundary. Static energy calculations w ere performed using three different potentials. The first potential we used is the Basak potential. This potential has the functional form given by:
122 ) r r ( exp ) r r ( [exp D r C r exp A ) r ( V* 2 26 (5 1) In addition to the Buckingham part described in Chapter 2, this potential has a Morse term (second term on the right hand side). This term represents covalence in the material. However, it does not follow the strict definition of covalent character in the sense that it is non -directional. This potential also includes partial c harges of +2.4 for uranium and 1.2 for oxygen as opposed to the formal charges of +4 and 2. In a review of empirical potentials by Govers et al.,55 this potential reasonably described the bulk and defect properties of UO2. This potential has also been us ed by Geng et al. to examine Xe (111) planar clustering. The Xe interaction to the uranium and oxygen atoms has been fit to DFT+U energies using a Lennard Jones potential for Xe -Xe and Xe O which is given as: 6 12rij rij ) r ( Vij ij (5 2) The Xe -U interaction is given by a Born Mayer type term: r exp ) r ( V (5 3) The second potential used in this work is the Grimes potential,14 which has the advantage of being parameterized for a variety o f fission product interactions with the UO2 host. This is a Buckingham type potential with formal charges for the cations and anions. This approach was used to perform cesium (Cs) incorporation and solution studies in UO2 as discussed in Chapter 4 .The thir d potential used is the Busker potential,93 which is similar in form to the Grimes potential except that it has no short range repulsion for the cations. One of the advantages of Busker model is the t ransferability of the potential for a variety of elements. I t also has parameters for U3+ and U5+, thereby allowing the effects of off -stoichiometry to be considered.
123 The Busker potential was originally developed to study defect properties in UO2 and ZrO2. The potentials with the Morse term have been shown to give better defect and lattice properties compared to formal charge Buckingham potentials that lead to overestimation of binding energies, as discussed in Chapters 3 and 4. This work is organized as fo llows: the Basak potential has been considered for the majority of this work and then the results are compared with the Grimes and Busker potentials in the discussion section. Results Structural Analysis The first step to validate the approach was to calc ulate the radial distribution function (RDF) for the amorphous grain boundary structure and compare it with bulk UO2. RDF is used in computational and statistical mechanics to describe how the density of matter changes with distance for a particular point. It is especially important for crystalline solids as it can describe the bonding environment for a particular atom. Figure 5 3 compares the RDF for the uranium uranium interactions for both bulk UO2 and the amorphous grain boundary structure. For the grai n boundary, the different parts are shown in different colors to illustrate the difference in bonding environments. For bulk UO2, which crystallizes in the fluorite structure, there are twelve nearest uranium neighbors at a distance of 3.88 and six secon d nearest neighbors (NN) at 5.49 For t he grain boundary structure, it was observe d that the bulk portion reproduces the RDF for crystalline UO2 reasonably. The NN peaks are at the same positions but are slightly more diffuse, especially the third NN pe ak where atoms in the bulk start interacting with those in the grain boundary. For the grain boundary part, the NN peaks are approximately at the same positions as the bulk part. This is surprising since the boundary has been termed amorphous so
124 far and hence should not have the spatial correlation of a crystalline solid. The NN peaks are more diffuse, as expected, since the bonding environment is slightly different from the bulk. Next, the uranium oxygen NN peaks were considered as shown in Figure 5-4. For crystalline UO2, there should be eight first NN oxygen atoms at a distance of 2.37 and twenty four second NN at a distance of 4.55 For t he grain boundary structure, it was observe d that the bulk part reproduces the crystalline character having the NN at the same distance. However, the first NN peak for the two grain boundaries is slightly shifted. For the first grain boundary (GB1), it is at 2.27 while for GB2, it is at 2.32 This suggests that the oxygen atoms surrounding the uranium atoms are slightly closer compared to the bulk. This kind of bonding analysis is important as it can provide atomic level explanations to the segregation energy trends. It is well known34 that the bonding environment around a site dictates how the segregation will take place to that site. It was expect ed that a similar RDF trend would be observed for the other two grain boundaries. Xe Segregation to the 5 Tilt Grain Boundary After this initial validation of approach, the segregation energies for Xe are presented. The initial test was on the structure shown in Figure 5 5 a). Here eight different uranium sites in the bulk as well as grain boundary were ea ch replaced with a Xe atom. The segregation energy trend calculated using GULP is plotted in Figure 5 5 b). Here the energy of each Xe site with respect to the lowest energy (normalized energy) is considered as a function of the direction perpendicular to the boundary plane. There appears to be no clear trend in the segregation energy behavior. The reason may be the small size of the supercell. The atoms in the grain boundary can influence the atoms in the bulk since the distance between the two sites is le ss than the cut -off of the potential (8 ).
125 The second initial test that was carried out to was to understand the effect of oxygen atoms on the segregation behavior. When any uranium atom is replaced by a Xe atom, the overall supercell is no longer charge neutral. The charge compensation in empirical potential calculations has to take place through the removal of oxygen vacancies such that a Schottky defect (one uranium and two oxygen atoms) is removed to preserve charge neutrality. The neglect of charge c ompensation in these calculations can have significant influence on the results since in ionic materials the ions can have long range effects. To test for this influence, some calculations were performed where a Schottky defect is removed from the cell and the Xe atom is placed at the removed uranium atom site. Here again, sites we re considered both in the bulk and grain boundary with different arrangements of oxygen atoms around the uranium site as seen in Figure 5 6 a). The position of the Xe is illustrated in green. The normalized energy (energy with respect to the most stable site) is plotted against position in Figure 5 6 b). As can be observed, the energies of the bulk sites are significantly higher than those of the grain boundary sites. The scatter f or sites in the boundary is also high (1 2 orders of magnitude). Even for the same uranium site, different oxygen atom configurations give different trends To understand the origin of this discrepancy, the optimized structures of both the configurations are illustrated in Figure 5 7. As the figure indicates, there is a significant amount of reconstruction of the Bulk_Schottky1 (see Figure 5 7) site and hence the energy of this site is lower. However, it is no longer a symmetrical tilt interface. Both thes e initial tests prove that larger supercells than those considered here are required for the segregation energy calculations. An additional complication is that it is well known that segregation can cause structural rearrangements at the grain boundary. T his is especially true for Xe. To alleviate both of these limitations, rigid body translations of the grain boundary structure were carried out These
126 consist of shifting one grain with respect to another in the plane of the boundary (x-y in this case) and minimizing the strain in the direction perpendicular to the boundary. In addition to the macroscopic degrees of freedom, there are microscopic degrees of freedom describing the relative translation between the two gr ains that must be minimized. T he surface, the map of the energy of the grain boundary as a function of the relative translation of the two grains, for both boundaries was constructed by shifting one grain with respect to another in the plane of the boundary and minimizing the strain in the direction perpendicular to the boundary. The lowest energy structure is shown in Figure 5 8. The lowest energy structures were then used in the segregation energy simulations. It is to be noted that these structures were still predicted to be the most stable ones even after the addition of Xe. Moreover, as discussed in Chapter 1, the grain boundary energy is an important measure of the sta bility of any boundary. T herefore the grain boundary energy for four different sizes of supercell for the minimized grain boundary structures was calculated The x and y dimensions are the same for all the four structures while the z direction that is perpendicular to the grain boundary plane is expanded. The normalized grain boundary energy is plotted against the num ber of atoms in Figure 5 9, which indicates that the energy is converged with system size at the lowest system size of 720 atoms. This kind of analysis gives an indication that the lowest system size may be large enough to carry out the segregation energy calculations. However, the segregation calculations were carried out for all the supercell sizes considered to observe the influence of cell size. After this initial analysis, the segregation of Xe to the 5 tilt boundary was considered. The relative seg regation energies for the four system sizes considered are plotted with respect to the direction perpendicular to the grain boundary plane in Figure 5 10. The relative energies are
127 plotted in blue while the average segregation energies are illustrated in red. As can be observed there is no symmetry in the segregation behavior on both sides of the grain boundary indicating no bulk behavior for the smallest structure as seen in Figure 5 10 a). Sites on the plane are the most stable sites as expected since th ere is enough space for the Xe atom due to the more open nature of the boundary. However, there are other sites that are e nergetically unfavorable. As progressively larger cells are approached there is a clear indication of bulk and grain boundary behavior and it becomes easier to define segregation energy (the difference in energy between an atom in the boundary and in the bulk). The structures consist of two boundaries: one in the middle and one at both ends due to periodic boundary conditions. The relat ive energy incr eases as the bulk was approached plateaus and then decreases as the boundary was approached although bulk-like behavior was never seen for this boundary This trend becomes more pronounced for larger cell sizes. The average segregation ene rgies as a function of system size are shown in Table 5 1. As can be observed, the average segregation energy is calculated to be ~ 2.3 eV. In other words, sites in the grain boundary on an average are more stable than those in the bulk by 2.3 eV. Moreover the average segregation energy is converged as a function of system size, providing further validation for the system sizes used in this work. The relative uranium atom energy is plotted against the first nearest neighbor distance in Figure 5 11. Xe Segr egation to the 5 Twist Grain Boundary The segregation of Xe to the 5 twist boundary is illustrated in Figure 5 12. Here again, the relative energies are shown in blue and the average energies shown in red. The segregation energy and grain boundary energy convergenc e was achieved for this particular size of the boundary. There are two symmetrical twist boundaries due to periodic boundary conditions. The relative energies do not change significantly with position and hence the average segregation
128 energy is ~ 0.4 eV. T his can be understood on the basis of the relative orientation of the grains. A twist boundary is created when the two grains are rotated relative to each other. Therefore, the relative positions atoms do not change significantly. To understand the origin of th e segregation energy trends, the relative energy of each uranium atom was calculated The relative uranium atom energy is plotted against its average first nearest neighbor distance in Figure 5 13. The relative energy depends on the surroundings of e ach atom and will be different for different structures. As can be observed, there is an inverse correlation between the two values. In other words, as the relative uranium atom energy increases, the average first nearest neighbor distance decreases. Major ity of the points are concentrated around 2.37 which is the 1st nearest neighbor distance of oxygen atoms in crystalline UO2. The results indicate that if the distance of the oxygen atoms is less than that in crystalline UO2, the relative energy of the uranium atom is high because it is constrained by the oxygen atoms. However, a similar trend between the relative segregation energy and the first oxygen nearest neighbor distance was not observed. This might be due to the movement of the xenon atom from t he original lattice position. Xe Segregation to the Amorphous Grain Boundary For the third case, the Xe segregation calculations for the amorphous grain boundary were carried out. Here again, every uranium site in the structure with Xe atoms was replaced These calculations are computationally expensive as the system size had to be significantly larger to observe amorphous behavior of the boundaries The relative segregation energies are plotted against the z coordinate of the structure in Figure 5 14. The relative energies are clearly distinguished by the bulk and grain boundary regi on. The energies get lower as the boundary was approached and the lowest energy site is at the boundary. These energies increase for the bulk region between the two boundaries and this trend is symmetrical. However, there are sites in the
129 boundary where the segregation is energetically unfavorable. In fact, the most energetically unfavorable site (~8 eV) is in the grain boundary. To understand the origin of these effects, various analyses were carried out One of them is illustrates the relative uranium atom energy versus the first nearest neighbor distance plotted in Figure 5 15. The majority of the data is concentrated around the nearest neighbor distance for a perfect UO2 cry stal (2.37 ) again suggesting that the bonding in this amorphous boundary is similar. Although there is some scatter in the data, there is a clear trend with a negative slope. In other words, as the deviation from the ideal nearest neighbor distance incre ases, the energy of the uranium atom decreases. This is physically reasonable because it means that the local arrangement of oxygen atoms surrounding the uranium atom is more open and hence the energy of that atom is low. When the average distance is less than ideal, the uranium atom is constrained and therefore its energy is high. The hypothesis behind this correlation is that if the energy of a uranium atom is high, it will be harder to incorporate a Xe atom at that site since its local structure is const rained and this process becomes energetically unfavorable. However, when the xenon segregation energy is plotted against the relative uranium atom energy (see Figure 5 16), there seems to be no apparent trend. There is significant scatter of the data although there is a bulk -like region for xenon segregation energy of about 5 7 eV. This result suggests that there is no direct correlation between the two energies and this aspect has to be investigated further. An important indicator of the local structure a round an atom is the stress state that particular atom is in. It will be energetically unfavorable for Xe to segregate to an atom that is under compressive stress due to size constrai nts. To test this assumption, the virial stress on each atom in the struc ture that is a measure of the mechanical stress at the atomic scale was calculated For calculation of the virial stress of the simulation cell, the total calculated stress on
130 the structure has to be divided by the simulation cell volume. In order to calcu late the stress/atom, the stress has to be divided by the volume occupied by the atom. As a first approximation, this is taken as the simulation cell volume divided by the number of atoms. The virial stress versus distance is plotted in Figure 5 17. The te nsile stress is defined as positive while the compressive stress is negative. As can be observed from the plot, half of the structure is under tensile stress while the other half is under compressive stress. This is expected since the structure has been ge nerated initially by allowing the two grains to move towards each other using MD. However, it is also to be noted that the part under compressive stress exhibits significant variation in the stress state. After this initial analysis, the relative uranium atom energy was plotted against the virial stress in Figure 5 18 for all the uranium atoms in the structure. For the uranium atoms under tensile stress, there is no apparent trend between the stress and the relative energy. This suggests that if the urani um atom is under tensile stress, xenon has no preference based on the position of the uranium atom. On the other hand, if the uranium atom is under compressive stress, there is direct correlation between the stress and the relative energy. This, in additio n to the conclusions based on Figure 5 15, suggests that if the oxygen atoms surrounding the uranium atom are far apart corresponds to a relatively low energy for that uranium atom and should correspond to low segregation energy of xenon to that particular uranium site. However when the virial stress is plotted against the relative xenon segregation energy in Figure 5 19, no direct relation is observed between the two quantities. There is also significant scatter of the data, especially for the uranium atom s under compressive stress. The above analysis suggests that to understand the origin of the segregation energy trends, further investigations are required and will form the basis of future work.
131 Discussion Comparison of Potentials The segregation behavi or at the surface can be understood on the basis of the number of bonds broken for the solute to segregate to a particular site. Accordingly, more open structures have lower segregation energy. The situation at grain boundaries is a bit more complex. The s egregation energy depends strongly on the orientation of the two planes. Theoretical predictions of segregation depend on the type of empirical potential used to describe the interaction between atoms. First comparing the Basak and the Busker potential for segregation of Xe to the 5 symmetric tilt boundary, both potentials show a symmetric bulk like behavior on either side of the grain boundary plane. Some sites in the grain boundary are predicted to be energetically favorable for segregation as compared to others. However, the quantitative trend is different between the two potentials. The relative segregation energies calculated using the Busker potential are shown in Figure 5 20. The plot shows identical segregation behavior for the bulk region across t he grain boundary plane. It is interesting to observe that the most stable segregation site and the site with the highest relative segregation energy (~ 12 eV) are both in the grain boundary. As compared to the boundary, some sites in the bulk have lower s egregation energies while some others are energetically unfavorable. This can be understood on the basis of the local structure around the segregating atom. Depending on the local structure, the chemical bonding will change and this leads to different segr egation energies. This phenomenon can be observed as illustrated by Briant in Figure 5 21.126 Here the hypothetical energy for a solute in a certain configuration at a boundary is depicted. The polyhedra or the structural arrangement of atoms around the solute is different for each site and the segregation energy correspondingly varies. If the energy of the trap site is lower than that of the bulk, that site will be enriched with that solute. This type of behavior can affect the diffusion properties of Xe in the grain boundary
132 as well as the bulk. For sites with low segregation energies or high stability, it will be harder for the Xe atom to escape these trap sites as the migration barrier would be high. Therefore, even though the migration of atoms is easier in the grain boundary as compared to the bulk, this may not always be the case. The relative segregation energies of the Busker pote ntial are significantly higher than those predicted by the Basak potential. Moreover, the Basak potential predicts a site in the bulk to the most stable solution site, as opposed to the grain boundary site predicted by the Busker potential. These two resul ts clearly show that the interatomic description matters when predicting quantitative trends. The Busker potential is a Buckingham type potential with formal charges for the two ions and this has been shown to lead to an overestimation of the binding energ ies. The Basak potential, on the other hand, has partial ionic charges with a Morse term to represent covalence. It is more long range in its description of ions and tends to give good lattice and defect properties. It is also probable that this long ran ge nature is not captured to a full extent by the current system size. However, both potentials give different local arrangements of oxygen atoms after the rigid body translation as illustrated in Figure 5 8, which seems to be the dominant effect in the di fferent segregatio n behavior predicted by the two potentials. As observed from these Xe segregation calculations and previous literature, the local structure or the oxygen atoms surrounding the Xe atom have an impact on the segregation energies. To unders tand this effect, the local structure of oxygen atoms around the segregating atom was analyzed as a continuation of the last discussion. Based on the arrangement of oxygen atoms, there are three possible sites when only taking the first nearest neighbors i nto account. The results are shown in Figure 5 22. Here the energy with respect to the lowest solution site (normalized energy) is plotted against the position of the Schottky defect. As can be observed,
133 one site (GB_Schottky1) is more stable as compared t o the other sites. However, for the same uranium site, different arrangements of oxygen atoms result in significantly different energies. This result can be seen as a validation of the effect of local arrangements on segregation energies. Moreover, this pl ot shows that some sites in the bulk are energetically favorable (bulk_Schottky2) as compared to the sites in the grain boundary, a trend previously observed for Xe segregation to uranium vacancies. This is an important point as it suggests that the negle ct of charge compensation (the cell not being charge neutral when segregation is considered only to uranium vacancies) does not change the qualitative picture of segregation but may affect the absolute segregation energy values. It is to be noted that the differences in energy are much lower as compared to similar Schottky calculations for smaller grain boundary sizes as previously shown in Figure 5 6 b). This suggests that the present grain boundary size is stable even after the addition of Xe and does not undergo reconstruction. Comparison of Different Boundaries When comparing the three boundaries, the qualitative trends appear to be similar. The relative energy increases from the grain boundary to the bulk region. All the three boundaries show bulk like behavior although this trend is more apparent in the amorphous boundary. There are some energetically favorable and unfavorable sites in both boundaries. This may affect the diffusion properties of Xe during grain boundary migration. However, it is clearl y seen that the relative segregation energies depend on the orientation. The average segregation energies and the grain boundary energies have been calculated in Table 5 2. As can be observed, when comparing the twist boundary to the other two boundaries, there is a direct correlation between the two energies. In other words, lower grain boundary energy corresponds to lower segregation energy. This is expected since a small value of grain boundary energy translates to smaller difference
134 between the bulk and the boundary region. However, when the tilt and amorphous boundaries are considered, the trend reverses. This suggests a strong coupling between the local structure around a segregating atom and the amount of segregation.
135 Table 5 1 Average segregation energy for the 5 tilt boundary as a function of the supercell size. 5 Tilt (310)/ grain boundary # of Xe atoms Avg. Segregation Energy (eV) 480 2.2 720 2.21 960 2.4
136 Table 5 2 Average segregation energy and the grain boundary energy for the three different types of boundaries considered. Type Avg. Segregation Energy (eV) GB Energy (eV) Twist 5 0.37 0.073 Tilt 5 2.2 0.099 Amorphous 1.95 0.116
137 Figure 5 1 A schematic of the grain boundary energy versus tilt angle for low energy ori entations. [Adapted from reference 34 (Figure 7).]
138 a) b) c) Figure 5 2 The three types of grain boundaries used in this work a) 5 (310)/ (001) symmetric tilt boundary, b) 5 twist boundary and c) amorphous grain boundary.
139 0 2 4 6 8 10 12 14 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 7.5Distance (A)g(r) Distance () a) 0 0.5 1 1.5 2 2.5 3 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 7.5 Distance (A)g(r) GB2 GB1 bulk Distance () b) Figure 5 3 Radial distribution functions for the uranium uranium interaction, a) perfect UO2 and b) amorphous boundary structure. The individual par ts of the amorphous grain boundary structure are shown in different colors.
140 0 5 10 15 20 25 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5Distance (A)g(r) Distance () a) 0 0.5 1 1.5 2 2.5 3 3.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7Distance (A)g(r) GB2 GB1 bulk Distance () b) Figure 5 4 Radial distribution functions for the uranium -oxygen interaction a) perfect UO2 and b) amorphous boundary structure. The individual parts of the amorphous gr ain boundary structure are shown in different colors.
141 a) 0 1 2 3 4 5 6 7 8 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 X Coordinate Normalized Energy (eV) Basak -Bulk Basak -GB Grimes -Bulk Grimes -GB b) Figure 5 5 Segregation of Xe a) Illustration of the 5 symmetric tilt boundary and b) normalized energy versus distance along the x coordinate.
142 a) 0 10 20 30 40 50 60 70 80 bulk_Schottky1 bulk_Schottky2 bulk_Schottky3 GB_Schottky1 GB_Schottky2 GB_Schottky3 Position Normalized energy (eV) b) Figure 5 6 Effect of oxygen vacancies on segregation a) Illustration of the position of the uranium atoms where Xe was placed and b) plot of normalized energy versus position.
143 a) b) Figure 5 7 Effect of reconstruction a) Bulk_Schottky1 and b) Bulk_Schottky_3.
144 a) b) Figure 5 8 Illustration of the optimized structures after performing rigid body translations a) Busker potential, b) Basak potential.
145 0 0.05 0.1 0.15 0.2 0.25 0 50 100 150 200 250 300 350 400 450 500 550 600 650 700 750 800 # of UO2 formula Units GB energy (eV/2) Busker Tilt Basak Tilt Figure 5 9 Grain boundary energy as a function of UO2 formula units in the supercell using two different potential for the 5 tilt boundary
146 Tilt5_240Xe 0 1 2 3 4 5 6 7 8 0 5 10 15 20 25 30 35 Distance (A)Energy (eV) Distance ( ) a) Tilt5_480Xe 0 0.5 1 1.5 2 2.5 3 3.5 4 0 5 10 15 20 25 30 35 40 45 50 55 60 65 70 Distance(A)Energy (eV) Relative Average Distance ( ) b) Figure 5 10 Relative s egregation energy as a function of distance perpendicular to the grain boundary plane using the Basak potential for different number of Xe atoms as illustrated in the caption of each figure.
147 Tilt5_720Xe 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 10 20 30 40 50 60 70 80 90 100Energy (eV) Relative Average Distance ( ) c) Tilt5_960Xe 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 20 40 60 80 100 120 140Distance (A)Energy (eV) Relative Average Distance ( ) d) Figure 5 10 Continued
148 Tilt5 2.2 2.25 2.3 2.35 2.4 2.45 2.5 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Relative U atom energy (eV) Avg. 1st NN distance Figure 5 11 Average first nea rest neighbors for the uranium atoms as a function of the relative energy of each uranium atom for the 5 tilt grain boundary.
149 Twist5 0 0.5 1 1.5 2 2.5 3 3.5 0 10 20 30 40 50 60 70 80 90 Energy (eV) Relative Average Distance ( ) Figure 5 12 Relative energy of xenon as a function of distance for the 5 twist grain boundary ; relative energies are illustrated in blue while average energies are shown in red.
150 Twist5 2.22 2.24 2.26 2.28 2.3 2.32 2.34 2.36 2.38 2.4 2.42 2.44 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Relative U atom energy (eV) Avg. 1st NN distance Figure 5 13 Average first nearest neighbors for the uranium atoms as a function of the relative energy of each atom for the 5 twist grain boundary.
151 Amorphous 0 1 2 3 4 5 6 7 8 9 0 10 20 30 40 50 60 70 80 Distance along z (A) Relative Energy (eV) Relative Segregation Energy Average Segregation Energy Figure 5 14 Plot of segreg ation energy as a function of distance for the amorphous grain boundary calculated using the Basak potential.
152 Amorphous Boundary 2.15 2.2 2.25 2.3 2.35 2.4 2.45 2.5 2.55 2.6 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 Relative U atom Energy (eV) Avg. 1st NN distance Figure 5 15 Average first nearest neighbors for the uranium atoms as a function of the relative energy of each atom for the amorphous grain boundary.
153 2 3 4 5 6 7 8 9 10 0 0.5 1 1.5 2 2.5 3 3.5 4 Relative U atom Energy (eV) Relative Xe Segregation Energy (eV) Figure 5 16 Relative segregation energy of Xe atoms as a function of the relative energy of each atom for the amorphous grain boundary.
154 -120 -100 -80 -60 -40 -20 0 20 40 60 0 10 20 30 40 50 60 70 80 Distance (A) Virial Stress (GPa) Distance ( ) Figure 5 17 Virial stress on each atom as a function of the distance along the direction perpendicula r to the grain boundary plane for the amorphous grain boundary.
155 -120 -100 -80 -60 -40 -20 0 20 40 60 0 0.5 1 1.5 2 2.5 3 Relative U atom energy (eV) Virial Stress (GPa) Figure 5 18 Virial stress on each atom as a function relative energy of each uranium atom for the amorphous grain boundary.
156 -120 -100 -80 -60 -40 -20 0 20 40 60 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 7.5 8 8.5 9 Xe Segregation Energy (eV) Virial Stress (GPa) Figure 5 19 Virial stress on each atom as a function of the relative Xe segregation energy for the amorphous grain boundary.
157 0 2 4 6 8 10 12 14 16 0 5 10 15 20 25 30 35 Distance along z (A) Relative Segregation Energies (eV) Figure 5 20 Relative segregation energy as a function of distance perpendicular to the grain boundary plane using the Busker potential for the 5 tilt boundary
158 Figure 5 21 Schematic of the variation in segregation energy with change in the position of solute according to the local arrangement of atoms around the s olute [Adapted from reference 113 (Figure 10). ]
159 0 1 2 3 4 5 6 7 bulk_Schottky1 bulk_Schottky2 bulk_Schottky3 GB_Schottky1 GB_Schottky2 GB_Schottky3 Position Normalized energy (eV) Figure 5 22 Effect of oxy gen vacancies on segregation: plot of normalized energy versus position. The energy is normalized with respect to the lowest energy site. The labels are the same as for Figure 5 6.
160 CHAPTER 6 CONCLUSIONS The aim of this work was to study the p oint defects, fission products and their interactions with grain boundaries in UO2. The first step was to study the common point defects and explore the effect of charge on the relative stability of these defects. In this work, a combined electronic struct ure DFT and thermodynamic approach was used to predict the stability of intrinsic point defects in UO2. These calculations were supplemented by empirical potential -calculations using larger supercells. The electronic structure and equilibrium structure of UO2 we re well reproduced. The formation energi es of individual point defects we re seen to depend on the exchange -correlation functional chosen within the DFT calculations, whether the defect containing supercell was allowed to relax or not, and on the chos en reference state chemical potential. The DFT and empirical potential calculated formation energies for neutral defect complexes match ed both qualitative and quantative trends reported in the literature in that the oxygen Frenkel pair is predicted to be t he most stable defect overall. The stability of charged defects was also analyzed by taking into account the various charge states they are known to possess. The +2 charged oxygen vacancy was predicted to be dominant when the Fermi level is near the valenc e band, but the 4 charged uranium vacancy was predicted to become the energetically favored defect as the Fermi level moves towards the conduction band. The stability of the anti -Frenkel and Schottky defect complexes depend to a significant degree on whet her the individual point defects that made them up were charged or neutral. The formation energies of the neutral complexes agree d better with the experimental values than the energies of the charged complexes. The influence of temperature and partial pres sure was also analyzed. In general, the oxygen interstitial remain ed dominant over the entire range of oxygen par tial pressures; however, it became increasingly difficult to form these defects as higher temperatures and
161 reducing conditions were approached These predictions are consistent with experimental observations and measurements. Next, the stability of selected fission products and their respective oxides in UO2 x were predicted using electronic structure and atomistic level simulations. The most st able solution site was found to depend on stoichiometry, the particular fission product and the exchange correlation functional used. Regarding the SP GGA+U calculations, for Xe, the most stable solution site was the bound Schottky for UO2 x and UO2 and th e uranium vacancy for UO2+x. For Cs, the preferred solution site was the bound Schottky for UO2x, the divacancy trap site for UO2 and the uranium vacancy for UO2+x corresponding to the progressive difficulty of oxygen vacancy formation. For Sr, the energe tically favorable solution site was the divacancy for UO2x and the uranium vacancy for UO2 and UO2+x. In general, higher charged defects were observed to be more soluble for all the three stoichiometries and based on this, the behavior of other fission pr oducts can be predicted. There were differences in the GGA and SP GGA+U methods for the solid fission products with regard to the most stable site. This was attributed to the limitation of GGA to capture the Coulombic effects completely in these charged f ission products. The scenarios where simplified expressions of trap formation energy could be used and where more complicated relations were required were also identified in this study. The solution energy model was then extended to include the binary oxides SrO and Cs2O. Here it was possible to make favorable comparisons with experimentally observed phenomena. It was shown that the increase in solubility of the fission products with stoichiometry was consistent with previous theoretical studies and experim ental observations. These studies can be extended to include other fission products and
162 microstructural features such as grain boundaries and can lead to a fundamental understanding of nuclear fuel phenomena. The work was then extended to study the segrega tion behavior of fission gas atoms at grain boundaries in UO2. In particular, Xe was chosen as the representative gas. The grain boundaries considered were symmetric 5 (310) tilt grain boundary and an amorphous boundary. By using radial distribution func tion as a tool, it was shown that the bonding in the amorphous grain boundary is similar to that of bulk crystalline UO2. Three different empirical potentials based on the interactions between ions and their charges were chosen. The effects of size of the grain boundary supercell and rigid body translations were first considered. It was shown for smaller supercells, that on addition of a Xe atom, the structure of the grain boundary changed significantly and the atom at the boundary influenced effects from t he bulk. Progressively larger supercells were then constructed till the grain boundary energy converged. Every atom in the bulk and the grain boundary was replaced by a Xe atom and the segregation energy calculated. It was found that the segregation behavior was dependent of the type of potential used. Both potentials predicted different sites to the most stable. Some sites in the grain boundary were observed to be more stable compared to others. This could affect the diffusion of Xe atoms in the grain boundary as the stable sites may serve as trap sites. The segregation was then studied to the amorphous boundary. It also showed similar behavior to the tilt boundary that there were some energetically favorable and unfavorable sites in both boundaries. This phenomenon was explained on the basis of local structure around the Xe atom. An inverse relationship was also found between the atom of each uranium atom in the supercell and the first oxygen nearest oxygen neighbors which could help explain the segregati on energy results. The segregation behavior was found to depend on the type of potential and the
163 orientation of the planes constituting the grain boundary. These trends agree with the published literature on the segregation studies.
164 LIST OF REFERENCES 1 IAEA, edited by IAEA, Austria, 2005. 2 C. K. Gupta, Materials in Nuclear Energy Applications (CRC Press, Boca Raton, Florida, 1989). 3 H. Collette, V. Deremincemathieu, Z. Gabelica, J. B. Nagy, E. G. Derouane, and J. J. Verbist, J. Chem. Soc., Faraday Trans. 2 83, 1263 (1987). 4 F. J. Farrell, T. G. Nevell, and D. J. Hucknall, J. Chem. Soc., Faraday Trans. 1 82, 3587 (1986). 5 P. Villars and L. D. Calvert, Pearson's Handbook of Crystallographic Dat a for Intermetallic Phases (ASM, Metals Park, Ohio, 1985). 6 R. D. Shannon, Acta Crystallogr. 32, 751 (1971). 7 L. Pauling, The Nature of the Chemical Bond (Cornell University Press, Ithaca, New York, 1960). 8 B. J. Lewis, W. T. Thompson, F. Akbari, D. M. Thompson, C. Thurgood, and J. Higgs, J. Nucl. Mater. 328, 180 (2004). 9 P. Y. Chevalier, Fischer,E., Cheynet, B., J. Nucl. Mater. 303, 1 (2002). 10 C. Gueneau, Baichi, M., Labroche,D., Chatillon, C.,Sundman,B., J. Nucl. Mater. 304 (2002). 11 M. W. Barsoum, Fundamentals of Ceramics (McGraw -Hill Companies, New York, 1997). 12 C. R. A. Catlow, Proc. R. Soc. London: Ser. A. 353, 533 (1977). 13 P. Kofstad, Non -Stoichiometry, Diffusion and Electrical Conductivity in Binary Metal Oxides, (Wiley Interscience, New York, 1972). 14 R. W. Grimes and C. R. A. Catlow, Philos. Trans. R. Soc. London, Ser. A 335, 609 (1991). 15 C. R. A. Catlow, J. Nucl. Mater. 79, 432 (1979). 16 S. Imoto, J. Nucl. Mater. 140, 19 (1986). 17 D. R. Olander, Fundamental Aspects of Nuclear Reactor Fuel Elements Springfield, Virginia, 1985). 18 F. T. Ewart, R. G. Taylor, J. M. Horspool, and G. James, J. Nucl. Mater. 61, 254 (1976). 19 J. K. Shultis and R. E. Faw, Fundamentals of Nuclear Science and Engineering (Marcel Dekker Inc., New York 2002).
165 20 D. R. Olander and P. Van Uffelen, J. Nucl. Mater. 288, 137 (2001). 21 R. G. J. Ball, Burns, W. G. Henshaw, J., Mignanelli, M.A., Potter, P.E., J. Nucl. Mater. 167, 191 (1989). 22 D. R. Oboyle, F. L. Brown, and J. E. Sanecki, J. Nucl. Mater. 29, 27 (1969). 23 H. Kleykamp, J. Nucl. Mater. 131, 221 (1985). 24 H. Matzke and H. Blank, J. Nucl. Mater. 166, 120 (1989). 25 M. G. Adamson and S. Vaidyanathan, Trans. Am. Nucl. Soc. 38, 289 (1981). 26 I. Sato, H. Furuya, T. Arima, K. Idemitsu, and K. Yamamoto, J. Nucl. Sci. Technol. 36, 775 (1999). 27 C. Sari, U. Benedict, and H. Blank, Thermodyn. Nucl. Mat. Proc. Symp. 3 587 (1968). 28 J. M. Howe, Interfaces in Materials (John Wiley & Sons, Inc, New York, 1997). 29 D. Wolf, Yip, S., Materials Interfaces: Atomic -level Structure and Properties (Springer, New York, 1992). 30 S. M. Allen, Thomas, E. L., The Structure of Materials :MIT Series in Materials Science anf Engineering (John Wiley & Sons, 1999). 31 A. Brokman, Balluffi, R. W., Acta Metallog r. 29, 1703 (1981). 32 W. Bollmann, Crystal Defects and Crystalline Interfaces (Springer Verlag, New York, 1970). 33 H. Miura, J. Mater. Sci. Lett. 13, 46 (1994). 34 P. A. Dowben, Surface Segregation phenomena (CRC Press, 1990). 35 H. Kleykamp, J. Nuc l. Mater. 344, 1 (2005). 36 R. M. Martin, Electronic Structure: Basic Theory and Practical Methods (Cambridge University Press, New York, 2004). 37 M. Payne, R. Teter, D. Allan, T. Arias, and J. Joannopulos, Rev. Modern Phys. 64, 1045 (1992). 38 P. Hohe nberg and W. Kohn, Phys. Rev. B 136, B864 (1964). 39 W. Kohn and L. J. Sham, Phys. Rev. 140, 1133 (1965). 40 A. E. Mattsson, P. A. Schultz, M. P. Desjarlais, T. R. Mattsson, and K. Leung, Model. Simul. Mater. Sci. Eng. 13, R1 (2005).
166 41 N. W. Ashcroft a nd N. D. Mermin, Solid State Physics (Holt Saunders, Philidelphia, Pensalvania, 1976). 42 H. J. Monkhorst and J. D. Pack, Phys. Rev. B 13, 5188 (1976). 43 G. Makov, and Payne, M.C., Phys. Rev. B 51, 4014 (1995). 44 J. Lento, J. L. Mozos, and R. M. Niemi nen, J. Phys. C 14, 2637 (2002). 45 N. Trouillier and J. L. Martins, Phys. Rev. B 42 (1991). 46 D. Vanderbilt, Phys. Rev. B 41, 7892 (1990). 47 P. E. Blochl, Phys. Rev. B 50, 17953 (1994). 48 Y. Baer and J. Schoenes, Solid St. Comm. 33, 885 (1980). 49 J. Faber and G. H. Lander, Phys. Rev. B 14, 1151 (1976). 50 P. J. Kelly and M. S. Brooks, J. Chem. Soc., Faraday Trans. 2 83 (1987). 51 A. I. Liechtenstein, V. I. Anisimov, and J. Zaanen, Phys. Rev. B 52, R5467 (1995). 52 V. I. Anisimov and O. Gunnar sson, Phys. Rev. B 43, 7570 (1991). 53 V. I. Anisimov, F. Aryasetiawan, and A. I. Lichtenstein, J. Phys. C 9 767 (1997). 54 S. L. Dudarev, D. N. Manh, and A. P. Sutton, Philos Mag. B 75, 613 (1997). 55 K. Govers, S. Lemehov, M. Hou, and M. Verwerft, J Nucl. Mater. 366, 161 (2007). 56 R. W. Grimes and C. R. A. Catlow, J. Am. Ceram. Soc. 73, 3251 (1990). 57 C. R. A. Catlow, Ann. Rev. Mater. Sci. 16, 517 (1986). 58 C. G. Van de Walle and J. Neugebauer, J. Appl. Phys. 95, 3851 (2004). 59 G. Kresse and J. Hafner, Phys. Rev. B 47, 558 (1993). 60 G. Kresse and J. Furthmuller, Phys. Rev. B 54, 11169 (1996). 61 J. D. Gale and A. L. Rohl, Mol. Simul. 29, 291 (2003). 62 J. D. Gale, J. Chem. Soc., Faraday Trans. 93, 629 (1997). 63 P. Nerikar, T. Watanabe, J. Tulenko, S. Phillpot, and S. Sinnott, J. Nucl. Mater. 384, 61 (2009). 64 H. Matzke, J. Chem. Soc., Faraday Trans. 2 83, 1121 (1987).
167 65 B. T. W. Willis, Proc. R. Soc. London Ser. A 274, 122 (1964). 66 S. Nicoll, J. Nucl. Mater. 226, 51 (1995). 67 M Stan, J. C. Ramirez, P. Cristea, S. Y. Hu, C. Deo, B. P. Uberuaga, S. Srivilliputhur, S. P. Rudin, and J. M. Wills, J. Alloys Compd. 444, 415 (2007). 68 P. Erhart, A. Klein, and K. Albe, Phys. Rev. B 72, 085213 (2005). 69 A. F. Kohan, G. Ceder, D. Morg an, and C. G. Van de Walle, Phys. Rev. B 61, 15019 (2000). 70 S. K. Xiang, H. C. Huang, and L. M. Hsiung, J. Nucl. Mater. 375, 113 (2008). 71 J. P. Crocombette, F. Jollet, L. N. Nga, and T. Petit, Phys. Rev. B 6410, 104107:1 (2001). 72 M. Freyss, T. Pet it, and J. P. Crocombette, J. Nucl. Mater. 347, 44 (2005). 73 T. Petit, C. Lemaignan, F. Jollet, B. Bigot, and A. Pasturel, Philos. Mag. B 77, 779 (1998). 74 M. Iwasawa, Y. Chen, Y. Kaneta, T. Ohnuma, H. Y. Geng, and M. Kinoshita, Mater. Trans. 47, 2651 (2006). 75 F. Gupta, G. Brillant, and A. Pasturel, Philos. Mag. 87, 2561 (2007). 76 G. Kresse and D. Joubert, Phys. Rev. B 59, 1758 (1999). 77 R. Laskowski, G. K. H. Madsen, P. Blaha, and K. Schwarz, Phys. Rev. B 69, 140408 (2004). 78 J. He, R. K. Beh era, M. W. Finnis, X. Li, E. C. Dickey, S. R. Phillpot, and S. B. Sinnott, Acta Mater. 55, 4325 (2007). 79 M. W. Finnis, A. Y. Lozovoi, and A. Alavi, Ann. Rev. Mater. Research 35, 167 (2005). 80 I. Barin, Thermochemical Data of Pure Substances, Parts II 1989). 81 K. Yamada, K. Kurosaki, M. Uno, and S. Yamanaka, J. Alloys Compd. 307, 10 (2000). 82 D. R. Olander, in Technical Information center Virginia, 1985). 83 G. E. Murch and C. R. A. Catlow, J. Chem. Soc., Faraday Trans. 2 83, 1157 (1987). 84 C. Catlow, R., A., Philos. Trans. R. Soc. London, Ser. A 364, 473 (1978). 85 R. G. J. Ball and R. Grimes, J. Chem. Soc., Faraday Trans. 2 86, 1257 (1990). 86 H. Kleykamp, J. Nucl. Mater. 206, 82 (1993).
168 87 J. -P. Crocombette, J. Nucl. Mater. 305, 29 (2002). 88 M. Freyss, N. Vergnet, and T. Petit, J. Nucl. Mater. 352, 144 (2006). 89 T. Petit, M. Freyss, P. Garcia, P. Martin, M. Ripert, J. P. Crocombette, and F. Jollet, J. Nucl. Mater. 320, 133 (2003). 90 T. Petit, G. Jomard, C. Lemaignan, B. Bigot, and A. Pasturel, J. Nucl. Mater. 275, 119 (1999). 91 G. Brilliant and A. Pasturel, Phys. Rev. B 77, 184110 (2008). 92 F. Gupta, A. Pasturel, and G. Brillant, J. Nucl. Mater. Accepted Proof (2009). 93 M. Abramowski, R. W. Grimes, and S. Owens, J. Nucl. Mater. 2 75, 12 (1999). 94 A. B. Lidiard, J. Nucl. Mater. 19, 106 (1966). 95 L. E. Thomas, C. E. Beyer, and L. A. Charlot, J. Nucl. Mater. 188, 80 (1992). 96 D. Farkas, J. Phys. C 12, R497 (2000). 97 H. Sautter, H. Gleiter, and G. Baro, Acta Metallogr. 25, 467 (1977). 98 T. Watanabe, S. Kitamura, and S. Karashima, Acta Metallogr. 28, 455 (1980). 99 A. Roy, U. Erb, and H. Gleiter, Acta Metallogr. 30, 1847 (1982). 100 O. I. Noboru, J. Takagi, H. Kim, and K. Park, Nucl. Eng. Technol. 2 127 (1965). 101 T. Sono da, M. Kinoshita, I. L. F. Ray, T. Wiss, H. Thiele, D. Pellottiero, V. V. Rondinella, and H. Matzke, Nucl. Instru. Methods Phys. Res. B 191, 622 (2002). 102 T. Kubo, Ishimoto, S., Koyoma, T., Nucl. Eng. Technol. 30, 664 (1993). 103 A. P. Sutton and V. Vi tek, Acta Metal. 30, 2011 (1982). 104 Z. G. Mao, S. B. Sinnott, and E. C. Dickey, J. Am. Ceram. Soc. 85, 1594 (2002). 105 E. C. Dickey, X. D. Fan, and S. J. Pennycook, J. Am. Ceram. Soc. 84, 1361 (2001). 106 L. van Brutzel and E. Vincent -Aublant, J. Nuc l. Mater. 377, 522 (2008). 107 T. Watanabe, S. B. Sinnott, J. S. Tulenko, R. W. Grimes, P. K. Schelling, and S. R. Phillpot, J. Nucl. Mater. 375, 388 (2008). 108 T. Watanabe, B. Ni, S. R. Phillpot, P. K. Schelling, and P. Keblinski, J. Appl. Phys. 102 (2007).
169 109 C. R. Stanek, M. R. Bradford, and R. W. Grimes, J. Phys. C 16, S2699 (2004). 110 W. T. Geng, Phys. Rev. B 68 (2003). 111 C. B. Basak, A. K. Sengupta, and H. S. Kamath, J. Alloys Compd. 360, 210 (2003). 112 H. Ogawa, Mater. Trans. 47, 2706 (2006). 113 C. L. Briant, Acta Metal. 31, 257 (1983).
170 BIOGRAPHICAL SKETCH Pankaj Vilas Nerikar was born in Hyderabad, Andhra Pradesh, India He received his B. E. degree in chemical e ngineering from A. I. S. S. M. S. C. O. E., Pune, Maharasht ra, India in 2003 and M.S. in c hemical e ngineering from the University of Louisiana, Lafayette in 2005. He then began his PhD study in the Department of Materials Science and Engineering at the University of Florida with Prof. Susan B. Sinnott as his advi sor. He plans to pursue a post doctorate at the Los Alamos National Laboratory.