Zero Mach Number Diphasic Equations for the Simulation of Water ...

Viewer
Transcript

Zero Mach Number Diphasic Equations for the Simulation of Water-Vapor High Pressure Flows Stéphane DELLACHERIE1

23

and Alain VINCENT2 3

1

Commissariat à l’Énergie Atomique 91191 Gif sur Yvette, France CERCA, Astrophysics Group, 5160 bv. Décarie bureau 400, Montréal QC, H3X 2H9 Physics Department, University of Montréal C.P. 6128, Succ. Centre-Ville, Montréal QC, H3C 3J7 2

3

Email: [email protected]; [email protected] Proceedings of the 11th Conference of the CFD Society of Canada, Vancouver, BC May 28-30, 2003

A BSTRACT A new set of equations to model two-phase water-vapor high pressure flows in cooling pipes heat exchangers of nuclear powerplants at low mach number (between 0.001 and 0.1) is proposed. These equations – namely the Zero Mach Number Diphasic equations (ZMND equations) – are similar to the Zero Mach Number Combustion equations [8, 9]. Simulations of temperature relaxation phenomena at low Mach number for two immiscible perfect gases are proposed to validate in monodimensional geometry this model without phase change.

1

I NTRODUCTION

Realistic modelling of water-vapor flows under high pressure and at low Mach number (between 0.001 and 0.1) in the core or in the cooling pipes heat exchanger of a nuclear powerplant is an important research topic. Numerical simulations help engineers to compute heat transfer, a necessary information for safety studies. Up to now, because of the complexity of such flows (Figs. 1 and 2), most of studies concerns large scales simulated with average models based on turbulence closure laws in which the interfaces between water and vapor are not accounted [1]. For this reason, laboratory experiments [2] as well as direct numerical simulations at scales small enough to detail both the structure of the water-vapor interface and the heat and mass transfer across the interface are necessary. Direct numerical simulations ([4], [5], [6]) have been made for incompressible and high compressible diphasic flows without phase change. However, for low

Mach number flows with phase change, it is essential to take into account compressibility and heat transfer phenomena although the Mach number is between 0.001 and 0.1. Thus, as in [7], the diphasic compressible Navier-Stokes equations are needed. But, the compressible equations include the acoustic waves and this puts a strong constraint on the stability CFL criterion when the convective terms are taken explicit in time. Moreover, the artificial viscosity of explicit or implicit compressible solvers is very high when the Mach number is low and this makes necessary to use preconditioning to reduce the numerical diffusion. To overcome these numerical difficulties and to simplify the numerical solvers, we propose to filter out the acoustic waves in the diphasic compressible NavierStokes equations by extending the Zero Mach Number Combustion equations proposed in [8, 9] (see also [11] for recent numerical simulations) – which models infinitely thin flames front structure constituted of perfect gases having the same equations of state and diffusion properties – to low Mach number water-vapor flows where the equations of state and the diffusion properties strongly depend of each phase. We call the resulting set of equations (see [12]) the Zero Mach Number Diphasic equations (ZMND equations). Indeed, from a qualitative point of view, we can understand the low speed combustion phenomena – where the transition from unburnt gas to burnt gas occurs in a thin layer – as a “phase change” at low Mach number. The modelling of the combustion transition layer cannot be the same than a realistic modelling of the water-vapor transition layer (see [7] for example). However, both transition phenomena are strong molecular or chemical reorganization of molecules spatialy localized in a thin layer. As a first step, to study the potential of the ZMND equations for the modelling of low Mach number

water-vapor flows, we have considered simple heat exchanges without phase change between two immiscible fluids at low Mach number. We show that pressure variations and heat exchanges between the two immiscible fluids are directly function of the discontinuity of the equations of state near the two fluids interface. Our model does not include at the present time any phase change neither for an infinitely thin transition layer (phase change model which could be deduced from [8, 9]), nor for a finite thin transition layer ([7]). Let us note that the Low Mach number equations were firstly proposed in [3] to filter acoustic waves in monophasic natural convection problems and for perfect gases. The plan of this paper is the following: in the second part, we present our Zero Mach Number Diphasic equations and its basic properties in a bounded domain. In the third part, we write these equations in monodimensional lagrangian coordinates for two immiscible perfect gases. In the fourth part, we detail the numerical scheme. In the fifth part, we show numerical simulations and, at last, we discuss the work underway.

