Zirconium Redistribution in U-Pu-Zr

Binary constituent redistribution model

We begin by considering a binary substitutional alloy, U-Zr for definiteness. Let $c_0$ be the density of lattice sites, which for simplicity we will assume is independent of composition and temperature, and let $x_{\text{Zr}}$, $x_{\text{U}} \ge 0$, and $x_{\text{Zr}}+ x_{\text{U}}=1$ be the atomic fractions of zirconium and uranium. In the presence of thermal gradients, the Zr atom flux in a single-phase region is, $${{\mathbf{J}}_{\text{Zr}}}=-c_0 D(\nabla {x_{\text{Zr}}}+S\nabla T),$$(1) where $D$ is the interdiffusion coefficient of Zr in the phase. The term involving $\nabla T$ is the so-called Soret effect or thermodiffusion. The coefficient $S(x_{\text{Zr}},T)$ may be positive or negative. This term contributes an advective component to the Zr flux directed toward lower temperatures when positive and higher temperatures when negative. To ensure physical compositions, i.e. $0\le x_{\text{Zr}} \le 1$, the flux ${{\mathbf{J}}_{\text{Zr}}}$ must tend to 0 as ${x_{\text{Zr}}}$ tends to either 0 or 1. This leads us to take as a leading order approximation, $$S={{x}_{\text{Zr}}}(1-{{x}_{\text{Zr}}})\frac{{{Q}^{\star }}}{R{{T}^{2}}}.$$(2) where ${Q^{\star}}$ is the heat of transport of Zr in the phase and $R$ is the gas constant. In the dilute limit ${x_{\text{Zr}}}\ll 1$, we recover the usual form $S={{x}_{\text{Zr}}}{{Q}^{\star }}/R{{T}^{2}}$ used in previous models. The additional factor $1-{{x}_{\text{Zr}}}$ ensures that ${{x}_{\text{Zr}}}\le 1$. Thus in a single- phase region the Zr flux is: $${{\mathbf{J}}_{\text{Zr}}}=-{{c}_{0}}{\mkern 1mu} D\left( \nabla {{x}_{\text{Zr}}}+{{x}_{\text{Zr}}}(1-{{x}_{\text{Zr}}})\frac{{{Q}^{\star }}}{R{{T}^{2}}}\nabla T \right),$$(3) where $D$ and ${{Q}^{\star}}$ are coefficients associated with the phase. It is important to recognize that the usual form of $S$ proportional to ${{x}_{\text{Zr}}}$ alone is not a fundamental law, but is itself merely a leading order approximation that is only valid for small ${{x}_{\text{Zr}}}$. This can be seen in two ways. First, as just noted, it gives a non-zero flux when ${{x}_{\text{Zr}}}=1$ that leads to unphysical solutions with ${{x}_{\text{Zr}}}>1$. Second, $S$ and its approximations should exhibit a certain symmetry. The problem could just as well be posed in terms of ${{x}_{\text{U}}}$ with the U flux ${{\mathbf{J}}_{\text{U}}}$ given by an expression analogous to Equation Eq. (1). On the other hand ${{\mathbf{J}}_{\text{U}}}=-{{\mathbf{J}}_{\text{Zr}}}$. Both should lead to the same flux, and reconciling the two leads to the symmetry. Our approximation Equation Eq. (2) exhibits this symmetry, and the necessary asymptotic behavior, while the usual form does not.

Next consider the flux in a 2-phase region very near a solubility curve defined by ${{x}_{\text{Zr}}}=X(T)$. For ${{x}_{\text{Zr}}}\le X(T)$ suppose the alloy is single- phase, but for ${{x}_{\text{Zr}}}>X(T)$ a second Zr-rich precipitant phase appears. Following Shewmon, P. C. (1958) we will assume that a local equilibrium between the two phases is maintained at each point in a temperature gradient through rapid adjustment of the phase fractions via local diffusion processes. This means the composition of the major continuous phase is fixed, ${{x}_{\text{Zr}}}=X(T)$ and as a consequence, $$\nabla {{x}_{\text{Zr}}}={X}'(T)\nabla T.$$(4)

