Thermal transport study in actinide oxides with point defects

Nuclear Engineering and Technology 51 (2019) 1398e1405

Contents lists available at ScienceDirect

Nuclear Engineering and Technology journal homepage: www.elsevier.com/locate/net

Original Article

Thermal transport study in actinide oxides with point defects n, Tien Yee Alex Resnick 1, Katherine Mitchell 1, Jungkyu Park*, Eduardo B. Farfa Kennesaw State University, Department of Mechanical Engineering, Kennesaw, GA, 30144, USA

a r t i c l e i n f o

a b s t r a c t

Article history: Received 25 October 2018 Received in revised form 11 March 2019 Accepted 14 March 2019 Available online 19 March 2019

We use a molecular dynamics simulation to explore thermal transport in oxide nuclear fuels with point defects. The effect of vacancy and substitutional defects on the thermal conductivity of plutonium dioxide and uranium dioxide is investigated. It is found that the thermal conductivities of these fuels are reduced significantly by the presence of small amount of vacancy defects; 0.1% oxygen vacancy reduces the thermal conductivity of plutonium dioxide by more than 10%. The missing of larger atoms has a more detrimental impact on the thermal conductivity of actinide oxides. In uranium dioxide, for example, 0.1% uranium vacancies decrease the thermal conductivity by 24.6% while the same concentration of oxygen vacancies decreases the thermal conductivity by 19.4%. However, uranium substitution has a minimal effect on the thermal conductivity; 1.0% uranium substitution decreases the thermal conductivity of plutonium dioxide only by 1.5%. © 2019 Korean Nuclear Society, Published by Elsevier Korea LLC. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

Keywords: Nuclear fuel Thermal conductivity Defects Molecular dynamics

1. Introduction Effective control of operating temperatures within a nuclear reactor is a fundamental requirement of power plant operations and is strongly correlated to the thermal conductivity of the nuclear fuel being used within the reactor core. Among other thermophysical characteristics of nuclear fuels, thermal conductivity is a key component when managing the temperature safety margins during a power plant's operation. During operation, a large temperature gradient is induced with the fuel during reactor cycle and results in cracking because of the increased thermal stress. Fuel pellet swelling also occurs due to the limited dissipation of thermal energy as uneven heat distribution increases the release rate of fission gases. Moreover, the increased fission gas release is expected to cause embrittlement in grain boundaries, significantly shortening the creep life of nuclear fuel pellets. Therefore, ensuring efficient thermal transport in nuclear fuel pellets will stand as a major hurdle for the next generation reactors currently in development. Temperature, lattice structure, point defect, etc. are all thermophysical or chemo-physical material characteristic traits that can alter thermal conductivity. Point defects in actinide oxides are interruptions in the crystalline materials. Point defects include

multi-point vacancies, interstitial vacancies, and substitutional atoms. While these defects alter the lattice vibrations in the structures, many of the actinide oxides research focuses on formation energy [1], diffusion, and uranium doping [2,3]. Ma and Asok [4] used the all-electron full-potential linearized augmented plane wave plus local orbitals basis method while Petit et al. [5] used local-density approximation to calculate formation energies. Duriez et al. [6] experimented on low-Pu content MOX fuel with oxygen vacancy point defects to develop quantitative input for a simple phonon conduction model. In the present study, uranium dioxide (UO2), plutonium dioxide (PuO2), and thorium dioxide (ThO2) are selected to study the effect of point defects on the thermal transport properties in actinide oxides. The both oxygen point defects and cation (plutonium/uranium) defects are investigated. Although the majority of researchers reported the effect of anion point defects in the fluoritetype oxides using atomistic simulations, cation defects such as uranium vacancy and plutonium vacancy can affect thermal transport in oxide fuels since they alter the lattice structures of oxide fuels largely. Additionally, the phonon density of states (PDOS) for UO2 and PuO2 were examined to better understand the effects of the uranium substitution. To obtain the thermal conductivities of bulk structures, simulations are run for the simulation structures with several sample lengths, ranging from 10 to 100 unit cells in length.

* Corresponding author. E-mail address: [email protected] (J. Park). 1 These authors contributed equally to the study. https://doi.org/10.1016/j.net.2019.03.011 1738-5733/© 2019 Korean Nuclear Society, Published by Elsevier Korea LLC. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/ licenses/by-nc-nd/4.0/).

A. Resnick et al. / Nuclear Engineering and Technology 51 (2019) 1398e1405

1399

2. Simulation method In this research, reverse non-equilibrium molecular dynamics (RNEMD) [7] is employed to estimate the thermal conductivity of actinide oxides with various defect types and concentrations. In RNEMD, a heat flux is imposed across a simulation box by exchanging energies between the hot bath and the cold bath (Fig. 1). Once the heat flow in the structure reaches a steady state, the thermal conductivity (k), is calculated using the Fourier's heat conduction law, k ¼  ∕

where is averaged heat flux and

is temperature gradient. The potential field, or set of interatomic energy potentials, developed by Cooper et al. [8e10] is used together with LAMMPS simulator [11] for the all simulations in the present study. The potential field is constructed by combining many body interaction and pairwise interactions as shown in the following equation.

Ei ¼

  1X fab rij  Ga 2 j

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi X  ffi sb rij

(1)

j