w a te r

h o m o g e n e o u s w a te r-v a p o r m ix tu re

g

s e e th e z o o m o n th e fig . 2

v a p o r

g

σ Y1

q Y1

Fig. 2

T HE Z ERO M ACH N UMBER D IPHASIC E QUATIONS

Dt Y1 0

(a)

∂t ρ ∇ ρu 0

(b)

ρDt u ∇P ∇ σ Y1

ρg

(c)

ρDt h Dt P σ Y1 : ∇u ∇ q Y1

(d)

1 if x belongs to fluid 1 (3) 0 if x belongs to fluid 2

λ Y1 Y1 λ1 1 Y1 λ2 µ Y1 Y1 µ1 1 Y1 µ2 ρ Y1 T P Y1 ρ1 T P 1 Y1 ρ2 T P ε Y1 T P Y1 ε1 T P

1 Y1 ε2 T P

#"

x Σ t : T $ Σ1 % t & T $ Σ2 % t & λ1 ∇T $ Σ1 % t & n1' 2 λ2 ∇T $ Σ2 % t & n1 ' x ∂Ω :

The diphasic compressible Navier-Stokes equations solved in a bounded domain Ω are

∇ u I

(4) knowing that Y1 t x ! 0 1 (µk and λk are supposed to be positive constants in this paper). Functions ρk T P and εk T P are the equations of state of the fluid k, T being the thermodynamic temperature. The boundary conditions are defined with

Two-phase flow in a vertical pipe

2

2 3

Equation (1)(a) with the initial condition (3) imposes that for any t x , Y1 t x 0 1 . We also define Y2 1 Y1. The function Y1 t x can be regarded as a color function or as the volumetric fraction of the fluid 1 whose discontinuity surface Σ t localizes at any time the interface between fluid 1 and fluid 2 (notice that Σ t can be non connex as on figure 2). If the fluids 1 and 2 are respectively vapor and water, Y1 can be seen as a void fraction. The diphasic character of the flow is taken into account through

d ire c tio n o f th e tw o -p h a s e flo w

Fig. 1

Y1 t 0 x

p ip e b o u n d a ry

µ Y1 ∇u ∇u t λ Y1 ∇T

(a) (b) (2) and where Dt ∂t u ∇ is the lagrangian derivative. The constant vector g 9 81 zˆ m.s 2 is the gravity (ˆz is the unit vertical direction). The physical quantities ρ, h ε P ρ (ε is the internal energy), σ and q are respectively density, enthalpy, newtonian viscosity tensor and Fourier heat flux (µ and λ are the viscosity and conductivity: they are defined below). The initial condition for Y1 t x is given by

w a te r-v a p o r in te rfa c e

p ip e b o u n d a ry

d ire c tio n o f th e tw o -p h a s e flo w

where

(1)

2

∇T n 0 (5)

and" with x ∂Ω :

u 0

(6)

Vector n1 ' 2 is the" unit normal vector on the interface Σ t oriented towards the fluid 2; Σ1 t and Σ2 t are respectively the fluid 1 side and the fluid 2 side of Σ t . Vector n is normal to the surface ∂Ω. The first and the second conditions in (5) impose the continuity of the temperature and of the heat flux in Ω, despite the existence of the infinitely thin interface Σ t .

The Zero Mach Number Diphasic equations: When the Mach number is close to zero, it is possible to filter out the acoustic waves in (1) by performing an asymptotic expansion as in [8, 9] (see [12]) in terms of the Mach number M u c (u and c being representative of the value of the flow velocity and of the sound velocity in Ω). We obtain that Y1 T u P is solution of the two coupled systems

Dt Y1 0

(a)

(7)

ρC p Dt T αT P t ∇ λ Y1 ∇T

(b)

(P t refers to the time derivative of P t ) and

∇ u G t x

(a)

ρDt u

