© G. Ricciardi et al., Published by EDP Sciences, 2024

1. Introduction

The design of advanced experimental nuclear reactors consists of integrating safety and operational requirements as well as reaching targets in terms of thermal power and neutron spectrum. In order to meet these constraints, slender structures with little support and crossing the entire reactor vessel are implemented in the reactor and are subjected to an axial flow that generates flow-induced vibration (FIV). This vibration mechanism results from the combination of many complex physical phenomena, such as the turbulence in the surrounding flow which leads to unsteady hydraulic loadings of the structure, the dynamic response of the structure which is generally not linear for industrial applications, and additional fluid forces resulting from the flow response due to the structure displacement. From the industrial point of view, the mastering of the occurrence of FIV and its associated wear in case of contacts is requested. For this purpose, one relies, on one hand, on experimental tests in a single channel with laser diagnostics for featuring structure dynamics. On the other hand, a numerical solution with computational cost compatible with the pace of industrial projects is also requested as it allows us to identify potential design options and rank them. The stationary fluid forces applied on slender tubular structures in response to its motion are of primary interest since they can significantly affect the vibration amplitudes and even the stability of the system, especially in case of confined axial flow with high velocities (typically higher than 10 m s⁻¹).

Fig. 1. RJH core and fuel element.

The study of fluid forces acting on a cylinder under inclined flow started with [1]. They proposed experimental results of galvanized steel wires in a wind tunnel with an inclination, measuring lift and drag forces for different diameters and angles. Based on these observations, [2] developed an empirical model for fluid forces on swimming animals depending on the roughness. Lighthill [3] proposed a different expression for fluid forces based on inviscid theory introducing the virtual mass concept. From the work of [2–4] proposed a decomposition of the fluid forces into viscous and pressure forces to describe the dynamic of a cylinder evolving in axial flow and to predict the potential fluid-elastic instabilities.

Fig. 2. Experimental device.

Narrow annular flows show important effect of the fluid forces with large added mass effect. The first to propose an extensive model of the annular flow was Hobson [5–7]. After some simplifications, accounting for a non-deformable solid in motion, he showed that boundary conditions have important effects on the stability of the system. Similar effects have been observed on axisymmetric bodies [8, 9]. Based on the Hobson’s work, one-dimensional models have been developed [10–13]. Bélanger et al. [14] proposed to determine the fluid forces using numerical integration of the Navier-Stokes equations of an annular flow with a simplified motion of the structure. Mateescu and Païdoussis [15] proposed to estimate the viscous effect using Hobson’s approach but based on the flow given by the potential flow theory. The model was extended to the case of deformable structure and showed similar results when compared to the original potential flow theory [16]. More recently, [17] proposed to apply Hobson’s work on a deformable structure and compare the proposed model to potential flow theory and experimental results.

This study, proposed to compare two approaches to estimate the vibration induced by turbulent excitation of an industrial device encountered in a research nuclear reactor. The control-rod guide-tube mock-up of the Jules Horowitz Reactor, previously tested in a hydraulic channel at the Technical Center from Le Creusot in France, is retained for this benchmark.

In the first approach, the model developed by [17], as an extension of Hobson’s work, is applied. In the second one, a simplified fluid-structure model using added mass and drag coefficients previously calibrated by means of steady-state Computational Fluid Dynamics (CFD) simulations is proposed. It must be noted that the first model needs to solve a complex coupled system accounting for velocity and pressure of the fluid whereas the second one is simplified since all the fluid effects are projected on the structural degrees of freedom. In both approaches, the turbulent excitation forces are estimated by refined unsteady fluid dynamics analyses assuming that this fluctuating part of the fluid forces is uncoupled with the tube motion. Displacements obtained numerically are compared to each other and confronted with the experimental results.

2. Experimental apparatus

