# Numerical Investigation of Fluid Forces Acting on a Confined Cylinder with Obstacle Subjected to Axial Flow

Authors: Guillaume Ricciardi
Journal: Science and Technology of Nuclear Installations
 (2020)

---

## Abstract

This paper focuses on fluid forces acting on a confined cylinder subjected to axial flow in application to fuel assembly dynamic behavior. From the literature, it is difficult to estimate the damping induced by the flow. Therefore, it is proposed to study numerically the damping fluid forces on a cylinder for various parameters. It has been observed that it increases with the smaller confinement and with the presence of an obstacle and decreases when the Reynolds number increases. Larger values correspond to a greater contribution of pressure forces. Dynamic simulations are compared to the steady ones and give different values, but the order of magnitude and general trend remain the same. Therefore, steady simulations are suitable to have a rough estimation of drag coefficients in dynamics.


---

## Body

## 1. Introduction

Fluid-structure interaction issues can be encountered in many nuclear power plant components [1]. Fluid-elastic instabilities can occur in the steam generator [2] in case of high flow rate. Turbulence can also induce vibrations of the rods in fuel assemblies [3] generating fretting wear. Axial flow strongly modifies the dynamic behavior of the whole fuel assembly [4, 5] by added stiffness, damping, and mass. Therefore, it is important to account for the fluid-structure phenomenon to model the dynamical behavior of a core in case of an external excitation. The model developed by Ricciardi [6] simulating the behavior of fuel assemblies under an external excitation accounts for fluid forces. These fluid forces are mostly responsible for the damping of fuel assembly. Therefore, knowing the effect of fuel assembly geometry on these forces is important. The study of fluid forces acting on a cylinder under inclined flow started with Relf and Powell [7]. They proposed experimental results of galvanized steel wires in a wind tunnel with an inclination of i, measuring lift and drag forces for different diameters and angles, and showed that the lift force was proportional to sin2i. Based on this observation, Taylor [8] developed an empirical model for fluid forces on swimming animals depending on the roughness. Lighthill [9] proposed a different expression for fluid forces based on inviscid theory introducing the virtual mass concept. From the work of Taylor [8] and Lighthill [9]; Païdoussis [10] decomposed the fluid forces into viscous and pressure forces to describe the dynamic of a cylinder and estimate instabilities. The theory is then expanded to a tube bundle [11]. The viscous force, responsible for the damping of the structure, has a component on the axis of the cylinder FL=1/2ρDV2Ct and on the normal axis FN=1/2ρDV2Cnsini. According to theory, drag coefficients Cn and Ct should be equal as proposed in Païdoussis [12], Païdoussis [11], and Triantafyllou and Chryssostomidis [13]. On the other hand, many studies find different values for these coefficients. Ortloff and Ives [14] give boundaries for the ratio 0.5<Cn/Ct<2, and the lower value is for rough cylinders and the upper bound for smooth ones. Chen and Wambsganss [15] and Lee and Kennedy [16] exceed this range with a ratio larger than 10, whereas Dowling [17] and Jamal et al. [18] proposed a ratio smaller than 0.25. The range for Cn in the literature is quite large 0.003<Cn<0.200 according to Païdoussis [19]; Moussou et al. [20] even show values up to 0.55. Divaret et al. [21] performed steady and dynamic tests and showed that the normal force responsible for the damping was due to not only viscosity but also a pressure distribution. Ersdal and Faltinsen [22] showed that the lift coefficient depends closely on the Reynolds number. De Ridder et al. [23] reached the same conclusions as Ersdal and Faltinsen [22] and Divaret et al. [21] numerically. From the literature, it is difficult to estimate the damping of a structure under axial flow as the range of value for Cn and the dependency on various parameters seems wide. Therefore, it is proposed in what follows to study numerically the damping fluid forces on a cylinder for various parameters. First, steady simulations are made with inclined walls to observe the influence of the Reynolds number, inclination angle, and confinement. Then, simulations are made with a small obstacle to account for the spacer grid effect. And finally, dynamic simulations are performed with the motion imposed of the cylinder.


## 2. Numerical Model