Assuming that the flux of Zr (or counter flux of U) occurs only through the continuous phase, we then have, $${{\mathbf{J}}_{\text{Zr}}}=-{{c}_{0}}{\mkern 1mu} D\left( {X}'(T)+X(T)(1-X(T))\frac{{{Q}^{\star }}}{R{{T}^{2}}} \right)\nabla T,$$(5) where $D$ and ${{Q}^{\star}}$ are the coefficients associated with the continuous phase. For simplicity we consider solubility curves $X(T)$ satisfying, $$\frac{dX}{dT}=X(1-X)\frac{\Delta H}{R{{T}^{2}}},\quad X({{T}_{0}})={{X}_{0}},$$(6) where $\Delta H$, ${{X}_{0}}$, and ${{T}_{0}}$ are constant model parameters. Note that this equation uniquely determines the solubility curve. With this choice Equation Eq. (5) becomes, $${{\mathbf{J}}_{\text{Zr}}}=-{{c}_{0}}{\mkern 1mu} D{\mkern 1mu} X(T)(1-X(T))\frac{\Delta H+{{Q}^{\star }}}{R{{T}^{2}}}\nabla T.$$(7) This approach mimics Shewmon, P. C. (1958) who took, in the dilute limit ${{x}_{\text{Zr}}}\ll 1$, $$\frac{dX}{dT}=X\frac{\Delta H}{R{{T}^{2}}}.$$(8) The models developed in Karahan (2009), Kim et al. (2004), and Kim et al. (2006) also adopted this approach, but the expressions they use for $\Delta H$ are inconsistent with their phase diagrams, and have the wrong sign in several cases. If one desires to use arbitrary solubility curves, then equation Equation Eq. (5) must be used for the flux.

The mirror situation where the alloy is single phase for ${{x}_{\text{Zr}}}\ge X(T)$and two-phase for ${{x}_{\text{Zr}}}<X(T)$ is precisely the same. The Zr flux in the two- phase region very near the solubility curve $X(T)$ is also given by Equation Eq. (7) when $X(T)$ is of the form Equation Eq. (6), with the coefficients $D$ and ${{Q}^{\star}}$ associated with the continuous phase.

Constituent migration model for U-Pu-Zr

For U-Pu-Zr fuel it has been argued that a pseudo-binary treatment, in which the Pu fraction is assumed fixed, is justified both theoretically and on the basis of experimental data that show the Pu is largely immobile (Kim et al., 2004; Kim et al., 2006). Thus we consider a ternary U-Pu-Zr alloy where the Pu fraction is a fixed, spatially uniform value ${{x}_{\text{Pu}}}$. As in a true binary alloy, the remaining U and Zr constituents flow counter to each other, and if we set ${{\tilde{x}}_{\text{Zr}}}={{x}_{\text{Zr}}}/(1-{{x}_{\text{Pu}}})$ we may apply the preceding binary model to the relative Zr atom fraction ${{\tilde{x}}_{\text{Zr}}}$. Expressing the result in terms of the original ternary atom fraction ${{x}_{\text{Zr}}}$ we obtain, $${{\mathbf{J}}_{\text{Zr}}}=-{{c}_{0}}{\mkern 1mu} D\left( \nabla {{x}_{\text{Zr}}}+\frac {{{x}_{\text{Zr}}}(1-{{x}_{\text{Pu}}}-{{x}_{\text{Zr}}})}{1-{{x}_{\text{Pu}}}}\frac{{{Q}^ {\star }}}{R{{T}^{2}}}\nabla T \right),$$(9) for the Zr flux in a single-phase region, where the coefficients $D$ and ${{Q}^{\star}}$are associated with the phase, and, $${{\mathbf{J}}_{\text{Zr}}}=-{{c}_{0}}{\mkern 1mu} D{\mkern 1mu} \left( \frac{X(T)(1-{{x}_{\text{Pu}}}-X(T))}{1-{{x}_{\text{Pu}}}} \right)\left( \frac{\Delta H+{{Q}^{\star }}}{R{{T}^{2}}} \right)\nabla T,$$(10) for the flux in a 2-phase region very near the solubility curve ${{x}_{\text{Zr}}}=X(T)$, where the coefficients $D$ and ${{Q}^{\star}}$ are those associated with the continuous phase. The solubility curves are chosen to satisfy, $$\frac{dX}{dT}=\frac{X(1-{{x}_{\text{Pu}}}-X)}{1-{{x}_{\text{Pu}}}}\frac{\Delta H}{R{{T}^{2}}},\quad X({{T}_{0}})={{X}_{0}},$$(11) for some choice of parameters $\Delta H$, ${{X}_{0}}$, and ${{T}_{0}}$, which has the solution: $$X(T)=(1-{{x}_{\text{Pu}}}){{\left\{ 1+\frac{1-{{x}_{\text{Pu}}}-{{X}_{0}}}{{{X}_{0}}}\exp \left[ \frac{\Delta H}{R}\left( \frac{1}{T}-\frac{1}{{{T}_{0}}} \right) \right] \right\}}^{-1}}.$$(12)

We use the simple pseudo-binary U-19Pu-Zr phase diagram shown in Figure 2 which is based on the phase diagrams used in Karahan (2009), Kim et al. (2004), and Kim et al. (2006). There are two 2-phase regions, $\alpha + \delta$ and $\beta + \gamma$, depending on the temperature range. The nominal transition temperatures are ${{T}_{\alpha}}=595\text{ C}$ and ${{T}_{\beta}}=650\text{ C}$, but these two values were varied in the simulations presented in this manual. The parameters defining the solubility curves are given in Table 1.

The 2-phase flux Equation Eq. (10) only applies very near a solubility curve and supposes that all the flux occurs through the major continuous phase, ignoring the precipitant phase. This is a reasonable assumption near a solubility limit, but we need an expression for the flux that spans the entire 2-phase region from the Zr-poor phase solubility limit to the Zr-rich phase solubility limit. For this we use a simple weighted average of the two fluxes that apply near the solubility limits. For the $\beta + \gamma$ region we take $$\begin{aligned} {{\mathbf{J}}_{\text{Zr}}}=& -w\ {{c}_{0}}{{D}_{\beta }}\left( \frac{{{X}_{\beta }}(T)(1-{{x}_{\text{Pu}}}-{{X}_{\beta }}(T))}{1-{{x}_{\text{Pu}}}} \right)\left( \frac{\Delta {{H}_{\beta }}+Q_{\beta }^{*}}{R{{T}^{2}}} \right)\nabla T \\ & -(1-w)\ {{c}_{0}}{{D}_{\gamma }}\left( \frac{{{X}_{\gamma }}(T)(1-{{x}_{\text{Pu}}}-{{X}_{\gamma }}(T))}{1-{{x}_{\text{Pu}}}} \right)\left( \frac{\Delta {{H}_{\gamma }}+Q_{\gamma }^{*}}{R{{T}^{2}}} \right)\nabla T \end{aligned}$$(13) The weight factor $w$ is set equal to the $\beta$ phase fraction $f=({{X}_{\gamma }}(T)-{{x}_{\text{Zr}}})/({{X}_{\gamma }}(T)-{{X}_{\beta }}(T))$. The flux in the $\alpha + \beta$ region is defined analogously. Some computational experiments were done using the more sophisticated weighting $w=\exp{(a)} f^b/[ \exp{(a)} f^b + (1-f)^b]$ for parameters $a$ and $b>0$, however in the end we opted for the simple linear weighting ($a=0$, $b=1$).

Numerical Implementation

The preceding pseudo-binary model for U-19Pu-10Zr fuel was implemented into BISON with two essential numerical modifications: addition of an artificial diffusion term in the 2-phase region to stabilize the algorithm, and smoothing of the model coefficients near phase diagram curves. The full model with these modifications is described in ZirconiumDiffusion.