Phase-Field Model of Hydride Blister Growth Kinetics on Zirconium Surface

The precipitation of zirconium hydride blisters on the surface of zirconium alloy cladding destroys the integrity of the material and leads to material failure, which has serious potential safety hazards. To study the corrosion kinetics of zirconium hydride blisters, we established a phase-field model coupled with anisotropic elastic, which can reflect the microstructure evolution process of zirconium hydride blisters. The model studies the natural growth of hydride blister and the structural evolution process after applying radial stress and discusses the corrosion characteristics, stress distribution, and displacement changes. Zirconium hydride blisters tend to be semi-ellipsoidal in natural growth. Applying radial outward stress will promote the growth of blisters and aggravate the harm of corrosion. The stress state changes with applying stress, which affects the elastic driving force, resulting in the change of blister morphology. This work is helpful to understand the corrosion kinetic mechanism of hydride blisters.

Introduction

Zirconium alloy has a low neutron absorption rate and good high-temperature corrosion resistance and is widely used in structural materials in fuel cladding of water-cooled reactors (Krishnan and Asundi, 1981; Bair et al., 2015; Motta et al., 2015). However, severe waterside corrosion reactions occur during use and storage, and its performance will also be limited (Hong et al., 2002; Han et al., 2019). During the corrosion process on the waterside, the outer layer of the zirconium alloy cladding is oxidized, releasing hydrogen atoms, and the hydrogen atoms migrate into the zirconium matrix. It moves along the hydrostatic stress gradient and moves against the thermal gradient (Kammenzind and Berquist, 2000; Reheman and Ståhle, 2018). In other words, hydrogen is attracted to high-stress areas and cold areas, or areas where the two interact together. When the hydrogen concentration exceeds the solid solution limit in the zirconium matrix, zirconium and solid solute hydrogen transform into the zirconium hydride phase (Ghosal et al., 2002; Bair et al., 2015; Reheman and Ståhle, 2018; Han et al., 2019). The formation of brittle hydride will not only cause the material to expand and increase its volume greatly, but also it is easy to initiate crack growth and penetrate the matrix (Kim et al., 2007; Puls, 2009; McRae et al., 2010). The hydride formed on the metal surface appears as a small bump, usually called a blister (Cheong et al., 2001; Singh et al., 2002). The growth of zirconium hydride blister is an important form of pitting corrosion. It is the main initiation point of fuel cladding failure in nuclear power plants, which will pose a serious threat to nuclear power safety. Therefore, the basic understanding of hydride blister growth behavior is of great significance for predicting material performance degradation and failure.

The size of the hydride blister is much larger than that of the single hydride, reaching the micrometer level (Cheong et al., 2001; Singh et al., 2002; Long et al., 2017; Kim et al., 2018). The formation of the hydride blister requires the existence of a cold spot in the surrounding environment, and the effect of cold spot aggregation on the growth of hydride blisters is also different (Singh et al., 2002). In the experiment of (Domizzi et al., 1996), the growth of blisters on the surface of the Zr-2.5% Nb alloy tube was studied. It divides the observed phenomenon into three areas (I) the main part of the zirconium hydride blister, whose growth is similar to that of a single hydride, and the anisotropic growth is ellipsoidal, (II) the radial hydride accumulation area, and (III) circumferential hydride accumulation area. The phase-field method has been widely used to study the microstructure changes such as hydride precipitation in Zr cladding materials (Bair et al., 2015; Bair et al., 2016; Bair and Zaeem, 2017; Bair et al., 2017; Han et al., 2019; Toghraee et al., 2021). Based on region (I) (Reheman and Ståhle, 2018), studies the growth and cracking of hydride blisters under elastic-plastic action by using the phase-field method, but its model is established based on isotropy, which can not represent the anisotropic growth of zirconium hydride blisters. In this work, we study the morphological evolution of zirconium hydride blisters by using the phase-field method coupled with anisotropic elastic strain energy. In addition, some studies (Cheong et al., 2001; Motta and Chen, 2012) show that during reactor operation, the fission gas and the helium filling gas in the cladding tube create a combined pressure that exposes the nuclear fuel cladding to constant stress. Even the swelling of nuclear fuel also produces a certain pressure on the zirconium alloy cladding tube (Rozhnov et al., 2011). Therefore, it is important to understand the evolution mechanism of zirconium hydride blisters under specific pressure. The evolution of hydride blisters under radial stress is considered in this work.

Phase-Field Model

Combining the zirconium hydrogen phase diagram (Yu et al., 2020) and experimental research (Santisteban and Domizzi, 2009; Vicente Alvarez et al., 2011; Lin et al., 2016; Motta et al., 2019), it can be found that the hydrides in the zirconium hydride blister are all δ-hydrides, and the ratio varies approximately in the range of 60–85%, and the rest is the Zr matrix and impurities. Since the chemical energy parameters of the zirconium hydride blister are not clear yet, we simulate the zirconium hydride blister by using the chemical energy (Bair et al., 2016; Bair and Zaeem, 2017; Bair et al., 2017; Toghraee et al., 2021) and elastic parameters (Bair et al., 2017; Han et al., 2019) of the δ-hydride. The size of the cladding tube in the model is similar to some studies (Hong et al., 2002; Hellouin de Menibus et al., 2014), and 11.25-degree section is simulated. Since the longitudinal plane of the blister is usually observed in the experiment, we convert the data in the z-y direction of the coordinate axis to x-y. The whole calculation process is carried out in two dimensions. Figure 1 shows the setting of the hydride blister growth model and boundary conditions. For the displacement field $u$, the top boundary is set as free boundary conditions, the left and right are periodic boundary conditions, and the bottom is fixed boundary conditions. The boundary settings of phase-field η and concentration field c are the same. All boundaries are set as zero flux boundary conditions, and hydride and the matrix phase are distinguished by the diffuse interface.

