Coupled Thermomechanical Responses of Zirconium Alloy System Claddings under Neutron Irradiation

: Zirconium (Zr) alloy is a promising fuel cladding material used widely in nuclear reactors. Usually, it is in service for a long time under the effects of neutron radiation with high temperature and high pressure, which results in thermomechanical coupling behavior during the service process. Focusing on the UO 2 /Zr fuel elements, the macroscopic thermomechanical coupling responses of pure Zr, Zr-Sn, and Zr-Nb binary system alloys, as well as Zr-Sn-Nb ternary system alloy as cladding materials, were studied under neutron irradiation. As a heat source, the thermal conductivity and thermal expansion coefﬁcient models of the UO 2 core were established, and an irradiation growth model of a pure Zr and Zr alloy multisystem was built. Based on the user material subroutine (UMAT) with ABAQUS, the current theoretical model was implemented into the ﬁnite element framework, and the consequent thermomechanical coupling behavior under irradiation was calculated. The distribution of temperature, the stress ﬁeld of the fuel cladding, and their evolution over time were analyzed. It was found that the stress and displacement of the cladding were sensitive to alloying elements due to irradiated growth.


Introduction
With sources of fossil fuels becoming increasingly exhausted and the environmental challenges associated with their use, efficient and clean nuclear energy will become a promising energy source in the future. The most significant issue with nuclear power is its safety and controllability in the long term. Cladding is one of the most important components in a nuclear reactor. It has been used in extreme environments with high radiation levels, high temperatures, and high pressure levels for a long time [1]. Zirconium (Zr) alloys have been applied as a cladding material because of their high temperature mechanical properties, small neutron absorption cross section, strong radiation resistance, and good irradiation stability [2]. The properties of the Zr alloy clad and UO 2 fuel are critical issues to be resolved for nuclear power [3,4]. However, the effects that extreme environments-e.g., with high irradiation levels, high temperatures, and high pressure levels-have on the mechanical and thermal properties of Zr alloy cannot be determined through conventional experimental methods [5,6]. Experimental data on optimizing component structure properties is lacking.
Numerical simulations provide an effective method for performance analyses of nuclear materials in extreme environments. Bison, a three-dimensional finite element analysis program developed by the National Laboratory of Idaho, USA, has been used 2 of 10 to simulate the thermal mechanical behavior of fuel rods in the reactor [7]. Alcyone, a three-dimensional finite element analysis program developed by the CEA of France, has been used to study the high burnup performance of fuel [8]. The user subroutine UMAT of ABAQUS finite element software was employed to study the radiation-thermal-mechanical coupling behavior of fuel rods by Tang et al. [9]. The effects of a nonhomogeneous irradiation environment were studied using the UMAT of ABAQUS. The three-dimensional deformation of a typical 17 × 17 fuel assembly was determined by numerical simulation [10]. The thermomechanical behavior simulation, coupled with the hydrostatic-pressure-dependent, grain-scale fission gas swelling calculation for a monolithic UMo fuel plate, was studied under heterogeneous neutron irradiation [11]. The thermomechanical behavior of a Zr-4 fuel rod was studied by experiments and a simulation. It was found that the creep behavior was not captured by a simple law [12]. The deformation behavior of Zr alloys was also investigated at the microscale level. A physically-based crystal plasticity model was proposed to study the contact force and gap opening of fuel rod assembly by Patra and Tome [13]. The effects of fast neutron irradiation on the tensile deformation mechanism of Zr-2.5 pressure tube samples was studied by in situ neutron diffraction [14]. It was found that the starting stress of base plane <a> slip and cylinder <a> slip increased due to neutron diffraction, while the influence of cone < c + a > slip was small and became the dominant deformation mechanism. After the model was established, the SRIM/TRIM software could be used to conduct an ion energy analysis on the thickness of the model to evaluate its safety and stability, and provide a scientific basis for the subsequent behavior model in the fuel rod reactor [15].
Irradiation-induced hardening was confirmed for the Zr alloy cladding. An increase of Nb was shown to be able to retard the formation of irradiation defects [16]. The thermomechanical behavior of Zr alloys with UO 2 /Zr fuel element is key for nuclear reactors. At the same time, the alloying elements in the Zr alloy for the nuclear reactor can also affect the properties of Zr alloys. For example, tin can improve the strength of Zr, reduce its plasticity, and improve its creep and corrosion. Niobium has a high strengthening effect on Zr; the thermal neutron absorption cross-section of niobium is very small, which can eliminate the detrimental effects of impurities such as C, Ti, and Al on the corrosion resistance of Zr, and reduce the hydrogen absorption capacity of Zr. Copper has also a significant strengthening effect on Zr. A small amount of copper can reduce the impact toughness of Zr, and improve the corrosion resistance of Zr-Nb alloy at high temperature [17]. The refractory metal Zr may therefore be a promising cladding material with excellent high-temperature properties [18,19]. The effect of alloy elements on the irradiation behavior of cladding is important in the development of cladding materials. Therefore, pure Zr, Zr-Sn, and Zr-Nb binary systems, and the Zr-Sn-Nb ternary system were considered as cladding materials in this study of the thermomechanical behavior of the UO 2 /Zr fuel element. The effect of neutron irradiation on the mechanical and thermal properties of Zr and Zr alloy, as well as the temperature field, displacement field, and stress field of Zr alloy cladding under neutron irradiation, were studied using the ABAQUS distribution and evolution mechanism.

