c
TU Delft, The Netherlands, 2006
A RELAXATION SCHEME TO APPROXIMATE THE
10-MOMENT GAUSSIAN CLOSURE
Christophe Berthon∗ ∗MAB, Universit´e Bordeaux 1
351, cours de la lib´eration, 33405 Talence Cedex, France e-mail: Christophe.Berthon@math.u-bordeaux.fr
Key words: hyperbolic system of conservation laws, Gaussian moment closure, relax-ation scheme, discrete entropy inequalities.
Abstract. A relaxation scheme is derived to approximate the weak solutions of the 10-moment gaussian closure system. This model is governed by an hyperbolic system of conservation laws supplemented by entropy inequalities. The proposed relaxation scheme is proved to satisfy both needed positiveness property and all discrete entropy inequalities. Numerical experiments illustrate the interest of the method.
1 INTRODUCTION
The relaxation schemes have recently been introduced by Jin-Xin [10] (see also Suliciu [15, 16]) to approximate the weak solutions of hyperbolic systems of conservation laws. In the sequel, this under consideration hyperbolic system will be called: equilibrium system. Motivated by the work of Liu [13] or Chen-Levermore-Liu [8], the weak solutions of the equilibrium system are approximated by the solutions of a relevant system of partial differential equations with singular perturbation: the so-called relaxation system. The relaxation model aims to restore the equilibrium model completed by the entropy in-equalities within the limit of an infinite relaxation parameter. The reader is referred to Aregba-Natalini [1] or Natalini [14] where a precise analysis is proposed into the context of the scalar conservation laws.
or entropy inequalities. In [3], the approximate solutions of the standard 3 × 3 Euler equations are proved to satisfy stability properties as soon as a relaxation scheme is used. In the present work, these results are extended to the 10-moment Gaussian closure model. The paper is organized as follows. In the next section, a brief description of the 10-moment Gaussian closure model is given. The third section is devoted to the relaxation model while in the section 4, we derive the relaxation scheme. The last section concerns the numerical experiments.
2 The 10-moment model
The usual local thermodynamic equilibrium cannot be assumed in the present work devoted to the simulations of low pressure rarefied gas flows. As a consequence, the stan-dard Euler equations cannot be considered and alternative models must be proposed. In a recent work, Levermore [11] derived a hierarchy of moment closure systems. The simplest model of this hierarchy coincides with the Euler equations. The second derived model, and investigated by Levermore-Morokoff [12], admits ten equations: the 10-moment Gaussian closure model. In one space dimension, this model is governed by the following system:
∂tρ + ∂x(ρu1) = 0, t > 0, x∈ R, ∂t(ρu1) + ∂x(ρu21+ p11) = 0, ∂t(ρu2) + ∂x(ρu1u2+ p12) = 0, ∂tE11+ ∂x((E11+ p11)u1) = 0, ∂tE22+ ∂x(E22u1+ p12u2) = 0, ∂tE12+ ∂x(E12u1+ 1 2(p11u2+ p12u1)) = 0, (1)
where the state laws are given by Eij =
1
2ρuiuj + 1
2pij, 1≤ i ≤ j ≤ 2.
With clear notations, the following condensed form is introduced for the sake of simplicity: ∂tw + ∂xf (w) = 0,
w =t(ρ, ρu1, ρu2, E11, E22, E12).
where w∈ Ω with the following definition of the space of the admissible states: Ω = w∈ R6; ρ > 0, (u
1, u2)∈ R2, p11> 0, p11p22− p212> 0
.
The system under consideration is easily seen to be hyperbolic over Ω. In addition, the weak solutions of (1) satisfy the following entropy inequalities:
where the functionsF, G : R → R verify F0(y) < 0, F00(y) F0(y) < 1 3, ∀y ∈ R, G0(y) < 0, G 00(y) G0(y) < 1 4, ∀y ∈ R, (2)
in order to enforce the functions w→ ρF(ln s(w)) and w → ρG(ln σ(w)) to be convex.
3 A relaxation model
The aim of the present work is to propose a numerical scheme which is robust, when preserving the numerical invariance of Ω, and stable with the satisfaction of all the discrete entropy inequalities. To access such an issue, a relaxation model is introduced. In this model, the equilibrium pressures p11 and p12 are approximated by two new unknowns Π11
and Π12. In the limit of an infinite parameter, the variables Π11 and Π12 relax to the
equilibrium pressures. The relaxation model under consideration reads as follows: ∂tρ + ∂x(ρu1) = 0, t > 0, x∈ R, ∂t(ρu1) + ∂x(ρu21+ π11) = 0, ∂t(ρu2) + ∂x(ρu1u2+ π12) = 0, ∂tE11+ ∂x((E11+ π11)u1) = 0, ∂tE22+ ∂x(E22u1+ π12u2) = 0, ∂tE12+ ∂x(E12u1+ 1 2(π11u2 + π12u1)) = 0, ∂tρπ11+ ∂x(ρπ11u1+ a2u1) = µρ(p11− π11), ∂tρπ12+ ∂x(ρπ12u1+ a2u2) = µρ(p12− π12), (3)
where the relaxation parameter a > 0 must satisfy the well-known Whitham sub-characte-ristic condition [19]:
a2
ρ > 3p11. (4)
Involving formal arguments, in the limit of µ to infinity, the system (3) aims to restore the equilibrium system (1). Indeed, as soon as µ is large enough, we formally obtain Π11= p11
and Π12 = p12. As a consequence, the equations of the momentum (ρu1 and ρu2) and
the conservation laws of the energies (E11, E22 and E12) in the equilibrium model (1) are
recovered from the momentum and energies equations of the relaxation model (3). With clear notations, we set
(
∂tW + ∂xFa(W) = µR(W),
W =t(ρ, ρu1, ρu2, E11, E22, E12, ρπ11, ρπ12),
where W belongs to V = {W ∈ R8; ρ > 0}.
The derivation of the relaxation model (3) is dictated by a linear degeneracy of all the fields which makes the Riemann problem easily solvable as stated in the following result: Lemma 1 Assume a > 0 and µ = 0. The first order extracted system from (3) is hyperbolic and it admits the eigenvalues λ1 = u1 − a/ρ, λ2 = u1 and λ3 = u1+ a/ρ. All
the fields are linearly degenerated. In addition, assume that the parameter a is such that (u1)L− a ρL < u?1 < (u1)R+ a ρR , u?1 = (u1)L+ (u1)R 2 + (π11)L− (π11)R 2a , (6)
where WL and WR are constant states in V to define the initial data of a Riemann
problem:
W0(x) =
WL if x < 0,
WR if x > 0. (7)
Then the solution of (3)-(7) with µ = 0 is given by
W(x, t) = WL if xt < λ1(WL), W1 if λ1(WL) < xt < λ2(W1), W2 if λ2(W2) < xt < λ3(WR), WR if λ3(WR) < xt, (8)
where λ2(W1) = λ2(W2) = u?1. Let us set for all 1≤ i ≤ j ≤ 2
u?i = (ui)L+ (ui)R 2 + (π1i)L− (π1i)R 2a , π1i? = (π1i)L+ (π1i)R 2 + a 2((ui)L− (ui)R), 1 ρ1 = 1 ρL +u ? 1− (u1)L a , 1 ρ2 = 1 ρR + (u1)R− u ? 1 a e1ij = (Eij)L ρL − (ui)L(uj)L 2 − 1 2a2 (π1i)L(π1j)L− π ? 1iπ?1j , e2ij = (Eij)R ρR − (ui)R(uj)R 2 − 1 2a2 (π1i)R(π1j)R− π ? 1iπ?1j , Eij1 = ρ1u ? iu?j 2 + ρ 1e1 ij, Eij2 = ρ2 u? iu?j 2 + ρ 2e2 ij.
Then the constant states W1 and W2 in V are defined as follows:
4 A relaxation scheme
Involving the above relaxation model, a numerical scheme is now derived when consid-ering a structured uniform mesh with size ∆x. At time tn, we note wh the approximate
solution of (1). This approximation is piecewise constant and it is defined as follows: wh(x, tn) = win=t(ρni, (ρu1)ni, (ρu2)ni, (E11)ni, (E22)ni, (E12)ni) ,
for all x ∈ (xi−1 2, xi+
1
2). To evolve in time this approximation, a two step method is
adopted.
Evolution: tn → tn+1,−
With 0 < t < ∆t, the Cauchy problem for the relaxation system (3) with µ = 0 is solved when considering the following initial data:
Wh(x, tn) = Wni
=t(ρni, (ρu1)ni, (ρu2)ni, (E11)ni, (E22)ni, (E12)ni, (ρπ11)ni, (ρπ12)ni) , x∈ (xi−1
2, xi+12),
where (π11)ni = (p11)ni and (π12)ni = (p12)ni.
Under the CFL like condition ∆t
∆xmax |λi(W
h(x, tn))|≤ 1
2, (9)
where (λi)1≤i≤3 denote the eigenvalues introduced in Lemma 1, the solution Wh at time
tn+ ∆t is made of the juxtaposition of the non interacting Riemann problem solutions
set at the cell interface xi+1
2. Let us emphasize that the relaxation parameter ai+12, set at
each cell interface xi+1
2, is evaluated to satisfy both Whitham condition (4) and condition
(6) with WL= Wni and WR= Wi+1n . Then we set
Win+1,− = 1 ∆x Z xi+ 1 2 xi− 1 2 Wh(x, tn+ ∆t)dx. Relaxation: tn+1,− → tn+1
The vector Wn+1,−i is projected on the equilibrium manifold {W ∈ V; Π11 = p11, Π12 =
p12}. Put in other words, this projection is equivalent to solve
∂tW = µR(W),
with Wn+1,i − as initial data, and to consider the solution in the limit of µ to infinity. After computations, at the end of the relaxation step the state vector is updated as follows:
with (π11)n+1i = (p11)n+1i and (π12)n+1i = (p12)n+1i .
The relaxation scheme can be summarized as follows: win+1= wni − ∆t ∆x fi+n 1 2 − f n i−1 2 , fi+n 1 2 = Fai+ 12 Wh(xi+1 2, t n+ ∆t)| [ρ,ρu1,ρu2,E11,E22,E12]. (10)
Let us note that the scheme (10) is nothing but a standard 3-point finite volume method. Indeed, the definition of Wh implies that the vector Wh(x
i+1 2, t
n+ ∆t) solely depends on
Wn
i and Wi+1n . Since Wni coincides with an equilibrium sate, it depends on wni.
The following result establishes the main stability properties satisfied by the considered scheme:
Theorem 2 Consider the relaxation scheme (10) and assume the CFL condition (9). Let us assume both sub-characteristic Whitham condition (4) and restriction (6).
i) Assume wn
i ∈ Ω for all i ∈ Z then win+1∈ Ω for all i ∈ Z.
ii) For all functions F and G such that (2) is satisfied, the following discrete entropy inequalities hold: ρn+1i F(ln sn+1 i )− ρniF(ln sni) + ∆t ∆x {ρF(ln s)u1}ni+1 2 − {ρF(ln s)u1} n i−1 2 ≤ 0, ρn+1i G(ln σn+1 i )− ρniG(ln σin) + ∆t ∆x {ρG(ln σ)u1}ni+1 2 − {ρG(ln σ)u1} n i−1 2 ≤ 0, with sni = (p11) n i (ρn i)3 , σni = (p11) n i(p22)ni − ((p12)ni) 2 (ρn i)4 . The numerical entropy flux functions are given by
{ρF(ln s)u1}ni+1 2 = (fρ) n i+1 2 × ( F(ln sn i) if (fρ)ni+1 2 > 0, F(ln sn i+1) otherwise, {ρG(ln σ)u1}ni+1 2 = (fρ) n i+1 2 × ( G(ln σn i) if (fρ)ni+1 2 > 0, G(ln σn i+1) otherwise, where (fρ)ni+1
2 denotes the first component of the numerical flux function involved in (10).
5 numerical experiments
A numerical experiment is performed with the relaxation scheme (10). We propose to consider a Riemann problems over the interval (−0.5, 0.5), the initial discontinuity being located at x = 0. We use an uniform mesh made of 500 cells. The CFL number fixed to 0.5 according to the CFL condition (9). The result is compared with the exact Riemann solution and the approximate solution performed with the HLL scheme and the Lax-Friedrichs scheme (see [18] and references therein). The initial data of this Riemann problem is given by
ρ u1 u2 p11 p12 p22
left state 1 0 0 2 0.05 0.6
right state 0.125 0 0 0.2 0.1 0.2
The numerical results for t = 0.125 are displayed in figure 1. The numerical solutions obtained with the HLL scheme and the relaxation scheme have the same accuracy while the Lax-Friedrichs scheme gives a very diffusive approximation. We note a better accuracy for the relaxation scheme when approximating the contact wave.
-0,4 -0,2 0 0,2 0,4 0,2 0,4 0,6 0,8 1 -0,4 -0,2 0 0,2 0,4 0,2 0,4 0,6 0,8 1 HLLE Lax Friedrichs Relaxation density -0,4 -0,2 0 0,2 0,4 0 0,2 0,4 0,6 0,8 -0,4 -0,2 0 0,2 0,4 0 0,2 0,4 0,6 0,8 x-velocity -0,4 -0,2 0 0,2 0,4 -0,1 0 0,1 0,2 0,3 0,4 -0,4 -0,2 0 0,2 0,4 -0,1 0 0,1 0,2 0,3 0,4 y-velocity -0,4 -0,2 0 0,2 0,4 0,5 1 1,5 2 -0,4 -0,2 0 0,2 0,4 0,5 1 1,5 2 (x,x)-pressure -0,4 -0,2 0 0,2 0,4 0 0,05 0,1 0,15 0,2 0,25 -0,4 -0,2 0 0,2 0,4 0 0,05 0,1 0,15 0,2 0,25 (x,y)-pressure -0,4 -0,2 0 0,2 0,4 0,2 0,3 0,4 0,5 0,6 -0,4 -0,2 0 0,2 0,4 0,2 0,3 0,4 0,5 0,6 (y,y)-pressure
REFERENCES
[1] D. Aregba-Driollet, R. Natalini, Convergence of relaxation schemes for conservation laws, Appl. Anal., 61, 163–190 (1996).
[2] M. Baudin, C. Berthon, F. Coquel, R. Masson, Q. H. Tran, A relaxation method for two-phase flow models with hydrodynamic closure law, Numer. Math., 99, no. 3, pp. 411–440 (2005).
[3] F. Bouchut, Entropy satisfying flux vector splittings and kinetic BGK models, Nu-mer. Math., 94, pp. 623–672 (2003).
[4] F. Bouchut, Nonlinear stability of finite volume methods for hyperbolic conserva-tion laws, and well-balanced schemes for sources, Frontiers in Mathematics series, Birkh¨auser (2004).
[5] C. Berthon, Numerical approximations of the 10-moment Gaussian closure, Math. of Comp., accepted.
[6] C. Berthon, In´egalit´es d’entropie pour un sch´ema de relaxation, C. R., Math., Acad. Sci. Paris 340, No.1, pp. 63-68 (2005).
[7] C. Chalons, Bilans d’entropie discrets dans l’approximation numrique des chocs non classiques. Application aux quations de Navier-Stokes multi-pression 2D et quelques systmes visco-capillaires, PhD Thesis, Ecole polytechnique (Paris, France) (2002). [8] G.Q. Chen, C.D. Levermore, T.P. Liu, Hyperbolic Conservation Laws with Stiff
Relaxation Terms and Entropy, Comm. Pure Appl. Math., 47 , pp. 787–830 (1995). [9] F. Coquel, B. Perthame, Relaxation of Energy and Approximate Riemann Solvers for
General Pressure Laws in Fluid Dynamics, SIAM J. Numer. Anal., 35, pp. 2223–2249 (1998).
[10] S. Jin, Z. Xin, The Relaxation Scheme for Systems of Conservation Laws in Arbitrary Space Dimension, Comm. Pure Appl. Math., 45, pp. 235–276 (1995).
[11] C.D. Levermore, Moment closure hierarchies for kinetic theory, J. Statist.. Phys., 83, pp. 1021–1065 (1996).
[12] C.D. Levermore and W.J. Morokoff The Gaussian moment closure for gas dynamics, SIAM J. Appl. Math., 59, pp. 72–96 (1998).
[14] R. Natalini, Convergence to equilibrium for the relaxation approximation of conser-vation laws, Comm. Pure Appl. Math., 49, pp. 1–30 (1996).
[15] I. Suliciu, On modelling phase transitions by means of rate-type constitutive equa-tions, shock wave structure, Internat. J. Engrg. Sci., 28, pp. 829–841 (1990).
[16] I. Suliciu, Some stability-instability problems in phase transitions modelled by piece-wise linear elastic or viscoelastic constitutive equations, Internat. J. Engrg. Sci., 30, pp. 483–494 (1992).
[17] E. Tadmor, A minimum entropy principle in the gas dynamics equations, Appl. Nu-mer. Math., 2, pp. 211–219 (1986).
[18] E. F. Toro, Riemann solvers and numerical methods for fluid dynamics. A practical introduction, Second edition. Springer-Verlag, Berlin (1999).
