# A Homogenization Method for a FSI Problem: Application to a Tube Bundle Row

Gianluca Artini, Daniel Broc

► To cite this version:

Gianluca Artini, Daniel Broc. A Homogenization Method for a FSI Problem: Application to a Tube Bundle Row. ASME 2017 Pressure Vessels and Piping Conference, Jul 2017, Waikoloa (Hawaii), United States. pp.PVP2017-65737, 10.1115/PVP2017-65737. hal-02418141

HAL Id: hal-02418141

https://hal.science/hal-02418141

Submitted on 13 Feb 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A Homogenization Method for a FSI Problem:
## Application to a Tube Bundle Row

Gianluca Artini, Daniel Broc

► To cite this version:

Gianluca Artini, Daniel Broc. A Homogenization Method for a FSI Problem: Application to a Tube Bundle Row. ASME Pressure Vessels and Piping Division Conference (PVP2017), Jul 2017, Waikoloa, United States. 10.1115/PVP2017-65737. hal-02418141

HAL Id: hal-02418141

https://hal-cea.archives-ouvertes.fr/hal-02418141

Submitted on 23 Jan 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
# AN HOMOGENISATION METHOD FOR A FSI PROBLEM
APPLICATION TO A TUBE BUNDLE ROW

Gianluca Artini*

Den-Service d'études mécaniques et thermiques
(SEMT), CEA,
Université Paris-Saclay,
F-91191, Gif-sur-Yvette, France
Email: gianluca.artini@cea.fr

Daniel Broc

Den-Service d'études mécaniques et thermiques
(SEMT), CEA,
Université Paris-Saclay,
F-91191, Gif-sur-Yvette, France
Email: daniel.broc@cea.fr

## ABSTRACT

A research program is ongoing at CEA in the framework of ASTRID project, for developing and validating tools for the analysis of the mechanical dynamic behaviour of Fast Breeder Reactor (FBR) cores. In this context, Fluid Structure Interaction is a recurring problem, especially in case of tube bundles in nuclear reactor cores or steam generators. The large number of elements makes modelling difficult to perform and any direct numerical simulation at whole core scale unachievable. In order to overcome this limitation, we resort to homogenisation techniques which allow to describe the global dynamic behaviour of the whole fluid-structure system. The purpose of this study is to better understand the effects of the spatial variation of the tube displacement field on the whole bundle dynamics. Thanks to multiple scales, a complete development of the problem is conducted. An analysis of the effects of spatial variations of the long scale displacement field on the tube bundle's dynamics is led. As a first step we consider the homogenisation of Euler linear equations that gives us the possibility to make a thorough assessment of the problem with the extension to follow to Navier-Stokes equations.

## INTRODUCTION

A research program is ongoing at the French Commissariat à l'Énergie Atomique et aux Énergies Alternatives (CEA) for the

development of the new prototypical Sodium Fast Breeder Reactor (Na-FBR) ASTRID.

The ASTRID reactor core presents a large amount of fuel elements (~300) ranged in a very narrow gap configuration (~3 mm) due to neutronic reasons. They are surrounded by about 300 reflectors elements and 500 PNL (Lateral Neutronic Protection) shielding elements.

In FBR configuration, the simulation of whole core results to be necessary for a complete safety assessment of the power plant. For example, phenomena occurring during seismic events or during core flowering involve a strong coupling in the whole core scale. Then, for a complete design and safety analysis understanding the fluid structure interaction phenomena occurring at whole core scale is of crucial relevance.

It is well known that the presence of a fluid strongly changes the dynamic behaviour of the whole reactor core. The fluid presence implicates a strong coupling between reactor core elements in the whole core scale. In addition, the presence of a large number of tubes and the narrow gap configuration make the dynamic analysis very difficult to assess. In this perspective, the homogenisation techniques have a crucial role in the dynamic analysis of a large number degrees of freedom systems, as the FBR cores. This is one of the main reason through we resort to homogenisation techniques [1] or [2].

In this paper, the authors obtain an homogenised equation system at long scale for pressure and displacement field in an enclosed tube bundle. Starting from the Euler linearised equations, thanks to the introduction of two scales development for pressure

<footnote>*Address all correspondence to this author.</footnote>
and displacement, a series of elementary problems describing the system at the short scale level is obtained. They report numerical simulations for these elementary problems and they discuss the main dynamic characteristics of the tube bundle displacement field obtained by the proposed homogenised model.

## HOMOGENISATION METHODS: STATE OF THE ART

The literature on fluid structure interaction for a single tube immersed in an inviscid or viscid fluid is largely widespread, as in, for example, [3] and [4].

It is well known, the presence of fluid involves on oscillating structure two different effects: inertial and dissipative effects.

The fluid flow around an oscillating body may be described by Navier-Stokes equations or under some simplifying hypothesis by Euler or linear Euler equations.

For example, Fritz, in [5] analyses an inviscid fluid initially at rest enclosed in two oscillating cylinders. Thanks to Lagrange formulation, the author obtains an analytical expression of the fluid reaction force on the solid body in the assumption of small displacement amplitude. A proportionality between solid body acceleration and fluid force is shown. The proportional coefficient is the so called *added mass* coefficient.

On the other hand, in [6], Chen et al. analyse the problem of a cylindrical tube immersed in a viscid fluid. In the case of imposed displacement of inner tube's wall, the authors obtain an analytical expression of the fluid force, where two different terms are shown. The first one, representing the inertial effects due to the fluid presence, is in phase with the solid acceleration; the second one, in phase with the solid body velocity, is due to the fluid internal friction and it is representative of dissipative effects.

In a large number of industrial facilities, periodic structures immersed in a fluid are often encountered (nuclear reactors, steam generators, etc.). Especially when the tubes number largely increases, the coupling provoked by the fluid represents the main issue. Moreover, the large amount of degree of freedom makes any direct numerical simulations unachievable. For this reason, homogenisation method are developed.