Thermomechanical Model of UO 2 Fuel
The thermal properties of UO 2 depend on the temperature, porosity, and burnup. In the present study, UO 2 fuel was the only considered heat source. The thermal conductivity of fuel UO 2 was calculated using the Fink-Lucuta model, which is valid from 298 to 3120 K [20,21]: where k 95 is the thermal conductivity of unirradiated UO 2 at 95% theoretical density in W/(m·K). In the present study, the temperature dependence of annealing on irradiation defects was considered in the thermal conductivity of UO 2 .
where T is temperature in K, and F D and F P are the impact factors of dissolved and precipitated fission products on thermal conductivity, respectively.
where F x is the deviation from stoichiometry (1.0 for urania fuel), F R and F M are the impact factor of irradiated effect and the porosity, respectively, and B is fuel consumption in atom%, B = 50%.
where P is the porosity and s is the shape factor that is spherical in this work, s = 0.5. The coefficient of thermal expansion was calculated using Equation (7): where E D and K T are the thermodynamic parameters and Boltzmann constant, E D = 6.9 × 10 − 20 J and K T = 1.38 × 10 −23 J/K respectively. The values of K i are K 1 = 1.0 × 10 −5 K −1 , K 2 = 3.0 × 10 −3 K −1 and K 3 = 4.0 × 10 −2 K −1 respectively. The density and Poisson's ratio are 10.96 g/cm 3 and 0.316, respectively.

Thermomechanical Models of Zr Cladding
In the current work, pure Zr, Zr-4 alloy, Zr-2.5Nb-0.5Cu alloy, and N36 alloy were studied. The thermal conductivity of Zr alloy from room temperature to melting point was shown to be temperature dependent [22]: Irradiation growth and creep of cladding are important mechanical behaviors under neutron irradiation [23]. The irradiation growth of Zr and its alloys were calculated according to a power-law empirical model [24,25]. The axial strain is described as The strain increment of the irradiation growth is given as where ε z is the axial strain, induced by irradiation growth, Φ is the fast neutron fluence rate, and A and n are the material parameters. The model parameters of irradiation growth for different Zr alloys are shown in Table 1. A = 1.2 × 10 −21 N/cm 2 , n = 0.794. Due to the volume conserving of the irradiation growth of the cladding [26], the growth in the other two directions was assumed to be consistent; the radial strain ε r is given by

Finite Element Modeling
The fuel rod contained 40 pellets with a radius of 4.09 mm. The inner diameter and thickness of the cladding were 4.16 mm and 0.57 mm, respectively. The length of the simulated tube was 200 mm, as shown in Figure 1a. The initial gap distance was set at 0.07 mm, as shown in Figure 1b. The reduced integration FEs (C3D8R) were employed in ABAQUS, and the number of meshes was 10,100. Due to the axial symmetry of cladding cube, axisymmetric boundary conditions were applied, and the node along the mid-circumference was constrained at the axial end of the cladding to prevent rigid body rotation.
where εz is the axial strain, induced by irradiation growth, Φ is the fast neutron fluence rate, and A and n are the material parameters. The model parameters of irradiation growth for different Zr alloys are shown in Table 1. A = 1.2 × 10 −21 N/cm 2 , n = 0.794. Due to the volume conserving of the irradiation growth of the cladding [26], the growth in the other two directions was assumed to be consistent; the radial strain εr is given by

