Continuum model for hydrogen pickup in zirconium alloys of LWR fuel cladding

Zirconium (Zr) alloys have been widely used as cladding materials for nuclear fuels in light-water reactors (LWRs). During the operation of the Zr cladded fuels, hydrogen generated by water and Zr corrosion reactions permeates through the protective Zr oxide layers, and diffuses and accumulates in Zr metal, potentially reaching or surpassing the hydrogen solubility limit in the Zr alloy.¹ This process is called hydrogen (H) pickup and it can lead to the formation of brittle hydrides that significantly reduce the ductility of Zr alloys.² Therefore, H pickup is one of the major issues potentially limiting the reliability and durability of cladding materials, especially under high burnups and accident conditions such as loss-of-cooling accidents and reactivity-initiated accidents.³ 

It is usually proposed that H pickup proceeds in three steps.^4,5 First, H₂O molecules adsorbed at the oxide/water interface react with anion oxygen vacancies to leave protons (H⁺) on the oxide surface. Some of the protons are discharged by electrons migrating from the oxide/metal interface and become H adsorbates (H_ad). Subsequently, the H_ad atom either reacts to form H₂ to be released as hydrogen gas or H_ad is absorbed into the oxide. Second, the absorbed hydrogen atoms, possibly along with protons, migrate through the barrier oxide layer and reach the Zr metal surface. Third, due to the H concentration gradient, hydrogen diffuses into Zr metal and hydride precipitates out when the hydrogen concentration is high enough.

Many studies have focused on understanding the transport of hydrogen through the barrier oxide layer, as this transport is often regarded as the rate-limiting step for H pickup.^5–7 A number of factors, including oxide morphology, alloy additive elements, and local stress, play important roles in the hydrogen transport process. It is widely accepted that the Zr oxide scale develops a double-layered structure during the oxidation.⁸ The outer layer is formed by porous oxide with cracks and pores that provide fast ingress routes for hydrogen, while the inner layer consists of dense oxide and it is usually regarded as a diffusion barrier.⁸ Recent transmission electron microscope (TEM) analysis found a suboxide (Zr₃O) region existing at the metal/oxide interface in some Zr alloys.⁹ According to density functional theory (DFT) calculations, the hydrogen migration energy in the suboxide is higher than that in pure Zr, so the suboxide layer may also slow down the H diffusion and contribute to the diffusion barrier.¹⁰ Using in-situ nuclear reaction analysis, Une et al. measured the deuterium concentration depth profile in oxide layers of Zr alloys corroded in the D₂O steam.⁷ The result shows a nearly flat concentration profile in the outside layer followed by a steeply decreasing concentration in the inner layer, which agrees well with the anticipated higher diffusivity in the porous oxide and lower diffusivity in the dense oxide. After growing to certain thickness, the dense oxide layer typically becomes porous (a change referred to as the “transition”) and the corrosion rate is suddenly increased.^8,11 Usually, a new dense oxide layer will start growing after the transition, so the oxidation of most Zr alloys shows a periodic feature. Recent studies have discovered that the H pickup process also follows the oxidation periodicity.^12,13 Besides the oxide morphology, additive elements (e.g., Fe, Cr, Nb, Sn, and Ni) in different Zr alloys have substantial effects on the H pickup fraction. These elements may either behave like trapping sites and directly decrease the hydrogen diffusivity in the oxide^14,15 or they may form second phase precipitates that have been hypothesized to be a preferred path for hydrogen migration or a source for pores or cracks formation.^16–18 Due to the lattice mismatch between Zr oxide and metal, high compressive stress is generated in the oxide near the interface.^19,20 Raman spectroscopy measurement revealed that the stress varies cyclically and can be as large as several GPa.⁷ DFT calculations found that under 1 GPa compressive stress hydrogen diffusion coefficient in tetragonal ZrO₂ is only about 60% of the coefficient without stress at 600 K.¹⁵ It is worth noting that there is also literature arguing that the diffusion of hydrogen through the barrier oxide layer is not the rate-limiting step.^21,22 Evidence from the chemical exchange experiments¹⁷ and transmission electron microscope (TEM) analysis²³ suggests that micropores/cracks exist even in the dense oxide layer. It has been proposed that hydrogen can penetrate to the oxide/metal interface via these flaws and the hydrogen cathodic reaction at the Zr metal surface is the rate-controlling process for H pickup. One technical difficulty in evaluating this hypothesis is that the observed crystallite boundary cracks or pores could also be formed by the TEM sample preparation process.⁸ Therefore, it is still an open question which process (or processes) is (are) rate-limiting for H pickup.