The domain studied accounts for one cylinder of diameterD in a rectangular confinement. The distance between the cylinder and the wall is noted Jx in the ex direction and Jy in the ey direction at the midplane where the cylinder is centered relatively to the walls (Figure 1). At the inlet, a homogeneous flat profile of the velocity Vbulk is imposed. After an entry length Lin, the confinement walls are inclined of the angle b over the length Lincl. Finally, a straight length Lout is imposed before the outlet where a homogeneous pressure is imposed (Figure 2). The entry length Lin is sufficiently long to have a developed velocity profile at the beginning of the inclined region, and the length Lout is chosen sufficiently long so that the outlet boundary condition does not affect the flow in the inclined region. All dimensions are given in Table 1.Figure 1
Confinement geometry at midplaneLincl/2.Figure 2
Geometry of the inclined walls.Table 1
Geometrical dimension and fluid property used in the simulations.
L
incl (cm)
L
in (cm)
L
out (cm)
D (mm)
ρ (kg/s)
ν (m2/s)
35
20
5
9
1.103
1.10−6The purpose of these simulations is to obtain the fluid forces acting on a cylinder subjected to a slightly inclined axial flow. In order to have reliable data, the length of the inclined region has to be long enough to have fluid forces that do not depend on the axial coordinate. Since this paper considers narrow confinement, only small anglesb<0.5° are accounted for. Let us define the following dimensionless parameters:(1)gx=JxD,gy=JyD,ReD=VbulkDν,ReDh=VbulkDhν,Dh=4JxJy−πD22Jx+2Jy+πD,v=VVbulk,r=x−D/2Jx,where gx and gy are the dimensionless confinements, ReD and ReDh are, respectively, the Reynolds numbers based on the diameter of the cylinder and on the hydraulic diameter Dh.Actual geometries of fuel assemblies’ spacer grids are complex. To avoid meshing difficulties, in this study, the grid is modeled by an obstacle with a simple shape. This allows to keep a structured mesh and it gives qualitative results. Thus, the obstacle is an increase, notedE, of the cylinder radius (Figure 3) on a small length Lob centered in the inclined region (Figure 4). Let us define the following dimensionless parameters:(2)e=ED,l=LobD.Figure 3
Confinement geometry at midplaneLincl/2 with an obstacle.Figure 4
Geometry of the inclined walls with an obstacle.Simulations are made for several values of confinementgx and gy, angle b, and Reynolds number ReD, with and without obstacle and for different sizes of obstacle e and l, and the range of each parameter can be found in Table 2. Over 500 simulations are performed with different combinations of the parameters.Table 2
Parameter range used for the simulations.
g
x
g
y
b (°)
Re
D
Re
D
h
e
l
0.22⟷ 3.33
0.22⟷ 3.33
0.05⟷ 0.5
9.103 ⟷ 9.104
5.103 ⟷ 5.105
0.033⟷ 0.11
1.11⟷ 4.44Steady turbulent simulations are made with the open-source CFD softwareCode_Saturne using the k−ωSST turbulence model with a logarithmic wall law. A structured mesh is used with a dimensionless wall distance y+∈ [15 : 40] accounting for half of the geometry with symmetry conditions. A simulation accounting for the whole geometry was performed on the smallest mesh to insure that the flow was symmetric. Depending on the parameters, meshes needed between 2 and 3.5 million elements (Figure 5). The half diameter was divided into 120 cells, the axial length into 770 cells, and the gap between 27 and 40 cells depending on the geometry. A convergence study has been made on the two extreme cases (narrowest and largest confinement) to ensure that fluid forces extracted from the numerical results were meaningful.Figure 5
Example of meshes used for the simulations.Fluid forces per unit of lengthFz and Fx, respectively, in the ez and ex directions are extracted from the simulation by integrating the local fluid forces over a circle located at z and projected on the corresponding axis. These fluid forces can be expressed for small angles:(3)Fz=12ρDCtVbulk2,Fx=12ρDCnbVbulk2,where Ct is the drag force coefficient and Cn is the normal coefficient:(4)Ct=2FzρDVbulk2,Cn=2FxρDbVbulk2.The transition between straight and inclined domains shows a singularity that induces significant spatial fluctuations of the fluid force. To account only on the effect of inclination,Cn and Ct are calculated based on the mean values of Fx and Fz on a length where their variations are negligible, between five diameters downstream of the transition inlet/inclination and five diameters upstream of the transition inclination/outlet.