(b)

∇Π ∇ σ Y1

ρg

(8)

is not a perfect gas (this term is not taken into account in [3], [9] and [11]), which is the case for water-vapor mixtures. The two coupled systems (7)-(8) closed with (4) and (11), the initial condition (3) and the boundary conditions (5)-(6) define the Zero Mach Number Diphasic equations (ZMND equations) without phase change in a bounded domain Ω. Nature of the equations: The system (7) is a mixed hyperbolic-parabolic equation with a source term. Moreover, by applying the Hodge decomposition on the velocity u, we can say that there is an unique v ∇φ such that u v ∇φ and ∇ v 0. Then, we can deduce from (8) and (2)(a) that there is an other dynamic pressure Π t x such that v ∇φ is solution of the system

∂T ρ % Y1 T P & ρ % Y1 T P &

α Y1 T P and C p Y1 T P ∂T h Y1 T P are respectively the compressibility coefficient and the calorific capacity at constant pressure. c Y1 T P is the sound speed of the fluids 1 or 2 and

∆φ G t x ∇ v 0 ρDt v ∇Π µ Y1 ∆v ρg ρDt ∇φ

with the boundary conditions x ∂Ω :

G t x

P t

β∇ λ Y1 ∇T ρc2

(9)

with α Y1 T P

ρ Y1 T P C p Y1 T P

β Y1 T P

(10)

The dynamic pressure Π t x in (8)(b) is a new unknow and is the first spatial pressure perturbation up to a constant. Equation (8)(a) is the mass conservation equation (1)(b). The thermodynamic pressure P t is now a function of time, its space dependance being close to zero at low Mach number ([3], [8], [9], [12]). In the case of a diphasic mixture with the most general equations of state and in the case of a bounded domain Ω with conditions (5), P t is solution of the nonlinear integro-differential equation (see [12]) P t

Σ% t&

β Y1 T P Σ % t & λ Y1 ∇T n1 '

Ω ∂T β

2 dΣ

dx Ω ρ % Y1 T P & c % Y1 T P & 2

Y1 T P λ Y1 ∇T 2 dx

dx Ω ρ % Y1 T P & c % Y1 T P & 2

(11) with β Σ β $ Σ1 β $ Σ2 . The first term in the right hand side of (11) is directly due to the discontinuity of the equations of state across the interface Σ t and, thus, to the diphasic character of the flow; the second term is not equal to zero when, at least, one of the two fluids

(a) (b) (c) (12)

∇φ n 0 v n 0 v n ∇φ n 0

"

(a) (b) (c)

(13)

Equation (12)(a) with the boundary condition (13)(a) is an elliptic equation and admits, by the Neumann compatibility condition, a solution with ∇φ uniquely determined if and only if

Ω

G t x dx 0

(14)

Equation (11) was obtained such that (14) is verified. The equations (12)(b)(c) are non homogeneous incompressible Navier-Stokes equations. At last, let us note that the Zero Mach Number Diphasic equations verify dt Ω ρεdx 0 (indeed, the energy u2 2 is neglected in this model) and dt Ωk % t & ρdx 0 where Ωk t

x Ω:

Yk t x 1

On the thermodynamic fluid balance conditions: A rigorous proof that the Zero Mach Number Diphasic equations is a well-posed system could be studied by extending the technique of P. Embid written in [10] for the study of the Zero Mach Number Combustion equations. In [10], it is shown that the initial velocity has to verify the chemical fluid balance condition introduced by A. Majda in [8, 9]. In the present case, this condition imposes that, for given initial temperature and thermodynamic pressure, the initial velocity

field u t 0 x has to verify (8)(a) with u $ ∂Ω 0 and that the initial Mach number has to be close to zero. Moreover, we have observed that the initial temperature has to satisfy the boundary conditions (5) on the interface Σ t 0 to avoid a non physical high acceleration of the interface. In reference to A. Majda, we name these conditions the thermodynamic fluid balance conditions.

3

Let us note that, in lagrangian geometry, the interfaces Σa and Σb and, thus, the domains Ωk are constants because of (15)(a). We can now write that η Y1 η m , C p Y1 C p m , γ Y1 γ m and R Y1 R m . The boundary conditions are now given by