Previous studies have provided a number of insights into the detailed mechanism of H pickup and acquired a large body of data under various corrosion conditions. However, few theoretical models have been developed that take the advantage of the accumulated knowledge to describe the overall H pickup process quantitatively. In particular, the accuracy for predicting H pickup has not been assessed for even what might be considered as the simplest model, which assumes only a rate limiting process of diffusion in the oxide barrier layer. In this paper, we have therefore developed and assessed such a model. This type of a model is important as a baseline for more complex models that invoke additional phenomena, e.g., rate limiting surface reaction processes, second phase precipitation, crack and pore formation, strain effects, etc. Using up-to-date diffusivities and corrosion measurements, H pickup fractions (f_H) in six different alloys (Zry-2, GNF-Ziron, VB, Zry-4, ZIRLO, and Zr-2.5Nb) were calculated and compared with the experimentally measured f_H.^7,12 Here, f_H is defined as the fraction of hydrogen absorbed by the Zr metal to the total hydrogen generated during corrosion. f_H has been widely used for comparing resistances of various Zr alloys to the H uptake under different corrosion conditions (e.g., temperature, corrosion solution).^7,12 The alloys are chosen to include all alloys for which the necessary data on H diffusivity, time-dependent weight gain, and experimentally measured f_H needed by the model for comparing reason are simultaneously available. In this study, corrosions in 360 °C pure water (for Zry-4, ZIRLO, and Zr-2.5Nb), 400 °C steam, and 290 °C LiOH-containing water (for Zry-2, GNF-Ziron, and VB) were analyzed. In Ref. 12, the samples of Zry-4 and ZIRLO were processed in both sheet and tube forms in order to test whether the sample geometry could affect the H pickup process. These data are also included and compared in this work. Previous investigations have confirmed that the H pickup rate is significantly accelerated by the LiOH addition to water.²⁴ Further TEM analysis discovered that extended networks of degraded grain boundaries were formed from the oxide surface to near the metal/oxide interface, probably due to the preferential dissolution of zirconia in LiOH solution.⁷ In this case, the H pickup process is controlled by the dissociation reaction of H₂O at the front of the degraded grain boundaries, rather than the hydrogen diffusion process. Including the LiOH case aims at showing the limitation of the current model and indicating possible improvement for future work.

In our model, hydrogen diffusion through the dense barrier oxide layer is taken as the rate-limiting step for H pickup.^5–7 The hydrogen diffusion follows the equation

Here, $CHZrO2(x,t)$ is the hydrogen concentration in the barrier oxide layer at the distance x to the oxide surface and time t. As shown on the right side of the equation, the H flux contains two terms. The first term is the flux from the concentration gradient across the oxide film and the second term is due to the electric field generated by other migrating charge particles (e.g., oxygen ions, electrons) during the oxidation process. Previous studies suggest that at least part of hydrogen migrating through the oxide layer is charged,^12,25 so it is necessary to include the effect of electric filed from other charged particles on the H pick up process. Here, D_H is the hydrogen diffusion coefficient. As mentioned earlier, alloy additives have substantial effects on the H diffusivity in Zr oxides and the measured D_H for different Zr alloys can vary significantly. When multiple D_H values are available for a single Zr alloy, we chose the D_H value that is closest to the averaged D_H values among all the Zr alloys. The chosen D_H for the six Zr alloys being studied here are summarized in Table I.^14,26–29 The effect of different D_H values on the final calculated H pickup fraction is analyzed in the discussion part. The electric field-induced H flux is calculated based on the steady state Nernst-Planck equation. In Equation (1), q_H is the charge of H ion (+1 unit charge), k_B is the Boltzmann constant, and T is the corrosion temperature. E represents the electric filed across the Zr oxide film and is calculated by