Finite Element Modeling
The fuel rod contained 40 pellets with a radius of 4.09 mm. The inner diameter and thickness of the cladding were 4.16 mm and 0.57 mm, respectively. The length of the simulated tube was 200 mm, as shown in Figure 1a. The initial gap distance was set at 0.07 mm, as shown in Figure 1b. The reduced integration FEs (C3D8R) were employed in ABAQUS, and the number of meshes was 10,100. Due to the axial symmetry of cladding cube, axisymmetric boundary conditions were applied, and the node along the mid-circumference was constrained at the axial end of the cladding to prevent rigid body rotation. The cladding material parameters of fuel rods are shown in Tables 1-4.  The cladding material parameters of fuel rods are shown in Tables 1-4. In the simulation process, the left and right end faces of the cladding were fixed, and the heat source was the core area of the pellet. The pellet core was heated to 1200 • C by radial irradiation. The cooling coefficient of the outer surface of the fuel rod was set at 20,000 W/(m 2 ·K), which made the outer surface of the fuel rod gradually cool down. The coolant pressure was kept at a constant value of 15.5 MPa during the entire simulation process. The initial pressure of helium precharged in the pellet cladding gap was 2 MPa. The simulated power was calculated according to the average power of the nuclear power plant, which was 200 W/cm. A power increase of 1 d resulted in a rise from 0 W/cm to 200 W/cm, and the steady state operation of 200 W/cm took 1200 days. During steady state operation, the fast neutron fluence rate of the fuel rod remained constant, and its value was 9.5 × 10 17 n/(mm 2 ·s). In the simulation, the input parameters were obtained from the practical working conditions of the reactor pressure vessel.