∂t Y1 0

P t

∂m η Y1 ∂m T ∂t C p Y1 T R Y1 T P t

∂m u G t m

(a)

∂t u

(b)

(16)

∂m Π

∂m T $

γ 1 R1 γ 2 R2 1 Y1

γ1 1 γ2 1 γ Y1 Y1 γ1 1 Y1 γ2 R Y1 Y1 R1 1 Y1 R2

C p Y1 Y1

L

(b) (c)

(18)

Ωa2 Ω1 Ωb2 Ω2

a L Σ a b Σ Σ Σb L

Ωa2 Ωb2

(20)

(21)

0

and by u $ L 0 (Σk Σak Σbk is the fluid k side of Σ: here, the fluid 1 is on the right side of Σa and on the left side of Σb ). By noting that ∂T β 0 for a perfect gas, (11) allows to write that the pressure P t is solution of the non linear ordinary differential equation

P t Γ12 P t η∂m T $ Σb η∂m T $ Σa L 1 R m T t m

dm γ m

L

(22)

Γ12

γ1 1 γ2 1 γ1 γ2

(23)

Of course, due to (22), the function G t m given by (17) verifies the Neumann compatibility condition (14) which means that

L

G t m dm 0

(24)

In monodimensional geometry, the Hodge decomposition u v ∂m φ is trivial since ∂m v 0 and v $ L 0 imply that v 0: thus, u ∂m φ and φ Π is solution of

∂2m φ G t m

2 ∂m Π

(a)

with γk 1 and Rk 0 (we have ρk T P P Rk T

for a perfect gas). Here, we suppose that Ω L L with L 0 and that the initial condition (3) is such that Ω1 is connex and Ω2 is non connex, Σa and Σb being respectively the left interface and the right interface. And we note Ω Ωa2 Ω1 Ωb2 (19) with

L

with

R Y1 T P t

γ Y1 1 G t m ∂m η Y1 ∂m T 2 γ Y1 P t

γ Y1 P t

(17) The positive functions C p Y1 , γ Y1 and R Y1 are given by

η∂m T $ Σb

(b)

The function G t m is defined by

T $ Σb T $ Σb 1 2 η∂m T $ Σb η∂m T $ Σb 1 2

(15)

with η ρλ, and by

(a)

η∂m T $ Σa

EQUATIONS IN LAGRANGIAN COORDINATES FOR TWO IMMISCIBLE PERFECT GASES

T $ Σa1 T $ Σa2 η∂m T $ Σa1 η∂m T $ Σa2

T HE ZMND

The equations in lagrangian coordinates: Suppose that the geometry is monodimensional and let us define the lagrangian variable m such that dm ρdx. Then, by supposing that the two fluids are immiscible perfect gases, the ZMND equations without viscosity and gravity, and written in lagrangian coordinates, are given by

with ∂m φ $

∂t G t m

0

L

with ∂m Π $

(a) 0

(b) (25) We deduce in this very simple case that, because v 0, the velocity u t m and the dynamic pressure Π t m

are solved independently. Let us also remark that in lagrangian geometry, the system (15) is not coupled to the system (16) and this considerably simplifies the discretization of the equations. Nevertheless, the resolution of (16) and thus the velocity u t m and the dynamic pressure Π t m depend on the solutions T t m of (15) and P t of (22) by the mean of function G t m . L

Simplification of the equation (15)(b): Let us now define the variable ψ t m T t m P t

γm 1 γm

(26)

It is easy to show ([12]) that ψ t m is solution of the parabolic equation

∂t C p m ψ ∂m η m ∂m ψ

(27)

with the boundary condition

γ1 1 γ1

ψ $ Σa1 P t

P t

P t

γ1 1 γ1

P t

γ1 1 γ1

P t

γ1 1 γ1

η∂m ψ $ Σa 1

ψ $ Σb P t

1

η∂m ψ $ Σb 1

∂m ψ $

ψ $ Σa2

γ2 1 γ2

P t

γ2 1 γ2