Here, Ei represents the potential energy of the ith atom when compared to all other atoms, and rij is the distance between atoms. a and b represents the species of ith atom and jth atom respectively. The first term describes pairwise interaction between ith atom and jth atom, given by fab(rij), while the second term, sb(rij), incorporates many body perturbations into the pairwise interaction. Ga is the constant of proportionality. This force field has successfully predicted thermos-physical properties of actinide solids [8] and the thermal conductivity of UO2 [12]. Rahman et al. [13] showed the effectiveness of this potential model. The authors also successfully estimated the thermal properties of ThO2 using the same potential [14,15]. Simulation structures with different types of defects and defect concentrations are constructed using a custom MATLAB program. Fig. 2 shows the unit cell of stoichiometric sample (Fig. 2a) and the unit cells with three different types of vacancy defects (Fig. 2bed) and substitutional defect (Fig. 2e) studied in this research study. Vacancy defects occur when atoms are missing from their lattice sites. For example, oxygen vacancy defects appear when some of the oxygen atoms in the simulation structure are absent from their lattice sites while plutonium vacancy defects are found when plutonium atoms are missing from their lattice sites. Simulation structures with oxygen vacancy defects are constructed for the both of PuO2 and UO2 with different sample lengths (Fig. 2b). Plutonium vacancy defects and uranium vacancy defects are simulated for

Fig. 2. Unit cell of PuO2 or UO2 with various defects. (a) Stoichiometric unit cell. (b) Oxygen vacancy. (d) Combined vacancies of primary and oxygen atom. (e) Uranium substitutional defect.

PuO2 and UO2 respectively (Fig. 2c). Combined vacancy defects occur when both primary atoms and oxygen atoms are missing in the molecular structure (Fig. 2d). For the case of PuO2, uranium substitutional defects are explored as shown in Fig. 2e; PuO2 structures with uranium substitutional defects are obtained by replacing plutonium atoms with uranium atoms. To obtain a better understanding on the effect of uranium substitution, ThO2 with uranium substitution is also simulated in the present study. The point defects shown in Fig. 2 are not likely to maintain their original positions during molecular dynamics simulations; vacancy and/or interstitial diffusions of these point defects will rearrange cations and anions in the structure. Since the diffusion process is expected to be accelerated at high temperatures, molecular arrangement (or the positions of cations and anions) in a structure with high temperature will be different from that in a structure with low temperature. This difference between the arrangements of cations and anions can cause a difference in the structure's thermal transport property by altering its phonon scattering mechanism. Therefore, more simulation studies are required to carefully track the rearrangement of cations and anions during molecular dynamics simulations. The simulation utilizes actinide dioxides of the FCC fluorite-type (F3m3) with the notation of M1-yO2±x where M ¼ U or Pu and x or y ¼ 0, 0.001, 0.005, 0.01, 0.02, or 0.05 [16] for all vacancy defect types. An example of a simulation sample with a 1% oxygen vacancy defect would be U1O2-0.01. Similarly, a simulation sample of PuO2 with uranium substitution would be P1-xUxO2 in which x ¼ 0, 0.001, 0.005, 0.01, 0.02, or 0.05. The defect concentration is defined by using the formula shown below.

%¼

Number of defect sites Total number of atoms in the stoichiometric structure  100: (2)

Fig. 1. The schematic of RNEMD in actinide oxide structure with vacancy defects.

Defect concentrations are determined by dividing the number of defect sites by the total number of atoms in the stoichiometric structure. In the simulation structure with 7200 unit cells (6  6  200 supercell) of PuO2, for example, 1% oxygen vacancy defects indicate that there are 864 oxygen atoms missing in the

1400

A. Resnick et al. / Nuclear Engineering and Technology 51 (2019) 1398e1405

simulation structure when compared to its stoichiometric structure since there are 86,400 atoms (28,800 oxygen atoms þ 57,600 plutonium atoms) present in the stoichiometric structure. Likewise, the concentration of substitutional defects is determined by dividing the number of substitutional defect sites by the total number of atoms in the stoichiometric structure. For the case of combined vacancies, two different concentrations are selected; PuO2 with 1% Pu and 1% O vacancies, PuO2 with 5% Pu and 5% O vacancies, UO2 with 1% U and 1% O vacancies, and UO2 with 5% U and 5% O vacancies are simulated. To estimate the bulk thermal conductivities of actinide oxides with defects, the size of simulation structures is increased progressively in the direction of the thermal conductivity measurement (x direction in this study). Here, the number of unit cells (N) is used to represent a structure size instead of its physical length (L) to avoid the complexity that stems from the different sample lengths after the equilibration of the simulation structure (or after the simulation system reaches a state that does not show changes in its thermophysical properties such as internal energy and temperature gradient); a simulation structure with a fixed number of unit cells often exhibits different simulation box sizes depending on the defect concentrations once the simulation structure is equilibrated. Nx, Ny, and Nz indicate the number of unit cells of a simulation structure in x, y, and z directions, respectively. In all structures simulated, the number of unit cells in y (Ny) and z (Nz) directions are fixed to be 6. In their previous study [14], the authors observed that thermal transport property is not affected much by the boundary scatterings from x-y or x-z planes. This is mainly because the simulations conducted here are carefully designed to maintain 1D heat conduction. The same phenomenon has been observed repeatedly by other researchers [17e19]. The half of the total number of unit cells in x direction is defined as Nx. In RNEMD algorithm, the hot bath is located at the center of the simulation box while two cold baths are located at two ends of the simulation box. Therefore, the characteristic length for the thermal conductivity estimation becomes the half of the total simulation box size. Seven different lengths (Nx ¼ 10, 20, 30, 40, 50, 75, and 100) are used to estimate the thermal conductivities of PuO2 and UO2 with various defect types and concentrations. It should be noted that the charge imbalance caused by missing atoms and/or substitutions are ignored in this research study because of the relatively small amount of defect concentrations. In reality, charge compensation is expected to occur in actinide oxides with randomly distributed point defects. To accurately analyze the effect of charged defects, it is necessary to consider the corresponding charge compensating defects. For example, actual charge compensating defects can be included in the same simulation structure to balance the charge of the simulation structure. In this case, however, another undesired effect such as phonon scattering from these additional defect sites need to be carefully examined. All structures simulated in this study (including structures with defects) are prepared by the following equilibration procedure. To achieve well-equilibrated simulation structures, energyminimization is performed first to get a near 0 K structure. Energy-minimization is a process of finding an arrangement of atoms that gives the almost zero net interatomic force. Then, the temperature is elevated to 500 K using isenthalpic (NPH) ensemble with Langevin thermostat to accelerate the equilibration process. Thereafter, the temperatures of the simulation structures are decreased to 300 K with isothermal-isobaric (NPT) ensemble. 0.5 fs is used as the time-step during the entire equilibration steps. RNEMD is carried out on the equilibrated structures for 2 ns with a time-step of 1 fs A constant amount of swapped energy is observed for all simulations in the present study and temperature profiles reach a steady state after 0.2 ns; a total time duration of 2 ns in