The Jules Horowitz Reactor (RJH) is a material test reactor under construction at Cadarache in France [18]. The reactor core is made of cylindrical fuel elements (Fig. 1). Each fuel element has a free space to place a control rod or an experimental device embedded in a guide tube. Fuel elements are approximately 0.7 m long for an outer diameter of 100 mm and the guide tubes are 2.1 m long for an outer diameter of 40 mm. The guide tube is made of aluminum and has an inner diameter of 36 mm. Fuel elements are subjected to an axial water flow of 15 m s⁻¹, therefore the guide tube is also subjected to an external axial flow of 15 m s⁻¹. Guide tubes can be divided into three regions, for each part the guide tube is immersed in water with cylindrical confinement. The lower part with a large confinement corresponds to the lower plenum of the reactor, in that region, the axial velocity is low because the cross section is much more important than in the second region. Then, the middle part is where the fuel elements are located with a high axial velocity of 15 m s⁻¹ and a narrow annular confinement with tree stiffeners, the distance between the guide tube and the fuel element is 2 mm. The stiffeners are part of the fuel element, their purpose is to maintain the fuel’s curved plates. They are therefore located in the middle region of the guide tube and are present all over this region. The stiffeners create locally narrower paths for the fluid. Finally, the upper part of the tube, located at the upper plenum of the reactor, is inserted inside a protective head cap with the same narrow confinement as the middle part (2 mm), but without significant axial flow, the main flow being deviated. This paper only focuses on the behavior of the guide tube based on a preliminary design. Although the actual design is slightly different, the experimental results still remain of interest with the aim to understand the dynamic of the guide tube. The guide tube has also an inner flow of about 2 m s⁻¹ and ends with a restriction of the cross section through a diaphragm (Fig. 2). Experiments reproducing the guide tube and the fuel element have been performed at Le Creusot (France) by Framatome with appropriate confinement and axial flow. The displacement of the guide tube was measured at three locations. Different boundary conditions of the structure were tested (Fig. 3): cantilevered; one support at the diaphragm location; and two supports (one at the diaphragm and the second one close to the outlet of the fuel element). Tests are performed at ambient temperature (17°C) which gives a Reynolds number of about 55 000 (the expected Reynolds number in real conditions is 100 000), therefore the flow is turbulent in that region.

Fig. 3. Support conditions of the experimental device.

3. Fluid excitation modeling

The fluid excitation is assumed to be induced by turbulence which is produced at the vicinity of flow singularities with stiff variations of the free flow section. The main singularities are located at the inlet and outlet of the fuel element where the discharge velocity varies from 7 m s⁻¹ up to 15 m s⁻¹ as mentioned in Figure 2. Due to the large complexity of the actual geometry, no satisfactory pressure fluctuation spectrum model can be found in the literature, and one better relies on Large Eddy Simulation (LES) employed in the Computational Fluid Dynamics (CFD) code STAR-CCM+. Indeed, the approach consisting in filtering the smallest turbulent scales while solving the largest is well suited for modeling the time-fluctuations of the hydraulic loadings induced by turbulence.

LES are performed within an industrial frame. This implies that the flow near walls is not explicitly solved but a standard wall-of-the-law is used as wall treatment. Consequently, the subgrid-scale closure is related to the WALE model from [19]. Solved elements are trimmed cells that correspond to hexahedrons cut by the boundaries of the computational domain, and the spatial resolution is sufficient for reaching a Pope criterion, i.e. relative amount of explicitly-solved turbulent energy, larger than 80%.

From a methodological point of view, the test vein is not treated in a single computational domain of the LES, but one takes advantage of the large pressure drop induced by the fuel element which decouples the flow conditions on both sides of the fuel element. Therefore, two LES computational domains are considered. The first domain related to the bottom of the test vein includes the lower area and the first 150 mm of the middle area identified in Figure 2. The second one ranges from the last 150 mm of the middle area up to the test vein outlet. The boundary conditions for the LES, especially the flow rate distribution along the water channels of the fuel elements, are provided by a steady RANS simulation of the whole test vein which is performed prior to the LES.

Fig. 4. Section colored by the velocity magnitude at several instants with a focus on the outlet of the fuel element

Fig. 5. Section colored by the velocity magnitude at several instants with a focus on the inlet of the fuel element