## 3. Inclination and Reynolds Number Effect

This part focuses on the effect of inclination and Reynolds number on a small confinement. Figure6 shows the evolution of Cn as a function of the angle b and the Reynolds number ReD for a small confinement. One can see that Cn reaches an asymptotic value as the angle goes close to zero. It can also be seen that it increases with the angle b as a power law, and the same trend has been observed in Divaret et al. [21] and De Ridder et al. [23]. Having an asymptotic value for small angles can be of great interest for instabilities or response to external excitation as it gives a lower bound for the damping value and thus a conservative estimation of the response.Figure 6
Evolution ofCn as a function of the angle b (a) and as a function of the Reynolds number ReD (b), for a small confinement.
(a)
(b)On the range of angle studied the evolution versusReD is much important with a decrease of Cn as ReD increases which is also observed in Divaret et al. [21] and De Ridder et al. [23]. The evolution of Ct as a function of ReD shows similar trend than Cn and does not depend on the angle, and thus, the ratio Cn/Ct can be considered as constant over the studied range. If Cn can be difficult to measure experimentally, Ct is much more affordable, so it would be very convenient to have a known relationship between these two coefficients. As a consequence, Cn could be estimated from a Ct measurement or a CFD simulation. In this narrow confinement case, a power law in Reynolds number can be found by least square method: Ct=0.105/ReD0.23, and the ratio approximately gives Cn/Ct=15.


## 4. Confinement Effect without Obstacle

In this section, the effect of the confinement is studied. Figure7 shows that Cn decreases as the size of the confinement gets larger, and one can observe a large range of values. One can also observe that the slope drastically increases for the lower values of gx and gy. The ratio Cn/Ct shows the same trend, and one can observe that for large values of confinement gx and gy, the ratio seems to converge toward 1 which conforms to the original theory from Taylor [8]. Since there is no variation of Cn for higher values than gy>2, the cylinder could be considered as in an infinite environment. Nevertheless, one has to keep in mind that these values could change for a lower Reynolds number especially in laminar flow.Figure 7
Evolution ofCn (a) and Cn/Ct (b) as a function of gy.
(a)
(b)Increasing the size of the confinement changes the geometry. Thus, the Reynolds numberReD based on the cylinder diameter may not be relevant since flow profile and turbulence can be significantly different from one geometry to another. Nevertheless, even accounting for the Reynolds number ReDh based on the hydraulic diameter Dh, Ct can take different values depending on the geometry and it is the same for the ratio Cn/Ct that can go from 4 to 13 for ReDh = 100 000 (Figure 8). No particular trend rises from this figure except that as the Reynolds number increases, both coefficients decrease.Figure 8
Evolution ofCt (a) and Cn/Ct (b) as a function ReDh.
(a)
(b)Figure9 shows the ratio of pressure and viscous contribution to the lift force as a function of Cn/Ct. One can observe consistent results with a very limited spread and two asymptotic branches. On one side when Cn/Ct tends toward 1, the portion of the pressure force goes to zero which means that only the viscous drag force projected is responsible for the lift force, this statement is in total agreement with the Taylor [8] theory. On the other side, the participation of pressure on the lift force seems to converge toward 90% as the ratio Cn/Ct increases. It confirms the observation made by Divaret et al. [21] that lift force is due to a pressure distribution. The present study confirms that both theories are correct and depend on the confinement.Figure 9
Ratio of pressure on the lift force as a function ofCn/Ct.Figure10 shows different velocity profiles from each side of the cylinder for different confinements. For the small confinement, the inclination of the walls modifies the velocity profile, and this phenomenon disappears as the confinement increases and the profile becomes flat. The inclination generates an unbalanced flow that induces an inhomogeneous pressure distribution around the cylinder and thus an extra lift force. Based on these observations, one could state that the lift force due to pressure is observable when the spatial variation of the velocity profile is of the same order of magnitude than the size of the confinement. This observation has to be related to the angle, since here only small angles are considered, one could expect to have a significant modification of the flow profile for large angles and thus an increase of the lift coefficient as observed in Divaret et al. [21].Figure 10
Velocity profile from each side of the cylinder forgx = gy = 0.22 (a), for gx = gy = 0.33 (b), and for gx = gy = 0.77 (c).
(a)
(b)
(c)Since the ratioCn/Ct shows regular evolution as a function of gx and gy, an empirical formula is proposed here to have a rough estimation of the ratio:(5)CnCtapprox=1+cxgxnx+cygyny,with cx=2.17, cy=1.22, nx=0.301, and ny=1.22.The approximation gives reasonable estimation (Figure11). Based on this formula, it is possible to have an estimation of Cn from Ct which can be easily estimated from simulation or experiment.Figure 11
Approximation ofCn/Ct.