A variety of methods for homogenisation are found in literature. The main goals of homogenisation is to obtain an equivalent model maintaining the local characteristics of the medium [1] and [2]. Hereafter we recall the main homogenisation methods for FSI problems.

In [7], Planchard proposes a method for homogenising the added mass matrix as the bundle is composed of a large number of tubes, identical and regularly spaced. This technique consists in homogenising the added mass matrix only, and maintaining unchanged the tubes dynamic equation.

Shinohara and Shimogo in [8] and Jacquelin and al. in [9], establish a bi-dimensional analytical homogenised model for square [8,9] or hexagonal [8] section tubes. The model proposed

is limited to elements with plane faces; in fact, it is based on the hypothesis that the translation movement of elements leads to fluid velocity fluctuations in the direction parallel to solid faces. This hypothesis is valuable only in case of high confinement configurations (BWR or FBR cores). Expressions for frequencies of the bundle are obtained as functions of the main physical characteristics and parameters.

Hammami in [10] and Conca et al. after, in [11], utilise the method of multiple scales based on the asymptotic development of the main variables, pressure *p* and structure displacement *s*, in order to obtain an homogenised two phase model. The model proposed by the authors is based on the periodicity of tube bundle geometry. The system obtained is representative of the long scale variation of the main quantities maintaining, at the same time, the effects of short scale phenomena.

In [12] and in [13], the authors utilise a physical approach to the homogenisation. Thanks to the analysis of an elementary cell corresponding to one tube surrounded by fluid, the authors are able to describe the global behaviour in terms of mean fluid structure velocity (the *Darcy velocity*) and homogenised global pressure. In [12], Broc introduces the Rayleigh damping in homogenised equation system, in order to take into account the friction effects within the fluid using an Euler equation based model. Moreover, in [13], Broc and Desbonnets propose an homogenised model based on the Navier-Stokes equations, where the small amplitude assumption is not made. Via volume average of main quantities, equations for long scale variations are obtained in terms of Darcy velocity, long scale homogenised pressure and tube displacement.

This last two approaches, the mathematical one of Conca and Hammami and the physical one of Broc and Desbonnets, are completely equivalent. In fact, it may be demonstrated that the final equation systems obtained with this two techniques are equivalent to each others.

All the previous homogenised models are able to describe the global fuel bundle movement and, to the author knowledge, no homogenised model can describe the effects of the tube differential movements on the dynamic of the bundle itself. In this paper, the authors propose a homogenised model where it is possible to take into account these latter effects.

## PROBLEM POSITION

In this paper, we analyse the dynamic behaviour of a tube bundle row, immersed in an inviscid fluid. The tubes are enclosed in a large vessel with a characteristic length *L* and there is no fluid flow within the bundle itself. The tubes are disposed in a periodic array with a pitch length equal to *ε* (Fig. 1).

From a mechanical point of view, the tubes are identical to each other and they are represented as rigid bodies and their dynamics is described by the Fundamental Principle of the Dynamic, where the forces acting on the tubes are the recalling elas-
<figure><img src="image_3.png" /><figcaption>FIGURE 1. TUBE BUNDLE ROW: GEOMETRICAL CONFIGURATION AND MAIN PARAMETERS</figcaption></figure>

tic force ($k_s s_c$) and the fluid force ($F_{F,i}$). This means that each tube of the bundle is assimilated to an oscillator of one degree of freedom with a mass $m_s$ and a stiffness $k_s$ representative of elastic recalling and constrain forces with no structure damping

$$m_s \ddot{s}_{c,i} + k_s s_{c,i} = F_{F,i} \quad (1)$$

where $s_c$ the tube's axis displacement. In the hypothesis of inviscid fluid, the fluid reaction force is given only in terms of pressure corresponding to the inertial effects due to the fluid presence [5].

Some hypothesis on the geometry and on the mechanical and dynamical behaviour are assumed:

The dimension of the pitch $\epsilon$ is small compared to the main characteristic length of the tube bundle $L$,

$$\epsilon < L \quad (2)$$

The small amplitude hypothesis is assumed for the tube displacement, i.e.:

$$\frac{|s_c|}{\epsilon} \ll 1 \quad (3)$$

Moreover, we assumed the displacement field of the tube bundle only in the horizontal direction $x_1$; thus we can write for each tube:

$$s_{c,i} = s_{c,1} = s_c \quad (4)$$

The fluid, in the inviscid and incompressible hypothesis, is governed by the Euler equations. In our case, the hypothesis of

<figure><img src="image_4.png" /><figcaption>FIGURE 2. QUALITATIVE INTERPRETATION OF LONG AND SHORT SCALE VARIATIONS AND “CHAIN DERIVATION RULE”</figcaption></figure>

small displacement leads to Euler linearised equations