Figure 4 gives an illustration of the LES results computed at the fuel element outlet. The cross sections are colored by the velocity magnitude and are separated of 1 ms. Regarding the flow phenomenology, vortex shedding develops at the outlet of the fuel plates and downstream in the top nozzle where the free flow section is contracted. The inlet of the fuel element is investigated in Figure 5, turbulent structures are produced near the bottom nozzle acting as an obstacle to the flow.

Associated Power Spectral densities of transverse forces are displayed in Figure 6. Except for the PSD of the transverse force applied to the diaphragm, the shape is typical of a wide-brand turbulent excitation with a cut-off frequency of about 800 Hz. Some peaks can also be pointed, they are related to the vortex shedding downstream various flow obstacles inside the fuel nozzles. The PSD for the diaphragm is more peculiar with a clear peak at about 11 Hz. Further analysis reveals that this peak is related to a large recirculation zone developing inside the cavity downstream.

Fig. 6. PSD of turbulent transverse hydraulic forces: inlet fuel assembly (up left), outlet fuel assembly (up right), upper plenum (bottom left), and at the upper end of the tube near the diaphragm (bottom right).

The time evolutions of the transverse hydraulic loadings applied on the outer surface of the guide tube are recorded in the simulation. The fluid forces are then applied as an external nodal transverse force on the guide tube, this force is noted F_turb in the following.

After examining the spatial correlation of the flow turbulence, the intensity of the fluctuating fluid forces, and their distribution along the tube, the hydraulic loading can reasonably be reduced to four uncoupled nodal forces located at both fuel element nozzles (inlet and outlet), at the cross-section restriction inside the upper plenum and at the upper end of the tube near the diaphragm. The spatial correlation is neglected between these flow regions since their distances from each other are clearly larger than the correlation length.

4. Annular flow modeling

The model presented here was developed by [17] based on the original work of Hobson [5–7] and the amalgam proposed by [20]. The assumptions and equations are briefly recalled below.

Fig. 7. System notation for the annular flow.

Let us assume that the fluid-structure interaction effects are preponderant in the middle part inside the fuel element, and to a lesser extent, in the upper plenum. It will be then assumed that the effect of the flow in the lower region is negligible. Let us note h the distance between the tube and the confinement, u the mean axial fluid velocity and v the azimuthal component of the velocity (Fig. 7). Assumptions are that h is small compared to the radius of the tube R and that fluctuations of u and v are small compared to the bulk velocity U. To account for the stiffeners (Fig. 8) it will be assumed that they induce an increase of the friction by accelerating the azimuthal velocity associated with local pressure drop. This friction is assumed to be distributed over the circumference of the tube. It will be also assumed that the stiffeners have no effect on the axial flow. Assuming a planar displacement at θ = 0, one can decompose variables into a mean part and a fluctuating one noted with the prime ′:

(1)

(2)

(3)

(4)

where H and h′ are respectively the means and the fluctuating part of the distance between the guide tube and the confinement h, U and u′ are respectively the mean and the fluctuating part of the axial velocity u, v′ is the fluctuating part of the azimuthal velocity and P and p′ are respectively the mean and the fluctuating part of the fluid pressure p. The fluctuating parts described here refer to the fluctuations induced by the motion of the tube. They are not related to fluctuations induced by turbulence which have different time and space scales.

Fig. 8. Stiffeners.

Removing the steady state equation and linearizing gives the continuity and momentum equations of the fluid:

(5)

(6)

(7)

with the boundary conditions:

(8)

(9)

(10)

where C_f is the friction coefficient, R_m  = R + H/2 is the mean radius, ρ is the fluid density, H_s is the distance between the guide tube and the stiffeners at rest, θ_s is the angle made by a stiffener (Fig. 8) and η is the pressure drop induced by the abrupt restriction and given by [21]:

(11)

The dynamic equation of the structure is obtained considering an Euler-Bernoulli beam with tension induced by the axial friction fluid force and its projection on the displacement direction associated with pressure and azimuthal friction fluid forces integrated over the tube surface:

(12)

where m is the mass per unit of length of the tube, E is Young’s modulus, I is the quadratic moment, μ is the structural damping, g is the gravity and T₀ is the tension at the extremity of the tube. In the following, the model constituted of (5)–(12) will be referred to as the Leakage Flow model (LF).

