International audienceWe introduce here some possible way to deal with the unsteady interfacial coupling of two distinct two-phase flow codes. The left code relies on the standard two-fluid approach, and the right code provides approximations of solutions of the Homogeneous Relaxation Model (HRM). The basic idea of the coupling method is to introduce a father model, which corresponds to the two-fluid two-pressure approach, and to define convective fluxes through the interface, once a prolongation of initial conditions on both sides of the coupling interface has been performed. A sketch of the algorithm is provided. More details are given in the paper [HER 06]. A few numerical results show a rather good behaviour of the coupling approach

Hurisse, Olivier

Hérard, Jean-Marc

English

Archive Ouverte en Sciences de l'Information et de la Communication

HAL Id: hal-01580997
https://hal.archives-ouvertes.fr/hal-01580997
Submitted on 5 Dec 2019
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.
The numerical coupling of a two-fluid model with an
homogeneous relaxation model
Jean-Marc Hérard, Olivier Hurisse
To cite this version:
Jean-Marc Hérard, Olivier Hurisse. The numerical coupling of a two-fluid model with an homogeneous
relaxation model. Finite Volumes for Complex Applications V, Jun 2008, Aussois, France. pp.503-510.
￿hal-01580997￿
The Numerical Coupling of a Two-Fluid Model
with an Homogeneous Relaxation Model
Jean-Marc Hérard — Olivier Hurisse
EDF-R&D, Department of Fluid Dynamics, Power Generation and Environment
6 Quai Watier, F-78401 Chatou Cedex
jean-marc.herard,olivier.hurisse@edf.fr
ABSTRACT. We introduce here some possible way to deal with the unsteady interfacial coupling
of two distinct two-phase flow codes. The left code relies on the standard two-fluid approach,
and the right code provides approximations of solutions of the Homogeneous Relaxation Model
(HRM). The basic idea of the coupling method is to introduce a father model, which corresponds
to the two-fluid two-pressure approach, and to define convective fluxes through the interface,
once a prolongation of initial conditions on both sides of the coupling interface has been per-
formed. A sketch of the algorithm is provided. More details are given in the paper [HER 06]. A
few numerical results show a rather good behaviour of the coupling approach.
KEYWORDS: interfacial coupling, unsteady two-phase flows, finite volumes
1. Introduction
Different two-phase flow models are used by the nuclear community, in order to
investigate various situations in the cooling circuit of nuclear power plants. Among
these, we should at least mention the standard six-equation two-fluid model (noted
STFM afterwards) and the four-equation homogeneous relaxation model (noted HRM
in the following). The latter HRM is mainly used for nuclear reactor core applica-
tions, but also for standard computations of flow patterns in steam generators, with
THYC, FLICA or GENEPI software. The former STFM is actually the keystone of
the CATHARE code for instance, which stands for the reference code for French nu-
clear authorities, especially for simulations of the whole coolant circuit (see [GUE 07]
for a more detailed presentation of the whole framework). Some situations may re-
quire the interfacial coupling of these codes, either for steady or unsteady numerical
experiments. This has motivated the development of new methods and algorithms to
cope with these coupled cases. The present work aims at providing a sketch of the
whole coupling strategy which has been proposed for such a purpose.
Some recent work on the interfacial coupling of hyperbolic systems of conser-
vation laws has emerged quite recently. Some references ([AMB 08a], etc.) can be
found for instance on the website [AMB 03], which basically focus on the interfa-
cial state coupling (see [GOD 96] among others), or alternatively on the flux coupling
techniques. Roughly speaking, the former favours the continuity of the conserva-
tive variable, whereas the latter aims at ensuring the intrinsic conservation laws. The
reader is for instance referred to [GOD 08], which provides an extensive description
of recent work on that topic. The method that is proposed in this paper relies on the
use of a coupling model, from which both models on each side of the coupling inter-
face may be at least formally recovered. We first briefly present both models to be
coupled through the interface in section 2. Then we show which coupling model is
used in the present situation, and we provide a sketch of the coupling algorithm in the
one-dimensional framework. Since all models involved here are invariant under frame
rotation, the extension to the 2D or 3D case is achieved in a straightforward manner
by using a projection along the normal to the coupling interface. A few numerical
results are eventually presented in a 2D framework. More useful details on the whole
approach can be found in [HER 06]. The reference [HER 07a] also provides some
recent ongoing work on the interfacial coupling of existing codes. Since the method
should apply for pre-existing codes, we do not detail the algorithms that are used in
sole codes, and focus on the treatment of boundary conditions through the coupling
interface in sections 3 and 4.
2. Models to be coupled
2.1. The left standard two-fluid model: STFM
The derivation of the standard two-fluid approach may be found in [ISH 75] for
instance. The governing equations of the STFM are the following:
∂mk
∂t
+
∂mkUk
∂x
= 0 ;
∂mkUk
∂t
+
∂mkU
2
k
∂x
+ αk
∂P
∂x
= Ik(W ) ;
∂αkEk
∂t
+
∂αkUk(Ek + P )
∂x
+ P ∗I
∂αk
∂t
= V ∗I Ik(W ) ;
[1]
where the state variable is: W tSTFM = (m1,m2,m1U1,m2U2, α1E1, α2E2). Here,
the void fractions αk ∈ [0, 1] comply with α1 + α2 = 1. Moreover, Uk, ρk, mk =
αkρk respectively stand for the velocity, density and partial mass within phase k.
Index 1 refers to the vapour phase, while the water phase is labeled by index 2. The
equilibrium pressure P is the same for both phases. A closure law for the interfacial
momentum transfer term Ik, which takes drag effects into account, must be given. We
thus introduce the classical contribution:
Ik = (−1)k m1m2
τU (m1 +m2)
(U1 − U2) [2]
where τU is the velocity relaxation time, which can be deduced from classical formu-
las of drag coefficients. A closure law relates P ∗I to WSTFM , we consider here the
rough closure: P ∗I = P , but other forms such asP ∗I = P+β(W )(m1+m2)(U2−U1)2
might be considered. The interface velocity V ∗I is usually taken as the velocity of the
dilute phase (which corresponds to the vapour phase in nuclear applications). The to-
tal energy within phase k equals Ek = ρkU2k/2 + ρkek(P, ρk). An equation of state
(EOS) is of course needed, which prescribes the internal energy ek(P, ρk) within each
phase. Hence the system is complete. We recall that the STFM is not hyperbolic,
unless U1 = U2, when using the closure P ∗I = P .
2.2. The right homogeneous relaxation model: HRM
The HRM is probably the most widespread model in nuclear two-phase thermohy-
draulics. The derivation of the model is straightforward when starting from the two-
fluid approach ([ISH 75]). If we assume a vanishing relative velocity Ur = U2−U1 =
0, the governing equations of the main variable W tHRM = (ρ, ρC, ρU,E) that de-
scribes the mixture are the following:
∂ρ
∂t
+
∂ρU
∂x
= 0 ;
∂ρC
∂t
+
∂ρUC
∂x
= 0 ;
∂ρU
∂t
+
∂ρU2
∂x
+
∂P
∂x
= 0 ;
∂E
∂t
+
∂U(E + P )
∂x
= 0 .
[3]
Of course, all interfacial transfer terms have vanished when summed up. Once more,
we use classical notations, and ρ stands for the density of the mixture (in other words
ρ = m1 +m2 = α1ρ1 + α2ρ2), while U,P,C respectively stand for the velocity of
the mixture, the mean pressure and the vapour concentration (that is U = (m1U1 +
m2U2)/(m1 +m2), P = α1P1 +α2P2 and C = m1/(m1 +m2) ). The total energy
of the mixture is E = ρU2/2 + ρe(ρ, P,C), where the closure law for e(ρ, P,C) is
given by the right code. Other closure laws for the relative velocity Ur can be found in
the literature. If we consider the closure law provided in [GUI 07] for instance, it must
be emphasized that the latter does not modify the convective form of the homogeneous
model [3]. System [3] is hyperbolic if:
(P/ρ− ρ∂e
∂ρ
|ρ,C)/( ∂e
∂P
|ρ,C) > 0.
A classic entropy inequality holds here for the HRM.
3. The coupling model
The coupling model (or father model) enables us to retrieve the left STFM model
by enforcing a pressure equilibrium, while it also leads to the right HRM model by
imposing a pressure-velocity-temperature equilibrium. The state variable of the cou-
pling model is: W t = (α1,m1,m2,m1U1,m2U2, α1E1, α2E2). The void fractions
αk ∈ [0, 1] still comply with α1 + α2 = 1. Moreover, Uk, Pk, ρk, mk = αkρk
respectively stand for the velocity, pressure, density and partial mass within phase k.
With some abuse of notation, we still note Ek = ρkU2k/2 + ρkek(Pk, ρk) the total
energy of phase k. We will thus focus herein on the following two-fluid model (see
[BAE 86, KAP 97] for instance):
∂α1
∂t
+ VI
∂α1
∂x
= φ1(W ) ;
∂mk
∂t
+
∂mkUk
∂x
= 0 ;
∂mkUk
∂t
+
∂mkU
2
k
∂x
+ αk
∂Pk
∂x
+ (Pk − PI)∂αk
∂x
= Ik(W ) ;
∂αkEk
∂t
+
∂αkUk(Ek + Pk)
∂x
+ PI
∂αk
∂t
= VIIk(W ) .
[4]
A closure law for the source terms φk must be added, and these should comply with:
φ1+φ2 = 0. We will prescribe: φk = (−1)kα1α2(P2−P1)/(τPΠ0), where τP stands
for the pressure relaxation time, while Π0 refers to some pressure constant given by
user. In practice, nuclear applications in pressurised water reactors obviously lead
to the choice (VI , PI) = (U1, P2), where the index 1 refers to the dilute (vapour)
phase, while the index 2 refers to the water phase. With this definition of (VI , PI), the
father model [4] is hyperbolic; moreover, shock solutions are uniquely defined though
the system contains some non-conservative contributions, and a meaningful entropy
inequality governs the motion of smooth solutions (see [COQ 02, GAL 04]). In the
following, we will note Z = (α1, ρ1, U1, P1, ρ2, U2, P2). The entropy inequality
which governs regular solutions of the seven-equation model is the exact counterpart
of the entropy inequality associated with the standard two-fluid model. If we define
the specific entropy sk(Pk, ρk) within each phase, which compels with: c2k
∂sk
∂Pk
|ρk +
∂sk
∂ρk
|Pk = 0, where: c2k = (
∂ek
∂Pk
|ρk)−1(Pk/ρ2k −
∂ek
∂ρk
|Pk) , we can introduce the
entropy η = m1s1 + m2s2, and the entropy flux fη = m1U1s1 + m2U2s2; thus,
regular solutions of both STFM and coupling model obey the inequality:
∂η
∂t
+
∂fη
∂x
≥ 0.
4. Main ideas of the coupling algorithm
Boundary conditions are imposed on each side of the coupling interface that sepa-
rates the STFM and the HRM finite volume codes, by solving the 1D Riemann prob-
lem associated with the LHS of [4] together with initial conditions ZL, ZR as defined
below. Eventually, the interface state Z∗ = ZRiemann(ZL, ZR, (x − xint)/∆t = 0)
is accounted for by both codes as a boundary condition through a flux prescription.
1) Defining initial conditions on both sides of the coupling interface
- Starting from the STFM-cell value on the left side of the interface:
Z−STFM = (α
−
1 , ρ
−
1 , U
−
1 , ρ
−
2 , U
−
2 , P
−)
we define ZL on the left side of the coupling interface as follows:
ZL = (α−1 , ρ
−
1 , U
−
1 , P
−
1 = P
−, ρ−2 , U
−
2 , P
−
2 = P
−)
- Starting from the HRM-cell value on the right side of the interface:
Z+HRM = (ρ
+, ρ+C+, ρ+U+, E+)
we define ZR on the right side of the coupling interface as:
ZR = (α+1 , ρ
+
1 , U
+
1 , P
+
1 = P
+, ρ+2 , U
+
2 , P
+
2 = P
+)
This is achieved by enforcing a pressure, velocity and temperature equilibrium be-
tween phases. This implies U+1 = U
+
2 = U
+
, noting: 1/Tk = (
∂sk
∂Pk
|ρk)/(
∂ek
∂Pk
|ρk):
P+1 = P
+
2 = P
+
T1(P+1 , ρ
+
1 ) = T2(P
+
2 , ρ
+
2 ) = T
+
α+1 ρ
+
1 = ρ
+C+
(1− α+1 )ρ+2 = ρ+(1− C+)
α+1 ρ
+
1 e1(P
+
1 , ρ
+
1 ) + (1− α+1 )ρ+2 e2(P+2 , ρ+2 ) = (E − ρU2/2)+
[5]
If EOS are chosen as stiffened gas EOS within each phase: (γk − 1)ρkek = Pk +
γkP∞,k, an explicit formula may be obtained. In that case, a direct elimination leads
to a non-linear equation wrt P+, which admits a unique solution (see [HER 06]) in
agreement with 0 ≤ P+ + P∞,φ. Hence, we obtain physically relevant values:
T+ = α
+
1 (P
++P∞,1)
ρ+C+CV,1(γ1−1)
α+1 =
1
1+ψ(P+)
ψ(P+) = (1−C
+)(CV,2)(γ2−1)(P++P∞,1)
C+(CV,1)(γ1−1)(P++P∞,2) .
[6]
The specific heat with constant volume has been noted CV,k within phase k.
2) Solving the Riemann problem at the coupling interface
The Riemann problem at the coupling interface associated with [4], ZL and ZR is
solved using either the Rusanov scheme or an approximate Godunov scheme such
as the one described in [GAL 04]. Once the intermediate state Z∗ defined above is
known, we may calculate the flux on the right side of the interface (HRM side) by
summing up both flux components associated with the total mass and momentum,
that is:
F+ρ = (α1ρ1U1 + α2ρ2U2)
∗
F+ρU = (α1(ρ1U
2
1 + P1) + α2(ρ2U
2
2 + P2))
∗
and the counterpart corresponding to the total energy flux:
F+E = (α1U1(E1 + P1) + α2U2(E2 + P2))
∗.
The mass flow rate for the vapour is F+ρC = (α1ρ1U1)∗. The implementation is
straightforward on the left (STFM) side of the coupling interface, since we only need
to skip the first component corresponding to the void fraction in [4].
5. Numerical results
In this first academic test case, the one-dimensional computational domain con-
tains 500 regular cells and the CFL number is CFL = 0.5. The coupling interface be-
tween the STFM domain (left side) and the HRM domain (right side) is located at x =
0. A Rusanov scheme has been used here to calculate convective fluxes through the
coupling interface. The initial condition is such that: P (x < x0, t = 0) = 150× 105,
P (x > x0, t = 0) = 155 × 105, where x0 = −0.15. The flow is at rest at the begin-
ning of the computation (U1(x, t = 0) = U2(x, t = 0) = 0), and the void fraction
is uniform (α1(x, t = 0) = 0.995). The pressure disequilibrium inside the STFM re-
gion generates a rarefaction wave which moves towards the right side. The grey-green
line in Figure 1 indicates the result obtained with the coupled simulation, whereas
the black-red line corresponds to the sole calculation using the STFM model on both
sides of the coupling interface x = 0. We have plotted here the state variable of the
HRM model (liquid mass fraction C, density ρ, velocity U , and pressure P ) in order
to obtain a relevant comparison. Obviously, we can notice a difference in the mass
fraction profile, which is essentially due to the fact that no drift velocity is accounted
for in the right HRM domain. Many other experiments can be found in [HER 06].
The second industrial test case corresponds to a rough representation of a transient
flow in the core of the reactor, after a sudden increase of the released power through
the fuel assembly. The central part of the core is governed by the HRM model, while
the flow in the coolant circuit in the three inlet and outlet pipes is represented by the
STFM model. The three coupling interfaces have been located inside the pipes, in
order to minimize transverse perturbations (see [HER 07b]). Though the interfacial
coupling procedure appears to be stable, the global analysis is difficult, and suffers
from the absence of any reference solution. More details can be found in [HER 07a].
Acknowledgments: this work has been achieved in the framework of the NEP-
TUNE project, with financial support from CEA (Commissariat à l’Énergie Atom-
ique), EDF (Électricité de France), IRSN (Institut de Radioprotection et de Sûreté
Nucléaire) and AREVA-NP. O. Hurisse benefited from a CIFRE/EDF contract during
his PhD thesis. This work has benefited from discussions with the group [AMB 03].
-1 -0.5 0 0.5 1
0.99951
0.99951
0.99952
0.99952
0.99952
0.99952
0.99952
liquid mass fraction
-1 -0.5 0 0.5 1
700
720
740
760
780
800
mean density
-1 -0.5 0 0.5 1
0
0.5
1
1.5
2
2.5
3
mean velocity
-1 -0.5 0 0.5 1
1.5e+07
1.52e+07
1.54e+07
1.56e+07
1.58e+07
1.6e+07
pressure
Figure 1. Simulation of waves passing through the coupling interface (at x = 0). The
left (resp. right) code corresponds to the two-fluid (resp. HRM) model. The green line
corresponds to the coupled simulation, the red one to the sole STFM computation.
6. References
[AMB 03] AMBROSO A., CHALONS C. , COQUEL F., GODLEWSKI E. , LAGOUTIÈRE
F., RAVIART P. A., SEGUIN N., Working group on the interfacial coupling of models,
http://www.ann.jussieu.fr/groupes/cea, 2003.
[AMB 08a] AMBROSO A., CHALONS C. , COQUEL F., GODLEWSKI E. , LAGOUTIÈRE F.,
RAVIART P. A., SEGUIN N., Coupling of general Lagrangian systems, Mathematics of
Computation, vol. 77, 2008, pp.909-941.
[AMB 08b] AMBROSO A., CHALONS C. , COQUEL F., GALIÉ T., Relaxation and numerical
approximation of a two-fluid two-pressure diphasic model, in these proceedings: Finite
Volumes for Complex Applications V, Aussois, June 8-13, 2008.
[BAE 86] BAER M.R., NUNZIATO J.W., A two-phase mixture theory for the deflagration to
detonation transition (DDT) in reactive granular materials, Int. J. Multiphase Flow, vol. 12-
6, 1986, pp. 861–889.
[COQ 02] COQUEL F., GALLOUËT T., HÉRARD J.M., SEGUIN N., Closure laws for a two-
fluid two-pressure model, C. R. Acad. Sci. Paris, vol. I-332, 2002, pp. 927–932.
[GAL 04] GALLOUËT T., HÉRARD J.-M., SEGUIN N., Numerical modelling of two phase
Figure 2. Pressure distribution (wrt a pressure reference P0 = 150.105Pa)
flows using the two-fluid two-pressure approach, Math. Mod. Meth. in Appl. Sci., vol. 14,
n◦5, 2004, pp. 663-700.
[GOD 08] GODLEWSKI E., Coupling of fluid models- Exploring some features of interfacial
coupling, invited lecture, in these proceedings: Finite Volumes for Complex Applications
V, Aussois, June 8-13, 2008.
[GOD 96] GODLEWSKI E., LE THANH K.C., RAVIART P.-A., The numerical interface cou-
pling of nonlinear hyperbolic systems of conservation laws: the case of systems, Math.
Mod. Numer. Anal., vol. 39, n◦4, 2005, pp. 649-692.
[GUE 07] GUELFI A., BESTION D., BOUCKER M., BOUDIER P., FILLON P., GRANDOTTO
M., HÉRARD J.M., HERVIEU E., PÉTURAUD P. , NEPTUNE : a new software platform
for advanced nuclear thermal hydraulics, Nucl. Sci. and Eng., vol. 156, 2007, pp. 281-324.
[GUI 07] GUILLARD H., DUVAL F., A Darcy law for the drift velocity in a two-phase flow
model, J. of Comput. Physics, vol. 224, 2007, pp. 288 –313.
[HER 06] HÉRARD J.-M., HURISSE O., Couplage interfacial d’un modèle homogène et d’un
modèle bifluide, EDF report H-I81-2006-04691-FR, in French, unpublished, 2006.
[HER 07a] HÉRARD J.-M., HURISSE O., Couplage d’un code basé sur un modèle HRM et du
code NEPTUNE-CFD, EDF report H-I81-2007-0863-FR, in French, unpublished, 2007.
[HER 07b] HÉRARD J.-M., HURISSE O., Coupling two and one-dimensional unsteady Euler
equations through a thin interface, Computers and Fluids, vol. 36(4), 2007, pp. 651-666.
[ISH 75] ISHII M., Thermo-fluid dynamic theory of two-phase flow, Collection de la Direction
des Etudes et Recherches d’Electricité de France, 1975.
[KAP 97] KAPILA A.K., SON S.F., BDZIL J.B., MENIKOFF R., STEWART D.S., Two-phase
modeling of a DDT: structure of the velocity relaxation zone, Phys. of Fluids, vol. 9-12,
1997, pp. 3885–3897.


The numerical coupling of a two-fluid model with an homogeneous relaxation model

https://hal.archives-ouvertes.fr/hal-01580997/document

The numerical coupling of a two-fluid model with an homogeneous relaxation model

Abstract

Similar works

Full text

Available Versions

Archive Ouverte en Sciences de l'Information et de la Communication