$$\left\{ \begin{array}{l} \displaystyle \frac{\partial w_i}{\partial x_i} = 0 \qquad \text{in } \Omega_F \\ \\ \displaystyle \frac{\partial w_i}{\partial t} = -\frac{1}{\rho} \frac{\partial p}{\partial x_i} \qquad \text{in } \Omega_F \end{array} \right. \quad (5a) \quad (5b)$$

where $w_i$ is the $i$-th scalar component of the fluid vector velocity, $\rho$ is the fluid density and $p$ is the pressure fluid field. Moreover, the small displacement assumption (Eq. 3) allows to declare the slip boundary conditions in the non-displaced geometrical configuration:

$$w_i n_i = \dot{s}_{c,i} n_i \quad \text{at each } \gamma_c \quad (6)$$

where $n_i$ is the scalar component of the normal vector $n$ at each surface $\gamma_c$.

Coupling together mass balance equation (Eq. 5a) and momentum equation (Eq. 5b), we obtain the equation governing the fluid flow in terms of pressure field only with appropriate boundary conditions:

$$\Delta p = 0 \qquad \text{in } \Omega_F \quad (7a)$$

$$\frac{1}{\rho} \frac{\partial p}{\partial n} = -\ddot{s}_{c,i} n_i \quad \text{at each } \gamma_c \quad (7b)$$

$$\frac{\partial p}{\partial n} = 0 \qquad \text{at } \partial \Omega_F \quad (7c)$$

where $\frac{\partial p}{\partial n}$ is the normal derivative of the pressure field to the surface $\gamma_c$ of each tubes.

## HOMOGENISATION PROCESS

The equations governing the problem are written for the whole system, where, due to the geometrical configuration of
the bundle, the main unknowns (pressure field $p$ and tube displacement $s$) strongly vary in the scale of an elementary cell. Therefore, the variations of the unknowns may be considered depending on two different scales: $(x,y)$. The “long scale” $x$ corresponds to a variation in the whole bundle, while the “short scale” $y$ is representative of variations in a single elementary cell $Y$ (Fig. 5).

The dimensionless short scale variable $y$ is defined by the following formula:

$$y_1 = \frac{x_1}{\epsilon}, \quad y_2 = \frac{x_2}{\epsilon} \qquad (8)$$

The introduction of two scales for the variation of main quantities leads to the “chain derivation rule” (Fig. 2) as seen in [11] or in [10]:

$$\frac{d}{dx_i} = \frac{\partial}{\partial x_i} + \frac{1}{\epsilon} \frac{\partial}{\partial y_i} \qquad (9)$$

The “chain derivation rule” leads to a new formulation for the Laplace operator and directional derivatives, that are now written in the following form:

$$\begin{aligned} \Delta &= \frac{\partial^2}{\partial x_1^2} + \frac{2}{\epsilon} \left( \frac{\partial^2}{\partial x_1 \partial y_1} \right) + \frac{1}{\epsilon^2} \left( \frac{\partial^2}{\partial y_1^2} + \frac{\partial^2}{\partial y_2^2} \right) \\ \frac{d}{dn} &= n_1 \frac{\partial}{\partial x_1} + \frac{1}{\epsilon} \left( n_1 \frac{\partial}{\partial y_1} + n_2 \frac{\partial}{\partial y_2} \right) \end{aligned} \qquad (10)$$

## Displacement Field

The solution for the displacement field $s_c$, in the perspective of homogenisation, is searched in term of a continuous function $s(x_1)$ in the long scale variable $x_1$ only. The function must satisfy the following relationship with the discrete vector field $s_c$ (Fig. 3):

$$s_c = s(x_c) \quad \forall x_c \qquad (11)$$

However, we need a piece-wise function for the successive short scale analysis. We can achieve this purpose thanks to the following elementary Taylor development in the two scale $x,y$ for the displacement field:

$$\begin{aligned} s^{\epsilon}(x_1, y_1) = & s(x_1) - \epsilon y_1 \frac{\partial s(x_1)}{\partial x_1} + \epsilon^2 \frac{y_1^2}{2!} \frac{\partial^2 s(x_1)}{\partial x_1^2} - \epsilon^3 \frac{y_1^3}{3!} \frac{\partial^3 s(x_1)}{\partial x_1^3} + \\ & \mathcal{O}(\epsilon^4) \end{aligned} \qquad (12)$$

<figure><img src="image_5.png" /><figcaption>**FIGURE 3.** SCHEMATIC REPRESENTATION OF THE DISPLACEMENT FIELD DEVELOPMENT. A MULTIPLICITY OF SOLUTIONS SATISFIES THE CONDITIONS ASSUMED ON THE FUNCTION $s(x_1)$ (SOLID AND DASHED LINES ARE TWO POSITIVE SOLUTION FOR THE DISPLACEMENT FIELD)</figcaption></figure>

where the function $s^\epsilon(x_1, y_1)$ is a piece-wise constant function in each cell (solid black line in Fig. 3). The development proposed is made of by a function in long scale variable $x_1$ $s(x_1)$ continuous on the whole tube bundle domain and its successive derivatives. The chosen development represents the link between the discrete vector displacement field $s_c$ of each tubes and the continuous function $s(x_1)$. In fact, independently from the form of $s(x_1)$, the development represents a piece-wise functions in each elementary cell. We can check it thanks to the “chain derivation rule”:

$$\frac{ds^{\epsilon}}{dx_1}(x_1, y_1) = \frac{\partial s^{\epsilon}}{\partial x_1} + \frac{1}{\epsilon} \frac{\partial s^{\epsilon}}{\partial y_1} = 0.$$

and so the function $s^\epsilon(x_1, y_1)$ is constant in each elementary cell $Y$. Nevertheless, the condition imposed by Eq. 11 to $s(x_1)$ leads to an infinity of possible solutions for the long scale function that may satisfy the proposed development of Eq. 12. An analysis of the uniqueness of the solution for $s(x_1)$ will be discussed hereinafter.

## Pressure Field

As made for the displacement field, we search for the pressure field a solution in the form of a development in long and short scales $x,y$. We choose the development of Eq. 13 that, due to the physics of the problem, describes a pressure field as a linear combination between shape functions and long scale tube acceleration and its successive derivatives. The development represents a pressure field where its variations are provoked only by
<figure><img src="image_6.png" /><figcaption>FIGURE 4. SCHEMATIC REPRESENTATION OF THE PRESSURE FIELD DEVELOPMENT WHERE THE LONG AND SHORT SCALES ARE POINTED OUT</figcaption></figure>

structure oscillations since there is no global fluid flow within the bundle.

$$p^{\epsilon}(x_1, y) = b_0 \rho \int_{0}^{x_1} \ddot{s}(x) dx + \epsilon \rho \ddot{s}(x_1) (b_1 + f_1(y_1, y_2)) + \\ \epsilon^2 \rho \frac{\partial \ddot{s}(x_1)}{\partial x_1} (b_2 + f_2(y_1, y_2)) + \\ \epsilon^3 \rho \frac{\partial \ddot{s}(x_1)}{\partial x_1} (b_3 + f_3(y_1, y_2)) + \mathcal{O}(\epsilon^4) \quad (13)$$

The proposed development for pressure field $p$ is a sum of variations in the two scales: the long scale behaviour is given by $b_i$ constant coefficients while the fast variations at the level of an elementary cell are described via the functions $f_i(y_1, y_2)$ defined in the unitary elementary cell fluid region $Y^*$ (see Fig. 5).

In order to ensure the continuity of the pressure field between adjacent cells, condition of periodicity are imposed on the functions $f_i$. Moreover, the functions $f_i$ have got a zero average value at the boundary $\partial Y$ of the unitary elementary cell $Y$ for ensure the uniqueness of the function itself.

A specific change of variable is introduced for the functions $f_i$:

$$\begin{align}
f_1(y) &= \chi_1(y) - b_0 y_1 \nonumber \\
f_2(y) &= \chi_2(y) + b_0 \frac{y_1^2}{2} - b_1 y_1 + \frac{b_0}{8} - y_1 \chi_1 \nonumber \\
f_3(y) &= \chi_3(y) - b_0 \frac{y_1^3}{6} + b_1 \frac{y_1^2}{2} - b_2 y_1 + \frac{b_1}{8} - y_1 \chi_2 + \frac{y_1^2}{2} \chi_1 - \frac{b_0}{8} y_1 \nonumber \\
&\dots \tag{14}
\end{align}$$

<figure><img src="image_7.png" /><figcaption>FIGURE 5. ELEMENTARY UNITARY ELEMENT $Y$ SCALED BY A FACTOR $\epsilon$. THE ELEMENTARY CELL REPRESENTS A SINGLE TUBE SURROUNDED BY FLUID. IT IS THE BASIC REPETITIVE ELEMENT CONSTITUTING THE TUBE BUNDLE. THE FIGURE SHOWS THE MAIN GEOMETRICAL PARAMETERS OF THE UNITARY CELL</figcaption></figure>

This specific change of variables find its justification in the following; in fact, thanks to it, we are able to declare an easier problems series for short scale analysis in terms of function $\chi_i(y)$ rather than in terms of $f_i(y)$.

Moreover, the change of variables of Eq. 14 with Eq. 12 leads to a new development for pressure field $p$ determined up to a constant value:

$$\begin{align}
p^{\epsilon}(x_c, y) = \epsilon \rho \ddot{s}(x_c) \chi_1(y_1, y_2) + \epsilon^2 \rho \frac{\partial \ddot{s}(x_c)}{\partial x_1} \chi_2(y_1, y_2) + \nonumber \\
\epsilon^3 \rho \frac{\partial^2 \ddot{s}(x_c)}{\partial x_1^2} \chi_3(y_1, y_2) + \mathcal{O}(\epsilon^3) \tag{15}
\end{align}$$

The latter development of pressure field $p$, as the previous one, represents a pressure field where tubes accelerations and its successive derivatives are the cause of pressure fluctuations and there is no global fluid flow across the bundle.

The pressure field, due to the physics of the problem, is assumed as a linear combination of structure acceleration $\ddot{s}(x_1)$ and shape function $\chi_i(y)$. However, in this development, the solid acceleration and its successive derivatives are evaluated at the tube's axis ($\ddot{s}(x_1) = \ddot{s}(x_c)$). The $\chi_i$ functions are functions of the short scale variable $y$ only and they are defined in $Y^*$, fluid region of the unitary elementary cell $Y$ (Fig. 5).
## Construction of the short scale problem series

We now proceed to the homogeneisation of the fluid structure interaction problem. Firstly the fluid problem is analysed from which we can obtain the description of the fluid force acting on the structures in terms of long scale variables only.

The introduction of Eq. 13 in the Laplace operator for p of Eq. 7a with the application of the “chain derivation rule” (Eq. 9) gives rise to a series of problems in short scale variable y in the fluid region $Y^*$ of the unitary elementary cell Y.

The expression of the problem series for the $f_i$ functions results to be complex where each problem recalls the previous ones in the series making its resolution difficult to achieve. However, the change of variable between the $f_i$ and $\chi_i$ functions (Eq. 14), makes the problem series largely simpler. In fact, as $\ddot{s}(x_1)$ act as a parameter, we obtain the following expression for each $\chi_i(y)$ problems, after regrouping them in terms of $\epsilon$ powers:

$$ 
\begin{aligned}
\epsilon^0) \quad & \Delta_y \chi_1 (y_1, y_2) = 0 \\
\epsilon^1) \quad & \Delta_y \chi_2 (y_1, y_2) = 0 \\
\epsilon^2) \quad & \Delta_y \chi_3 (y_1, y_2) = 0 \\
\dots
\end{aligned} 
\qquad (16) 
$$