The system of equations is solved with a finite element method for the spatial discretization coded in Python. Unknowns of the fluid (velocity and pressure) and structure (displacement) have the same spatial discretization. An explicit time scheme is used and boundary conditions are handled with Lagrange multipliers. All the fluid and structure unknown are solved at the same time with a monolithic approach.

5. Potential flow modeling

In this section, a simplified model that does not take into account the fluid variable is proposed. The model proposed by [22] is based on the decomposition of fluid forces as inviscid and viscous effects. Inviscid terms are given by potential flow and viscous effects are accounted for by the drag force of which the projection on displacement direction gives rise to damping and stiffness terms. The axial effect of the drag force and the gravity effect are neglected here knowing that previous sensitivity analyses had shown these terms are not predominant and can be omitted for JHR guide-tube application.

By separating the effect of viscous and inviscid effect, this model does not take into account the effect of fluid friction forces on the fluid dynamic which is the main difference with the Leakage Flow model. This model will be referred to as Potential Flow model (PF):

(13)

where C_D1 and C_D2 are drag coefficients and C_H1, C_H2, and C_H3 are added mass coefficients. In the original theory no distinction is made between the drag coefficients C_D1 and C_D2 and between the added mass coefficient C_H1, C_H2, and C_H3. To compensate for the simplicity of this model we choose to identify each coefficient separately based on CFD simulations.

Fig. 9. RANS computations with prescribed motions.

Fig. 10. U-RANS computation with pure translational motion of the guide tube.

Three kinds of U-RANS computations (Fig. 9) are needed to successively get the five required coefficient values.

As a first step (Fig. 10), a pure translational motion is imposed on the tube. A harmonic lateral displacement is enforced, and the CFD analysis performed with STAR-CCM+ gives in return a phase-shifted fluid force on the tube. This temporal evolution of the hydraulic force is used to calibrate the two first coefficients C_H1 and C_D1 of the FSI model. An overview of the CFD results is available in Figure 10, one can point out that the transverse flow establishes in the opposite direction of the guide tube displacement. The friction of the transverse flow on the walls allows for dissipating the kinetic energy of the guide tube.

As a second step, the same methodology is applied by considering a static deformed shape with a given beam curvature is considered in order to get the second set of parameters C_H2 and C_D2.

The last step is dedicated to a harmonic displacement on mode 1 in order to calibrate the last missing coefficient C_H3 related to the Coriolis force.

The coefficient values obtained in this way are provided in Table 1 for both confined flow regions, the first one being inside the fuel element and the second one inside the protective head cap. Outside these two narrow annular spaces, precisely in the lower plenum of the reactor, the flow is unconfined, therefore the added mass coefficients C_H are taken equal to 1 and the drag force effects are neglected (C_D = 0). At this stage, the following comments can be made. Inside the protective head cap, at the upper part of the guide tube, the three added mass coefficients C_H1, C_H2, and C_H3 have almost the same value which is very close to the theoretical estimation using [23] formulae applicable to the case of two coaxial cylinders separated by an inviscid annular fluid at rest:

(14)

Inside the fuel element, higher added mass coefficients are found even though the internal and external diameters of the fluid space are the same as the upper section ones. Actually, the presence of the fuel element stiffeners creates a locally narrow annular flow passage with a small gap (around 0.5 mm) producing an increase of the fluid confinement effect. Another important finding is that the drag coefficients C_D strongly exceeds the usual values found in the literature. According to [20], the range for C_D is quite large but rarely exceeds 0.2 except for very confined flow. Basically, the high values obtained here are the direct consequence of the very narrow annular geometry such as R ≫ H. This result is consistent with the leakage flow model developed by [17] and described in Section 4 which introduces drag forces depending on (h′) and ∂h′/∂x, proportional to R ²/H ². The larger C_D values are obtained inside the fuel element where the stiffeners provide an additional annular pressure drop by dividing the flow outside the guide tube into three sectors.