RNEMD is sufficient for thermal conductivity estimation. Table 1 lists the lattice constants for the various structures simulated in this study after the equilibration. As can be seen in Table 1, the lattice constants of simulation structures with defects do not exhibit a notable difference from the lattice constant of their stoichiometric structures because of the small amount of defect concentrations simulated in this study. Fig. 3 shows the representative temperature profiles after 2 ns of RNEMD simulations. The distinct irregularities in the simulation structures with lower defect concentrations make thermal conductivity estimation more uncertain when compared to the case of higher defect concentrations, thus resulting in bigger standard deviations in the structures with lower defect concentrations. This indicates that the thermal transport in the simulation structures with high purity is ballistic and phonons in these structures encounter more disruptive phonon scatterings at defect locations than in the structures with many defects. As the defect concentration increases, the thermal transport in the structures becomes diffusive, resulting in a smooth linear temperature profile. Temperature gradient is obtained from the linear region of temperature profile (red line shown in Fig. 3a) using a linear regression function in MATLAB. It should be noted that Pu-239 has been selected as a constituent for PuO2 instead of the most stable isotope of plutonium, i.e. Pu-244 in the present study since the composition of Pu-239 is the largest among the isotopes of Pu in MOX fuels, therefore most commonly used in nuclear reactors. To examine the effect of atomic mass in Pu, the authors simulated the 100 unit cells of PuO2 with Pu-244 and identified that the thermal conductivity of PuO2 with Pu-239 is 5% higher than the thermal conductivity of PuO2 with Pu244. The authors will provide a more systematic study on the effect of different isotopic compositions in MOX fuels in the future.

3. Results and discussion The thermal conductivities of UO2 and PuO2 with different vacancies of various concentrations are plotted together with the thermal conductivities of their stoichiometric samples in Fig. 4 as a function of the number of unit cells in x direction, Nx. The highest conductivities of stoichiometric PuO2 and UO2 are found to be 17.15 ± 0.29 W/m-K and 14.99 ± 0.42 W/m-K respectively for the simulation structures with the longest sample length, i.e. Nx ¼ 100. For the case of stoichiometric UO2 and PuO2, the only source of phonon scattering is intrinsic phonon-phonon scattering and the phonon scattering from the simulation boxes. Consequently, the thermal conductivities of stoichiometric UO2 and PuO2 increase with increasing the sample length because of the ballistic thermal transport when the sample length is smaller than their intrinsic phonon mean free path. When the sample length is smaller than the phonon mean free path, phonons are scattered mainly by the simulation boundaries, and the length of the simulation structure

Table 1 Representative lattice constants after equilibration. Lattice constant, a (Å) % Defect

0% 0.1% 0.5% 1% 2% 5%

UO2 (Nx ¼ 100)

PuO2 (Nx ¼ 100)

O vacancy

U vacancy

O vacancy

U substitution

5.469 5.468 5.466 5.466 5.462 5.445

5.469 5.468 5.469 5.468 5.469 5.468

5.395 5.394 5.392 5.390 5.385 5.362

5.395 5.394 5.395 5.396 5.399 5.405

A. Resnick et al. / Nuclear Engineering and Technology 51 (2019) 1398e1405

1401

Fig. 3. Representative temperature profiles during RNEMD.

