#### Abstract

In a fast spectrum reactor, the fuel rod bundle is mainly positioned radially by the wire which can make contact with the adjacent fuel rods, and then it is inevitable that flow-induced vibration (FIV) will cause fretting wear and vibration fatigue of the fuel cladding at the contact position. Therefore, the computational model of fretting wear and fatigue life about the fuel rod bundle caused by FIV will be studied in this paper. Based on the random vibration model of the fuel rod bundle, Hertz contact theory, and Archard wear theory, the fretting wear life computational model and the fatigue life computational model of the wire-to-adjacent fuel rod (WAFR) contact have been established. Finally, taking CEFR design parameters as an example, the fretting wear life and vibration fatigue life of the cladding are calculated, and it is found that fatigue affects the service life of the fuel rod more seriously than fretting wear. The calculation model and method lay a foundation for further accurate prediction and analysis of the fuel rod life.

#### 1. Introduction

In a pressurized water reactor (PWR), the FIV causes grid-to-rod fretting wear, and the fretting wear is responsible for 70% of fuel cladding leakage accidents [1]. Similarly, there are some typical accidents caused by flow-induced vibration in fast spectrum reactors [2, 3]. In general, the fuel rods in fast reactors such as CEFR in China, FFTF in the United States, JOYO in Japan, and China initial Accelerator Driven System (CiADS) use the wire wrapping for radial positioning. With the change of fuel burnup under the reactor operation, pellet swelling and thermal expansion of the cladding will lead to the wire-to-adjacent fuel rod (WAFR) contact. Turbulent flow-induced vibration occurs when the coolant flows through the fuel rod bundle. Therefore, the turbulence-induced vibration is an important factor that causes the fretting wear and vibration fatigue of the cladding [4, 5]. In order to study the wear life and fatigue life of the contact between the wire and adjacent fuel rods, a theoretical computational model must be established. Bakewell, Clinch, and Chen and Wambsganss obtain pressure data and semiempirical pressure expressions [6–8]. According to the empirical formula of the pressure spectrum, the displacement power spectrum of the rod bundle can be obtained and the stress power spectrum can be calculated, which lays a foundation for the establishment of fretting wear life calculation and random vibration fatigue life calculation models. At present, there is little research on the complex life analysis of the contact between the wire and the fuel rod. In order to ensure the safe operation of the reactor, the theoretical computational model of the WAFR contact fretting wear and vibration fatigue under the random excitation of turbulence pulsation in the China Experimental Fast Reactor (CEFR) is studied in this paper.

Firstly, the FIV dynamic model of the fuel rod bundle is established. According to the relationship between the pressure power spectrum and the displacement power spectrum, the displacement power spectrum of the fuel rod is obtained. Based on the Hertz contact theory, the WAFR contact stress model is established. Furthermore, the WAFR contact wear model is established by using Archard wear theory, and the fatigue life model of the WAFR contact is established, which is based on the random vibration fatigue theory. Finally, the theoretical computational model is used to predict the life of the contact danger points of fuel rods in CEFR.

#### 2. Problem Description

In normal fast reactor component design, the initial installation gap is reserved between one fuel rod and another, as shown in Figure 1. During reactor operation, the fuel cladding will undergo radial deformation because of thermal expansion and other factors, such as the error of installation clearance, irradiation swelling, and irradiation creep. About the effect of the error of installation clearance, we will use different equivalent spring stiffness to analyze the effect in Section 3.2. In general, the natural deformation which is caused by other factors for the cladding is mainly considered plastic after some operation time [9], and its contribution to the contact stress is very little to be considered so that the point contact between the wire and adjacent fuel rods occurs. Because of FIV of the turbulent pulsating pressure, the fretting wear and vibration fatigue will appear in the WAFR contact position. For the safe operation of the reactor, the fretting wear and vibration fatigue at the contact position must be calculated and studied. Table 1 shows the main parameters of the cold state design of fuel rods for the China Experimental Fast Reactor (CEFR). The following sections will take these parameters as examples for analysis.

**(a)**

**(b)**

In the fast reactor, 316 stainless steel is regarded as the fuel cladding material (316SS). The thermal expansion coefficient is [12]where *T* is the cladding average temperature (K).

For the sodium-cooled fast reactor CEFR, the cladding average temperature is about 500°C and the installation clearance *δ*_{0} is 0–0.05 mm. Utilizing equation (1) to compute the radial thermal deformation of the cladding, Δ, it is found that Δ = 0.06 mm > *δ*_{0}. As shown in Figure 1, the installation WAFR gap will close and make contact. Then, the stress concentration will occur at the contact point, which is the design danger point. The fretting wear and vibration fatigue of random vibration at dangerous points will be analyzed below.