## 5. Obstacle Effect

In this section, the effect of an obstacle is assessed by slightly increasing the radius of the cylinder over a small length. The purpose of this obstacle is to generate a singularity as encountered in the grid design of a fuel assembly. It induces spatial fluctuations of the fluid force and thus significantly changes the force induced by regular inclination. Simulations show that the presence of a small obstacle (l = 1.1 and e = 0.033) on a small confinement induces a large increase (more than 50%) of the coefficient Cn (Figure 12). Cn follows the same trend as b or ReD increases with or without obstacle.Figure 12
Evolution ofCn as a function of the angle b for a small confinement, comparison with and without small obstacle (l = 1.1 and e = 0.033).As discussed in the previous section, the effect of the obstacle is due to the modification of the flow profile, and thus, for a large confinement, the obstacle has no effect onCn and the ratio Cn/Ct goes back to 1 like in the situation without obstacle (Figure 13).Figure 13
Evolution ofCn (a) and Cn/Ct (b) as a function of a square confinement gx=gy, comparison with and without small obstacle (l = 1.1 and e = 0.033).
(a)
(b)For small values ofgx and gy, the obstacle locally reduces the cross section, and thus, the lift coefficient Cn increases as the length l and the thickness of the obstacle increase (Figure 14) to reach the value of 0.35 which is in the range [0.3–0.55] encountered in Moussou et al. [20]. One can conclude that a structure will be more damped with the combination of a small confinement and obstacle and that an increase of the pressure drop, which is related to the drag force, will induce an increase of the damping. This also explains why a fuel assembly [20], which presents many obstacles and a confined flow, shows high values for Cn.Figure 14
Evolution ofCn with an obstacle for a small confinement gx = gy = 0.22 as a function of the length of the obstacle l.


## 6. Dynamic Simulations