The Thermomechanical Behavior
The temperature field distribution in°C of the fuel elements after irradiation for 1200 days is shown in Figure 2, which are pure Zr, N36, Zr-4, and Zr-2.5Nb-0.5Cu. The temperature distribution in the radial direction at the middle part of the cladding and the axial direction at the ends of the cladding are presented, respectively. It was observed from the simulation results that the temperature of the inner surface was the highest under irradiation, and that this decreased gradually from the inner to outer surface along the radial direction due to the action of the coolant. The distribution of temperature showed a gradient increasing trend at the end of the cladding in the axial direction. In the simulation, the temperature of the four cladding materials was maintained at 300 to 320 • C under the action of the coolant, which was consistent with the value and distribution of the cladding temperature of the LWR fuel rod [27]. The distribution of temperature was the same for pure Zr and Zr alloy. The results show that the composition of alloying elements has little effect on the temperature distribution and magnitude of the cladding tube under neutron irradiation.
In the present study, the thermomechanical behavior of the cladding was studied under normal operating conditions. In the initial stage, the pressure on the cladding was due to the compressive stress from the outer surface to the inner surface, induced by the coolant pressure (15.5 MPa) and the helium pressure (2 MPa) in the gap. Figure 3 shows the Mises stress field distributions in the MPa of cladding for different alloy systems during stable operation for 1200 days under irradiation conditions; these data can be employed to evaluate the cladding safety. It may be observed from the simulation results that the stress field distribution of the cladding was similar to the temperature field distribution; namely, the stress gradually decreased from the inner to the outer surface of the cladding in the radial direction. The distribution of stress along the radial direction was the same for the four alloys, which is consistent with the SiC/SiC composite cladding at the end of its life [28]. The stress was related to the temperature of the cladding produced by fuel burnup. Therefore, the gradient distribution of stress along the radial direction was shown to be independent of the alloys of cladding materials.
However, the axial stress distribution of the cladding was found to be dependent on the alloy element. In the pure Zr and Zr-4 alloy, the stress level was low at the end. The stress gradually increased from the end to center, along the axial direction, as shown in Figure 3a,c. The maximum magnitudes of both pure Zr and Zr-4 alloys were about 30.2 MPa on the 1200th day. In contrast, the stress was high at the end for the N36 alloy. The stress distribution of the N36 alloy showed that the stress decreased from the end to the center along the axial direction, as shown in Figure 3b. Uniform stress distribution was observed along the axial direction for the Zr-2.5Nb-0.5Cu alloy, as shown in Figure 3d. The maximum magnitudes of the N36 and Zr-2.5Nb-0.5Cu alloys were about 78.9 MPa and 44.3 MPa on the 1200th day, respectively. This axial stress distribution led to different axial strains for different cladding materials. Therefore, axial deformation should be considered as an important parameter for selecting a cladding material.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 6 of 11 the radial direction due to the action of the coolant. The distribution of temperature showed a gradient increasing trend at the end of the cladding in the axial direction. In the simulation, the temperature of the four cladding materials was maintained at 300 to 320 °C under the action of the coolant, which was consistent with the value and distribution of the cladding temperature of the LWR fuel rod [27]. The distribution of temperature was the same for pure Zr and Zr alloy. The results show that the composition of alloying elements has little effect on the temperature distribution and magnitude of the cladding tube under neutron irradiation.
(a) pure Zr (c) Zr-4 (d) Zr -2.5Nb-0.5Cu  In the present study, the thermomechanical behavior of the cladding was studied under normal operating conditions. In the initial stage, the pressure on the cladding was due to the compressive stress from the outer surface to the inner surface, induced by the coolant pressure (15.5 MPa) and the helium pressure (2 MPa) in the gap. Figure 3 shows the Mises stress field distributions in the MPa of cladding for different alloy systems during stable operation for 1200 days under irradiation conditions; these data can be employed to evaluate the cladding safety. It may be observed from the simulation results that the stress field distribution of the cladding was similar to the temperature field distribution; namely, the stress gradually decreased from the inner to the outer surface of the cladding in the radial direction. The distribution of stress along the radial direction was the same for the four alloys, which is consistent with the SiC/SiC composite cladding at the end of its life [28]. The stress was related to the temperature of the cladding produced by fuel burnup. Therefore, the gradient distribution of stress along the radial direction was shown to be independent of the alloys of cladding materials.
However, the axial stress distribution of the cladding was found to be dependent on the alloy element. In the pure Zr and Zr-4 alloy, the stress level was low at the end. The stress gradually increased from the end to center, along the axial direction, as shown in Figure 3a,c. The maximum magnitudes of both pure Zr and Zr-4 alloys were about 30.2 MPa on the 1200th day. In contrast, the stress was high at the end for the N36 alloy. The stress distribution of the N36 alloy showed that the stress decreased from the end to the center along the axial direction, as shown in Figure 3b. Uniform stress distribution was observed along the axial direction for the Zr-2.5Nb-0.5Cu alloy, as shown in Figure 3d. The maximum magnitudes of the N36 and Zr-2.5Nb-0.5Cu alloys were about 78.9 MPa and 44.3 MPa on the 1200th day, respectively. This axial stress distribution led to different axial strains for different cladding materials. Therefore, axial deformation should be considered as an important parameter for selecting a cladding material.
(a) pure Zr

The Change of Gap Distance
The displacement fields, in mm distribution of cladding, of different alloy systems irradiated under stable operation for 1200 days, are shown in Figure 4. The results show that the displacement field of the four cladding materials presented a gradient distribution along the radial direction. It was found that the outer surface of the cladding had a large displacement deformation due to the irradiation of pure Zr, N36, and Zr -2.5Nb-0.5Cu alloys. In contrast, a large displacement deformation was observed on the inner surface for the Zr-4 alloy, as shown in Figure 4c. The change of the displacement deformation was affected by the alloy system due to thermal expansion and irradiation-induced growth. Normally, there is no contact between the inner surface of the cladding and the outer surface of the fuel pellet for safety reasons; however, this may occur due to displacement deformation along radial direction.