FIGURE 1

FIGURE 1. Schematic diagram of zirconium-zirconium hydride blister model and boundary condition settings.

Governing Equations and Bulk Free Energy

Interfacial reaction, hydrogen transport driven by diffusion potential, and internal stress caused by precipitation and growth of zirconium hydride blister. Multi-physical field coupling based on dynamic equations can simulate these effects and reflect the dynamic characteristics (Shi et al., 2021; Liu et al., 2022; Shi et al., 2022). The phase field, hydrogen concentration field, and elastic strain energy field were coupled to establish the zirconium hydride blister growth model in this study. These fields are necessary to track the movement of the diffusion interface, calculate the diffusion/redistribution of hydrogen (Yang et al., 2021), and reveal the effect of radial stress on hydride growth. The advantage of the phase-field model in this study is that the shape of the grown precipitates can be calculated (Gao et al., 2021; Gao et al., 2022; Xu et al., 2022). The phase transition from metal to hydride is given by a phase variable η that continuously covers a range of values, where the ends of the range represent pure metal or pure hydride, respectively. The analysis is based on totally free energy, including strain energy, chemical potential energy, and gradient energy describing the diffusion interface. In the phase-field model, both the diffusion of hydrogen atoms and the precipitation of hydrides are considered. The time evolution of the concentration c is determined by the Cahn-Hilliard diffusion equation (Kim et al., 2008), and the evolution of the microstructural parameter η is governed by the time-dependent Allen-Cahn equation (Allen and Cahn, 1979). In order to solve the evolution equation, the finite element method is adopted, and the mesh adaptive method is used to speed up the calculation.

$\frac{\mathrm{\partial }\mathit{\eta }}{\mathrm{\partial }\mathit{t}}=-\mathit{L}\frac{\mathit{\delta }{\mathit{F}}_{\mathbf{total}}}{\mathit{\delta \eta }}(1)$

$\frac{\partial \mathit{c}}{\partial \mathit{t}}=\nabla \cdot \mathit{M}\nabla \frac{\mathit{\delta }{\mathit{F}}_{\mathbf{total}}}{\mathit{\delta c}}(2)$

The driving force of the zirconium-hydride blister phase-field model comes from the minimization of the total free energy ($F_{\text{total}}$) of the system, which includes the uniform bulk chemical free energy ($F_{\text{bulk}}$) the interface gradient energy ($F_{\text{int}}$) and the elastic energy ($F_{\text{el}}$). For simplicity, it is assumed that the diffusion coefficients of hydrogen in the matrix and precipitate are the same.

${\mathit{F}}_{\mathbf{total}}={\mathit{F}}_{\mathbf{bulk}}+{\mathit{F}}_{\mathbf{int}}+{\mathit{F}}_{\mathbf{el}}(3)$

Here, all the parameters used in this simulation are listed in Table. 1 The bulk chemical free energy $F_{\text{bulk}}$ is a function of the concentration $c$ and the order parameter $\eta$ at T = 600 k, while $F_{\text{int}}$ is the gradient energy due to the inhomogeneity of the interface and $F_{\text{el}}$ is the elastic strain energy during the Zr matrix phase to the hydride phase transition.

${\mathit{F}}_{\mathbf{bulk}}={\displaystyle \underset{\mathit{V}}{\int }{\mathit{f}}_{\mathbf{bulk}}\mathit{dV}}(4)$

TABLE 1

TABLE 1. Model parameters and conditions used in the simulations.

The Kim-Kim-Suzuki (KKS) binary alloy model (Kim et al., 1999) is used to determine the expression form of bulk chemical free energy density. The KKS model assumes that each material point is considered to be a mixture of two or more phases with different concentrations, and the diffusion potentials of the same concentration between different phases are equal. The assumption can be expressed as follows:

$\mathit{c}=\mathit{h}\left(\mathit{\eta }\right){\mathit{c}}_{\mathit{\alpha } }+\left[1-\mathit{h}\left(\mathit{\eta }\right)\right]{\mathit{c}}_{\mathit{\delta }}(5)$

$\frac{\mathrm{\partial }{\mathit{f}}_{\mathit{\alpha }}\left({\mathit{c}}_{\mathit{\alpha }}\right)}{\mathrm{\partial }{\mathit{c}}_{\mathit{\alpha }}}=\frac{\mathrm{\partial }{\mathit{f}}_{\mathit{\delta }}\left({\mathit{c}}_{\mathit{\delta }}\right)}{\mathrm{\partial }{\mathit{c}}_{\mathit{\delta }}}(6)$

${\mathit{f}}_{\mathbf{bulk}}\left(\mathit{c},\mathit{\eta }\right)=\mathit{h}\left(\mathit{\eta }\right){\mathit{f}}_{\mathit{\alpha }}\left({\mathit{c}}_{\mathit{\alpha }}\right)+\left[\mathrm{1}-\mathit{h}\left(\mathit{\eta }\right)\right]{\mathit{f}}_{\mathit{\delta }}\left({\mathit{c}}_{\mathit{\delta }}\right)+\mathit{wg}\left(\mathit{\eta }\right)(7)$

where $h\left(\eta \right)=3{\eta }^2-2{\eta }^3$ is the interpolation function, the double potential well function $g\left(\eta \right)=\frac{1}{2}\left({\eta }^2\times {\left(1-\eta \right)}^2\right)$ , and $w$ is the potential barrier. The free energy density in this study is constructed according to the quadratic fitting method adopted by (Bair et al., 2016; Bair and Zaeem, 2017; Bair et al., 2017; Toghraee et al., 2021). Taking advantage of the equal slope at the two-phase equilibrium point, the two-phase free energy curve can be constructed by using the artificial parameters $A_{\alpha }$ and $c_{\alpha }$.