Here, ρ is the electrical resistivity of the oxide (in MΩ cm) and j is the oxidation current density (in A/cm²). According to previous studies, the values of ρ vary as oxidation proceeds. Different Zr oxides can also have quite different ρ values with a range of about 2–138 MΩ cm.^30–32 In order to estimate the maximum possible effects of the electric field on the H pickup process, the maximum ρ (138 MΩ cm) among all the reported values in the literature is chosen in our calculation. We assume that the oxidation current density j is proportional to the oxidation rate. In Equation (2), wg(t) is the time-dependent weight gain of the Zr specimen. If the contribution of absorbed H to the weight gain is neglected, then dwg(t)/dt is the oxidation rate (in mg/dm² s). q_o is the charge of oxygen ion (+2 unit charge), N_A is Avogadro's constant, and Mo is the molar mass of the oxygen ion (15.9994 g/mol). Our calculations show that the contribution of the electric filed to the H pickup is negligible when compared to the contribution of concentration gradient, since the oxidation current density j decreases rapidly as the oxidation proceeds. The detailed analysis will be presented in Section III.

In Equation (1), L_b(t) is the time-dependent thickness of the barrier oxide layer. In order to solve the equation, the value of L_b(t) must be accessible at all times as an analytical function since Eq. (1) must be evaluated at different time steps for the numerical solution. Therefore, we need to fit an analytical form for L_b(t). It is difficult to directly measure the oxide thickness during the corrosion as the alloy specimen must be periodically taken out of the autoclave and analyzed, typically by electron microscopy. Instead, weight gains of the specimen as a function of corrosion time are usually reported in the literature.^12,21,22 The weight gain can be related to the oxide thickness based on the overall corrosion reaction $Zr+2H2O →ZrO2+2H2$⁠. If assuming the weight gain arises only from the added oxygen and the oxide has a uniform thickness, the time-dependent oxide thickness L(t) (including the protective barrier oxide layer L_b(t) and the non-protective porous layer) can be calculated by

where M_X is the molar mass of chemical X and $ρZrO2$ is the zirconia density (5.68 g/cm³). As similar to Equation (2), wg(t) is the time-dependent weight gain of the Zr specimen. For diffusion controlled growth, the oxidation kinetics follows a power law yielding^33,34

By fitting a series of weight gains measured at different times during the corrosion, the K and q can be obtained for all six alloys being studied here under steam or water corrosion. More specifically, the weight gain data for fitting K and q values of Zry-2, GNF-Ziron, and VB are from Ref. 7 and the data for fitting K and q values of Zry-4, ZIRLO, and Zr-2.5Nb are from Ref. 12. For the corrosion in LiOH-containing water, as the surface reaction of H₂O and Zr is the controlling process, the wg vs. time follows a simple linear relationship, so q is set equal to one and only K is fitted.⁷ The fitted coefficients are summarized in Table I. K has the same unit as wg (mg/dm²) and q is a numerical factor corresponding to corrosion time in days. As mentioned earlier, the entire oxide layer (L(t)) contains both the dense protective layer (L_b(t)) and the porous non-protective layer. To determine L_b, the nuclear reaction analysis was used to measure the deuterium concentration profiles in the oxide layer of Zry-2, GNF-Ziron, and VB corroded in 400 °C D₂O steam before the transition.^7,14 In the concentration profile, the region of a flat deuterium concentration is regarded as corresponding to the porous layer and the region of a decreasing concentration is regarded as being due to the dense protective oxide layer.^7,14 The measurement has found that the thickness of the dense oxide layer is about 0.53–0.60 of the entire oxide layer. A large number of theoretical and experimental analyses suggest that the dense oxide undergoes a transition when it reaches its maximum thickness L_bm.^7–9,35,36 Therefore in our calculation, it is assumed that before the transition, the barrier oxide layer grows with a constant thickness fraction (η) of the entire oxide layer, but after the transition, the dense oxide layer becomes porous and is no longer a barrier to H diffusion. We also assume that after the transition, a new dense oxide layer starts growing, following the same kinetics as before the transition. Our assumption is consistent with the periodic feature of Zr oxidation and the fact that a thicker oxide layer is usually associated with a superior resistance to oxidation and H pickup.³⁷ With all these assumptions, we get L_b(t) = ηL(t) for t < t_transition. Here, t_transition is the time when the oxidation transition happens, which is indicated by a sudden increase of the oxidation rate shown in the measured weight gain-time curve.^7,12 The values of η and L_bm for Zry-2, GNF-Ziron, and VB from Ref. 7 are listed in Table I. For Zry-4, ZIRLO, and Zr-2.5Nb, the values of L_bm and η are not reported and only the entire oxide thickness right before the transition was calculated in Ref. 12. For these three alloys, we assume that the fraction η is equal to the average value of η (0.57) of the other three Zr alloys in Ref. 7. For corrosion in LiOH, the dense oxide layer is very thin or possibly does not exist at all.^7,14 However, in order to compare with the water corrosion case, we treat the entire oxide layer thickness L(t) as a barrier layer for the H diffusion calculation. The failure of our model to explain the LiOH data both supports the model by showing it fails for systems where it does not include the correct physics and supports the hypothesis that little or no dense oxide exists in the LiOH system.