Contrary to the theory from [20], one can see that the drag force coefficients C_D can be slightly different depending on whether the tube moves laterally or is inclined. The same observation is made for the different added mass coefficients, related to the purely inertial effect C_H1, the centrifugal (C_H2), or the gyroscopic (C_H3) fluid forces, especially inside the fuel element.

Sensitivity analyses have shown that the set of coefficients given in Table 1 depends neither on the fluid velocity, nor the tube motion frequency provided that the latter is greater than 2 or 3 Hz (below this value, the fluid viscosity tends to increase the added mass effects). This results show that the simplified model described in this section, based on the [22] work (1981), can be reasonably used, at least for the main parts of the flow outside the geometrical singularities.

6. Comparison with experimental results

The guide-tube response to the fluid excitation is calculated by means of non-linear transient analyses on an Euler-Bernoulli beam Finite Element model. Only the guide tube is meshed. The fuel element and the protective head cap in the upper plenum are assumed motionless, they are not introduced in the Finite Element-model. Both fluid-structure interaction models, the Leakage Flow (LF) and the Potential Flow (PF) models described respectively in Sections 4 and 5, are used separately and coupled with the guide-tube beam model. The initial damping when the fluid is at rest (i.e. without flow) is set to 5% for the first dynamic modes of interest, which matches with the available measurements coming from modal analyses performed on the tube inside the hydraulic single-channel mock-up at Le Creusot Technical Center.

Different boundary conditions of the structure were investigated, cantilevered beam; one support at the diaphragm location; and two supports, one at the diaphragm and the second one close to the core outlet. In all cases, the lower end of the tube is clamped. The same turbulent excitation loads, derived from the unsteady CFD analysis discussed in Section 3, are applied for all the support conditions since the mechanical devices used for ensuring the tube lateral support function do not modify significantly the flow.

The frictional contact against the fuel element stiffeners and the head cap is taken into account at the smallest annular clearance locations. The common penalty method is used for the contact algorithm. During the impact durations, the Coulomb model is used to describe the adhesion phases by friction and the sliding phases. The contact normal stiffness K_c , is taken equal to twice the tube ovalization stiffness given by the empirical formula , where D and e are respectively the external diameters and the thickness of the section, and E is Young’s modulus. The friction coefficient used in the Coulomb model is set to 0.3 when the contact is bonded and 0.2 when sliding occurs.

Fig. 11. Comparison of experimental and numerical displacements with no support (left), 1 support (middle) and 2 supports (right) at the H (top), M (middle), and B (bottom) positions.

Figure 11 compares the spectral displacement of experiments and both models for all the support conditions at the three levels of measurements, and Table 2 gives the RMS values.

For the cantilevered case without support at the top of the guide tube, both models give a reasonable estimation of vibration obtained by the experiment. However, both models overestimate the displacement at low frequencies, especially the PF model and at the top position H, whereas the opposite trend is found at higher frequencies. An important result must be highlighted here, the dynamic response of the tube at 13 Hz observed in the experimental spectra is not reflected by the PF model when the LF model succeeded in reproducing it (in air the first two natural frequencies are 2.9 Hz and 21.3 Hz). A possible explanation is that the PF model used is valid only for the main part of the flow, outside the geometrical singularities like the sudden expansion of the cross-section at the fuel element exit. The fluid-elastic forces computed with the simplified PF model are locally erroneous near these particular zones where the flow conditions change drastically. Thanks to the convection effect accounted for with the axial fluctuation u′ terms in (5)–(7), and the fluid boundary conditions given in (8)–(10), the LF model is most likely able to provide a better prediction of the fluid-elastic forces at the ends of the confined axial flow.

The same observations can be made for the second test when a lateral support is added at the upper end of the tube. But in this case, the comparison with the experiment is better since the main frequencies at 20 Hz are well reproduced by both simulations, especially at the lower sensor in position B.

For the last test, two additional supports are implemented, one at the upper end and the second one near the fuel element outlet. The oscillations become small, and the number of impacts on the fuel element stiffeners is limited now. So, the tube response is quite linear in this case and that is most likely why both models’ results are really satisfactory here and fit well with the experimental data, except at the high frequencies.