Two-phase molar free energy $G_{\alpha }$ and $G_{\delta }$ at T = 600 K can be expressed as:

${\mathit{G}}_{\mathit{\alpha }}={\mathit{A}}_{\mathit{\alpha }}{\left({\mathit{c}}_{\mathit{\alpha }}-{\mathit{c}}_{\mathit{\delta }}^{\prime }\right)}^2(8)$

${\mathit{G}}_{\mathit{\delta }}=\Delta {\mathit{G}}_{\mathit{\delta }}^0\left({\mathit{A}}_{\mathit{\delta }}{\left({\mathit{c}}_{\mathit{\delta }}-\left({\mathit{c}}_{\mathit{\delta }}^{\prime }-0.05\right)\right)}^2-1\right)(9)$

$\Delta {\mathit{G}}_{\mathit{\delta }}^0=46.37-0.0414\mathit{T}\left[\frac{\mathit{KJ}}{\mathit{mol}}\right](10)$

${\mathit{c}}_{\mathit{\alpha }}^{\prime }=3.27\times {10}^4\mathbf{exp}\left(-\frac{25024}{\mathit{RT}}\right)(11)$

Here, $c_{\alpha }$ and $c_{\delta }$ are the molar fractions of atomic H in the matrix and blister. $c_{\alpha }^{\prime }$ is the concentration of maximum hydrogen solubility in the matrix before precipitation occurs (Une et al., 2009); $c_{\delta }^{\prime }$ is the hydrogen concentration in the hydride blister; $\Delta G_{\delta }^0$ is the molar Gibbs free energy of hydride blister formation. The constant $A_{\alpha }$ and $A_{\delta }$ are the contribution of the free energy of the two phases to the interface, which are determined by constructing a common tangent between the molar Gibbs free energy curves of the two phases (Yang et al., 2020; Yang et al., 2021; Sheng et al., 2022), from which an approximate expression of the free energy can be obtained. $f_{\alpha }$ and $f_{\delta }$ are the free energy densities of the zirconium matrix and hydride blister, respectively, as shown in Figure 2 The blue dashed line is the common tangent when the two phases are in equilibrium.

${\mathit{f}}_{\mathit{\alpha }}=\frac{{\mathit{G}}_{\mathit{\alpha }}}{{\mathit{V}}_{\mathit{m}}}=\mathrm{14}{\mathit{c}}_{\mathit{\alpha }}^{\mathrm{2}}-\mathrm{2}\mathrm{.80}{\mathit{c}}_{\mathit{\alpha }}\mathrm{+}0\mathrm{.14}(12)$

${\mathit{f}}_{\mathit{\delta }}=\frac{{\mathit{G}}_{\mathit{\delta }}}{{\mathit{V}}_{\mathit{m}}}=\mathrm{36}\mathrm{.9683}{\mathit{c}}_{\mathit{\delta }}^{\mathrm{2}}-\mathrm{46}\mathrm{.5801}{\mathit{c}}_{\mathit{\delta }}\mathrm{+}13\mathrm{.1349}(13)$

${\mathit{F}}_{\mathbf{int}}={\displaystyle \underset{\mathit{V}}{\int }{\mathit{f}}_{\mathbf{int}}\mathit{dr}}(14)$

FIGURE 2

FIGURE 2. Two-phase chemical free energy density curves of zirconium matrix and zirconium hydride blister.

The interfacial energy density can generally be written as a function of the gradient of the phase-field variable,

${\mathit{f}}_{\mathbf{int}}=\frac{\mathrm{1}}{\mathrm{2}}{\mathit{\kappa }}_{\mathit{\eta }}{\left(\nabla \mathit{\eta }\right)}^{\mathrm{2}}(15)$

where ${\kappa }_{\eta }$ is the gradient energy coefficient related to the phase-field variable $\eta$.

Elastic Strain Energy and Applied Stress

When hydride blister precipitate and grow in the matrix, their elastic strain energy comes from the structural difference between hydride blister and matrix. The elastic strain energy $F_{e\text{l}}$ is expressed as follows:

${\mathit{F}}_{\mathbf{el}}={\displaystyle \underset{\mathit{V}}{\int }{\mathit{f}}_{\mathbf{el}}\mathit{dr}}=\frac{\mathrm{1}}{\mathrm{2}}{\displaystyle \underset{\mathit{V}}{\int }{\mathit{\sigma }}_{\mathit{ij}}{\mathit{\epsilon }}_{\mathit{ij}}^{\mathit{el}}\mathit{dr}}(16)$

where $f_{\text{el}}$ is the mismatch elastic strain energy density, ${\sigma }_{ij}$ and ${\mathit{ϵ}}_{ij}^{el}$ are the stress and elastic strain, respectively. Considering that the cladding tube will be subjected to the radial outward stress from the swelling of nuclear fuel and so on, we add a new cylindrical coordinate system and add the radial outward stress ${\sigma }_{app}$ to the model to study the evolution difference of the blister under different stresses. The radial stress ${\sigma }_{app}\mathit{ }$ is added through the stress coordinate conversion in elasticity.

${\mathit{\sigma }}_{\mathit{ij}}={\mathit{\sigma }}_{\mathit{app}}+{\mathit{C}}_{\mathit{ijkl}}\left(\mathit{\eta }\right){\mathit{\epsilon }}_{\mathit{ij}}^{\mathit{el}}(17)$