η∂m ψ $ Σa 2

(28)

ψ $ Σb

2 γ2 1 γ2

P t

η∂m ψ $ Σb 2

0

L

γ2 1 γ2

And, the pressure P t is now solution of

P t Γ12 P t

2γ1 1 γ1

η∂m ψ $ Σb

1

R m ψ t m P t

γ m

L

1 γm 1 γm

L

η∂m ψ $ Σa

dm

1

Remark on the final equilibrium: We can easily show that LL γRm% m& & T1 dm E and LL RP% mt& & T dm L are % % positive constants fixed by the initial conditions. This proves that sup P t

0 when ∞ and that inf P t

P t 0

0. And by supposing that P t converges to a pressure P∞ and that ψ t m converges to ψ1 in Ω1 and to P∞Γ12 ψ1 in Ω2 when t goes to ∞ (which remains to be proved), we obtain that the temperature T t m converges also to a positive constant T∞ in Ω. And, the equilibrium T∞ P∞ is given by

T∞

L

P∞

E L

4

E

R % m& dm L γ % m& 1

L

(29) Thus, the parabolic equation (27) and the nonlinear ordinary differential equation (29) are coupled through the time dependent boundary conditions (28) on the fluid 1 / fluid 2 interface Σ. Let us remark that if γ1 γ2 (in that case, Γ12 0), we have P t P t 0 and (27)-(28) become a classical parabolic equation.

we deduce the temperature with (26), we compute G t m with (17) and we obtain the velocity with (16)(a). The scheme that we use to solve (27) is the simple first order in time explicit scheme for the monodimensional heat equation and thus imposes a classical stability criterion on the time step ∆t. To verify a discretized version of the Neumann compatibility condition (24), it is easy to see that the resolution of (29) has also to be a first order in time explicit scheme. Thus, we obtain a second stability criterion which makes positive the pressure P t at each time step. An important point of the previous scheme is the way we take into account the boundary conditions (28) at the discretized level on the interfaces Σa and Σb . Let us define the subscripts I b and I b 1 of the mesh points mI b and mI b 1 which are respectively on the left side and on the right side of the interface Σb whose mesh point is noted mI b 1 2 . Let us now define the value of ψ and the discretized gradients of ψ on the left side and on the right side of Σb with

(30)

L R m dm L R % m& dm L & γ m 1 %

T HE

NUMERICAL SCHEME IN LAGRANGIAN COORDINATES

To obtain the discretized solution of (15)-(16)-(17)(21)-(22), we firstly solve (27)-(28)-(29). Secondly,

ψ $ Σb 1

η1 ψI b η2 P t Γ12 ψI b η1 η2

ψ $ Σb 2

(31)

∆m 2

and with

1

ψ $ Σb ψI b 1

∂m ψ $ Σb 1

Γ12 ψ

η1 P t

∂m ψ $ Σb 2

ψI b

I b η2 ψ I b η1 η2

1

ψ $ Σb

∆m 2

2

1

(32)

Here, we suppose that the mesh size ∆m is constant in L L (we recall that the fluid 1 is on the left side of Σb and that Γ12 is defined by (23)). Such discretization for the gradients is needed to define the numerical scheme near the interface Σb and, of course, near the interface Σa by using similar definitions. Good properties of (31) and (32) are that they verify (28) and that the scheme preserves now the equilibrium when it is reached. Numerical results show that (31) and (32) give good numerical results without any new stability criterion, even when Γ12 0, η1 η2 and when the number of mesh points is very low. At last, let us recall that the initial conditions u t 0 m , T t 0 m and P t 0 have to verify the thermodynamic fluid balance conditions at the discretized level. In the case of perfect gases and when ηk Ckste , our numerical results show that for a given initial temperature, the higher the initial pressure, the lower the initial Mach number. At last, we compute Gi in each

mesh point mi of L L and we solve (16)(a) with ui

1 2

ui

1 2

∆mGi

Second test case compression of the “bubble” Ω1 : The initial conditions are now

knowing that u1 2 u $ L 0 and that the velocity is defined at the interfaces mi 1 2 of the mesh.

5