The second support added in the 3rd test, at the middle part, is achieved by means of 3 screws. Two of these screws were equipped with force sensors and it is possible to compare the recorded data with the numerical predictions.

Fig. 12. Reaction forces at the middle support.

The spectral distribution of the force fluctuations is provided in Figure 12. Both models are able to satisfactorily reproduce the experimental data. Since the force measurement is very close to the turbulent fluid force application point on the beam model, this comparison supports qualitatively the robustness of the LES method used for the CFD simulation.

It should be mentioned that the PF model seems to give better results here. Actually, the middle support condition was modeled with a rigid tangential stiffness for the LF computation whereas a more realistic tangential stiffness taking into account the bending of the lateral screws was implemented when the PF model was used. This distinction in the middle support modeling is probably the reason why the PF model provides more accurate results for the reaction forces. Nevertheless, this will not significantly affect the vibration results in both situations the contact remains stiff compared to the flexibility of the guide tube.

7. Conclusion

A numerical method using CFD and mechanical models was developed in support of the hydraulic qualification of the JHR internals. Two models were proposed, one based on the derivation of leakage flow theory and the other one based on potential flow theory with adjusted coefficients given by CFD simulations. The hydraulic tests performed on a single-channel mockup equipped with a control rod guide tube at Le Creusot Technical Center were numerically reproduced for validation. The flow induced vibration amplitude is consistent with the experimental data. Also, the calculation and experiment provide similar trends when the boundary conditions are changed. However, the spectral distribution of the vibration at lower frequencies is often overestimated, especially with the PF model. Moreover, some frequencies observed in tests and preponderant for the tube displacement are sometimes not reflected by the PF model whereas the LF model succeeds in restoring them. For these reasons, the LF model seems to be slightly more efficient to predict FIV of JHR internals but its use is more difficult in an industrial environment since it needs to solve a complex coupled system accounting for the fluid velocity and pressure. Nonetheless, the use PF model is also complex as it requires to make series of URANS simulations to identify the coefficients.

In spite of the flaws of the models, the modeling strategy could be extended to the other internals to have a reasonable estimation of vibration amplitudes. Parametric studies could be done using the proposed models in this paper to transpose the experimental data which validates the design to the actual reactor conditions taking into account the complex 3D global flow, the interaction between neighboring internals, and the potential component geometrical imperfections due to the manufacturing, the assembly, and the operating conditions.

Funding

This research did not receive any specific funding.

Conicts of interest

The authors declare that they have no competing interests to report.

Data availability statement

Data associated with this article cannot be disclosed due to legal reason.

Author contribution statement

Guillaume Ricciardi: writing – review and editing, formal analysis. Christophe Collignon: writing – review and editing, formal analysis. Jérémy Galpin: writing – review and editing, formal analysis. Mamadou Cissé: writing – review and editing, investigation. Philippe Piteau: writing – review and editing, formal analysis. Domnin Erard: writing – review and editing, supervision. Patrick Blanc: writing – review and editing, supervision.

References

All Tables

All Figures

	Fig. 1. RJH core and fuel element.
In the text

	Fig. 2. Experimental device.
In the text

	Fig. 3. Support conditions of the experimental device.
In the text

	Fig. 4. Section colored by the velocity magnitude at several instants with a focus on the outlet of the fuel element
In the text

	Fig. 5. Section colored by the velocity magnitude at several instants with a focus on the inlet of the fuel element
In the text

	Fig. 6. PSD of turbulent transverse hydraulic forces: inlet fuel assembly (up left), outlet fuel assembly (up right), upper plenum (bottom left), and at the upper end of the tube near the diaphragm (bottom right).
In the text

	Fig. 7. System notation for the annular flow.
In the text

	Fig. 8. Stiffeners.
In the text

	Fig. 9. RANS computations with prescribed motions.
In the text

	Fig. 10. U-RANS computation with pure translational motion of the guide tube.
In the text

	Fig. 11. Comparison of experimental and numerical displacements with no support (left), 1 support (middle) and 2 supports (right) at the H (top), M (middle), and B (bottom) positions.
In the text

	Fig. 12. Reaction forces at the middle support.
In the text