limits the thermal conductivity. The authors estimated the phonon mean free path of ThO2 to be 7e8.5 nm at 300K in their previous study [14], and they believe the phonon mean free paths of UO2 and PuO2 are not deviated much from that of ThO2. This length dependent thermal conductivity of nanostructured solids is constantly identified by researchers [18,20,21] who used atomistic simulations. The thermal conductivity ceases to increase linearly once the structure size becomes larger than the phonon mean free path, and eventually reach its bulk thermal conductivity. The simulation method used in this research study has already been validated in the authors’ previous articles [14,15]. The thermal conductivities of stoichiometric UO2 and PuO2 obtained in the present research study are supported by other researchers who reported the thermal conductivities of UO2 and PuO2 using atomistic simulations. For example, Arima et al. [3] reported 15.2 W/m-K as the thermal conductivity of PuO2 using equilibrium molecular dynamics. Qin et al. [12] reported 14.6 W/m-K as the thermal conductivity of UO2 using molecular dynamics simulations. Cooper et al. [22] used the same RNEMD and estimated the thermal conductivities of UO2 and PuO2 to be ~17.5 W/m-K and ~22.5 W/m-K respectively. Fink [23] reported ~13.5 W/m-K as the thermal conductivity of UO2 by adding the contribution of a small polaron to molecular dynamics estimation. Wang et al. [24] reported ~12.5 W/ m-K and 7.5 W/m-K for the thermal conductivities of UO2 and PuO2 respectively from first-principles simulations. The similar dependency on sample length is observed in UO2 and PuO2 with the lowest defect concentrations, i.e. 0.1%. This indicates that thermal transport is still ballistic when a small number of defects is introduced in the simulation structure. This is well supported by the temperature profiles (Fig. 3); when the defect concentration is only 0.1% (Fig. 3a and c), slips at the locations of the hot and cold bath are notable unlike the simulation structures with the highest defect concentrations, i.e. 5% (Fig. 3b and d). However, the size effect on the thermal conductivity becomes negligible

when defect concentration is high (1%, 2%, and 5%) (Fig. 4). The diminishing slopes of thermal conductivity as defect concentrations increase demonstrate that thermal transport becomes more diffusive in structures with higher vacancy defect concentrations; the diffusive thermal transport in these structures can be identified in Fig. 3b and d where temperature slips are not identifiable at the locations of hot bath and cold bath. This diffusive thermal behavior in actinide oxides with large amount of defects indicate that the intrinsic phonon mean free path in defective actinide oxides becomes much shorter when compared to stoichiometric samples. Vacancy defects decrease the thermal conductivities of UO2 and PuO2 for all defect concentrations and lengths simulated in this study. It is interesting to note that oxygen vacancies degrade the thermal conductivities of these nuclear fuels significantly in spite of their small sizes. For example, the thermal conductivity of UO2 with Nx ¼ 100 is decreased by 19.4% when 0.1% oxygen vacancies are introduced into the simulation structure. For the case of PuO2 with Nx ¼ 100, the thermal conductivity decreases by 16.5% when 0.1% oxygen vacancy is added. This detrimental effect of oxygen vacancies have been identified in the authors’ previous research study [15] that estimated the thermal conductivity of ThO2 with defects. Since the thermal transport in actinide oxides is largely dependent on the lattice vibrations, any small perturbations to the lattice structure is expected to alter phonon transport notably. Park et al. [25]. also showed that a small hydrogen functionalization brings a detrimental impact on the thermal conductivity of carbon nanotubes, comparable to the impact caused by large phenyl groups. The thermal conductivities of UO2 and PuO2 decrease rapidly by increasing the concentration of vacancy defects. For example, the thermal conductivity of UO2 decreases by nearly 90% when 5% of oxygen vacancy is present in the sample (Fig. 4a). The complete thermal conductivity estimation results for UO2 and PuO2 with Nx ¼ 100 are listed in Table 2. It is apparent that a vacancy of primary element degrades thermal conductivity more

1402

A. Resnick et al. / Nuclear Engineering and Technology 51 (2019) 1398e1405

Fig. 4. The thermal conductivities of UO2 and PuO2 with various vacancies as a function of sample length.

Table 2 Thermal conductivities of UO2 and PuO2 with vacancy. Thermal conductivity, k (W/m-K) % Defect

0% 0.1% 0.5% 1% 2% 5%

UO2 (Nx ¼ 100)

PuO2 (Nx ¼ 100)

U vacancy

O vacancy

Pu vacancy

O vacancy

14.99 ± 0.42 11.30 ± 0.33 6.02 ± 0.18 3.92 ± 0.08 2.30 ± 0.05 1.22 ± 0.01

14.99 ± 0.42 12.08 ± 0.47 7.46 ± 0.13 5.12 ± 0.09 3.42 ± 0.05 2.03 ± 0.01

17.15 ± 0.29 12.82 ± 0.14 6.67 ± 0.20 4.38 ± 0.26 2.22 ± 0.03 1.32 ± 0.01

17.15 ± 0.29 14.32 ± 0.35 8.83 ± 0.13 5.82 ± 0.07 3.93 ± 0.12 2.16 ± 0.03

severely than the oxygen vacancy equivalent. For example, 0.1% uranium vacancies decrease the thermal conductivity of UO2 with Nx ¼ 100 by 24.6% while the same concentration of oxygen vacancies decreases the thermal conductivity by 19.4%. This can be easily explained by the difference between atomic sizes of oxygen atom and primary atom, i.e. uranium or plutonium; the absence of larger atoms is expected to alter the lattice vibrations in the simulation structure more critically when compared to the absence of smaller atoms. The same phenomenon has been observed in the authors’ previous research [15] that studied thorium vacancies and oxygen vacancies in ThO2, and also constantly observed in other nanostructured materials [25e27]. Using experiments, moreover, researchers repeatedly reported that the thermal transport properties of oxide nuclear fuels are degraded by the deviation from the stoichiometry [28,29]. However, the direct quantitative comparison