We evaluate the boundary conditions at the surface $\gamma$ (tube's surface $\gamma_c$ scaled by a factor $\epsilon$) in the unitary elementary cell Y (Fig. 5). For this reason, we introduce the development of Eq. 13 in boundary conditions of Eq. 7c and thanks to the relationship of Eq. 10 we obtain a series of boundary conditions for each function $f_i$. These boundary conditions are expressed, as for the definition of the functions $f_i$, in a complex form. Once again, thanks to the change of variable of Eq. 14, we are able to obtain a series of simple boundary conditions for the functions $\chi_i$ at the tube surface $\gamma$ in the elementary cell Y:

$$ 
\begin{aligned}
\epsilon^0) \quad & \frac{\partial \chi_1}{\partial n_y} (y_1, y_2) = n_1 (y_1, y_2) \\
\epsilon^1) \quad & \frac{\partial \chi_2}{\partial n_y} (y_1, y_2) = 0 \\
\epsilon^2) \quad & \frac{\partial \chi_3}{\partial n_y} (y_1, y_2) = 0 \\
\dots
\end{aligned} 
\qquad (17) 
$$

where $n_1$ is the horizontal component of the normal vector at the boundary $\gamma$ and $\frac{\partial \chi_i}{\partial n_y}$ are the normal derivatives with respect to y.