In order to solve Equation (1), boundary conditions at the H₂O/oxide interface (x = 0) and oxide/metal interface (x = L_b(t)) must be set. As the hydrogen diffusion is assumed to be the rate-limiting step, the hydrogen chemical potentials (μ_H) on both sides of each interface are treated as equal. Under the equilibrium condition, the boundary conditions are

where $CH,1atmZrO2$ is the hydrogen solubility in zirconia at the standard atmospheric pressure and $pH2$ is the effective H₂ partial pressure just outside the oxide. All concentrations in this work are given as mole fractions of the host unless stated otherwise. Equation (5) simply follows the classic Sievert's law, and both Equations (5) and (6) assume ideal-mixing behavior of the dissolved hydrogen in ZrO₂ and Zr, respectively. According to Ref. 38, $CH,1atmZrO2$ (unit: mol H/mol ZrO₂) equals to 2.78 × 10⁻⁴ at 400 °C, 3.80 × 10⁻⁴ at 360 °C, and 7.29 × 10⁻⁴ at 290 °C. It is worth mentioning that here $pH2$ is only an effective pressure that represents the activity of H for entering the oxide. It is not meant to represent a real gas pressure at the water/oxide interface or the overall activity of H in the surrounding water. The $pH2$ represents the H activity that comes from the detailed reactions and H generation occurring right at the surface of the oxide. Currently, no experimental data about the time-dependent $pH2$ or H activity at the oxide surface are available, so we hypothesize that $pH2$ is a single constant value during the entire corrosion process for all Zr alloys corroded by steam or pure water. The value of $pH2$ is fitted by minimizing the sum of the squares of the calculation error (defined as difference between the calculated f_H and the measured f_H). The fitted $pH2$ is 3.35 × 10⁶atm. Again, we note that this is an effective value representing the local H activity at the interface and potentially weakly related or unrelated to the H activity of the surrounding water. Furthermore, our calculations show that increasing or decreasing the $pH2$ value only makes the total calculation errors larger, but does not change the trend of the calculated f_H on which our conclusions are based. Similar fitting has been performed for corrosion in LiOH and the fitted $pH2$ equals to 6.10 × 10⁷atm. In Equation (6), $CHZr$ is the time-dependent hydrogen concentration in Zr alloy and can be calculated by