between their experimental results and the current simulation results is not possible since it is unattainable to control the amount of vacancies independently during experiments. For example, thermogravimetry method allows the control of oxygen-to-metal ratio of oxide fuels only and does not provide information about the actual amount of oxygen vacancies; with the same amount of deficit in oxygen, different amounts of oxygen vacancies are possible depending on the oxygen interstitials present in the fuel. Fig. 5 compares the effects of vacancy defects and uranium substitutional defects on the thermal conductivities of PuO2 and ThO2. Both PuO2 and ThO2 show that the impact of uranium substitution on thermal conductivity is nominal compared to any other vacancy defect. For the case of PuO2 with Nx ¼ 100, 0.5% oxygen vacancy causes a decrease of 45.6% in thermal conductivity when compared to stoichiometric PuO2, while 5% uranium substitution causes a thermal conductivity decrease of only 12.3%. The similar phenomenon is observed for the case of ThO2; the effect of uranium substitution on the thermal conductivity of ThO2 is minimal when compared to the effect of vacancy defects. However, the impact of the uranium substitution on the thermal conductivity of ThO2 is much larger when compared to its effect on the thermal conductivity of PuO2; 5% uranium substitution decreases the thermal conductivity of ThO2 by 33% while the uranium substitution with the same concentration decreases the thermal conductivity of PuO2 by 12.3% only. This difference can be attributed to the difference between the atomic sizes and ionic radii; the atomic mass of plutonium (239 amu) is nearly the same as the atomic mass of uranium (238 amu) while the atomic mass of thorium (232 amu) is slightly different from the atomic mass of uranium (238 amu). The larger atomic mass interferes with the propagation of phonons more severely, when compared to the smaller atomic mass. It is also important to note the difference in ionic radii of each element, based on its lattice constant, the difference between ThO2 (5.5974 Å [30]) and UO2 (5.471 Å [31]) lattice constant is about 12 p.m., while the difference in PuO2 (5.3975 Å [32]) and UO2 lattice constant is 8 p.m. The effect of the difference in ionic radii is consistent with other studies conducted involving thermal conductivity and vacancy defect relationships as observed in Tada et al. [33], wherein it is examined that the greater difference in ionic radius results in a more pronounced change in thermal conductivity. As observed in Klemen's perturbation theory [34], the perturbation energy in

Fig. 5. Comparison between the effects of oxygen vacancy and uranium substitutional on the thermal conductivity of PuO2 and ThO2. The thermal conductivities are normalized by the thermal conductivities of their stoichiometric samples. Kmax is the thermal conductivity of the stoichiometric sample.

A. Resnick et al. / Nuclear Engineering and Technology 51 (2019) 1398e1405

phonon scattering by point defects is proportional to the difference between the mass of the substituted atom and average atomic mass. The difference between the mass sizes of thorium and plutonium compared to uranium supports the data in that the effect of the uranium substitution is less substantial in PuO2. Important mechanisms that influence thermal transport in actinide oxides with defects include magnitude of phonon frequencies; phonon group velocities; and phonon relaxation times. Point defects reduce thermal conductivity by limiting intrinsic phonon mean free path, resulting in reduced phonon relaxation times as observed in previous sections; the thermal conductivities of actinide oxides with high defect concentrations reach their bulk thermal conductivities when the sample length is increased while the thermal conductivities of stoichiometric actinide oxides continue to increase with an increase in their sample lengths. However, the effect of point defects on magnitude of phonon frequencies and phono group velocities still need to be examined. Although phonon dispersion is a powerful tool that is used to examine these parameters, the production of phonon dispersion relations is prohibited for defective structures with randomly distributed point defects. In light of this, we compare the phonon density of states (PDOS) in actinide oxides with different types of defects that are obtained from the Fourier transform of the autocorrelation function of atom velocities. The velocity vectors are obtained at 300K using molecular dynamics simulations together with the force field used in RNEMD in the present study. Fig. 6a and b shows the PDOS of UO2 with oxygen and uranium vacancies. The shaded area represents the PDOS of stoichiometric UO2. Fig. 6c, d, and e illustrate the PDOS of PuO2 with oxygen vacancies, plutonium vacancies, and uranium substitutional defects. The shaded area represents the PDOS of

Fig. 6. Phonon density of states (PDOS) as a function of phonon energy for (a) UO2 with oxygen vacancy, (b) UO2 with uranium vacancy, (c) PuO2 with oxygen vacancy, (d) PuO2 with plutonium vacancy, and (e) PuO2 with uranium substitutional defects. The area plot indicates the PDOS of stoichiometric samples as a function of the energy of the phonon.

1403