To properly close the problem, repetitive boundary conditions are obtained for left and right boundaries of the elementary cell $\partial Y$. We declare the boundary conditions in order to assure

the continuity of pressure field p and fluid flow exiting and entering the unitary cell for the $\chi_i$ functions. For the functions $\chi_i$ we obtain the following relationship:

$$ 
\begin{aligned}
\epsilon^0) \quad & \chi_1^R = \chi_1^L \\
\epsilon^1) \quad & \chi_2^R - \hat{Y} \chi_1^R = \chi_2^L + \hat{Y} \chi_1^L \\
\epsilon^2) \quad & \chi_3^R - \hat{Y} \chi_2^R + \frac{\hat{Y}^2}{2!} \chi_1^R = \chi_3^L + \hat{Y} \chi_2^L + \frac{\hat{Y}^2}{2!} \chi_1^L \\
\dots
\end{aligned} 
\qquad (18) 
$$

and for the derivative with respect to $y_1$ we obtain:

$$ 
\begin{aligned}
\epsilon^0) \quad & \frac{\partial \chi_1^R}{\partial y_1} = \frac{\partial \chi_1^L}{\partial y_1} \\
\epsilon^1) \quad & \frac{\partial \chi_2^R}{\partial y_1} - \hat{Y} \frac{\partial \chi_1^R}{\partial y_1} = \frac{\partial \chi_2^L}{\partial y_1} + \hat{Y} \frac{\partial \chi_1^L}{\partial y_1} \\
\epsilon^2) \quad & \frac{\partial \chi_3^R}{\partial y_1} - \hat{Y} \frac{\partial \chi_2^R}{\partial y_1} + \frac{\hat{Y}^2}{2!} \frac{\partial \chi_1^R}{\partial y_1} = \frac{\partial \chi_3^L}{\partial y_1} + \hat{Y} \frac{\partial \chi_2^L}{\partial y_1} + \frac{\hat{Y}^2}{2!} \frac{\partial \chi_1^L}{\partial y_1} \\
\dots
\end{aligned} 
\qquad (19) 
$$

where $\chi_i^L$ and $\chi_i^R$ are the values of $\chi_i$ functions at left and right boundaries of the unitary elementary cell Y, i.e. :

$$ \chi_i^L = \chi_i(-\hat{Y}, y_2) \quad \chi_i^R = \chi_i(+\hat{Y}, y_2) \quad \text{where } \hat{Y} = 0.5 \qquad (20) $$

## Numerical Simulations of Short Scale Problems

We are going to show the numerical simulations for the elementary cell problem series of Eq. 16 with the appropriate boundary conditions. These simulations are obtained thanks to the Finite Elements code Cast3M [14] developed at the CEA of Saclay.

In Fig. 6, we show results for the first three problems of the series and in Fig. 7 we present fluid accelerations fields relative to each elementary problem. The numerical simulations show that the pressure field of the first problem only is due to the tube's movement whereas the others elementary pressure field are due to pressure (odd problems) or to fluxes (even problems) at the boundary of the elementary cell $\partial Y$ (see boundary conditions imposed by Eq 17 Eq. 19 and Eq. 18).

Figure 7 shows that the structure movement in the first elementary problem $\chi_1$ generates in the fluid domain $Y^*$ a mean acceleration proportional to the fluid volume displaced by the tube itself $|\gamma|$ (the surface of the tube in the elementary unitary cell Y). Instead, for the successive elementary problems, where the structure it is assumed fixed, the fluid presents a zero mean acceleration field: it maintains a stagnant global behaviour.
<figure><img src="image_8.png" /><figcaption>FIGURE 6. ISOVALUES OF THE FUNCTIONS $\chi_i$ FOR THE FIRST 3 ELEMENTARY PROBLEMS</figcaption></figure>

The boundary conditions imposed for $\chi_i$ even problems lead to an elementary pressure field provoked by opposite fluid fluxes imposed at the boundary $\partial Y$. The obtained pressure field generates a zero resultant force on the tube wall. On the other hand, the boundary conditions for the odd $\chi_i$ problems impose to the fluid a pressure drop across the elementary cell. In this case, the fluid reaction force is proportional to this pressure drop.

The convergence of the $\chi_i$ problems is also illustrated by the

<figure><img src="image_9.png" /><figcaption>FIGURE 7. FLUID ACCELERATIONS FOR THE FIRST 3 ELEMENTARY PROBLEMS. ONLY THE $\chi_1$ PROBLEM PRESENTS A NON ZERO MEAN ACCELERATIONS FIELD</figcaption></figure>

values of the coefficients $a_i$ (see below) corresponding to the fluid reaction force of the pressure field for the odd problems (Fig. 8). From Fig. 8 we can see that the coefficient's series presents a trend that can be approximated by the convergent series:

$$a_{2i} \approx \sum_{i=0}^{\infty} \frac{1}{(2i)!} \qquad (21)$$
<figure><img src="image_10.png" /><figcaption>FIGURE 8. COEFFICIENTS $a_{2i}$ FOR THE FIRST 18 PROBLEMS OF THE SERIES</figcaption></figure>

## Fluid Force

We introduce the development of the pressure field of Eq. 13 in the fluid force expression:

$$F_{F,i} = \varepsilon \int_{\gamma} \left( \varepsilon \rho \ddot{s}(x_c) \chi_1(y) + \varepsilon^2 \rho \frac{\partial \ddot{s}(x_c)}{\partial x_1} \chi_2(y) + \varepsilon^3 \rho \frac{\partial^2 \ddot{s}(x_c)}{\partial x_1^2} \chi_3(y) + \mathcal{O}(\varepsilon^4) \right) n_i(y) d\gamma \quad (22)$$

where we use the relationship between the integral at the surface $\gamma_c$ and the structure surface $\gamma$ in the unitary elementary cell Y ($\int_{\gamma_c}(\cdot) d\gamma_c = \varepsilon \int_{\gamma}(\cdot) d\gamma$).