Here, d is the thickness of the Zr cladding, which is 600 μm for Zry-2, GNF-Ziron, and VB,⁷ and 800 μm for the remaining three alloys according to the real sample size.¹² The first term in Equation (7) is the intrinsic initial hydrogen concentration in the alloys and the second term represents the accumulated hydrogen due to the H flux from oxide into metal. Equation (7) assumes that no hydride precipitation occurs, which is consistent with the concentrations according to our model calculations. Based on Refs. 12 and 14, $CH,0Zr$ is 9 weight ppm for both Zry-2 and VB, 6 weight ppm for GNF-Ziron, and about 12.5 ppm for the remaining three alloys. For initial conditions, we take $CHZrO2=CH,1atmZrO2pH2$ for x = 0 and $CHZrO2=ρZrO2ρZrCH,0Zr$ for $0<x≤Lb$⁠, in which equations fulfill the equilibrium condition for hydrogen chemical potentials at the water/oxide and oxide/alloy interface, respectively. With the boundary and initial conditions, the time evolution of the hydrogen concentration profile in zirconia is solved using the standard finite difference method implemented by us in Matlab. Based on the obtained hydrogen concentration in Zr, the H pickup fraction is calculated by

Here, wg_exp(t) is the experimentally measured weight gain at time t. Note that the H pickup fraction in Equation (8) is always calculated with respect to an experimentally measured weight gain. The power law expression in Equation (4) for the weight gain is not used in Equation (8) and is only used to estimate the oxidation current density in Equation (2) and the barrier layer thickness in Equation (3). The time-dependent f_H calculated by our model is then compared with the experimentally measured f_H values at the same time. In the experiments, the H concentrations in the Zr alloys are measured at a certain corrosion time by either the vacuum hot extraction method or the cold neutron prompt gamma activation analysis. Subsequently, the H concentrations are converted into the H pickup fractions with the weight gain measured at the same time.^7,12

Figure 1 shows the calculated time-dependent H pickup fractions of all Zr alloys corroded in 400 °C steam or 360 °C water. In general, the pickup fractions increase with time, which trend agrees with the experiments. One exception is that for Zry-4 and ZIRLO, the f_H has a slight but a sharp decrease at the transition time (90 days for Zry-4 sheet, 120 days for ZIRLO sheet, and 135 days for both Zry-4 tube and ZIRLO tube samples) and f_H values for these alloys resume increasing afterwards. The sudden decrease is due to the fact that the oxidation rate, as well as the H generation rate, is suddenly accelerated at the oxidation transition. The same effect of the oxidation transition has also been observed in the experimentally measured H pickup fractions for these alloys, although the decrease is not as obvious as in the model and only a plateau of f_H appears before the oxide transition time.¹² For Zry-2, GNF-Ziron, VB, and Zr-2.5Nb, the increase in the oxidation rate at the transition is not that pronounced, so no similar f_H decrease is shown in the calculations, and no decrease or plateau of f_H was observed experimentally.^7,12

To directly compare the model calculations to the experimental measurements at the same corrosion time and conditions, the calculated fraction $fHcal$ vs. the experimentally measured fraction $fH exp$ are plotted in Figure 2(a). The first letter in the symbol represents different alloys and the number is the corrosion time in days. The letter with prime symbol means that the sample is in the tube form, and otherwise, it is in the sheet form. Since the data points of Zry-4, ZIRLO, and Zr-2.5Nb are concentrated in the low f_H corner, that region (marked by dotted line) is magnified in Figure 2(b). According to Figures 2(a) and 2(b), most of the data points are relatively close to the $fHcal$ = $fH exp$ dashed line, which would correspond to a perfect agreement between the modeling predictions and experimental values. The values of the $fHcal$ and $fH exp$ shown in Figure 2 are summarized in Table II. According to the data, the average absolute error, which is defined as the average of absolute difference between $fHcal$ and $fH exp$ for all Zr alloys and at all different time, is 4.9%, and the average relative error, which is defined as the average of absolute difference between $fHcal$ and $fH exp$ divided by $fH exp$ for all Zr alloys and at all different time, is 61.0%. Given the uncertainties in the experimental measurements and the input parameters (e.g., H diffusivity and dense oxide layer thickness), as well as the various simplifying assumptions in the model, the agreement between $fHcal$ and $fH exp$ is relatively good. This agreement shows that the dense oxide layer plays a significant role in the H pickup process, and suggests that hydrogen diffusion through the dense oxide barrier layer may be the rate-limiting process for H pickup in many situations.