#### 3. Mathematical Model and Analysis Method

In general, the displacement response of the fuel rod is caused by turbulence excitation at the subcritical flow velocity (defined as the flow velocity before buckling and divergence). Therefore, in order to calculate the fretting model and vibration fatigue of the WAFR contact, the random vibration response model of the rod bundle caused by the turbulence of the fuel rod is firstly determined, and then the stress response spectrum is calculated through the contact stress and displacement response spectrum so as to establish the fretting wear life model to predict the failure life of the fuel rod. Finally, the fatigue life computational model of the WAFR contact is established through the stress response spectrum.

##### 3.1. Vibration Response Model of Bundle

The WAFR contact stress is composed of static stress and dynamic stress, and turbulence pressure pulsation is the main factor causing dynamic stress change. In order to obtain the dynamic stress power spectral density, the dynamic equation governing the vibration motion of the fuel rod bundle can be expressed in the matrix form as [13, 14]where **M**, **C**, and **K** are, respectively, the mass, damping, and stiffness matrices and **F** is a generalized force matrix caused by turbulent pressure fluctuation. For a given system, the elements of the matrices are constant and the dimension is 2*KN*, where *K* is the number of cylinders and *N* is the Galerkin discrete order. Through the Fourier transform and frequency response function, the power spectral density (PSD) of fuel rod bundle displacement in equation (2) is expressed as [14]where *Sη*_{i}*η*_{j}(*ξ*, *ξ′*, *ω*) is defined as the cross-spectral density (CSD), **H**(*ω*) is the frequency response function and equal to [−*ω***M** + *iω***C** + **K**]^{−1}, and *ϕ*(*ξ*) is the beam comparison function; here, the simply supported condition is applied. *N* is the discrete order, *K* is the number of rods, EI is the flexural rigidity, *R* is the rod radius, *L* is the length, *a*(*p*) = {largest integer ≤ (2*K* + *p* *−* 1)/2*K*}, *b*(*p*) = *p* *−* 2*K*[*a*(*p*) *−* 1], and *β*(*n*, *j*) = *2K*(*n* *−* 1) *+* *j*. It is difficult and expensive to measure the pressure field of a wire-wrapped rod bundle. Therefore, numerical simulation can be used instead of the experiment to obtain the turbulent pressure spectrum. Abbasian et al. have done 43 fuel rods of CANDU reactor wall pressure test experiments and CFD simulation based on LES, DES, and RSM turbulence models. The comparative analysis shows that the numerical results are in good agreement with the experimental pressure spectrum, but the numerical simulation is time expensive [15]. And thus in this paper, the more accurate pressure spectrum formula, which is based on Curling and Païdoussis measured data of the bundle pressure field, is applied, and it is written as [16, 17]where Ω is the circular frequency, *U* is the average axial flow velocity, and parameters *K*_{1}, *K*_{2}, and *K*_{3} are 0.25, 2.0 × 10^{−6}, and 7.94 × 10^{−6}. According to the relation *S*_{pp}(*ω*) = [*EI*/(*m* + *ρA*)/*L*^{^4}]^{0.5}Φ_{PP}(Ω), the displacement power spectrum of the arbitrary rod *i* can be calculated by substituting equation 4 into equation 3. In the fast reactor, the rods of the fuel assembly are arranged in a triangular arrangement, and for the convenience of calculation, only a three-cylinder system is selected for analysis, as shown in Figure 2. And the in-plane direction means the selected motion direction is in the same plane. The out-plane direction means the selected motion direction is in the parallel plane and perpendicular to the in-plane direction.

Because of the very small amplitude, it is reasonable to assume that the change of contact depth between adjacent rods is caused by the displacement of one of the rods. The stress power spectrum can be obtained from the displacement power spectrum *S*_{x}(*ω*). As shown in Figure 3, the theoretical PSD is in agreement with the experimental PSD, where the system parameters are the same as those in the experiment [18].

##### 3.2. Contact Stress Model of WAFR Contact

The WAFR contact will inevitably produce contact stress. In order to obtain the stress power spectrum, it is necessary to establish the contact stress model. According to Hertz contact theory [19], for upper and lower surfaces, the contact area is generally an elliptical region, as shown in Figure 4.

The relation between contact pressure and stress along the pressure direction is expressed aswhere *a* is the long half-axis of the ellipse, *b* is the short half-axis of the ellipse, eccentricity of the ellipse *e* = (1 − (*b*/*a*)^{2})^{1/2}, and *T* = ((*b*^{2} + *y*^{2})/(*a*^{2} + *y*^{2}))^{1/2}. *p*_{0} is the maximum contact pressure (MPa). At *y* = 0, the maximum contact pressure is the maximum contact stress. The maximum contact pressure and maximum indentation depth, *δ*, are expressed as follows [19]:from Figure 4.4 in Reference [19], we can get the correction factors *F*_{1}(*e*) and *F*_{2}(*e*). and *R*_{e} are the equivalent Young’s modulus and the equivalent radius, respectively. They satisfy the relationswhere *E*_{1}, *E*_{2} and *ν*_{1}, *ν*_{2} are Young’s modulus and Poisson’s ratio of the wire wrap and cladding. The relationship between positive real constants *A* and *B* iswhere *R*_{1}, and *R*_{2}, are the maximum and minimum curvature radius of the wire wrap and fuel rods, respectively, and *α* is the slope angle between the wire wrap and the fuel rod axis (Figure 3). Curvature radius is expressed aswhere *H* is the helical pitch, *R* is the rod radius, and *r*_{w} is the wire radius. The wire and cladding materials are 316L stainless steel. Table 2 shows the mechanical parameters of 316L. = 109.89 GPa is obtained from equation (8).