The Change of Gap Distance
The displacement fields, in mm distribution of cladding, of different alloy systems irradiated under stable operation for 1200 days, are shown in Figure 4. The results show that the displacement field of the four cladding materials presented a gradient distribution along the radial direction. It was found that the outer surface of the cladding had a large displacement deformation due to the irradiation of pure Zr, N36, and Zr -2.5Nb-0.5Cu alloys. In contrast, a large displacement deformation was observed on the inner surface for the Zr-4 alloy, as shown in Figure 4c. The change of the displacement deformation was affected by the alloy system due to thermal expansion and irradiation-induced growth. Normally, there is no contact between the inner surface of the cladding and the outer surface of the fuel pellet for safety reasons; however, this may occur due to displacement deformation along radial direction.
There is typically a gap between pellet and cladding (i.e., the difference between the radial distance on the inner surface of the cladding and the radial distance on the outer surface of the pellet). In order to study the change of gap distance between the pellet and the cladding, the radial distance between the outer surface of the pellet and the inner surface of the cladding is also an important index. The radial distances of cladding and pellet are presented in Figure 5. It was shown that the radial distance increased with running time. At the beginning of the simulation, the radial distance of the fuel increased rapidly, which was ascribed to the heating and expansion of the fuel. The gap distance gradually decreased due to irradiation growth from the simulation results. The initiation of Pellet Cladding Mechanical Interaction (PCMI) can take place at the end due to burnups of pellets, which will induce additional cladding growth from pellet expansion [29]. The change rule of the gap is consistent with the reactor behavior of PWR SiC cladding fuel rods [30]. This is mainly due to the rapid increase of the local thickness of the pellet under irradiation. As the irradiation time increases, the gas distance decreases, which leads to a Appl. Sci. 2021, 11, 1308 8 of 10 short service period of the four cladding materials. The performance of irradiation growth resistance for Zr-based alloy cladding could be improved by using a multi-alloy system. irradiated under stable operation for 1200 days, are shown in Figure 4. The results show that the displacement field of the four cladding materials presented a gradient distribution along the radial direction. It was found that the outer surface of the cladding had a large displacement deformation due to the irradiation of pure Zr, N36, and Zr -2.5Nb-0.5Cu alloys. In contrast, a large displacement deformation was observed on the inner surface for the Zr-4 alloy, as shown in Figure 4c. The change of the displacement deformation was affected by the alloy system due to thermal expansion and irradiation-induced growth. Normally, there is no contact between the inner surface of the cladding and the outer surface of the fuel pellet for safety reasons; however, this may occur due to displacement deformation along radial direction.
(a) pure Zr (b) N36 (c) Zr-4 (d) Zr -2.5Nb-0.5Cu There is typically a gap between pellet and cladding (i.e., the difference between the radial distance on the inner surface of the cladding and the radial distance on the outer surface of the pellet). In order to study the change of gap distance between the pellet and the cladding, the radial distance between the outer surface of the pellet and the inner surface of the cladding is also an important index. The radial distances of cladding and pellet are presented in Figure 5. It was shown that the radial distance increased with running time. At the beginning of the simulation, the radial distance of the fuel increased rapidly, which was ascribed to the heating and expansion of the fuel. The gap distance gradually decreased due to irradiation growth from the simulation results. The initiation of Pellet Cladding Mechanical Interaction (PCMI) can take place at the end due to burnups of pellets, which will induce additional cladding growth from pellet expansion [29]. The change rule of the gap is consistent with the reactor behavior of PWR SiC cladding fuel rods [30]. This is mainly due to the rapid increase of the local thickness of the pellet under irradiation. As the irradiation time increases, the gas distance decreases, which leads to a short service period of the four cladding materials. The performance of irradiation growth resistance for Zr-based alloy cladding could be improved by using a multi-alloy system. (c) Zr-4 (d) Zr-2.5Nb-0.5Cu

Conclusion
In this study, radiation growth models of Zr alloys of different types, i.e., pure Zr, Zr-Sn, a Zr-Nb binary system, and a Zr-Sn-NB ternary system, were established for UO2/Zr fuel elements. The thermal mechanical coupling behavior of different system alloys under irradiation was simulated and calculated. The distribution and evolution laws of temperature field, displacement field, and stress field were analyzed:

Conclusions
In this study, radiation growth models of Zr alloys of different types, i.e., pure Zr, Zr-Sn, a Zr-Nb binary system, and a Zr-Sn-NB ternary system, were established for UO 2 /Zr