As an example of how the model can provide insights on what mechanisms are dominant for H pickup, our calculations show that the effect of electric filed on the H pickup is negligible when compared to the effect of the H concentration gradient across the oxide layer. To demonstrate this point, we calculated the H pickup fractions without considering the electric filed-induced H flux (i.e., the second term on the right side of Equation (1)). The H pickup fractions with the effect of the electric filed (i.e., $fH, eleccal$⁠) and without the effect (⁠$fH, non−eleccal$⁠) are also summarized in Table II. The average relative difference between $fH, eleccal$ and $fH, non−eleccal$⁠, which is defined as the average of the absolute differences between $fH, eleccal$ and $fH, non−eleccal$ divided by $fH, eleccal$ for all six types of Zr alloys and at all different corrosion time, is only 6.9%. It is worth noting that in order to estimate the maximum possible effect of the electric field, the maximum oxide resistivity we found measured in the literature (138 MΩ cm) has been used in the calculation. If the average value of the oxide resistivity (70 MΩ cm) is used, the average difference between $fH, eleccal$ and $fH, non−eleccal$ will be only 3.5%. The effect of the electric field on the H pickup is significant at early times, e.g., within the first hour or so (in this paper, all times are relative to the start of oxidation), but becomes very small when considered over the more relevant longer time scales of many days. The effect of electric filed on the H pickup process is negligible over long times because the oxidation current density j, as well as the associated electric filed E, decreases rapidly to less values as the oxidation proceeds, since j is proportional to the time derivative of weight gain, which evolves with a negative exponent (q − 1) of time. For example, j in Zry-4 oxide is 5.44 × 10⁻⁵ A/cm² at 30 min, 3.32 × 10⁻⁵ A/cm² at 1 h, and just 3.48 × 10⁻⁶ A/cm² at 1 day. Correspondingly, the ratio of electric filed-induced H flux over the concentration gradient-induced H flux is 17.7 at 30 min, 0.16 at 1 h, and only about 0.02 at 1 day. This ratio remains lower than 2% in the rest of the corrosion process.

A careful examination of the distribution of the data points can reveal a few systematic discrepancies between the experimental and calculated f_H. First, as shown in Figure 2, the calculated H pickup fractions are higher than the experimental values for most data points of Zry-4 and ZIRLO (both sheet and tube samples), whereas lower for Zr-2.5Nb, GNF-Ziron, and VB. Second, as shown in Figure 2(a), the $fHcal$ of Zry-2 is lower than $fH exp$ for shorter corrosion time but higher than $fH exp$ for longer corrosion time. This means that the slopes of f_H with time are steeper in the model calculations than in the experiments. In fact, a similar error in the slope of f_H with time also exists for GNF-Ziron and VB. Third, some experiments showed that the speed of the H pickup suddenly increases just before the end of the oxidation transition.^12,39 However, such a phenomenon is not obviously observed in our calculated H pickup fractions.