As the structure acceleration field appears here as a function of the long scale variable $x_1$ only, the integral in the short scale y of the forces can be decomposed as follows:

$$\int_{\gamma} \chi_i n_j d\gamma = \int_{Y^*} \frac{\partial \chi_i}{\partial y_j} dy - \int_{\partial Y} \chi_i n_j d\partial Y \quad (23)$$

thanks to the Green theorem, where $i = 1, 2, 3, ...$ and $j = 1, 2$.

The previous decomposition allows us to express the fluid force acting on the tube's wall by constant values deriving from short scale analysis in the elementary cell Y.

For the first elementary problem $\chi_1(y)$ we obtain:

$$\int_{\gamma} \chi_1 n_1 d\gamma = |\bar{\gamma}| - a_0 \quad (24)$$

where $|\bar{\gamma}|$ is the structure's surface in the elementary unitary cell and $a_0$ is a constant value (see below).

The structure motion generates in the fluid an acceleration that puts on movement the fluid itself. This fluid acceleration is

proportional to the structure surface ($|\bar{\gamma}|$ in the recursive elementary element Fig. 5). The other terms of Eq. 24 ($a_0$) represents the mean value of the pressure drop across the elementary cell.

For the successive problems ($i = 2, 3, ...$), the fluid reaction force is expressed proportionally to the pressure drop in the elementary cell Y only:

$$\int_{\gamma} \chi_i n_1 d\gamma = -a_{i-1} \quad \text{for } i = 2, 3, ... \quad (25)$$

where the $a_{i-1}$ are constant coefficients.

In conclusion, the component of the fluid force in the horizontal direction is expressed by the following formulation:

$$F_{F,1}(x_1) = \varepsilon^2 \rho \ddot{s}(x_1) (|\bar{\gamma}| - a_0) - \varepsilon^4 \rho \frac{\partial^2 \ddot{s}(x_1)}{\partial x_1^2} a_2 - \varepsilon^6 \rho \frac{\partial^4 \ddot{s}(x_1)}{\partial x_1^4} a_4 + \mathcal{O}(\varepsilon^8) \quad (26)$$

where $|\bar{\gamma}|$ is the surface on the elementary unitary cell, the $a_i$ are constant values obtained by the numerical simulations and $\ddot{s}(x_1)$ is the continuous acceleration field. Moreover, in Eq. 26 we resort to the result that the mean pressure drop is zero for each even problems: in fact the odd $a_i$ are equal to zero.

As expected, the vertical component of the fluid force is a zero resultant force ($F_{F,2} = 0$).

## Dynamic Equation for the Tube Bundle

Thanks to the analysis at the level of a single elementary cell, we are able to express the fluid force in terms of the long scale variable $x_1$ only.

The hydrodynamic force of Eq. 26 is composed by a first term proportional to the tube acceleration field $\ddot{s}(x_1)$ representing the inertial effects caused on the oscillating structure by the presence of a dense fluid [5]. In addition, the expression of the fluid force presents a series of successive terms involving the spatial derivatives of the continuous acceleration field $\ddot{s}(x_1)$.

We introduce the following quantities,

$$m_s = m_0 \varepsilon^2 \quad k_s = k_0 \varepsilon^2 \quad (27)$$

valuable due to the proportionality of the tube mass and stiffness to the tube surface and the expression of the hydrodynamic force (Eq. 26) in the dynamic equation (Eq. 1). Thus, we obtain the following expression for the dynamic equation of the whole tube
bundle:

$$m_0\ddot{s}(x_1) + k_0s(x_1) = -\ddot{s}(x_1)M_{add} - \varepsilon^2\rho \underbrace{\frac{\partial^2\ddot{s}(x_1)}{\partial x_1^2}a_2 - \varepsilon^4\rho a_4 \frac{\partial^4\ddot{s}(x_1)}{\partial x_1^4} + \mathcal{O}(\varepsilon^6)}_{\text{Fluid Reaction Force}} \quad (28)$$

where:

$$M_{add} = -\rho (|\tilde{\gamma}| - a_0) \quad (29)$$

is a constant positive value proportional to the solid acceleration hence the name of *Added Mass* coefficient.

## DISCUSSION: MODAL ANALYSIS OF THE DYNAMIC EQUATION

As already said, the function $s(x_1)$ is not unequivocally defined: in fact, there are an infinity of solutions that can satisfy Eq. 28 and Eq. 11 (see section “Displacement Field”).

We illustrate the effects of spatial variations of the continuous displacement field $s(x_1)$ on bundle's water frequency. In order to achieve this purpose, a modal analysis of Eq. 28 is performed by replacing time and spatial derivatives with the formula:

$$\frac{\partial}{\partial t} = i\omega, \quad \frac{\partial}{\partial x_1} = i\frac{\lambda}{\varepsilon} \quad (30)$$

where $\omega$ is the tube angular frequency, $\varepsilon/\lambda$ the spatial period of $s(x_1)$ and $i$ is the imaginary unit. The assumptions made on the spatial and time derivatives correspond in searching a solution for $s(x_1)$ in the form of:

$$s(x_1, t) = \sin\left(\frac{\lambda}{\varepsilon}x_1\right) e^{i\omega t} \quad (31)$$

The introduction of the formula for space and time derivatives in Eq. 28 leads to:

$$-\omega^2 m_0 s(x_1) + k_0 s(x_1) = \omega^2 M_{add} s(x_1) + \omega^2 \rho \sum_{n}^{\infty} (-1)^n \lambda^{2n} a_{2n} s(x_1) \quad (32)$$

In order to have a solution of Eq. 32, whatever values are assumed by the function $s(x_1)$, the following statement must be

<figure><img src="image_11.png" /><figcaption>FIGURE 9. WATER TUBE BUNDLE FREQUENCY VS $\lambda$ FOR $m_0 = 10.\,\text{kg}\,\text{m}^{-2}$ AND $f_{air} = 10.\,\text{Hz}$ IN THE APPROXIMATION FOR THE COEFFICIENTS $a_i$ OF Eq. 21</figcaption></figure>

