## 1. Introduction

In Korea, since the first operation of a commercial nuclear power plant in 1978, pressurized water reactors (PWRs) have produced more than 10,000 tons of spent nuclear fuel (SNF) and spent fuel pools for storage of the SNF are expected to be saturated soon [1]. In the present circumstances, dry storage of SNF is the most promising candidate to overcome the situation. Integrity evaluation of SNF is a legal requirement [2-4]. The integrity of SNF must be maintained intact under normal conditions and implemented under hypothetical accident conditions to evaluate the safety of the SNF like a criticality. The US Department of Energy (DOE) carried out a damage evaluation of SNF under various conditions for low burn-up SNF [5]. They conducted computer simulations and specified three possible failure modes for SNF cladding. Modes I and II are transverse tearing caused by a combination of axial tension and bending, and Mode III is longitudinal tearing due to the pinch force. Following the DOE study, the Electrical Power Research Institute (EPRI) performed an integrity evaluation for high burnup SNF [6]. EPRI’s method statistically assessed the damage probability using the strain energy density (SE) of the cladding considering the hydride structure for only Mode III, not for Modes I and II. A more advanced method was developed to define a conceptual theory by analogizing the hydride structure of the cladding with a composite to evaluate the normal conditions and accident conditions separately, and to assess probabilistically by introducing a reliability method. That improved approach was the first evaluation for Mode III [7].

In this study, we conducted an integrity evaluation related to Modes I and II for 14 × 14 SNF using an advanced method. We used the critical strain energy density (CSE) as the failure criteria, which is affected by material degradation, for example, the circumferential hydride concentration. SE is used as an assessment measurement determined by the stress and strain. SE is influenced by loading and geometric conditions. The purpose of this study was to generate the probability distribution of the CSE and the SE then calculate the damage probability of 14 × 14 SNF cladding under drop events. A Monte Carlo method was used to produce the probability distributions, and the probability of transverse tearing was estimated by the reliability evaluation method.

## 2. Analysis Regime

### 2.1 Mixed Composite Model Setup With Cladding Tube and Hydride

The mechanical behavior of SNF cladding is affected by hydrides through irradiation and dry storage. When we define characteristics of SNF, we regard the cladding mixed with hydrides as a composite structure [7]. Hydrides subjected to perpendicular stress can be regarded as an isostress model, and the equivalent elastic modulus (*M _{eq}*) of that model is described as follows:

where *M _{m}* and

*M*are the elastic modulus or Poisson’s ratio of a zircaloy matrix and hydride, respectively, and

_{h}*V*and

_{m}*V*are the volume fraction of the matrix and hydride. Generally, the equivalent elastic modulus and Poisson’s ratio of an iso-stress model are lower than those of an isostrain model. In case of transverse tearing, axial and bending loads act along the direction perpendicular to the cladding and hydride, so the mechanical properties deteriorate. The general system of stress-strain relationships combining each phase is represented as follows:

_{h}

where *D ^{EP}* is the elastoplastic constitutive tensor in each phase and

*ijkl*, each index has 1, 2, 3.

Stress calculations for each phase are performed using a radial-return method with the following assumptions.

where *ε _{m}*,

*ε*,

_{h}*ε*are the strains of the matrix, hydride and composite, and

_{c}*σ*,

_{m}*σ*,

_{h}*σ*are the stresses of the matrix, hydride and composite. Equation (3) indicates that the strains are the same in the plane of the composite. Similarly, equation (4) denotes that stresses are the same normal to the plane of the composite [7]. These theoretical backgrounds are employed in specifying the material properties and establishing the damage evaluation procedures as described in the next sections.

_{c}### 2.2 Structural Integrity Evaluation Approach

The procedure for transverse damage evaluation is shown in Fig. 1. The first step is to define the cladding damage mode from the mechanical response analysis under the drop event during the SNF transportation. The three potential cladding damage modes are shown in Fig. 2 [8]. This study focuses on the damage evaluation of Modes I and II. Mode I damage may occur when the strain exceeds the material's ductility limit. Mode II damage is a complete fracture of the entire section. Damage for Modes I and II are the main route of transverse tearing. Unlike Mode III, radial hydride has little effect on Modes I and II because radial hydride platelets are parallel to axial loading, which causes the transverse failure.