stoichiometric PuO2. The peak locations of the PDOS is similarly found in other works, such as Antropov et al. [35], where the PDOS for the solid uranium examined at 300K had a peak around 10 meV. That is similar to the peak locations of UO2 in Fig. 6a and b. The PDOS for plutonium as observed by Manley et al. [36] has a peak location also around 10 meV. This observation is consistent with the low-energy peak for PuO2 found in this study. It is noted that vacancy defects reduce mode degeneracies identified as peaks, especially at energies around 20 meV. However, these reduced peaks are not likely to affect the thermal conductivities of these nuclear fuels much as long as phonon group velocities remain unchanged. In the both UO2 and PuO2, the high frequency (energy) peaks (~75 meV) are shifted in all cases although they are not expected to contribute to the reduction in thermal conductivity because high energy phonons generally have low group velocities. More important peak shift we have to focus on are observed in low frequencies (energies) below 20 meV. Peak shift (or frequency reduction) is more noticeable for the vacancies of larger atoms (UO2 with uranium vacancy and PuO2 with plutonium vacancy). This reduction in phonon mode frequencies may result from reduced phonon group velocities that are caused by increased structural compliance with vacancy defects. The phonons from low frequency (energy) modes contribute more to thermal transport when compared to higher frequency (energy) phonon modes, thus it is expected that the uranium and plutonium vacancies contribute to the previously observed reduction in thermal conductivity for both fuel types. When observing the PDOS for UO2 and PuO2 with oxygen vacancies, the same differences between the graphs can be seen, however, the degree to which the graphs vary is significantly smaller compared to PDOS for UO2 and PuO2 with vacancies of larger atoms, i.e. uranium or plutonium vacancies. This implies that group velocities and phonon frequencies are not largely affected by oxygen vacancies, and reduced phonon relaxation times are assumed to be the dominant mechanism that limits the thermal conductivity of UO2 and PuO2 with vacancy defects. The PDOS in PuO2 with uranium substitutional defects do not show the same extent of variations in the low energy region as vacancy defects did. Although there is a slight peak shift in the high energy region in between 70 meV and 80 meV, it is unlikely to have an effect on the phonon transport of PuO2 since high-frequency (energy) phonons carry vibrational energy inefficiently when compared to low-frequency (energy) phonons because of their low group velocities; phonon dispersion relations in PuO2 repeatedly showed that acoustic phonon modes (the most important heat carrier) have high group velocities when their frequencies are low [37,38]. This minimal effect of substitutional defects on phonon spectrum is well supported by our thermal conductivity calculations that showed a small change in the thermal conductivity of PuO2 by uranium substitutional defects observed in Fig. 5. To provide a better understanding on phonon scattering in defective actinide oxides, we repeat the RNEMD calculations at sample temperatures varying from 100 K to 1200 K. The characteristic length during the thermal conductivity calculation is fixed to be Nx ¼ 50. (Note that Table 2 lists the thermal conductivities for the structures with Nx ¼ 100). Defect concentrations are controlled to be 0.1%, 1%, and 5% for uranium vacancies and oxygen vacancies. In stoichiometric samples, as shown in Fig. 7a, high temperature degrades the thermal conductivity of UO2. This is explained by Umklapp phonon-phonon scattering. High temperatures increase the population of phonons with large momenta which often lose their energy by being collided with other high momentum phonons. Therefore, this Umklapp phonon-phonon scattering reduces the phonon mean free path of the system, resulting in decreased thermal conductivities. Kim et al. [39] also demonstrated that the acoustic contribution to the thermal conductivity of UO2 has a

1404

A. Resnick et al. / Nuclear Engineering and Technology 51 (2019) 1398e1405

conductivity of structures with a large number of defects cannot be increased much by decreasing the temperature. In Fig. 7b, thermal conductivity calculations for UO2 and PuO2 by other researchers are plotted together with the results in this study. The results in the current study show a good match with the results by other researchers. However, our thermal conductivity calculations for PuO2 exhibit a deviation from the calculations by Cooper et al. [22] although the same interatomic potential has been used. This can be explained by different simulation factors we have chosen in this study such as simulation box size, equilibration methods, and amount of energy swap during RNEMD. The authors will provide more systematic study on the effect of different simulation settings on the thermal conductivity calculation by RNEMD in the near future. 4. Conclusions

Fig. 7. (a) Temperature effect on the thermal conductivity of UO2 with defects. Data are curve-fitted using a correlation, k ¼ k100 K/(A þ BT) where k100 K is the thermal conductivity at 100 K, and A and B are fitting parameters. (b) Comparison with other experimental (Fukushima et al. [41], Gibby [42]) and simulation results (Arima et al. [3], Cooper et al. [22]).

strong T1 dependence and dominates the lattice thermal conductivity of UO2. The fitting curve for the stoichiometric sample in Fig. 7a also demonstrates this T1 dependence. This is explained more precisely by considering dominant phonon wavelength. Dominant phonon wavelength, ld can be evaluated from ħud ¼ hn/ ld z kBT [40], where ħ is the modified Planck constant, kB is the Boltzmann constant, v is the phonon group velocity, ud is the dominant phonon frequency, and T is the temperature. As shown in this relation, the dominant phonon wavelength increases with decreasing the temperature. Therefore, at low temperatures, long wavelength phonons are expected to dominate the thermal transport in the system and they are not likely to participate in Umklapp scattering because of their small wavenumbers. At high temperatures, however, the dominant phonons are short wavelength phonons that are usually engaged in Umklapp scattering, and thus the intrinsic phonon mean free path is decreased with an increase in temperature. Unlike the stoichiometric sample, the thermal conductivities of samples with high defect concentrations do not exhibit a clear dependency on temperature (Fig. 7a). The parameter “B” which is the fitting parameter related to the temperature dependence is only 0.00016 for the UO2 with 5% oxygen vacancy while the parameter “B” is 0.00476 for the stoichiometric UO2. In samples with high defect concentrations, the activation of long wavelength phonons is physically prohibited. Therefore, the thermal