Several possible reasons may contribute to the above discrepancies. First, the uncertainty in the measured H diffusion coefficients in various Zr alloy oxides is relatively large, even in the same Zr alloy. For example, according to Ref. 29, D_H = 1.81 × 10⁻¹⁹ m²/s for Zr−2.5Nb at 360 °C, whereas in Ref. 40 D_H is as large as 1.13 × 10⁻¹⁷ m²/s. This large difference may be caused by different techniques for measuring H concentration, details of sample preparations, and methods for deriving the diffusion coefficient. The large uncertainties of D_H can substantially influence the calculated H pickup fractions. A detailed estimation of the uncertainties of D_H and other input parameters, as well as their effects on the calculated H pickup fractions, is presented in the next paragraph. In addition, in our model, a fixed D_H is applied for calculating f_H in the entire H pickup process. This fixed D_H is usually derived by fitting the H concentration profiles measured at certain period of time after the corrosion starts.^29,41,42 Therefore, the obtained D_H is actually “averaged effective” diffusivity, while the real diffusivity in different regions of the oxide barrier layer or at different times may not be necessarily the same. For example, as mentioned in the introduction, the compressive stress in the oxide can decrease the H diffusivity. It is possible that during the initial oxidation, the stress in the thin oxide layer is still small so that the real H diffusivity is actually larger than the effective diffusivity, which would lead to a higher H pickup rate at short corrosion time. Second, an accurate model for the oxide growth is missing in our calculation. Here, we assume that the oxidation rate is proportional to the speed of weight gain and the weight gain follows the simple power law in Equation (4). However, the growth of the oxide layer may not necessarily follow this simple kinetics, especially during the oxidation transition. Third, due to the difficulty in measuring H activity at the corroding surface during the corrosion test, a constant effective H₂ partial pressure $pH2$ is used for all the f_H calculations. Here, the partial pressure $pH2$ is simply representing the activity of H to enter the oxide. As the rate of corrosion reaction changes during the H pickup process, it is quite possible that $pH2$ also varies with time. Finally, as mentioned in Section I, some other physical mechanisms that are not included in the current model may also substantially affect the H pickup process. For example, the formation of pores/cracks may provide a fast ingress route for H, and the water splitting and H cathodic reactions on the Zr surface may also be a rate-limiting step for the H pickup.^17,23 Furthermore, second phase precipitates can act as short-circuit paths for H through the dense oxide layer.¹¹ Given all these uncertainties in parameters and assumptions in oxidation kinetics, we do not claim that the current model can give quantitatively accurate predictions of H pickup fractions. However, this model provides useful qualitative guidance and acts as a baseline for more complex models with more accurate diffusion parameters, oxidation kinetics, and physical mechanisms.

Here, we analyzed the uncertainty of the model input parameters, including H diffusion coefficient, dense oxide layer thickness, and H₂ partial pressure on the oxide surface, and their effects on the calculated H pickup fractions. For the H diffusion coefficient, as summarized in Table I, the D_H values chosen in our calculations are between 1.0 × 10⁻¹⁸ m²/s and 1.0 × 10⁻¹⁹ m²/s (at T = 360 °C). However, the total range reported for D_H is somewhat larger (a summary of all H diffusion coefficients reported by the previous literature for various Zr alloy oxides can be found in Table III). Based on the values in Table III, we can calculate the standard deviation σ of ln(D_H) and the mean λ of ln(D_H). If we take D_H, Max = exp(λ + σ) and D_H, Min = exp(λ − σ), the range of D_H is (3.14 × 10⁻²⁰ m²/s, 7.06 × 10⁻¹⁸ m²/s) at T = 360 °C and (6.27 × 10⁻²⁰ m²/s, 1.13 × 10⁻¹⁷ m²/s) at T = 400 °C. This large uncertainty of D_H can substantially affect the calculation results. If D_H, Max is used, the $fHcal$ for all Zr alloys will quickly increase to nearly 100% within about 40 days. If D_H, Min is used, the $fHcal$ for all Zr alloys will remain nearly 0% during the entire corrosion process. Under these extreme cases, our model cannot give a reasonable prediction of the H pickup fractions. Although we chose the D_H values in Table I following a logical approach (described in Section II), the large uncertainty of D_H and its substantial effects on $fHcal$ suggest that we should regard the D_H as a partially fitted parameter. For the dense oxide layer thickness, previous studies show that the ratio of the dense oxide layer to the total oxide layer ranges from η = 0.53 to 0.60.⁷ By using this range, we can estimate that the potential uncertainty in the thickness of the dense oxide layer is within 0.2 μm for all kinds of Zr alloys studied here. This variance in the oxide layer thickness has relatively small effects on the $fHcal$⁠. More specially, if the maximum oxide thickness is used (based on η = 0.60), the average absolute error between the $fHcal$ and $fH exp$ is 4.6% and the average relative error is 56.2%. If the minimum oxide thickness is used (based on η = 0.53), the average absolute error is 5.3% and the average relative error is 68.4%. Both errors in each case are very similar to the original calculation results, where the absolute error is 4.9% and the relative error is 61.0%. The calculated H pickup fractions using different oxide thicknesses are listed in Table IV. For the effective H₂ partial pressure $pH2$⁠, a ten-fold cross validation was performed to evaluate its uncertainty range and the effects on the H pickup calculations. The original dataset was partitioned into ten groups. Each time, nine groups were selected as the training set to get the fitted $pH2$ using the same procedures as described in Section II, and then the fitted $pH2$ was applied to calculate the H pickup fractions in the final single group, which is called the validation set. This process was repeated ten times so each of the ten groups had been used as the validation set once. The detailed calculation results are summarized in Table V. The range of the fitted $pH2$ is between 2.6 × 10⁶atm and 3.7 × 10⁶atm, which are quite close to the original value of $pH2$ = 3.35 × 10⁶atm. The average absolute error of the validation set is 5.3% and the average relative error is 65.1%, which are also close to the corresponding values in the original calculations. We also plotted the $fHcal$ vs. $fH exp$ of the validation set in Figure 4. Most of the data points are relatively close to the $fHcal$ = $fH exp$ dashed line, which is also similar to Figure 2(a). Therefore, the cross validation shows that the $pH2$ is well constrained by the data and the likely errors from fitting have limited influence on the H pickup fractions. In summary, the H diffusion coefficient seems to affect the calculations most significantly. The large uncertainty range of D_H suggests that more high-quality diffusion data are necessary for fully assessing the accuracy of our model.