The SE is used as an assessment measure. We determined the circumferential hydride concentration from the oxide layer thickness and the hydrogen content database. The circumferential hydride is the parameter for defining the CSE. The SE was obtained by drop analysis using a global model [9] and the detailed fuel rod model, which had a significant effect on transverse damage. Finally, the failure percentage was acquired by the reliability evaluation with probability distributions of the CSE and SE.

### 2.3 Finite Element Modeling and Analysis

The detailed cladding finite element (FE) model was developed with ABAQUS. Design and in-reactor combustion features of 14 × 14 SNF are reflected on the model. The FE model used in the analysis was a 1/2 fuel rod crosssectional symmetry model, as shown in Fig. 3. Table 1 lists the parameters for cladding material properties [7]. The properties of the Zircaloy material and hydrides were defined using the MATPRO [10]. The cladding-to-pellet gap was assumed to be very small to simulate the end-of-life conditions of the SNF. In order to apply the axial force and bending moment, axial displacement and bending rotation were controlled at the reference node in the fuel rod cross section, as shown in Fig. 3.

Fig. 4 shows the fracture behavior of the cross section and the axial stress distribution of cladding as functions of displacement and rotation. As shown in Fig. 4(a), the fracture was initiated at the outer diameter of the cladding, so the axial stress at the outer diameter began to decrease. Fig. 4(b) shows that the fracture progressed from the outer diameter to the inner diameter at the completion of Mode I. It can be seen that fracture advanced from 0° to 90° in the clockwise direction, as shown in Fig. 4(c). Fig. 4(c) shows that full rupture did not occur and the analysis reached equilibrium state.

Fig. 5 shows the SE analysis results for the cladding. Along with the stress as shown in Fig. 4, the SE shows the highest value at the outer diameter of the cladding. It can also be seen that the SE rises continuously until the termination of the analysis.

Fig. 6 shows the SE results of various axial forces and bending moments. The SE is the maximum value of the integration points of the finite element at the outer diameter of the cladding. In the figure, it can be seen that the SE rises as the axial force and bending moment increase.

### 2.4 Regression Analysis

By considering the SE as the dependent variable of the two independent variables axial force and bending moment, a third-degree polynomial regression model with two independent variables can be presented as follows:

where *β* is the regression variable, *x*_{1} and *x*_{2} are the independent variables for the bending moment and axial load, *y* is the dependent variable for SE, and *ϵ* is the error term. Equation (5) is a third-degree polynomial regression model, but can be represented as follows by a multi-linear regression model with 9 independent variables.

By replacing ${x}_{1}^{2}={x}_{3},{x}_{2}^{2}={x}_{4},{x}_{1}{x}_{2}={x}_{5},{x}_{1}^{3}={x}_{6},{x}_{2}^{3}={x}_{7},{x}_{1}^{2}{x}^{2}={x}_{8},{x}_{1}{x}_{2}^{2}={x}_{9}$, the polynomial can be represented as a linear regression function. As shown in Fig. 7, if the axial force and moment is “0”, the SE is “0”, so *β*_{0} = 0. Regression coefficients are determined using EXCEL 2019. The selected values are shown in Table 2.

Fig. 7 shows SE as functions of the axial force and bending moment, respectively using Equation (5) and the regression results of Table 2. As the axial force and bending moment escalate, the SE increases.

In order to evaluate the integrity of cladding, we compared the SE with the CSE. In normal conditions, the elastic region of the CSE was used as the failure criterion and the total CSE was used as the failure criterion under accident conditions [7]. The CSE is determined by a simple FE model with hydride distribution. The hydride concentration increases from the ID to the OD. The amount of hydride at the ID was fixed as 50 ppm, and those at the OD were varied from 0 to 4,000 ppm. Mechanical properties of Zircaloy and hydride in MATPRO were used [10]. As shown in Fig. 8(a), the SE for the elastic region at the OD decreases as the amount of circumferential hydride increases, while the SE for the elastic region at the ID remains the same due to the fixed amount of hydride. In case of CSE, the presence of hydride has a large effect. As shown in Fig. 8(b), the CSE, the resistance to transverse failure, declines as the circumferential hydride rises.