N UMERICAL

m Ωa2 : T t 0 m 1 5 m Ω1 : T t 0 m 2 5 m Ωb2 : T t 0 m 2

Again, we compute the equilibrium with (30) and we find T∞ 2 28 P∞ 237 84

RESULTS

The two perfect gases are defined by the physical constants γ1 1 2 R1 0 2 η1 1 and

γ2 3 R2 0 7 η2 2

(we recall that η ρλ). The interfaces Σa and Σb (see (19) and (20)) are such that the mass in Ωa2 is equal to the mass in Ωb2 which simply means that Σa L L Σb since the coordinates are lagrangian. Here, we take L 1 Σa 1 2 Σb 1 2 At last, the initial thermodynamic pressure is equal to P t 0 200 and the number of mesh points is equal to 100 with ∆m 1 50. First test case dilatation of the “bubble” Ω1 : The initial conditions are

m Ωa2 : T t 0 m 1 5 m Ω1 : T t 0 m 1 m Ωb2 : T t 0 m 2

We regularize the initial temperature T t 0 m on some mesh points (5 in the present test cases) and we compute the initial velocity u t 0 m such that the thermodynamic fluid balance conditions is fulfilled. Nevertheless, it remains to verify that the initial Mach number is close to zero which is the case when P t 0

1 (it is inferior to 0.05 in the present case; notice that when P t 0 1, the initial Mach number is superior to 7 near the interface Σb ). We then deduce that T∞ 1 21 P∞ 153 90

by using (30).

after the regularization of the initial temperature. We also verify that the initial Mach number is close to zero (it is again inferior to 0.05). Simulations: To illustrate the dilatation or the compression of the “bubble” Ω1 , we have to show the results in the eulerian space t x since the interfaces are fixed in the lagrangian space t m . Thus, the numerical results are projected in the eulerian mesh which is computed at each time step: since we know the velocity ui 1 2 tn at the time tn and at the interface mi 1 2 of the lagrangian mesh, it is possible to compute the interface xi 1 2 tn 1 of the eulerian mesh with

x t x t

∆t u t Let us note that x x $ C and that x x $ C since u $ 0: here, we take x 0, x x being imposed by the physical constants and by the initial conditions. Moreover, to simi 1 2 n 1

L

L

1 2

ste

imax 1 2

i 1 2 n

i 1 2 n

ste

L

imax 1 2 1 2

1 2

plify the visualization of the results, we normalize the eulerian domain 0 ximax 1 2 to 0 1 . The first and second test cases show respectively a dilatation and a compression of the “bubble” Ω1 (Figs. 3 and 4) which is advected towards the middle of Ω: indeed, in the first test case, the temperature in Ω1 increases and, in the second test case, it decreases. Figures [5-10] show that the temperature T t x and the thermodynamic pressure P t go to the equilibrium T∞ and P∞ , and that the velocity u t x goes to zero when t goes to ∞.

6

C ONCLUSION

We have proposed in this paper a new set of equations to model water-vapor high pressure flows. These equations – namely the Zero Mach Number Diphasic equations (ZMND equations) – are similar to the Zero Mach Number Combustion equations proposed in [8, 9] to model infinitely thin flames front structure. Numerical results show that these equations can simulate dilatation and compression of a monodimensional bubble at low Mach number. A next step would be to propose a numerical scheme for monodimensional eulerian geometry (see [12]) using a level set technique to

capture the interface ([4]). A further step would be to extend the previous monodimensional eulerian scheme to bidimensional and tridimensional simulations. The introduction of a phase change model derived from [9] for an infinitely thin interface or from [7] for a finite thin interface will have to be carefully studied. In this last case, the introduction of adaptative mesh refinement techniques to optimize the resolution of the finite thin interface in bidimensional and tridimensional simulations should be also considered as it is done in [11] for combustion problems. At last, another interesting question concerns the wellposedness of our ZMND equations: the ideas exposed in [10] should be a good start point to give answers to that theoretical question. Acknowledgments: This research is supported by Commissariat à l’Énergie Atomique CEA, France. Computations are done with the machines of the Réseau Québécois de Calcul Haute Performance RQCHP.