To further validate the model, we demonstrate that it fails where it is expected not to work. In Figure 3, we plot the calculated the H pick up fraction $fHcal$ as a function of the experimentally measured fraction $fH exp$ for LiOH solution corrosion experiments. The detailed calculation results are listed in Table VI. Here, the model is expected to show poor correlation to experiment, as no dense oxide layer is believed to form on the Zr surface in the LiOH solution.⁷ As expected, most calculated H pickup fractions are far from the experimental values, which demonstrates that when diffusion through the dense oxide is not a rate limiting process, our model does not predict reasonable f_H values. More quantitatively, the average absolute error is 23.1% and the average relative error is 138.5% in the LiOH solution case. Both errors are much larger than the corresponding errors in the water/steam case. The failure of the model when applied to the LiOH solution case (where other physical phenomena are expected to affect H transport) supports the assertion that the model captures real physics when it successfully matches experimental data in the steam/water corroded materials. In addition, the success of the model for the steam/water corroded systems further supports the notion that the dense oxide plays a significant role in controlling f_H in those measurements.

A continuum model for calculating the time-dependent H pickup fractions in different Zr alloys has been developed in this paper. To the best of our knowledge, all available experimental data that simultaneously measure time-dependent H pickup fractions and weight gain for Zr alloys of available H diffusivities in the oxide layer have been collected and compared with the f_H calculated by our model. The model's predictions match qualitatively with the experimental values of steam/water corrosion, which supports the hypothesis that hydrogen diffusion through the dense oxide barrier layer plays a significant role in affecting the H pickup process. The model breaks down when attempting to predict f_H in LiOH containing water corrosion, providing further evidence for the absence of dense oxide in the LiOH environment and supporting the model through demonstrating that it fails when expected. This model offers a primary framework for developing more sophisticated models in the future when more accurate parameters (e.g., H diffusivity and solubility, oxide barrier layer thickness, and effective H₂ partial pressure) are available, and for incorporating more physical mechanisms that may play an important role in the H transport (e.g., pores/cracks formation, water dissociation on oxide surface, second phase precipitation, and strain effects).

This research was primarily supported by the Consortium for Advanced Simulation of Light Water Reactors (http://www.casl.gov) and an Energy Innovation Hub (http://www.energy.gov/hubs) for Modeling and Simulation of Nuclear Reactors under U.S. Department of Energy Contract No. DE-AC05-00OR22725. The support for Xing Wang to collect additional data, perform analysis, refine the codes, and write parts of the manuscript was provided by the U.S. Department of Energy, Office of Basic Energy Sciences (Grant No. DE-FG02-08ER46493). We are very grateful to Adrien Couet for helpful discussions.