Through the application of RNEMD simulations, the thermal conductivities for three actinide oxides (UO2, PuO2, and ThO2) with vacancy and uranium substitutional defects were found. The effects of these defects were observed to be significant across all three nuclear fuels, and even the smallest percentage of defect has a noteworthy impact on thermal conductivity. Among the three defect types examined, the vacancies of primary atoms (uranium, plutonium, and thorium atoms) decreased the thermal conductivities of actinide oxide fuels the most, compared to the oxygen vacancy and the uranium substitutional defect. In the case of UO2 with 100 unit cell in length, for example, a 0.5% uranium vacancy decreased the thermal conductivity by 85.4% while the oxygen vacancy of the same concentration decreased the thermal conductivity by 67.1%. It can be concluded that although oxygen vacancies have a substantial negative effect on the thermal conductivity of actinide oxides, the vacancies of primary atoms degrade the thermal conductivity of actinide oxides more severely. Further conclusions can be drawn by investigating the uranium substitutional defects on PuO2 and ThO2. A 1.0% uranium substitutional defect in PuO2 (Nx ¼ 100) results in a thermal conductivity of 16.89 W/m-K; this is a decrease of 1.5% compared to the stoichiometric PuO2 of the same length (17.15 W/m-K). ThO2 of the same length and the same defect concentration has a thermal conductivity of 14.02 W/m-K while the thermal conductivity of its stoichiometric sample is 15.96 w/m-K. The resulting difference of 12.9% reveals that uranium substitution has a greater effect on ThO2 than PuO2. Acknowledgements This study was part of collaborative research effort among members of the Nuclear Energy, Science, and Engineering Laboratory at Kennesaw State University. This work made use of the High Performance Computing Resource at Kennesaw State University. Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.net.2019.03.011. References [1] Y. Lu, Y. Yang, P. Zhang, Thermodynamic properties and structural stability of thorium dioxide, J. Phys. Condens. Matter 24 (22) (2012) 225801. [2] P. Martin, D.J. Cooke, R. Cywinski, A molecular dynamics study of the thermal properties of thorium oxide, J. Appl. Phys. 112 (7) (2012), p. 073507. [3] T. Arima, S. Yamasaki, Y. Inagaki, K. Idemitsu, Evaluation of thermal properties of UO2 and PuO2 by equilibrium molecular dynamics simulations from 300 to 2000 K, J. Alloy. Comp. 400 (1e2) (2005) 43e50.

A. Resnick et al. / Nuclear Engineering and Technology 51 (2019) 1398e1405 [4] L. Ma, A.K. Ray, Formation energies and swelling of uranium dioxide by point defects, Phys. Lett. 376 (17) (2012) 1499e1505. [5] T. Petit, C. Lemaignan, F. Jollet, B. Bigot, A. Pasturel, Point defects in uranium dioxide, Phil. Mag. B 77 (3) (1998) 779e786. [6] C. Duriez, J.-P. Alessandri, T. Gervais, Y. Philipponneau, Thermal conductivity of hypostoichiometric low Pu content (U, Pu) O2 x mixed oxide, J. Nucl. Mater. 277 (2e3) (2000) 143e158. [7] F. Müller-Plathe, A simple nonequilibrium molecular dynamics method for calculating the thermal conductivity, J. Chem. Phys. 106 (14) (1997) 6082e6085. [8] M. Cooper, M. Rushton, R. Grimes, A many-body potential approach to modelling the thermomechanical properties of actinide oxides, J. Phys. Condens. Matter 26 (10) (2014) 105401. [9] M.W. Cooper, S.T. Murphy, P.C. Fossati, M.J. Rushton, R.W. Grimes, Thermophysical and anion diffusion properties of (Ux, Th1x) O2, in: Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, vol. 470, The Royal Society, 2014, p. 20140427, no. 2171. [10] M. Cooper, S. Murphy, M. Rushton, R. Grimes, Thermophysical properties and oxygen transport in the (U x, Pu 1x) O 2 lattice, J. Nucl. Mater. 461 (2015) 206e214. [11] S. Plimpton, Fast parallel algorithms for short-range molecular dynamics, J. Comput. Phys. 117 (1) (1995) 1e19. [12] M. Qin, et al., Thermal conductivity and energetic recoils in UO2 using a many-body potential model, J. Phys. Condens. Matter 26 (49) (2014) 495401. [13] M. Rahman, B. Szpunar, J. Szpunar, The induced anisotropy in thermal conductivity of thorium dioxide and cerium dioxide, Mater. Res. Express 4 (7) (2017), p. 075512. [14] J. Park, E.B. Farf an, C. Enriquez, Thermal transport in thorium dioxide, Nuclear Engineering and Technology 50 (2018) 731e737. n, K. Mitchell, A. Resnick, C. Enriquez, T. Yee, Sensitivity of [15] J. Park, E.B. Farfa thermal transport in thorium dioxide to defects, J. Nucl. Mater. 504 (2018) 198e205. [16] B. Willis, Structures of UO2, UO2þ x andU4O9 by neutron diffraction, J. Phys. 25 (5) (1964) 431e439. [17] P.K. Schelling, S.R. Phillpot, P. Keblinski, Comparison of atomic-level simulation methods for computing thermal conductivity, Phys. Rev. B 65 (14) (2002) 144306. [18] T. Watanabe, S.B. Sinnott, J.S. Tulenko, R.W. Grimes, P.K. Schelling, S.R. Phillpot, Thermal transport properties of uranium dioxide by molecular dynamics simulations, J. Nucl. Mater. 375 (3) (2008) 388e396. [19] M. Cooper, S. Middleburgh, R. Grimes, Modelling the thermal conductivity of (U x Th 1x) O 2 and (U x Pu 1 x) O 2, J. Nucl. Mater. 466 (2015) 29e35. [20] J. Park, V. Prakash, Phonon scattering and thermal conductivity of pillared graphene structures with carbon nanotube-graphene intramolecular junctions, J. Appl. Phys. 116 (1) (2014) 014303. [21] J. Park, V. Prakash, Thermal transport in 3D pillared SWCNTegraphene nanostructures, J. Mater. Res. 28 (7) (2013) 940e951. [22] M. Cooper, S. Middleburgh, R. Grimes, Modelling the thermal conductivity of (UxTh1x) O2 and (UxPu1x) O2, J. Nucl. Mater. 466 (2015) 29e35. [23] J. Fink, Thermophysical properties of uranium dioxide, J. Nucl. Mater. 279 (1) (2000) 1e18.