2.1

1.9

1.7

1.5

1.3

1.1

0.9 0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1.0

Fig. 5: temperature profile T x (test case 1) at time iterations 0, 450 and 10000

2.6

1100

2.4 1000 900

2.2

800

2.0

700 600

1.8

500 400

1.6 300 200

1.4 0

100 0.0 0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1.0

1.0

.

Fig. 6: temperature profile T x (test case 2) at time iterations 0, 450 and 10000

Fig. 3: density profile ρ x (test case 1) at time iterations 0, 450 and 10000

25e−3 600

21e−3 17e−3 500

13e−3 9e−3

400

5e−3 1e−3

300

−3e−3 200

−7e−3 −11e−3

100 0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

Fig. 4: density profile ρ x (test case 2) at time iterations 0, 450 and 10000

1.0

−15e−3 0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

Fig. 7: velocity profile u x (test case 1) at time iterations 0, 20, 450 and 10000

1.0

R EFERENCES [1] Clerc S. and Lebaigue O. – Simulation numérique locale des écoulements diphasiques liquide-vapeur - Scientific report of the Nuclear Reactor Division, Commissariat à l’Énergie Atomique, France, p. 135-143, 1998.

32e−3 28e−3

[2] Le Corre J.M., Hervieu E., Ishii M. and Delhaye J.M. – Benchmarking and improvements of measurements techniques for local time-averaged two-phase flow parameters – Fourth International Conference on Multiphase Flows (ICMF 2001), New-Orleans, USA, 2001.

24e−3 20e−3 16e−3 12e−3 8e−3

[3] Paolucci S. – On the filtering of sound from the NavierStokes equations – Sandia National Laboratories report (Livermore), SAND82-8253, 1982.

4e−3 0 −4e−3 −8e−3 0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1.0

[4] Osher S., Sussman M. and Smereka P. – A level set approach for computing solutions to incompressible twophase flow – J. Comput. Physics, 114, p. 146-159, 1994.

Fig. 8: velocity profile u x (test case 2) at time iterations 0, 20, 450 and 10000

[5] Lagoutière F. – Modélisation mathématique et résolution numérique de problèmes de fluides compressibles à plusieurs constituants – Ph.D. Thesis of Université Pierre et Marie Curie (Paris VI), 2000.

200

[6] Allaire G., Clerc S. and Kokh S. – A five-equation model for the simulation of interfaces between compressible fluids – J. Comput. Physics, 181, p. 577-616, 2002.

190

180

[7] Jamet D., Lebaigue O., Coutris N. and Delhaye J.M. – The second gradient method for the direct numerical simulation of liquid-vapor flows with phase change – J. Comput. Physics, 169, p. 624-651, 2001.

170

160

150 0

1e3

2e3

3e3

4e3

5e3

6e3

7e3

8e3

9e3

10e3

Fig. 9: pressure profile P t (test case 1) at time iterations 1 10000

[8] Majda A. – Equations for low mach number combustion – Center of Pure and Applied Mathematics, University of California at Berkeley, report no. 112, 1982. [9] Majda A. and Sethian J. – The derivation and numerical solution of the equations for zero Mach number combustion – Combust. Sci. and Tech., 42, p. 185-205, 1985. [10] Embid P. – Well-posedness of the nonlinear equations for zero Mach number combustion – Comm. Partial Differential Equations, 12(11), p. 1227-1283, 1987.

240 236 232 228

[11] Core X. – Méthode adaptative de raffinement local multi-niveaux pour le calcul d’écoulements réactifs à faible nombre de Mach – Ph.D. Thesis of Université de Provence (France), 2002.

224 220 216 212 208 204 200 0

1e3

2e3

3e3

4e3

5e3

6e3

7e3

8e3

9e3

10e3

Fig. 10: pressure profile P t (test case 2) at time iterations 1 10000

[12] Dellacherie S. – Dérivation des équations diphasiques à bas nombre de Mach. Simulation numérique en géométrie monodimensionnelle – CERCA technical report. In preparation.

Simulation of Water Pollution by Finite Difference Method