The equivalent spring stiffness at the contact point *K*_{stiff} *=* 1.42 N/mm [21], initial displacement *h*∼10 *μ*m, and then the normal initial load is expressed as

Considering turbulent excitation, the normal load is a function of rod displacement; then,where *x*(*t*) is the rod displacement caused by turbulent excitation.

Substituting equation 13 into equation (6), and carrying out the first-order Taylor expansion at *x*(*t*) *=* 0, the formula is written as

Therefore, the stress power spectral density corresponding to equation (14) is

The ratio of the principal relative radii of the curvature and the ratio of the long half-axis and short half-axis are defined as [19]

According to the relationship between load and pressure distribution on the semiellipsoid, the maximum pressure can be expressed as

Combining equations (17) and (18), the long half-axis *a* and the short half-axis *b* can be determined. When the maximum static stress exceeds the yield strength of the material, as shown in Table 2, plastic deformation occurs at the contact position until the contact area can resist external forces.

In order to verify the correctness of the contact stress model in this paper, we use the contact stress model to compute the contact problem as in Reference [22], as shown in Figure 5, and it is found that both results are in good agreement. This example verifies the correctness of the model in this paper.

##### 3.3. Wear Life Computational Model

In the previous section, the calculation of the WAFR contact stress has been analyzed. However, turbulence-induced vibration will cause the WAFR contact wear. In engineering, the maximum wear depth of the fuel rod cladding has to be less than 10% cladding thickness [23]. Therefore, the fretting wear analysis of the cladding is very important for reactor safety. Based on the Archard wear formula, written as equation (19), the fretting wear life calculation model of the cladding will be established in this section:where *A* is the WAFR wear contact area, *A*_{0} is the reference contact area, and *l* is the total sliding distance. According to the experiment, the wear coefficient of 316L stainless steel for a sliding speed of 0.06 m/s at a frequency of 2 Hz is *K* (mm^{3}/N/m), i.e., the slope of the fitting straight line, as shown in Figure 6. *F*_{N} is the normal load on the contact surface, *R*_{d} is the roughness of the abrasive paper, and the diameter of the reference contact surface *D*_{0} = 8 mm in Reference [24]. The WAFR wear contact area *A* is approximately a semiellipsoid surface, and the formula is written aswhere *a*, *b*, and *c* are the half-axis lengths of the contact semiellipsoid, respectively.

For the turbulent random vibration, the full-contact sliding friction distance of the wire on the cladding surface per second, *s*, is defined. According to the displacement power spectral density, the mean square root *dr*_{rms} of the wire displacement amplitude in the compression direction can be obtained.

As shown in Figure 7, from equation (13), the normal load of the WAFR contact is expressed as

The contact trajectory is written as

By deriving formula (22) with respect to time, the relative friction rate between the wire and the cladding can be obtained as follows:

From equation (19), we can obtain the wear volume of *T*_{0}/4 time, and it is written aswhere the average period *T*_{0} = 1/*f*_{ave}. The average frequency *f*_{ave}, that is, the average number of times the positive slope crosses the mean value in unit time, is consistent with the definition in Section 3.4:

After *N*_{c} cycles, as shown in Figure 7, the wear morphology can be determined according to the sliding path of the wire on the cladding surface. The wear volume can be divided into three parts: two quarter ellipsoids and a body whose cross section is a semiellipse rotating at a certain angle 2arcsin (*s*_{rms}/R). And thus, the total wear volume is written as

Then, the maximum wear depth *c* is expressed as

The longest allowable life of cladding wear can be calculated by equation (27), as follows:

Finally, according to equation (28), the wear life of the cladding can be determined in units of hours.

In order to verify the correctness of the wear model in this paper, we use Wang et al.’s example [9] which calculates the contact wear of the fuel rod to the grid. As shown in Figure 8, the selected parameters are in Reference [9]; from the comparison of both results, there is good agreement between the established model in this paper and the Wang computational model.

##### 3.4. Random Vibration Fatigue Life Computational Model

The turbulence excitation will generate random contact stress at the contact danger points, which will lead to the cladding fatigue. The following will establish the vibration fatigue model caused by random contact stress based on the statistical average and probability function method. The random stress PSD *W*(*f*) can be calculated through equation (15). The peak probability density function of the stationary stochastic process can be written as [25]where erf(*x*) is the error function and *σ* is the standard deviation, which means the root mean square value of stress *σ*_{rms}, and the random signal irregular factor *α* is expressed as

For narrowband random signals, *α* = 1, the probability density function of stress amplitude obeys Rayleigh distribution. However, for broadband random signals, the amplitude probability density function is the combination of Rayleigh distribution and Gauss distribution, and especially *α* = 0; it obeys Gauss distribution.

The mean square values of the continuous spectrum and the corresponding discrete spectrum are expressed as follows:

For the discrete spectrum, the zero-crossing expectation of stress *S*_{σ}(*t*) is expressed as equation (33), that is, the average number of times during unit time passing through the mean value at a positive slope:

For the discrete spectrum, the expected peak value of stress *S*_{σ}(*t*), i.e., the average number of peak occurrences per unit time, is written as

Then, according to equation (30), the irregular factor of random signals isand the bandwidth coefficient is defined as *ε* = (1 − *α*^{2})^{0.5}. Two parameters, *α* and *ε*, characterize the irregularity and bandwidth of stochastic processes, respectively. When *α* = 1, i.e., *ε* = 0, it is narrowband. It is generally believed in engineering that when *ε* ≤ 0.35, the random process can be approximated as a narrowband treatment [26].

According to Miner’s linear cumulative damage theory, the cumulative damage of the material can be expressed as

From equation (36), the following structural failure criteria can be established:where the subscript *i* represents the *i*th operating condition, *n*_{i} denotes the number of true cycles when the peak stress value is *s*, and *N*_{i} denotes the number of cycles for material failure when the peak stress is *s*.

The S-N curve describing the material is generally modeled by power functions as follows:where *C* and *m* are parameters related to the material, stress ratio, loading, and so on.

The Basquin model of the stress control conditions can be written as [27]

In order to consider the average stress effect, the SWT parameter method introduces equivalent stress as follows [28]:

Combining equations (39) and (40), the fatigue life under any type *i* working condition is written as

Let *m* = −1/b, *C* = 1/2()^{−1/b}, and *s*_{i} = . That is, equation (41) is written in the form of equation (39), where the strength coefficient of 316 stainless steel is 2100 MPa and the fatigue strength index *b* is −0.1 [29], and then substituting equation (41) into equation (37) leads to

Assuming *N*_{t} is the total number of cycles, the number of cycles in the *s* to *s* + *ds* interval is

By combining equation (42) with equation (43), the total cycles of material failure are written as

For the narrowband process, from equation (19), the peak stress probability density function is

Substituting equation (45) into equation (44), the gross number of stress cycles is obtained aswhere Γ(*x*) is the gamma function.

In practical engineering, most of them are broadband random vibration. For broadband stochastic processes, Wirsching and Light introduce a correction factor *λ* as [30]

The values of *a* and *b* in equation (47) are related to the material parameter *m*, which can be expressed as

Thus, a life estimation formula of the broadband stochastic process is obtained as

According to the value of the bandwidth coefficient *ε*, the random vibration fatigue life of the structure can be determined by using equation (46) or (49).

In order to verify the correctness of the random vibration fatigue life model, we use the model to compare with the experimental data [31], as shown in Figure 9, and it is found that most of the data are in the double error band and very few data are within the quadruple error band. The life prediction errors are acceptable in engineering, and thus, the applicability of the model is verified.

#### 4. Results and Discussion

In this section, the results of the fretting wear model and vibration fatigue model are analyzed, taking CEFR design parameters as an example. The physical parameters of the CEFR are *ρ*UO2 = 10 g/cm^{3}, *ρ*_{Na} = 0.8443 g/cm^{3} at 500°C, the equivalent mass *m* of the fuel rod per unit length of 0.2934 kg, and the average flow velocity of 4.51 m/s, and other parameters are shown in Table 1. When the parameters are dimensionless, the corresponding dimensionless parameters in the matrix of equation (2) are as follows: *β* = 0.0752, *ε* = 225, *γ* = −0.7505, *c*_{f} = 0.016, *u*_{0} = 0.2984, *α* = *c* = *c*_{b} = Π = 0, and *G*_{c} = *G*_{w} = 0.17. If the fuel rod is regarded as the Euler–Bernoulli beam model, the natural frequencies of the vibration of the pinned-pinned beam in vacuum are as follows: *f*_{v} = *π*/2/L^{2}(*EI*/*m*)^{1/2} = 4.37 Hz.

According to equation (3), PSD of rod No. 1 in the in-plane direction and out-plane direction can be obtained, as shown in Figure 10.

Figure 10 shows that the first-order frequency is 3.49 Hz, which is about 20% less than that in vacuum because of the effect of the added mass of the coolant. The corresponding frequency points of peak values in both directions are basically the same. The root mean square value of displacement is about 70 *μ*m. The value is used as amplitude *dr*_{rms}. As shown in Figure 7, *S*_{rms} = (*R*^{2} − (*R* − *dr*_{rms})^{2})^{0.5}. From equations (5)–(28), the calculation parameters of contact wear can be obtained, as shown in Table 3.

The equivalent spring stiffness *K*_{stiff} has a great influence on the wear and also varies with the operation of the reactor, as shown in Figure 11, and the maximum stress varies with the depth under different equivalent spring stiffness. The maximum depth is about 50 *μ*m, and the stress is 0. According to equation (28), the relationship between the wear life and the equivalent spring stiffness is calculated. As shown in Figure 12, it is found that, with the increase of equivalent spring stiffness, wear life decreases, maximum stress increases, and the effect on wear life becomes smaller.

Table 4 shows that *ε* ≤ 0.35, and then the stress spectrum corresponding to the pressure spectrum, i.e., equation (4), is a narrowband stochastic process. By substituting the equivalent stress amplitude of 214 MPa into equation (46), it is calculated that, after 81 h, the fatigue crack initiation of the cladding at the contact position begins. As shown in Figure 13, the S-N curve calculated from the random load model in this paper is compared with the S-N curve calculated from the SWT model. It is found that the random load has a more severe effect on fatigue life. By comparing the wear life with the fatigue life, it can be judged that fatigue failure is the main cause of fuel rod failure. But there are no corresponding theoretical and experimental results about the WAFR contact.

#### 5. Conclusion

Considering the dynamic behavior of rod bundles under turbulent axial flow and the environment in which the fuel rods are mainly subjected to high temperature and strong irradiation effect, it is found that the wire wrap and adjacent fuel rods will make contact by estimation. A theoretical model for life analysis of the WAFR contact is established. Taking the parameters of the CEFR as an example, the conclusion can be summarized as follows.

Based on the study of Païdoussis’ team, the dynamic equation of flow-induced vibration of rod bundles is established and the displacement power spectrum is obtained. The displacement amplitude is about 70 *μ*m; because of the effect of the additional mass of the fluid, the natural frequency in the coolant is about 20% less than that in vacuum.

If the contact equivalent stiffness is 1.42 N/mm, the wear life of the cladding is 1517 days, and it is found that the effect of the contact equivalent spring stiffness on the wear life decreases with the increase of the contact equivalent spring stiffness.

The random load irregularity factor *α* = 0.98 is calculated; that is, the bandwidth coefficient *ε* = 0.0196. Therefore, the random process is treated as a narrowband process, and the random vibration fatigue life of the contact position is calculated to be 81 h according to the random load model.

The results show that fretting wear and vibration fatigue are important problems in WAFR contact life analysis. The work preliminarily established the mechanical model of cladding failure about the WAFR contact for the first time. In the follow-up, the elastic-plastic problems in the WAFR contact will be further analyzed so as to establish a more accurate model to predict the life of fuel rods.

#### Data Availability

Reference for CEFR design is available at http://www.chinafra.com/newsitem/271381546. The data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

This work was supported by the Strategic Priority Research Program of the Chinese Academy of Sciences (Grant no. XDA21010202).