Dynamic simulations are performed; while the walls of the confinement remain still and are not inclined, the cylinder oscillates without inclination. All the points of the cylinder have the same displacementUx on the ex direction:(6)Uxt=aDcos2πfVbulkDt,f=FDVbulk,where a is the dimensionless amplitude, F is the imposed frequency, f is the dimensionless frequency, t represents the time, and Vbulk is the mean value of the velocity at the cylinder location.A uniform flow is imposed at the inlet (Figure15), the cylinder has the same diameter as in the steady simulations, and confinement values are taken among those used in static simulations. The same mesh parameters are used for the simulations taking an ALE approach to account for the displacement of the cylinder, and the time step used is 0.0011/F. Dimension of the domain can be found in Table 3. Around 100 simulations are run for two confinement values, various amplitudes, frequencies, and Reynolds number (Table 4).Figure 15
Test configuration for dynamic simulations.Table 3
Confinement values used for the simulations.
L
cylin (cm)
L
cyl (cm)
L
cylout (cm)
8
40
8Table 4
Confinement values used for the dynamics simulations.
a
f
Re
D
g
x
g
y
0.01⟷ 0.12
0.004⟷ 0.028
14.103 ⟷ 86.103
0.22, 1.11
0.22, 1.11Simulations are run for ten oscillations to obtain the harmonic response. The force is extracted over seven oscillations discarding the transient state. The lift forceFx is extracted from the simulation on a length sufficiently far from the extremities of the cylinder to avoid the effect of the singularities and to have a force that does not depend on the axial position. Assuming a decomposition of the force in terms of added damping, stiffness Kf and mass Mf and replacing the angle b=U˙x/Vbulk gives the following equation:(7)Fxt=12ρDCnVbulkU˙xt+MfU¨xt+KfUxt.The integral over a period of the work rate of the added stiffness and mass is equal to zero:(8)∫TMfU˙xtUxtdt=0,∫TKfUxtUxtdt=0,which finally leads to(9)Cn=2ρDVbulk∫TFxtUxtdt∫TU˙xtUxtdt.The angleb is defined as the relative inclination between the flow and the cylinder, since it remains straight, and for small angle, it can be simplified as(10)b=2πaf,and thus, one angle can be reached for several combination of dimensionless amplitude a and frequency f.Figure16 compares the ratio Cn/Ct given by steady and dynamic simulations, for two values of confinement. One can observe that, for narrow confinements, the steady simulations overestimate significantly the ratio, whereas for a larger confinement, it is the opposite and the difference is less important. If the variation with the amplitude a in dynamic is consistent with the steady results, the ratio Cn/Ct shows a much more important variation with the frequency, and this trend is even more observable for large confinement. This increase with frequency has been observed by Ricciardi and Boccaccio [24]. This could be due to a modification of the flow profile induced by the dynamic of the cylinder. Steady and dynamic simulations show a similar trend by decreasing when the Reynolds number increases (Figure 17).Figure 16
Evolution ofCn/Ct as a function of the angle b, comparison between steady and dynamic simulations.Figure 17
Evolution ofCn as a function of ReD, comparison between steady and dynamic simulations.Dynamic simulations with a small obstacle (l = 1.1 and e = 0.033) and a narrow confinement were performed. One can observe that the steady simulations overestimate the lift coefficient Cn. On the other hand, the comparison without obstacle confirms what has been observed in steady simulations, that the obstacle induces a significant increase of Cn.One dimensionless frequencyf can be reached for several combinations of frequency F and velocity Vbulk, and simulations show superimposed results for the ratio Cn/Ct (Figure 18) with and without obstacle. One could conclude that, for small amplitude a and frequency f ration, Cn/Ct only depends on the frequency.Figure 18
Evolution ofCn/Ct as a function of the dimensionless frequency f.


## 7. Conclusion

A parametric study has been conducted numerically to obtain liftCn and drag Ct coefficients varying the confinement, flow velocity, and inclination for steady simulations and frequency and amplitude for dynamic simulations. Moreover, the effect of an obstacle has been assessed. A wide range has been obtained for both coefficient 0.005<Ct<0.013 and 0.005<Cn<0.35 almost reaching the highest values found in the literature for Cn. It has been observed that Cn increases with the narrowest confinements and with the presence of an obstacle and decreases when the Reynolds number increases.The larger values of the ratioCn/Ct correspond to large values of Cn and it shows a range between 1<Cn/Ct<21. The lower bound corresponds to the case where only the viscous term participates to the lift force, and as the ratio increases, the share of pressure contribution in the lift force increases. The higher values of the ratio are obtained for the smaller confinements and show the largest modification of the velocity profile. Therefore, large values of Cn are obtained when the inclination of the cylinder locally modifies the flow either by a small confinement as shown in the present study or by large inclination angle as in Divaret et al. [21] and De Ridder et al. [23]. From this observation, 1 seems to be a theoretical lower bound for the ratio Cn/Ct, and thus, the values proposed by Ortloff and Ives [14], Dowling [17], and Jamal et al. [18] cannot be explained by the present study.Dynamic simulations give different values than the steady ones, but the order of magnitude and general trend remain the same. Therefore, steady simulations are suitable to have a rough estimation of drag coefficients in dynamics. The effect of amplitude oscillation is less important than expected as observed in Ersdal and Faltinsen [22]. On the other hand, the dimensionless frequency seems to have a much higher effect on the fluid forces.In a further study, simulations will be made to account for more representative fuel assembly geometry.



---

*Source: 5809175-2020-01-09.xml*