scheme. Namely, the time integration has used the three stage TVD Runge–Kutta discretization described in[11],
Qð1Þi ¼ Qni þ DtLðQniÞ ðaÞ
Qð2Þi ¼34Qni þ14Qð1Þi þ14DtLðQð1Þi Þ ðbÞ Qnþ1i ¼13Qni þ23Qð2Þi þ23DtLðQð2Þi Þ ðcÞ 8>
><
>>
:
ð57Þ
while the space discretization has used the MUSCL tech- nique described in Section 3.3, with f¼ 1=3 and the Spe- kreijse limiter. The results of this computation are shown inFig. 4. Although the results of this fine mesh computa-
tion are not totally identical to those ofFig. 3, one can note the close similarity between these results and the one ob- tained with the preconditioned upwind scheme.
4.2. Broken dam problem
Here we present a computation of the well-known bro- ken dam problem of Martin and Moyce [7]. Initially a water column with a¼ 0:06 m wide and g2a¼ 0:12 m high is a rest. Under the effect of the gravity g¼ 9:81 m s2, the column collapses. All the boundaries are solid walls. The
Fig. 2. Bubble ascension: isovalues of the volume fraction for the 100 100 mesh computation: classical upwind scheme.
Fig. 3. Bubble ascension: isovalues of the volume fraction for the 100 100 mesh computation: preconditioned scheme with b ¼ 0:1.
A. Murrone, H. Guillard / Computers & Fluids 37 (2008) 1209–1224 1217
mesh we have used for this test, is regular with Dx¼ Dz¼ 5 103m. The Mach number during the computa- tion is low and of the order of 1 101. The implicit scheme has been used with a CFL number equal to 2.5 in order to compute with a sufficient accuracy the unsteady pattern of the flow. The linear system is solved by an iter- ative method with a linear residual h¼ 1 102. The space discretization has used the MUSCL technique with f¼ 1=2
and the Van Albada–Van Leer limiter. We compute the solution with the standard and the preconditioned method.
For the preconditioned method, the parameter of the Tur- kel’s matrix b is chosen equal to 0.1 and remains constant in space and time.
Figs. 5 and 6 show the isovalues of the volume frac- tion at the different dimensionless times t ffiffiffiffiffiffiffiffiffiffi
p2g=a
¼ 0; 1:19; 1:98; 2:97; 4:02; 5:09 corresponding to the physical
Fig. 4. Bubble ascension: isovalues of the volume fraction for the 400 400 fine mesh explicit computation.
Fig. 5. Broken dam problem: isovalues of the volume fraction: classical upwind scheme.
Fig. 6. Broken dam problem: isovalues of the volume fraction: precon- ditioned scheme with b¼ 0:1.
1218 A. Murrone, H. Guillard / Computers & Fluids 37 (2008) 1209–1224
times t = 0, 0.066 s, 0.109 s, 0.164 s, 0.222 s, 0.281 s, for the standard (Fig. 5) and the preconditioned method (Fig. 6).
It is clear that the upwind preconditioned scheme predicts a faster development of the flow and for instance the front position at time 0.281 s is clearly in advance with respect to the results obtained with the standard scheme.
To quantify the difference between the two schemes, we compare inFig. 7, the two solutions with the experimental results of[7], for the front position x=a¼ F1 g2; t ffiffiffiffiffiffiffiffiffiffi
p2g=a
and the height of the column z=ðg2aÞ ¼ F2 g2; t ffiffiffiffiffiffiffiffi pg=a
. It is clear that the preconditioned method is more accurate than the standard one. For example at time t ffiffiffiffiffiffiffiffiffiffi
p2g=a
¼ 2:97 the error compared to the experimental data for the front position is the order of 1% for the preconditioned method while it is the order of 10% for the classical upwind scheme.
4.3. Two-phase flow in a nozzle
Finally we present a sequence of computations of two- phase flows in a symmetric nozzle where the Mach number
tends to zero. These computations are similar to the ones presented in[4]in single phase situation. The implicit scheme has been used with a CFL number equal to the inverse of the nonlinear residual of the mixture density CFL¼ 1=ResðqÞ and at maximum equal to 106. The discrete solutions pre- sented are the one obtained at convergence, i.e., for a resid- ual equal to 109. For the preconditioned method in this case, the parameter b is not taken as a constant but locally computed at each interface of the mesh and at each time step.
In these computations, an air–water two phase mixture defined by ða11 ¼ 0:5; q11 ¼ 1 kg m3;q12 ¼ 1000 kg m3Þ is injected in the nozzle with an horizontal imposed velocity equal to u1¼ 1 m s1. The law states(56)for air and water are the same than previously and a representative Mach number of the flow is defined by
Ma2¼ðu1Þ2 ð^a1Þ2¼ 1
ð^a1Þ2 ð58Þ
Then the inlet pressure is taken as solution of equation Ma2 1=ð^a1Þ2¼ 0 which can be rewritten using the equa- tions of state of the two fluids:
Ma2 ða11 q11 þ ð1 a11Þq12 Þ
a11
c1ðp1þ p1Þþ 1 a11 c2ðp1þ p2Þ
¼ 0 ð59Þ
We are interested in the situation where Matends to zero.
In this case, the asymptotic analysis given in Section2ap- plies and shows that the equations governing the flow are the two-phase incompressible Euler equations (14). Thus if we take at time t¼ 0; a1ðx; 0Þ ¼ 0:5 for all x, we will get a1ðx; tÞ ¼ 0:5 for all x and t > 0 and (1) is simply the incompressible Euler equation with a constant density gi- ven byðq1þ q2Þ=2. We then expect that the limit solution of (1)will be given by an incompressible potential flow of density ðq1þ q2Þ=2. In particular, the solution has to be symmetric with respect to the axis of the nozzle. To test the preconditioned scheme with respect to the classic non-preconditioned one, we realize three computations, respectively, at Ma¼ 0:1; Ma¼ 0:01 and Ma¼ 0:001.
Fig. 8 shows the isovalues of the normalized pressure p pmin=pmax pmin for the discrete stationary solutions obtained.Fig. 9shows the profile of pressure in the upper and lower boundaries of the nozzle. We present from left to right the results obtained with the classical and the precon- ditioned scheme.
Fig. 9shows clearly that the solution given by the clas- sical discretization is not symmetric and consequently could not be a reasonable approximation of the incom- pressible solution. In addition, one can notice that the pres- sure fluctuations with the classic scheme (of order Ma) are larger than which obtained with the preconditioned one (of order Ma2). To illustrate this difference in the behaviour of the pressure fluctuations, we plot in the first and second column of Fig. 9 the result for the classic and precondi- tioned method with the same pressure unit. We could note
0 1 2 3 4 5 6 7 8 9 10
0 1 2 3 4 5 6
Dimensionless front position x/a
Dimensionless time t*sqrt(2g/a)
Experimental data Standard Method Preconditioned Method B=0.1
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0 0.5 1 1.5 2 2.5 3 3.5
Dimensionless height z/(2a)
Dimensionless time t*sqrt(g/a)
Experimental data Standard Method Preconditioned Method B=0.1
Fig. 7. Comparison between numerical solutions of the classical and preconditioned scheme and experimental results for the broken dam problem. Front position (top) and height of the column (bottom).
A. Murrone, H. Guillard / Computers & Fluids 37 (2008) 1209–1224 1219