1405

[24] B.-T. Wang, J.-J. Zheng, X. Qu, W.-D. Li, P. Zhang, Thermal conductivity of UO2 and PuO2 from first-principles, J. Alloy. Comp. 628 (2015) 267e271. [25] J. Park, M.F. Bifano, V. Prakash, Sensitivity of thermal conductivity of carbon nanotubes to defect concentrations and heat-treatment, J. Appl. Phys. 113 (3) (2013), p. 034312. [26] H.-p. Li, R.-q. Zhang, Vacancy-defecteinduced diminution of thermal conductivity in silicene, EPL (Europhysics Letters) 99 (3) (2012) 36001. [27] N. Wei, Y. Chen, K. Cai, J. Zhao, H.-Q. Wang, J.-C. Zheng, Thermal conductivity of graphene kirigami: ultralow and strain robustness, Carbon 104 (2016) 203e213. [28] R. Kavazauri, S. Pokrovskiy, V. Baranov, A. Tenishev, Thermal properties of nonstoichiometry uranium dioxide, in: IOP Conference Series: Materials Science and Engineering, vol. 130, IOP Publishing, 2016 no. 1, p. 012025. [29] T. Pavlov, et al., Measurement and interpretation of the thermo-physical properties of UO2 at high temperatures: the viral effect of oxygen defects, Acta Mater. 139 (2017) 138e154. [30] T. Yamashita, N. Nitani, T. Tsuji, H. Inagaki, Thermal expansions of NpO2 and some other actinide dioxides, J. Nucl. Mater. 245 (1) (1997) 72e78. [31] G. Leinders, T. Cardinaels, K. Binnemans, M. Verwerft, Accurate lattice parameter measurements of stoichiometric uranium dioxide, J. Nucl. Mater. 459 (2015) 135e142. [32] J. Haschke, T.H. Allen, L.A. Morales, Surface and corrosion chemistry of plutonium, Los Alamos Sci. 26 (2) (2000) 252e273. [33] M. Tada, M. Yoshiya, H. Yasuda, Effect of ionic radius and resultant twodimensionality of phonons on thermal conductivity in M x CoO 2 (M¼ Li, Na, K) by perturbed molecular dynamics, J. Electron. Mater. 39 (9) (2010) 1439e1445. [34] P. Klemens, The scattering of low-frequency lattice waves by static imperfections, Proc. Phys. Soc. 68 (12) (1955) 1113. [35] A. Antropov, K. Fidanyan, V. Stegailov, Phonon density of states for solid uranium: accuracy of the embedded atom model classical interatomic potential, in: Journal of Physics: Conference Series, vol. 946, IOP Publishing, 2018 no. 1, p. 012094. [36] M. Manley, et al., Phonon density of states of a-and d-plutonium by inelastic x-ray scattering, Phys. Rev. B 79 (5) (2009), p. 052301. [37] M. Manley, et al., Measurement of the phonon density of states of PuO 2 (þ 2% Ga): a critical test of theory, Phys. Rev. B 85 (13) (2012) 132301. [38] P. Zhang, B.-T. Wang, X.-G. Zhao, Ground-state properties and high-pressure behavior of plutonium dioxide: density functional theory calculations, Phys. Rev. B 82 (14) (2010) 144110. [39] H. Kim, M.H. Kim, M. Kaviany, Lattice thermal conductivity of UO2 using abinitio and classical molecular dynamics, J. Appl. Phys. 115 (12) (2014) 123510. [40] R. Prasher, T. Tong, A. Majumdar, An acoustic and dimensional mismatch model for thermal boundary conductance between a vertical mesoscopic nanowire/nanotube and a bulk substrate, J. Appl. Phys. 102 (10) (2007), pp. 104312e10. [41] S. Fukushima, T. Ohmichi, A. Maeda, H. Watanabe, The effect of yttrium content on the thermal conductivity of near-stoichiometric (u,y) o2 solid solutions, J. Nucl. Mater. 102 (1e2) (1981) 30e39. [42] R. Gibby, The effect of plutonium content on the thermal conductivity of (U, Pu) O2 solid solutions, J. Nucl. Mater. 38 (2) (1971) 163e177.