${\mathit{\sigma }}_{\mathit{app}}=\{\begin{array}{c}{\mathit{\sigma }}_{\mathit{x}}=\frac{{\mathit{\sigma }}_{\mathit{r}}+{\mathit{\sigma }}_{\mathit{\phi }}}{2}+\frac{{\mathit{\sigma }}_{\mathit{r}}-{\mathit{\sigma }}_{\mathit{\phi }}}{2}\mathit{cos}\mathrm{2}\phi -{\mathit{\tau }}_{\mathit{r\phi }}\mathit{sin}2\phi \\ {\mathit{\sigma }}_{\mathit{y}}=\frac{{\mathit{\sigma }}_{\mathit{r}}+{\mathit{\sigma }}_{\mathit{\phi }}}{2}-\frac{{\mathit{\sigma }}_{\mathit{r}}-{\mathit{\sigma }}_{\mathit{\phi }}}{2}\mathit{cos}\mathrm{2}\phi +{\mathit{\tau }}_{\mathit{r\phi }}\mathit{sin}2\phi \\ {\mathit{\tau }}_{\mathit{xy}}=\frac{{\mathit{\sigma }}_{\mathit{r}}-{\mathit{\sigma }}_{\mathit{\phi }}}{2}sin\mathrm{2}\phi +{\mathit{\tau }}_{\mathit{r\phi }}\mathit{sin}\mathrm{2}\phi \end{array}(18)$

${\sigma }_r$ and ${\sigma }_{\phi }$ are the radial stress and circumferential stress applied in the cylindrical coordinate system respectively. $\phi$ is the angle value in the cylindrical coordinate system. ${\sigma }_x$, ${\sigma }_y$, ${\tau }_{xy}$ are the stress value in the rectangular coordinate system.

${\mathit{\epsilon }}_{\mathit{ij}}^{\mathit{el}}={\mathit{\epsilon }}_{\mathit{ij}}\left(\mathit{r}\right)-{\mathit{\epsilon }}_{\mathit{ij}}^{00}\left(\mathit{r}\right)(19)$

${\mathit{\epsilon }}_{\mathit{ij}}\left(\mathit{r}\right)=\frac{1}{2}\left(\frac{\partial {\mathit{u}}_{\mathit{i}}\left(\mathit{r}\right)}{\partial {\mathit{r}}_{\mathit{j}}}+\frac{\partial {\mathit{u}}_{\mathit{j}}\left(\mathit{r}\right)}{\partial {\mathit{r}}_{\mathit{i}}}\right)(20)$

${\mathit{\epsilon }}_{\mathit{ij}}^{00}\left(\mathit{r}\right)=\mathit{h}\left(\mathit{\eta }\right){\mathit{\epsilon }}_{\mathit{ij}}^{00}\left(\mathit{r}\right){\mathit{\delta }}_{\mathit{ij}}=\left[\begin{array}{ccc}\mathit{h}\left(\mathit{\eta }\right){\mathit{\epsilon }}_{11}^{00}& 0& 0\\ 0& \mathit{h}\left(\mathit{\eta }\right){\mathit{\epsilon }}_{22}^{00}& 0\\ 0& 0& \mathit{h}\left(\mathit{\eta }\right){\mathit{\epsilon }}_{33}^{00}\end{array}\right](21)$

The elastic strain ${\epsilon }_{ij}^{el}$ can be obtained from Equation 19, where ${\epsilon }_{ij}\left(r\right)$ is the total strain represented by the components $u_i$ and $u_j$ of the displacement field variable $u$ obtained by solving the mechanical equilibrium Equation 23 and ${\epsilon }_{ij}^{00}\left(r\right)$ is related to the lattice inherent strain ${\epsilon }_{ij}^{00}$. Tummala (Tummala et al., 2018) and Singh (Singh et al., 2007) give the eigenstrain of hydride at different temperatures. The eigenstrain of the hydride blister at T = 600 k can be obtained.

${\mathit{\epsilon }}_{\mathit{ij}}^{00}=\{\begin{array}{c}{\mathit{\epsilon }}_{11}^{00}={\mathit{\epsilon }}_{22}^{00}=0.03888+2.315\times {10}^{-5}T\\ {\mathit{\epsilon }}_{33}^{00}=0.06646+1.9348\times {10}^{-5}T\end{array}(22)$

$\frac{\partial {\mathit{\sigma }}_{\mathit{ij}}^{\mathit{el}}}{\partial {\mathit{r}}_{\mathit{j}}}=0(23)$

${\mathit{C}}_{\mathit{ijkl}}\left(\mathit{\eta }\right)=\mathit{h}\left(\mathit{\eta }\right){\mathit{C}}_{\mathit{ijkl}}^{\mathit{\delta }}+\left[1-\mathit{h}\left(\mathit{\eta }\right)\right]{\mathit{C}}_{\mathit{ijkl}}^{\mathit{\alpha }}={\mathit{C}}_{\mathit{ijkl}}^{\mathit{\delta }}+\mathit{h}\left(\mathit{\eta }\right)\mathrm{\Delta }{\mathit{C}}_{\mathit{ijkl}}(24)$

where $C_{ijkl}\left(\eta \right)$ is the elastic constant, a function of the order parameter $\eta$, $C_{ijkl}^{\alpha }$ and $C_{ijkl}^{\delta }$ are the elastic constants of the matrix and precipitate, respectively, and $\Delta C_{ijkl}$ is the difference between the elastic constants of $C_{ijkl}^{\alpha }$ and $C_{ijkl}^{\delta }$.

Results and Discussion

The effect of anisotropic elastic energy on the growth morphology of hydride blisters is explored, as shown in Figures 3A,C, which are the evolution results of hydride blisters after initial state t = 0 s and final state t = 1500 s, respectively. Due to the effect of elastic energy, the corrosion width of hydride blister is greater than the corrosion depth, the growth characteristics of anisotropy appear, and the morphology is semi ellipsoid. In addition, there is a slight bulge at the top of the blister. The growth of the hydride blister is closely related to the intrinsic strain, and to minimize the elastic properties, different growth rates will appear in different directions. Figures 3B,D show that the concentration in the hydride blister is relatively high, and there is a diffusion area at the interface between the matrix and the hydride blister. Figures 3E,F are the stress distribution results in the radial and circumferential directions of the coupled anisotropic elastic energy, respectively. The hydride blisters in both directions are in a compressed state, and there are stress concentration area near the phase interface and the tips on both sides. The stress concentration will redirect the δ-hydrides around the hydride blisters, form more harmful radial hydride clusters, and even cause the initiation and propagation of microcracks.

FIGURE 3

FIGURE 3. The evolution results of the coupled anisotropic elastic energy zirconium hydride blister, where (A) and (C) are the morphology evolution results (order parameter field) of the hydride blister growth at t = 0 s and t = 1,500 s, (B), and (D) is the distribution of hydrogen concentration field after t = 0 s and t = 1,500 s, (E) and (F) are the stress distribution of ${\sigma }_{rd}$ (Radial stress) and ${\sigma }_{\phi }$ (Circumferential stress).

Zirconium alloy cladding tube will be subjected to a certain radial outward stress during operation, which will have a great impact on the formation and growth of hydride. We applied stress of 1, 2, and 3 GPa along the radial direction in this model. Because this model is established under the ideal condition of no defects, the applied stress is much greater than that of the actual cladding tube. We can qualitatively study the effect of radial stress on the growth morphology of hydride blisters. The evolution results after t = 1,500 s are shown in Figure 4. Figures 4A–C shows that the width and depth of hydride blister increase with the increase of stress under the application of different radial stress. Figure 4D shows the contour comparison diagram of hydride blister under t = 1,500 s applied stress and no applied stress. The results in Figure 4 confirm that applying radial stress will change the morphology and dynamic evolution of blisters, promote the growth of blisters and enhance the process of pitting corrosion.

FIGURE 4

FIGURE 4. At t = 1,500 s, (A–C) are the growth morphologies of hydride blisters under 1, 2, 3 GPa radial stress, respectively (The direction of stress application is indicated by blue arrows). (D) is the contours of hydride blisters.

Figure 5A shows the change of corrosion depth of the hydride blister with time. It can be seen intuitively that the growth of hydrogenated blister is affected by different radial stresses, and the corrosion depth of hydride blister will increase with the increase of time and radial outward stress. Figure 5B shows the changes of the area of hydride blister and the average corrosion rate along the pitting depth under different radial stresses. It can be found that the greater the applied radial stress is, the more the corrosion area will be increased, which also shows that the radial stress will increase the influence of pitting corrosion. The corrosion rate changed from fast to slow, which may be due to the continuous consumption of hydrogen in the surrounding environment during the evolution process. The reduction of hydrogen content has become one of the key factors restricting its growth. After the stress is applied, the initial growth speed becomes larger, and the larger the value is, the more obvious the speed increase is, and then it will tend to be consistent. The matrix has a compressive effect on the expansion and growth of hydride blisters. The expansion state of the hydride blister in the simulation area is realized by the bending of the grid (Sheng et al., 2022). After the external stress is applied, the external stress will affect the overall stress state of the blister, and then affect its growth morphology. The growth of hydride bubbles is controlled by these three effects. Figure 5C,D show the distribution of displacement $u_x$, $u_y$ and η along the red arrow line after applying different radial stresses at t = 1,500 s. On both sides of the hydride blister, $u_x$ has the same size and changes in the opposite direction. The application of stress will change the initial distribution of displacement, but will not affect its change trend. Combined with the change of order parameters, it can be found that the displacement changes rapidly at the matrix-blister interface and reaches the maximum at the tip of the blister. The application of stress will promote the movement of the interface, thus changing the change of displacement at the same grid point, showing a larger expansion trend on both sides. In the depth direction of the hydride blister, when no stress is applied, as shown in the black line in the figure, the displacement $u_y$ first starts from zero (Due to the displacement of the lower surface $u_y$ is constrained) and then the value starts to increase near the phase interface, which means that this part of the matrix point moves along the z-axis and reaches a maximum at the upper free surface. The protrusion of hydride blister can also be realized by this part. This part can realize the bulge of hydride blister, which is the reason for the bulge in Figure 3. When the radial stress is applied, the stress distribution of the hydride blister and the matrix changes, and the displacement will change from the fixed bottom along the direction of the red arrow. The positive displacement changes after the stress are applied. A general trend is a form of “check”, and the inflection points appear at the interface.

FIGURE 5

FIGURE 5. (A) The corrosion pit depth, (B) the average corrosion rate and the area of hydride blister as a function of time under different radial stresses. Variation of (C) displacement component $u_x$ (solid line) and order parameter η (dotted line), (D) displacement component $u_y$ (solid line) and order parameter η (dotted line) along the red arrow in the diagram under different radial stresses at t = 1,500 s (The displacement of the green dotted line and the order parameter of the light blue dotted line are zero).

As mentioned above, applying different of radial outward stress to the whole simulation area will have a great impact on the growth of the hydride blister. Figure 6A and Figure 6B are t = 1,500 s respectively, under the application of 2 GPa radial stress ${\sigma }_{rd}$ and ${\sigma }_{\phi }$ stress distribution. Figure 6C and Figure 6D are t = 1,500 s respectively under different stress conditions ${\sigma }_{rd}$, ${\sigma }_{\phi }$ and η with the change of the black arrow in the diagram. Combined with the evolution process of order parameters, it can be found that when hydride blisters grow normally, the matrix bears compressive stress. The stress value changes rapidly at the interface and reaches the maximum, and there is a stress concentration area around the blister. The radial stress ${\sigma }_{rd}$ reverses at the interface, and part of the matrix bears tensile stress. The tensile effect caused by tensile stress is very weak, which can not prevent the expansion and compression of hydride. The circumferential stress ${\sigma }_{\phi }$ has a large stress concentration at the tip of the blister, which may lead to the cracking of the blister. It can be found from the change of the displacement $u_y$ in Figure 5D. When we apply stress and change the initial stress distribution, the overall stress of blister and matrix will change from compression to tension. No matter how the applied force changes, there will be stress concentration areas at the tip and interface of the hydride blister. The side of the upper surface of the blister close to the oxide layer will bear outward and upward stresses, accelerating the cracking of the upper surface. The growth process of zirconium hydride blisters is the result of the combined effect of chemical energy driving force and elastic energy driving force (Simon et al., 2021). The applied stress affects the growth of hydride blisters by affecting the elastic driving force, and the radial stress reduces the effect of the elastic driving force hindering the growth of hydrides. This leads to the difference in the morphology of hydride blisters under different stress states.

FIGURE 6

FIGURE 6. The distribution of (A) ${\sigma }_{rd}$ (Radial stress) and (B) ${\sigma }_{\phi }$ (Circumferential stress) under 2 GPa radial stress at t = 1,500 s. Variation of (C) ${\sigma }_{rd}$ and η, (D) ${\sigma }_{\phi }$ and η along the black arrow in the figure under different radial stresses at t = 1,500 s.

Conclusion

In this work, the anisotropic elastic properties are introduced into Zr-Hydride blister system, and the growth and evolution dynamics of zirconium hydride blister are studied by using phase-field model. The simulation results show that the natural growth morphology of hydride blister tends to be semi-ellipsoidal and in a compressed state, while the matrix near the interface is in a tensile state. There is a stress concentration area on the outer contour of the blister, and the maximum value appears at the end points on both sides of the hydride. The external stress affects the growth morphology of zirconium hydride blister by affecting the elastic stress distribution. Applying radial outward stress will promote the growth of blister, and its corrosion width and depth will increase. The greater the applied stress, the more obvious the promotion effect and the faster the growth rate. This model also has certain application value in other fields of pitting corrosion.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.

Author Contributions

HH and YL conceived and guided the work. SW conducted simulation and data integration. The results were analyzed and discussed by all authors. JS gave guidance in the theoretical part. SW, JS, HH and YL jointly wrote the manuscript, with contributions from all the authors.

Funding

This work was supported by the Science Challenge Project (Grant No. TZZT2019-D1-03) and the Foundation of LCP.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

Authors would like to thank Changqing Guo, Yuechao Wang and Guomin Han for their help on the model and basic knowledge in this work.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Allen, S. M., and Cahn, J. W. (1979). A Microscopic Theory for Antiphase Boundary Motion and its Application to Antiphase Domain Coarsening. Acta Metall. 27 (6), 1085–1095. doi:10.1016/0001-6160(79)90196-2

CrossRef Full Text | Google Scholar

Bair, J., and Zaeem, M. A. (2017). Effects of Applied Strain on Formation, Shape Evolution, and Reorientation of Multiphase Zirconium Hydrides: A Multiphase Field Modeling Study. arXiv e-prints. 2017. doi:10.1016/j.commatsci.2021.110367

CrossRef Full Text | Google Scholar

Bair, J., Zaeem, M. A., and Tonks, M. (2016). A Phase-Field Model to Study the Effects of Temperature Change on Shape Evolution Ofγ-Hydrides in Zirconium. J. Phys. D Appl. Phys. 49 (40). 405302. doi:10.1088/0022-3727/49/40/405302

CrossRef Full Text | Google Scholar

Bair, J., Asle Zaeem, M., and Schwen, D. (2017). Formation Path of δ Hydrides in Zirconium by Multiphase Field Modeling. Acta Mater. 123, 235–244. doi:10.1016/j.actamat.2016.10.056

CrossRef Full Text | Google Scholar

Bair, J., Asle Zaeem, M., and Tonks, M. (2015). A Review on Hydride Precipitation in Zirconium Alloys. J. Nucl. Mater. 466, 12–20. doi:10.1016/j.jnucmat.2015.07.014

CrossRef Full Text | Google Scholar

Cheong, Y. M., Gong, U. S., and Chooea, K. N. (2001). Formation and Growth of Hydride Blisters in Zr-2.5Nb Pressure Tubes. Nucl. Eng. Technol. 33 (2), 192–200.

Google Scholar

Domizzi, G., Enrique, R. A., Ovejero-García, J., and Buscaglia, G. C. (1996). Blister Growth in Zirconium Alloys: Experimentation and Modeling. J. Nucl. Mater. 229, 36–47. doi:10.1016/0022-3115(95)00204-9

CrossRef Full Text | Google Scholar

Gao, R., Shi, X., Wang, J., and Huang, H. (2022). Understanding Electrocaloric Cooling of Ferroelectrics Guided by Phase‐field Modeling. Journal of the American Ceramic Society. 105 (6). 3689–3714. doi:10.1111/jace.18370

CrossRef Full Text | Google Scholar

Gao, R., Shi, X., Wang, J., Zhang, G., and Huang, H. (2021). Designed Giant Room‐Temperature Electrocaloric Effects in Metal‐Free Organic Perovskite [MDABCO](NH₄)I₃ by Phase-Field Simulations. Adv. Funct. Mat. 31 (38), 2104393. doi:10.1002/adfm.202104393

CrossRef Full Text | Google Scholar

Ghosal, S. K., Palit, G. C., and De, P. K. (2002). Corrosion of Zirconium Alloys in Nuclear Applications – A Review. Mineral Process. Extr. Metallurgy Rev. 22 (4-6), 519–546. doi:10.1080/08827500208547428

CrossRef Full Text | Google Scholar

Han, G. M., Zhao, Y. F., Zhou, C. B., Lin, D.-Y., Zhu, X. Y., Zhang, J., et al. (2019). Phase-Field Modeling of Stacking Structure Formation and Transition of δ-hydride Precipitates in Zirconium. Acta Mater. 165, 528–546. doi:10.1016/j.actamat.2018.12.009

CrossRef Full Text | Google Scholar

Hellouin de Menibus, A., Auzoux, Q., Dieye, O., Berger, P., Bosonnet, S., Foy, E., et al. (2014). Formation and Characterization of Hydride Blisters in Zircaloy-4 Cladding Tubes. J. Nucl. Mater. 449 (1-3), 132–147. doi:10.1016/j.jnucmat.2014.03.006

CrossRef Full Text | Google Scholar

Hong, S. I., Lee, K. W., and Kim, K. T. (2002). Effect of the Circumferential Hydrides on the Deformation and Fracture of Zircaloy Cladding Tubes. J. Nucl. Mater. 303 (2-3), 169–176. doi:10.1016/s0022-3115(02)00814-0

CrossRef Full Text | Google Scholar

Kammenzind, B. F., Berquist, B. M., Bajaj, R., Kreyns, P. H., and Franklin, D. G., R. B. (2000). The Long-Range Migration of Hydrogen through Zircaloy in Response to Tensile and Compressive Stress Gradients. Zirconium Nucl. Industry Twelfth Int. Symposium, 196–233. doi:10.1520/stp14301s

CrossRef Full Text | Google Scholar

Kim, H.-J., Han, J.-H., Kaiser, R., Oh, K. H., and Vlassak, J. J. (2008). High-throughput Analysis of Thin-Film Stresses Using Arrays of Micromachined Cantilever Beams. Rev. Sci. Instrum. 79 (4), 045112. doi:10.1063/1.2912826

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, S.-D., Kim, J.-S., and Yoon, J. (2018). Phase Analysis of Hydride Blister in Zirconium Alloy. J. Alloys Compd. 735, 2007–2011. doi:10.1016/j.jallcom.2017.11.359

CrossRef Full Text | Google Scholar

Kim, S. G., Kim, W. T., and Suzuki, T. (1999). Phase-field Model for Binary Alloys. Phys. Rev. E 60 (6), 7186–7197. doi:10.1103/physreve.60.7186

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, Y. S., Ahn, S. B., and Cheong, Y. M. (2007). Precipitation of Crack Tip Hydrides in Zirconium Alloys. J. Alloys Compd. 429 (1-2), 221–226. doi:10.1016/j.jallcom.2006.09.034

CrossRef Full Text | Google Scholar

Krishnan, R., and Asundi, M. K. (1981). Zirconium Alloys in Nuclear Technology. Proc. Indian Acad. Sci. Engg. Sci.) 4 (1), 41–56. doi:10.1007/bf02843474

CrossRef Full Text | Google Scholar

Lin, J.-l., Han, X., Heuser, B. J., and Almer, J. D. (2016). Study of the Mechanical Behavior of the Hydride Blister/rim Structure in Zircaloy-4 Using In-Situ Synchrotron X-Ray Diffraction. J. Nucl. Mater. 471, 299–307. doi:10.1016/j.jnucmat.2015.12.048

CrossRef Full Text | Google Scholar

Liu, D., Wang, J., Jafri, H. M., Wang, X., Shi, X., Liang, D., et al. (2022). Phase-field Simulations of Vortex Chirality Manipulation in Ferroelectric Thin Films. npj Quantum Materials 7, 1–8. doi:10.1038/s41535-022-00444-8

CrossRef Full Text | Google Scholar

Long, F., Kerr, D., Domizzi, G., Wang, Q., and Daymond, M. R. (2017). Microstructure Characterization of a Hydride Blister in Zircaloy-4 by EBSD and TEM. Acta Mater. 129, 450–461. doi:10.1016/j.actamat.2017.03.016

CrossRef Full Text | Google Scholar

McRae, G. A., Coleman, C. E., and Leitch, B. W. (2010). The First Step for Delayed Hydride Cracking in Zirconium Alloys. J. Nucl. Mater. 396 (1), 130–143. doi:10.1016/j.jnucmat.2009.08.019

CrossRef Full Text | Google Scholar

Motta, A. T., Capolungo, L., Chen, L.-Q., Cinbiz, M. N., Daymond, M. R., Koss, D. A., et al. (2019). Hydrogen in Zirconium Alloys: A Review. J. Nucl. Mater. 518, 440–460. doi:10.1016/j.jnucmat.2019.02.042

CrossRef Full Text | Google Scholar

Motta, A. T., and Chen, L.-Q. (2012). Hydride Formation in Zirconium Alloys. Jom 64 (12), 1403–1408. doi:10.1007/s11837-012-0479-x

CrossRef Full Text | Google Scholar

Motta, A. T., Couet, A., and Comstock, R. J. (2015). Corrosion of Zirconium Alloys Used for Nuclear Fuel Cladding. Annu. Rev. Mat. Res. 45 (1), 311–343. doi:10.1146/annurev-matsci-070214-020951

CrossRef Full Text | Google Scholar

Puls, M. P. (2009). Review of the Thermodynamic Basis for Models of Delayed Hydride Cracking Rate in Zirconium Alloys. J. Nucl. Mater. 393 (2), 350–367. doi:10.1016/j.jnucmat.2009.06.022

CrossRef Full Text | Google Scholar

Reheman, W., and Ståhle, P. (2018). Self-inflicted Fracture of Expanding Surface Precipitates. Fatigue Fract. Eng. Mater Struct. 41 (12), 2614–2628. doi:10.1111/ffe.12872

CrossRef Full Text | Google Scholar

Rozhnov, A. B., Belov, V. A., Nikulin, S. A., and Khanzhin, V. G. (2011). Stress Corrosion Cracking of Zirconium Cladding Tubes: I. Proximate Local SCC Testing Method. Russ. Metall. 2010 (10), 979–983. doi:10.1134/s0036029510100186

CrossRef Full Text | Google Scholar

Santisteban, J. R. S. A., and Domizzi, G. (2009). Determination of Residual Stresses Around Blisters in Zr-2.5% Nb Pressure Tubes. Powder Diffr. 24 (S1), S72–S76. doi:10.1154/1.3139056

CrossRef Full Text | Google Scholar

Sheng, J., Wang, Y. C., Liu, Y., Wu, S., Xu, K., Chen, Z. H., et al. (2022). A Phase-Field Model for Simulating Hydrogen-Induced Pitting Corrosion with Solid-Solid Phase Transformation in the Metal. doi:10.48550/arXiv.2203.07289

CrossRef Full Text | Google Scholar

Shi, X., Wang, J., Cheng, X., and Huang, H. (2021). Ultrafast Ferroelectric Domain Switching Induced by Nano‐Second Strain‐Pulse, Advanced Theory and Simulations. 5 (3). 2100345. doi:10.1002/adts.202100345

CrossRef Full Text | Google Scholar

Shi, X., Wang, J., Cheng, X., and Huang, H. (2022). Phase Diagram of Sub‐GHz Electric‐Field‐Induced Polarization Oscillation. Phys. Rapid Res. Ltrs 16, 2100416. doi:10.1002/pssr.202100416

CrossRef Full Text | Google Scholar

Simon, P. C. A., Aagesen, L. K., Jokisaari, A. M., Chen, L-Q., Daymond, M. R., Motta, A. T., et al. (2021). Investigation of δ Zirconium Hydride Morphology in a Single Crystal Using Quantitative Phase Field Simulations Supported by Experiments. J. Nucl. Mater. 557. 153303. doi:10.1016/j.jnucmat.2021.153303

CrossRef Full Text | Google Scholar

Singh, R. N., Kishore, R., and Sinha, T. K.. ea. (2002). Hydride Blister Formation in Zr–2.5wt%Nb Pressure Tube Alloy. J. Nucl. Mater. 301 (2-3), 153–164. doi:10.1016/s0022-3115(02)00708-0

CrossRef Full Text | Google Scholar

Singh, R. N., Ståhle, P., Massih, A. R., and Shmakov, A. A. (2007). Temperature Dependence of Misfit Strains of δ-hydrides of Zirconium. J. Alloys Compd. 436 (1-2), 150–154. doi:10.1016/j.jallcom.2006.07.049

CrossRef Full Text | Google Scholar

Toghraee, A., Bair, J., and Asle Zaeem, M. (2021). Effects of Applied Load on Formation and Reorientation of Zirconium Hydrides: A Multiphase Field Modeling Study. Comput. Mater. Sci. 192. 110367. doi:10.1016/j.commatsci.2021.110367

CrossRef Full Text | Google Scholar

Tummala, H., Capolungo, L., and Tomé, C. N. (2018). Quantifying the Stress State in the Vicinity of a δ-hydride in α-zirconium. J. Nucl. Mater. 511, 406–416. doi:10.1016/j.jnucmat.2018.08.050

CrossRef Full Text | Google Scholar

Une, K., Ishimoto, S., Etoh, Y., Ito, K., Ogata, K., Baba, T., et al. (2009). The Terminal Solid Solubility of Hydrogen in Irradiated Zircaloy-2 and Microscopic Modeling of Hydride Behavior. J. Nucl. Mater. 389 (1), 127–136. doi:10.1016/j.jnucmat.2009.01.017

CrossRef Full Text | Google Scholar

Vicente Alvarez, M. A., Santisteban, J. R., Domizzi, G., and Almer, J. (2011). Phase and Texture Analysis of a Hydride Blister in a Zr–2.5%Nb Tube by Synchrotron X-Ray Diffraction. Acta Mater. 59 (5), 2210–2220. doi:10.1016/j.actamat.2010.12.024

CrossRef Full Text | Google Scholar

Xu, S., Shi, X., Pan, H., Gao, R., Wang, J., Lin, Y., et al. (2022). Strain Engineering of Energy Storage Performance in Relaxor Ferroelectric Thin Film Capacitors. Advcd Theory Sims, 2100324. doi:10.1002/adts.202100324

CrossRef Full Text | Google Scholar

Yang, C., Wang, X., Wang, J., and Huang, H. (2020). Multiphase-field Approach with Parabolic Approximation Scheme. Comput. Mater. Sci. 172. 109322, doi:10.1016/j.commatsci.2019.109322

CrossRef Full Text | Google Scholar

Yang, C., Liu, Y., Huang, H., Wu, S., Sheng, J., Shi, X., et al. (2021). Hydride Corrosion Kinetics on Metallic Surface: a Multiphase-Field Modeling. Mat. Res. Express 8 (10), 106518. doi:10.1088/2053-1591/ac1c32

CrossRef Full Text | Google Scholar

Yu, Q., Reyes, M., Shah, N., and Marian, J. (2020). Kinetic Model of Incipient Hydride Formation in Zr Clad under Dynamic Oxide Growth Conditions. Mater. (Basel) 13 (5). 1088. doi:10.3390/ma13051088

PubMed Abstract | CrossRef Full Text | Google Scholar