verified:

$$-\omega^2 \left( m_0 + M_{add} + \rho \sum_{n}^{\infty} (-1)^n \lambda^{2n} a_{2n} \right) + k_0 = 0 \quad (33)$$

Two different occurrences are possible corresponding to two distinct conditions on the parameter $\lambda$; in fact, we may have $\lambda = 0$ or $\lambda \neq 0$. When $\lambda = 0$, the bundle is oscillating as a whole and there are no spatial variations of the displacement field; instead when $\lambda \neq 0$ spatial variations of displacement field are important and not neglected.

In the first case, where $\lambda$ is equal to zero, the dynamic equation is reduced to a simple formulation and the resulting water frequency is given by:

$$f_{H_2O} = \frac{1}{2\pi} \sqrt{\frac{k_0}{m_0 + M_{add}}} \quad (34)$$

Thanks to this formula, the name of *Added Mass* for $M_{add}$ is well explained: the presence of the fluid acts on the structure as an added mass that makes the structure heavier leading to a decrease of bundle's water frequency (inertial effects [5]). This condition is represented by the gray dashed line of Fig. 9.

On the other hand, when spatial variations of the displacement field are not neglected ($\lambda \neq 0$), the hydrodynamic force explicitly depends on the spatial period $\lambda$. In this case, the bundle water frequency is expressed as follows:

$$f_{H_2O} = \frac{1}{2\pi} \sqrt{\frac{k_0}{m_0 + M_{add} + \rho \sum_{n}^{\infty} (-1)^n \lambda^{2n} a_{2n}}} \quad (35)$$

where the appearance of an extra contribution to the added mass
<figure><img src="image_12.png" /><figcaption>FIGURE 10. EXTRA ADDED MASS VS $\lambda$ IN THE APPROXIMATION OF THE COEFFICIENTS $a_i$ OF Eq. 21</figcaption></figure>

term, depending on the parameter $\lambda$ and the coefficients $a_{2n}$, is shown.

From the numerical results, the series of the coefficients $a_{2n}$ appearing in the extra added mass term of Eq. 35, can be approximated by the series $1/2n!$ (Eq. 21). Consequently, the contribution $\rho \sum_{n}^{\infty} (-1)^n \lambda^{2n} a_{2n}$ to the added mass due to the spatial variations of the tube bundle presents a periodic appearance as seen in Fig. 10. Consequently, periodicity of the extra added mass leads to periodicity in water frequency as shown in Fig. 9.

As Fig. 9 shows, when the parameter $\lambda$ increases, corresponding to a decreasing of the spatial period of $s(x_1)$, the frequency in water decreases until it reaches a minimum value. The minimum value is obtained for a value of $\lambda$ equal to $\pi$. This is the case of the black solid line of Fig. 11. In this configuration, the spatial period is equal to two times the size of an elementary cell. Tubes are oscillating in opposition provoking the higher added mass leading to the higher water frequency drop. The same configuration may be obtained with a spatial period smaller than the size of an elementary cell as shown by the grey solid line of Fig. 11. This configuration corresponds to the case where the parameter $\lambda$ assumes the value $3\pi$. From Fig. 9, we see that, when the parameter $\lambda$ is equal to $3\pi$, the water frequency presents a minimum of its values. This minimum value, due to periodicity condition, is the same of $\lambda = \pi$. More solutions of Eq. 28 are possible for a donor condition on $s(x_1)$, but the frequency obtained, due to the periodicity, is the same.

From Fig. 9, as the parameter $\lambda$ increases over $\pi$, water frequency increases, until it reaches a maximum value when $\lambda = 2\pi$. This is the case of black solid line of Fig. 12. In fact, the spatial period is equal to the size of an elementary cell thus each tubes present the same acceleration value and they are moving as an unique body. This condition is comparable to a constant value of tube displacement field, where no spatial variations are taken into account. Indeed, the frequency obtained is equal to the frequency of no spatial variations (grey dashed line of Fig. 9).

Moreover, in this latter case it may be demonstrated that the

<figure><img src="image_13.png" /><figcaption>FIGURE 11. UNIQUENESS OF THE SOLUTION $s(x_1)$. TWO DIFFERENT SOLUTION GIVE THE SAME DYNAMIC EFFECTS ON THE TUBE BUNDLE FREQUENCY</figcaption></figure>

problem describing the pressure field is the elementary problem $\chi_1(y)$ only.

In fact, for a sinusoidal movement with a spatial period of $\varepsilon/2\pi$ (see Fig. 12), we can introduce a new elementary problem $\zeta_1$, such that:

$$\zeta_1 = \chi_1(y) + \sum_{n}^{\infty} (-1)^n (2\pi)^{2n} \chi_{2n+1} \quad (36)$$

with its proper boundary conditions:

$$\begin{align}
\zeta_1(-0.5, y_2) &= \zeta_1(0.5, y_2) \quad \text{at } \partial Y \nonumber \\
\frac{\partial \zeta_1}{\partial n_y} &= n_i \quad \text{at } \gamma \tag{37}
\end{align}$$

The elementary problem defined by the function $\zeta_1$ corresponds precisely to the $\chi_1$ elementary problem. This results confirms what found with the previous analysis; in fact, as the bundle oscillates as a whole (for each tube the acceleration is the same) the pressure field correspond exactly to the elementary pressure field obtained by the function $\chi_1$ only.

In conclusion, we can deduce some fundamental results from the modal analysis conducted.

Firstly, the spatial variation of the tube displacement field leads to extra contribution on the hydrodynamic force, affecting directly the tube bundle frequency. The extra contribution gives rise to a decrease in tube bundle frequency corresponding to an extra added weight to the structure.

