5. Axial flow (example: fuel assemblies)#
The integration of this fluid-elastic excitation model in*Code_Aster* was discussed in the specification note [bib. 2]. Model concept note MEFISTEAU [bib. 6] is the theoretical reference documentation.
5.1. Description of the studied configuration#
We consider a bundle of
circular cylinders that are mobile in bending and subject to an incompressible flow of viscous fluid, limited by a rigid cylindrical chamber of circular or rectangular cross section [Figure 5.1-a].
Figure 5.1-a: Beam under axial flow
The cylinders are all parallel, directed along the axis of the chamber. They have a common length, noted
. To simplify the notations, it is subsequently assumed that
is the guiding axis. The steady flow is axial and assumed to be uniform in each section. The density of the fluid can be variable along the axis
(thermal gradients), the speed of the steady flow also depends on the variable
.
5.2. Calculation steps#
The first step concerns the determination of the air modal base of the beam. This operation is carried out by the operator CALC_MODES. This step is essential because fluid-elastic forces are projected onto this base.
The second step involves taking fluid-elastic forces into account by operator CALC_FLUI_STRU. This step is broken down into 7 sub-tasks:
5.2.1. Pre-treatments#
1°/ |
By means of the topology of the mesh, deduction of the coordinates of the centers of the cylinders in the bundle and then verification of the correct arrangement of the cylinders in relation to each other (in particular, it is verified that there is no overlap between two cylinders) and in relation to the rigid chamber. |
2°/ |
Determination of the fluid excitation length, common to all cylinders, as well as of an associated discretization along the direction axis. |
3°/ |
Composition of the tables giving the modal deformations in air of each cylinder of the bundle, for each of the modes taken into account for the fluid-structure coupling. To do this, the deformed ones are interpolated at the points of the discretization determined previously. |
5.2.2. Solving the modal problem under flow#
4°/ |
Solving the disturbed fluid problem. Determining the potential of disturbed velocities requires the reversal of high-order linear systems calling for the implementation of the Crout method. |
For each flow speed
5°/ |
Calculation of the added mass, damping and stiffness matrices giving the fluid-elastic force transfer matrix in the air modal base:
|
6°/ |
Solving the modal problem under flow; we solve the complete problem with vectors and eigenics |
We do not overlook the extra-diagonal terms of
. After reformulation, the resolution is carried out using the QR algorithm: obtaining the masses, frequencies and reduced modal damping under flow
, complex modal deformations
expressed in the air base; of these, only the real part is retained after minimizing the imaginary part (calculating a criterion on the imaginary part).
7°/ |
Restitution of deformations under flow in the physical base.
|
End of loop on flow velocities
Notes:
Knowledge of the coordinates of the centers of the cylinders (pre-treatment 1°) is necessary to solve the disturbed fluid problem (sub-task 4°). This resolution leads to the estimation of the terms of the fluid-elastic force transfer matrix (sub-task 5°), which involve pressure and speed disturbances.
The determination of a common excitation length and the creation of an associated discretization (2° pre-treatment) make it possible to define a domain of integration on the structures for the projection of fluid-elastic forces on a modal basis. The interpolation of the modal deformations at the same points is therefore required (3° pre-treatment) .
The dynamic behavior of the beam under flow can also be studied using a simplified representation of the beam (with equivalent tubes). The calculation steps for taking into account the fluid-structure coupling are then identical to those described above, the only differences appearing in the pre-treatments. This second approach is described more precisely in the note [bib. 2]. In step 1° of the pre-treatments, the coordinates of the centers of the cylinders in the bundle are then specified directly by the user, who also establishes the correspondence between the cylinders of the beam and the beams in the simplified representation given by the mesh. In step 3° of the pre-treatments, this correspondence makes it possible to assign to the cylinders of the beam, to the discretization points determined in step 2°, the modal deformations of the beams in the simplified representation.
5.3. Solving the unsteady fluid problem#
5.3.1. Simplifying hypothesis#
H1 |
The field of unsteady fluid velocities is determined analytically by assuming that the disturbed flow
: |
eq. 5.3.1-1
Such a speed field allows sliding on the walls of the cylinders which will make it possible to calculate the viscous stress by a law of friction.
H2 |
The movement of the cylinders only generates speed disturbances
|
H3 |
The pressure field is decomposed into stationary and disturbed parts according to
|
eq. 5.3.1-2
where |
, the hydraulic diameter of the bundle, and the roughness of the walls. This coefficient is deduced from Nikuradze’s law characterizing pipe flows; |
refers to the gravity field. Its action on the stationary pressure field depends on the inclination of the .. image:: images/Object_348.svg
beam. |
|
5.3.2. Determining the potential of disturbed speeds#
We are looking for an analytical solution for
in the form of a superposition of elementary singularities that can be written as:
eq. 5.3.2-1
at the center of each cylinder \(k\) and:
eq. 5.3.2-2
at the center of the rigid enclosure when the latter is circular where:
|
refers to the truncation order of the Laurent series , |
|
designate the polar coordinates in a plane perpendicular to the axis , centered at the center of the cylinder k, |
|
designate the polar coordinates in a plane perpendicular to the axis , centered at the center of the circular rigid enclosure. |
The coefficients
and
expressions [:ref:` éq. 5.3.2-1 < éq. 5.3.2-1 >`] and [:ref:` éq. 5.3.2‑2 < éq. 5.3.2‑2 >`] are determined by applying the condition to the non-penetration limits:
on the outline of each mobile cylinder
, this condition is written as:
where |
, |
, |
|
. |
|
on the outline of a circular rigid enclosure, it is written:
where
refers to the radius of the enclosure.
In the case of a rectangular rigid enclosure, this condition is taken into account by a method derived from the « images » method [bib. 6]; the fluid problem confined by the rectangular enclosure is made equivalent to the problem in an infinite environment by creating images of the cylinders moving in the beam with respect to the sides of the enclosure. This method leads to the introduction of new singularities of the form [:ref:` éq. 5.3.2-1 < éq. 5.3.2-1 >`], placed at the center of the « image » cylinders, in the expression of
. However, it does not add any unknown to the problem since the coefficients for these new singularities are derived from those of the moving cylinders of the beam through the interplay of images.
Finally, the potential of disturbed speeds is written as:
eq. 5.3.2-3
Where
refers to the number of moving cylinders in the bundle. The functions
are linear combinations of
and
whose coefficients are determined by the preceding boundary conditions. This requires the resolution of linear systems of high orders and full matrices. The inversions are carried out by implementing the Crout method.
5.3.3. Fluid force modeling#
First, the forces due to pressure disturbances are retained.
, related to the potential of speeds disturbed by:
eq. 5.3.3-1
The resultant of the disturbed pressure field around each moving cylinder is a linear force
Acting next
. This force depends linearly on
and
, thus generating added terms of mass, damping and stiffness.
The forces linked to the viscosity of the fluid are then taken into account.
In a quasistatic approach, we consider the action of the fluid velocity field
around a cylinder at the moment
: in the coordinate system linked to the cylinder, the flow, of speed
in order 0, has an incidence in relation to the cylinder that is a function of speed disturbances and the movement of the cylinder itself. This results in drag force and lift force. It is shown that the following components
Of the resultant linear force
Spell, for the cylinder
:
eq. 5.3.3-2
eq. 5.3.3-3
where |
refers to the zero incidence slope of the lift coefficient around a very slightly inclined cylinder ( .. image:: images/Object_397.svg
= 0.08). |
and .. image:: images/Object_399.svg
designate the mean of the speed disturbances along the axes .. image:: images/Object_400.svg
and .. image:: images/Object_401.svg
along the axes and around the cylinders, which depend linearly on the .. image:: images/Object_402.svg
and .. image:: images/Object_403.svg
(cf. [:ref:` éq. 5.3.2-3 < éq. 5.3.2-3 >`]). |
|
These forces generate additional damping and stiffness terms.
Finally, the action of the stationary pressure field on deformed mobile structures is taken into account. It is shown that the resultant linear force
On the cylinder
has as components, in order 1:
eq. 5.3.3-4
eq. 5.3.3-5
These forces only generate terms of added stiffness and no coupling between cylinders.
The expressions [:ref:` éq. 5.3.3-1 < éq. 5.3.3-1 >`], [:ref:` éq. 5.3.3-2 < éq. 5.3.3-2 >`], and [:ref:` éq. 5.3.3-3 < éq. 5.3.3-3 >`] highlight the need to solve the disturbed fluid problem before estimating fluid-elastic forces.
5.3.4. Expression of terms in the fluid-elastic force transfer matrix#
Summary of linear forces
For each cylinder
, fluid-elastic forces are written as follows
and
:
eq. 5.3.4-1
and are linear combinations of:
Breakdown of the movement on a modal basis
The movement of the cylinder bundle is broken down as follows:
air vibration modes. We note
(
and
) the following deformations
Of the cylinder
corresponding to the j th mode of the beam. The components of the displacement of the neutral fiber of the cylinder
At the abscissa
can then be written as:
eq. 5.3.4-2
eq. 5.3.4-3
where
is the vector of generalized movements.
Force projection on a modal basis
Note
the projection of fluid-elastic forces according to the i th mode of the beam.
eq. 5.3.4-4
is a linear combination of
Note
the vector of modal fluid-elastic forces:
Which is written:
eq. 5.3.4-5
Where |
refers to the matrix of terms of mass added by the fluid, |
refers to the matrix of damping terms added by the fluid, |
|
refers to the matrix of stiffness terms added by the fluid. |
|
These matrices are real square in order
and their terms are independent of the movement of the structures. The matrix
is symmetric; the matrices
and
are not necessarily a step.
The projection of the equations of motion on a modal basis provides:
eq. 5.3.4-6
where |
refer to the matrices of masses, damping and stiffness of structures in air; these matrices are of order .. image:: images/Object_440.svg
and diagonal. |
In Laplace’s domain, the [:ref:` éq. 5.3.4-6 < éq. 5.3.4-6 >`] relationship becomes:
eq. 5.3.4-7
The fluid-elastic force transfer matrix is then introduced
defined by:
eq. 5.3.4-8
And we find the relationship [:ref:` éq. 1.2-1 < éq. 1.2-1 >`] from paragraph [§ 1.2]:
5.4. Solving the modal problem under flow#
The modal problem under flow is formulated by the relationship [:ref:` éq. 5.3.4-7 < éq. 5.3.4-7 >`] from the previous paragraph.
This problem is solved after rewriting as a standard problem with vectors and eigenvalues of the type
.
The new wording is as follows:
eq. 5.4- 1
Notes:
The dimension of the problem is doubled compared to that of the initial problem.
The properties of the matries
and
allow for reversal.
The solution to this problem is done by means of the QR algorithm. The modules implemented by the CALC_FLUI_STRU operator are the same as those used by CALC_MODES.
The problem with its own elements that we solve is a complex problem. We therefore obtain an even number of complex eigenvalues conjugated two by two. We only keep those whose imaginary part is positive or zero.
Eigenvectors are complex, defined to within a complex multiplicative constant. As only real modes are taken into account, the first step is to determine, for each eigenvector, the constant that minimizes the imaginary part of the vector with respect to its real part, in the sense of the Euclidean norm. The eigenvectors are then redefined with respect to this norm. Given the normalization used, it is then possible to keep only the real part of the eigenvectors in the mode_meca concept. However, in file MESSAGE, indicators on the relationships between the imaginary part and the real part of the eigenvectors thus normalized are restored, so that the user can judge the bias introduced by not taking into account the imaginary part of the normalized vectors.
5.5. Taking into account the presence of the tube bundle grids#
The modeling described above, of the forces induced by an axial flow on a bundle of cylinders, does not take into account the presence of the grids of the bundle (for example, the grids for mixing and maintaining fuel assemblies). A comparison of this model with tests carried out on model CHAISE (in the configuration of a bundle of nine flexible tubes including a grid) is presented in a summary note [bib. 8]: we note that the fluid-elastic coupling between the grid and the axial flow is not negligible step and that it generates an increase in the reduced modal damping of the tubes. The purpose of this paragraph is to describe the additional effects due to grids and how they are taken into account in model MEFISTEAU.
5.5.1. Grid configuration description#
The study is restricted here to two types of grids:
the support grids which are located at the ends of the bundle,
the mixing grids which are distributed between the holding grids.
The grids are all positioned perpendicular to the cylinder bundle and are in the form of a prismatic network with a square base on each side.
And in height
(along the \(x\) axis of the cylinders). Grids of the same type are characterized by identical dimensions.
5.5.2. Additional calculation steps#
The first additional step concerns the specification of the type of configuration of the grids by the operator DEFI_FLUI_STRU, then the verification of the correct arrangement of the grids in relation to one another, and in relation to the ends of the beam.
The second step concerns the resolution of the modal problem under flow. In the flow velocity loop, the fluid-elastic force transfer matrix in the air modal base is completed by the calculation of an added damping matrix and an added stiffness matrix, linked to the grids.
5.5.3. Modeling the fluid forces exerted on the grids#
Pressure jump calculation
First of all, the presence of grids disrupts the stationary pressure field
; each grid is considered to be a singularity causing a pressure jump, whose expression is in the form:
eq. 5.5.3-1
where |
refers to the coefficient of pressure loss due to the grid, |
refers to the stationary speed of the flow at the level of the grid, |
|
refers to the density of the flow at the level of the grid, |
|
refers to the axial position of the middle of the grid along the beam. |
|
Density
is calculated by linear interpolation of the density profile
of the flow in the absence of a grate. Stationary speed
is calculated using the conservation of the mass flow rate, which results in the following equation:
Calculation of the point fluid forces exerted on each grid
Using the same quasistatic approach as that carried out in paragraph [§5.3.3], we show that the action of the fluid velocity field
around a grid involves a drag force and a lift force, depending on the incidence of the flow in relation to the grid. The y and z components of the resulting point force f
are therefore written, for each elementary cell
Of a
where |
and .. image:: images/Object_466.svg
respectively designate the density and stationary speed profile of the flow at the foot of the beam, |
designates the fluid section of the beam in the absence of a grid, |
|
designates the fluid section of the beam at the level of the grid: .. image:: images/Object_469.svg
with .. image:: images/Object_470.svg
the solid section of the grid. |
|
From this we deduce the expression:
The pressure loss coefficient
is calculated using the expression for the total hydrodynamic force that applies to the grid, and we get:
eq. 5.5.3-2
The 1st term (en
) comes from drag effort;
is the drag coefficient of the grid. The 2nd term (en
) is a term correcting the friction force applied to the beam alone at the altitude of the grid (
is the wet perimeter of the beam in the absence of a grid).
By introducing the expression [:ref:` éq. 5.5.3-2 < éq. 5.5.3-2 >`] into the relationship [:ref:` éq. 5.5.3-1 < éq. 5.5.3-1 >`], we therefore get the expression for the pressure jump
for each altitude grid
. This pressure jump is taken into account when calculating the stationary pressure field.
, in the following manner:
grill:
where |
refers to the zero incidence slope of the lift coefficient around a very weakly inclined grid. |
). |
|
These forces will therefore generate additional terms of added damping and stiffness, which are obtained after modal decomposition of the movement and projection of these forces on the modal basis.
5.6. Taking into account the amortization of fluid at rest#
Until now, the damping brought to a bundle of tubes by the presence of a fluid at rest has not been taken into account in the modeling. We therefore propose here a model of damping in fluid at rest, for which annex 1 of the summary note of tests CHAISE [bib8] constitutes the reference documentation.
5.6.1. Modeling the fluid force at rest exerted on a bundle of tubes#
The method for calculating depreciation in fluid at rest, which is implemented here, is a generalization of the method of CHEN [bib. 9].
The aim is to calculate the resultant force on each pipe from the stresses due to shear in the boundary layer. This is a non-linear problem because the fluid damping coefficient depends on the frequency. The following simplifications are therefore introduced:
the problem is written using the frequencies in water at rest calculated without taking into account fluid damping,
the coupling between modes is neglected.
The linear force
Exercising on the tube
subject to a harmonic movement of the beam according to the mode
At the frequency
is given by the following relationship:
eq. 5.6.1-1
where |
|
refers to the sliding speed between the tube
|
|
Denotes the drag coefficient of a cylinder with radius , subject to a harmonic flow with an amplitude of infinity , and is defined by:
|
|
The relationship obtained by replacing [:ref:` éq. 5.6.1-2 < éq. 5.6.1-2 >`] and [:ref:` éq. 5.6.1-3 < éq. 5.6.1-3 >`] in equation [:ref:` 5.6.1-1 < 5.6.1-1 >`] is linearized by a Fourier series expansion (of the term
) from which we only remember the first term; it comes from:
Projection on a modal basis
By projection on a modal basis and by neglecting the coupling between modes, we obtain the generalized force exerted on the tube bundle according to the mode
:
is therefore proportional to
and the associated modal force vector
takes the form:
where
refers to the damping matrix added by the fluid at rest.