Probabilistic assessment approaches were implemented to consider the variance of the damage related parameters. We first computed the probability density distributions of SE and CSE of 14 × 14 SNF using a Monte Carlo method under the drop analysis to evaluate the fuel rod transverse tearing damage probability under normal and hypothetical accident conditions, respectively. In order to generate the probability distribution of CSE, it was necessary to calculate the average circumferential hydride concentrations. The circumferential hydride concentration is determined by the hydrogen content due to corrosion of the cladding tube. The CSE was calculated using the defined average circumferential hydride concentrations. Fig. 9 shows the probability distributions of the hydride concentrations generated using the Monte Carlo method.

The damage probability was evaluated by applying a reliability assessment method with a stress-strength model, which was defined as follows [11]:

where *P _{f}* denotes the failure probability, and

_{f}and

*F*are the probability density function (PDF) and cumulative density function (CDF) for SE and CSE, respectively.

## 3. Results and Discussion

To calculate the failure reliability, the CDF terms of Equation (7) must be calculated using the PDFs for each SE and CSE. The axial forces and bending moments from the drop analysis are shown in Figs. 10 and 11. Fig. 12 shows the probability distribution of the SE using the regression analysis in Fig. 7. As shown in Fig. 9, the average circumferential hydride in the cladding is about 145 ppm, but we assumed 600 ppm of the average circumferential hydride for SE analysis as shown in Table 1. This assumption produced conservative SE results. In order to achieve more accurate results, we conducted additional SE analysis using various hydride conditions. Fig. 13 shows the probability distribution of the CSE using the hydride information of Fig. 9 and the CSE of Fig. 8.

The CDF for SE and CSE were calculated with Equation (7) under normal (0.3-m drop) and accident conditions (9.0-m drop) as shown in Fig. 14. Under the normal conditions, Fig. 14(a) reveals that the minimum CSE is greater than maximum SE, which means that no damage occurred. On the other hand, Fig. 14(b) shows the behavior of 9-m drop accident conditions. The integral of Equation (7) is approximated with a trapezoidal rule by dividing the SE values as sub-intervals. The finally calculated damage probabilities (*P _{f}*) are null for normal and 1.6 × 10

^{−2}for accident conditions. Table 3 summarizes evaluation results including the previous evaluation results of the DOE and EPRI. EPRI calculated the failure probability under accident conditions using load distribution instead of SE. They defined the failure probability as the number of load cases exceeding a specific axial force value among all cases [12]. All the results are very close to zero-damage probability under normal conditions. On the other hand, in accident conditions, the reason why the analysis result is higher than that of the DOE is mainly due to the cutting-edge material model difference of the SNF and the fuel feature. The DOE also carried out the finite element analysis and an extensive parametric study including the effects of temperature, irradiation, crack length, etc., except for the reflection of hydride effects. In this study, we applied a conservative amount of circumferential hydride to the SE analysis. If various circumferential hydrides and loading conditions are taken into account in the SE calculation, we expect a lower failure probability.

## 4. Conclusions

In this study, we evaluated Modes I and II of three possible failure modes for SNF fuel rods. The failure modes were due to the transverse deformation of SNF cladding caused by a combination of axial force and bending moment. We calculated the amount of damage from the degree of SE and also computed the CSE to build evaluation criteria. Through FE analysis, the SE and CSE were first calculated deterministically. Then, for the probabilistic approach based on the reliability assessment, we first yielded the probability distributions of the CSE and SE under axial and bending loadings for 14 × 14 SNF transportation. We evaluated the probability of damage by introducing the reliability evaluation scheme. The SE results were matched with axial loads and bending moments for regression analysis, and all probability distributions were executed with a Monte Carlo method, namely, the failure probability was calculated using the reliability assessment method of the stress-strength model. The results of this study on 14 × 14 commercial SNF revealed no damage possibility under the normal conditions, and the DOE’s was also very closed to null. Meanwhile, under the hypothetical accident conditions, the damage probability of this study was 1.6 × 10^{−2} and the DOE was 2.0 × 10^{−4}. These consequences were inferred from the difference of the sophisticated modeling technique in which an SNF cladding model composed of a matrix and imbedded hydrides made it possible to reflect more precise behavior of ductility degradation and SNF design features. Therefore, it is clear that the proposed approaches are practical and reliable for predicting behavior and responses of SNF cladding tube integrity, which reflect long term dry-storage degradation effects under axial and bending loads.