The contribution to the fluid force from spatial variations of $s(x_1)$ presents a periodic behaviour. The maximum frequency
<figure><img src="image_14.png" /><figcaption>FIGURE 12. TUBE BUNDLE OVERALL MOVEMENT: A TUBE MOVEMENT $s(x_1)$ DESCRIBED BY CONSTANT VALUE OR BY SPATIAL VARIATIONS GIVES RISE TO THE SAME EFFECTS IN TERM OF FLUID REACTION FORCE AND FREQUENCY</figcaption></figure>

drop is obtained when the spatial period of $s(x_1)$ is equal to two times the size of an elementary cell corresponding to tubes opposition movement. Instead, the minimum frequency drop is obtained when the bundle oscillates as a whole. However, due to periodicity of the bundle's frequency, the maximum and minimum values are obtained also when the displacement field spatial period is lower than the cell size. In this case, the water frequency assumes values already assumed in $\lambda \in [0, \pi]$.

Therefore, for an infinity solutions of $s(x_1)$ for a given configuration $s_c$, we can find one value of water frequency in the interval $[0 : \pi]$: i.e. the dynamic behaviour of the bundle is not influenced by the infinity of solutions for the continuous function $s(x_1)$.

## CONCLUSION AND PERSPECTIVES

In the present paper, the authors analysed, thanks to an homogenisation method, a fluid-structure interaction problem. The simple geometrical configuration allows to investigate and better understand the effects of spatial variations of tube bundle displacement field on its dynamics.

The authors developed the displacement and pressure field, unknown of the FSI problem, according to two different scales: long and short scale. These latter describe different kinds of variations of the main unknowns due to the particular geometrical configuration of the tube bundle. Thanks to the study at the local scale of an elementary unitary cell, the authors are able to introduce a series of problems, by means of which they describe fluid pressure and tube bundle displacement at long scale.

An equation characterising the homogenised dynamic behaviour of a tube bundle immersed in a stagnant fluid, where there is no global flow across the bundle, is obtained. The au-

thors point out the appearance in the expression of the fluid force of a term proportional to solid acceleration corresponding to the inertial effects. In addition, the expression obtained for the fluid reaction force presents a series of terms directly proportional to spatial variations of long scale tube bundle acceleration.

A modal analysis of the homogenised dynamic equation is conducted, where the lowering of the water frequency due to spatial variations is illustrate. Furthermore, thanks to the modal analysis, we highlighted the uniqueness of the solution for the dynamic equation in terms of frequency.

The analysis conduct allows the authors to assess the problem of fluid structure interaction via an homogenisation method and to identify the main dynamic characteristic described by an homogenised method. Moreover, the analysis conducted leads to a solid base for the further step in whole reactor core homogenisation: the application to Navier-Stokes equations in order to describe with an higher detail the fluid flow and introduce the dissipative effects naturally described with this later model.

## ACKNOWLEDGMENT

The present work has been supported by the CEA of Saclay in the frame of the RNR-Na (Sodium Fast Cooled Reactor) research and development program. The authors wish to thank Hélène Dumontet and Arnault Monavon for the precious support in the advancement of these works.

## REFERENCES

[1] Bensoussan, A., Lions, J. L., and Papanicolaou, G., 1978. *Asymptotic Analysis for Periodic Structure*. North-Holland Publishing Company.

[2] Bergman, D., Lions, J. L., Papanicolaou, G., Murat, F., Tartar, L., and Sanchez-Palencia, E., 1985. *Les méthodes de lèhomogénéisation: théorie et applications en physique*. Collection de la Direction des études et recherches d'Électricité de France. Eyrolles.

[3] Paidoussis, M. P. *Fluid-Structure Interactions: Slender Structures and Axial Flow*, Vol. 1. Academic Press.

[4] Axisa, F., and Antunes, J., 2006. *Modelling of Mechanical Systems: Fluid-Structure Interaction*. Modelling of mechanical systems. Elsevier Science.

[5] Fritz, R. J., 1972. "The Effect of Liquids on the Dynamic Motions of Immersed Solids". *Journal of Engineering for Industry, 94*(1), pp. 167-173.

[6] Chen, S. S., Wambsganss, M. W., and Jendrzejczyk, J. A., 1976. "Added Mass and Damping of a Vibrating Rod in Confined Viscous Fluids". *Journal of Applied Mechanics, 43*(2), pp. 325-329.

[7] Planchard, J., 1984. Homogénéisation des matrices de masse ajoutée dans le cas de grands faisceaux tubulaires. Tech. rep., Note Electricité de France.
[8] Shinohara, Y., and Shimogo, T., 1981. "Vibrations of Square and Hexagonal Cylinders in a Liquid". *Journal of Pressure Vessel Technology*, **103**(3), pp. 233-239.

[9] Jacquelin, E., Brochard, D., Trollat, C., and Jézéquel, L., 1996. "Homogenization of a Non-Linear Array of Confined Beams". *Nuclear Engineering and Design*, pp. 231-223.

[10] Hammami, L., 1990. "Contribution à l'Étude de l'Interaction Fluide Structure Dans Les Faisceaux de Tubes Par Une Technique d'Homogénéisation". PhD thesis, Université Paris VI.

[11] Conca, C., Planchard, J., Thomas, B., and Vanninathan, M., 1994. *Problèmes Mathématiques En Couplage Fluide-Structure Applications Aux Faisceaux Tubulaires*. Collection de la Direction des Etudes et Recherches d'Electricité de France.

[12] Broc, D., and Sigrist, J. F., 2010. "Fluid Structure Interaction for Tube Bundles, Physical Meaning of a Rayleigh Damping". *Proceedings of the ASME 2010 Pressre Vessel & Piping Conference*.

[13] Broc, D., and Desbonnets, Q., 2011. "Fluid Structure Interaction Modelling for the Vibration of Tube Bundles: Part II-Homogenization Method Based on the Navier Stokes Equations". *Proceedings of the ASME 2011 Pressre Vessel & Piping Conference*, pp. 33-40.

[14] Cast3M, 2016. "Cast3M Code de Calcul Aux Éléments Finis".