Group Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Title: 5.4.3 - A low-dissipation method for DNS of compressible turbulent multicomponent and multiphase flows with shocks
ALL VOLUMES CITATION THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00102023/00130
 Material Information
Title: 5.4.3 - A low-dissipation method for DNS of compressible turbulent multicomponent and multiphase flows with shocks Computational Techniques for Multiphase Flows
Series Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Physical Description: Conference Papers
Creator: Johnsen, E.
Larsson, J.
Publisher: International Conference on Multiphase Flow (ICMF)
Publication Date: June 4, 2010
 Subjects
Subject: DNS
multifluid algorithm
Richtmyer-Meshkov instability
 Notes
Abstract: High-speed turbulentmulticomponent andmultiphase flows with shocks are of importance in a number of applications, including supersonic and hypersonic propulsion, supernovae explosion and inertial confinement fusion. Numerical simulations can provide insights into these phenomena, but simulations of turbulent flows with shocks and interfaces are challenging due to the contradictory requirements to treat the different phenomena: while numerical dissipation is added to stabilize the solution near discontinuities, this dissipation tends to overwhelm the physical dissipation. To overcome this drawback, a hybrid WENO/central difference method is employed to solve the compressible multicomponent Navier-Stokes equations with mass and heat transfer. A non-dissipative central difference scheme is used in smooth regions while a weighted essentially non-oscillatory (WENO) scheme is employed near shocks; a shock sensor specifies regions in which each scheme is applied. In the present multifluid extension, an additional sensor based on the WENO weights of the mass fraction detects regions in which the fluid composition is discontinuous. In these regions, a fully conservative interface-capturing WENO scheme is employed, in which spurious pressure oscillations and temperature errors at interfaces are prevented. The resulting high-order accurate hybrid method can therefore simulate turbulence, shock waves and interfaces in a stable and accurate fashion and is used to study turbulent multi-material mixing in shock-accelerated flows. In particular, the set up of the experiments by Balakumar et al. (2008) is considered, in which a shock wave travelling in air interacts with a layer of SF6. Two- and three-dimensional calculations are discussed.
General Note: The International Conference on Multiphase Flow (ICMF) first was held in Tsukuba, Japan in 1991 and the second ICMF took place in Kyoto, Japan in 1995. During this conference, it was decided to establish an International Governing Board which oversees the major aspects of the conference and makes decisions about future conference locations. Due to the great importance of the field, it was furthermore decided to hold the conference every three years successively in Asia including Australia, Europe including Africa, Russia and the Near East and America. Hence, ICMF 1998 was held in Lyon, France, ICMF 2001 in New Orleans, USA, ICMF 2004 in Yokohama, Japan, and ICMF 2007 in Leipzig, Germany. ICMF-2010 is devoted to all aspects of Multiphase Flow. Researchers from all over the world gathered in order to introduce their recent advances in the field and thereby promote the exchange of new ideas, results and techniques. The conference is a key event in Multiphase Flow and supports the advancement of science in this very important field. The major research topics relevant for the conference are as follows: Bio-Fluid Dynamics; Boiling; Bubbly Flows; Cavitation; Colloidal and Suspension Dynamics; Collision, Agglomeration and Breakup; Computational Techniques for Multiphase Flows; Droplet Flows; Environmental and Geophysical Flows; Experimental Methods for Multiphase Flows; Fluidized and Circulating Fluidized Beds; Fluid Structure Interactions; Granular Media; Industrial Applications; Instabilities; Interfacial Flows; Micro and Nano-Scale Multiphase Flows; Microgravity in Two-Phase Flow; Multiphase Flows with Heat and Mass Transfer; Non-Newtonian Multiphase Flows; Particle-Laden Flows; Particle, Bubble and Drop Dynamics; Reactive Multiphase Flows
 Record Information
Bibliographic ID: UF00102023
Volume ID: VID00130
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: 543-Johnsen-ICMF2010.pdf

Full Text



7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


A low-dissipation method for DNS of compressible turbulent
multicomponent and multiphase flows with shocks


E. Johnsen* and J. Larssont

Department of Mechanical Engineering, University of Michigan, Ann Arbor, MI 48109, USA
t Center for Turbulence Research, Stanford University, Stanford, CA 93405, USA
ejohnsen @umich.edu and jola @stanford.edu
Keywords: DNS, multifluid algorithm, Richtmyer-Meshkov instability




Abstract

High-speed turbulent multicomponent and multiphase flows with shocks are of importance in a number of applications,
including supersonic and hypersonic propulsion, supernovae explosion and inertial confinement fusion. Numerical
simulations can provide insights into these phenomena, but simulations of turbulent flows with shocks and interfaces
are challenging due to the contradictory requirements to treat the different phenomena: while numerical dissipation is
added to stabilize the solution near discontinuities, this dissipation tends to overwhelm the physical dissipation.
To overcome this drawback, a hybrid WENO/central difference method is employed to solve the compressible
multicomponent Navier-Stokes equations with mass and heat transfer. A non-dissipative central difference scheme
is used in smooth regions while a weighted essentially non-oscillatory (WENO) scheme is employed near shocks;
a shock sensor specifies regions in which each scheme is applied. In the present multifluid extension, an additional
sensor based on the WENO weights of the mass fraction detects regions in which the fluid composition is discon-
tinuous. In these regions, a fully conservative interface-capturing WENO scheme is employed, in which spurious
pressure oscillations and temperature errors at interfaces are prevented. The resulting high-order accurate hybrid
method can therefore simulate turbulence, shock waves and interfaces in a stable and accurate fashion and is used
to study turbulent multi-material mixing in shock-accelerated flows. In particular, the set up of the experiments by
Balakumar et al. (2008) is considered, in which a shock wave travelling in air interacts with a layer of SF6. Two- and
three-dimensional calculations are discussed.


Nomenclature

Roman symbols
cp Specific heat at constant pressure
c, Specific heat at constant volume
D Mass diffusivity
c Internal energy per unit mass
E Total energy per unit volume
h Enthalpy per unit volume
j Mass flux
k Thermal conductivity
M Molecular mass
N Number of fluid components
p Pressure
q Heat flux
qd Enthalpy diffusion flux
R Gas constant
t Time
T Temperature


Ui Velocity vector
x Coordinate
z Mass fraction
Greek symbols
13 Bulk viscosity
6ij Kronecker delta
p Density
7 Ratio of specific heats
p Shear viscosity
Tij Viscous stress tensor
Subscripts
i, j, k Vector or tensor components
Superscipts
a", 3 Fluid component











Introduction

High-speed turbulent multicomponent and multiphase
flows with shocks waves occur in a variety of prob-
lems, such as in supersonic and hypersonic propulsion
systems, inertial confinement fusion, supernovae explo-
sion, and bubble clouds. In such problems, complex in-
teractions take place between shock waves, interfaces
separating different fluids and turbulence. One of the
key phenomena in such flows is the development of
the Richtmyer-Meshkov instability, which is a hydrody-
namic instability that occurs when a shock wave inter-
acts with a perturbed interface separating fluids of dif-
ferent densities. Because of the misalignment of the
pressure gradient (across the shock) and the density gra-
dient (across the interface), baroclinic vorticity is de-
posited along the interface, thus causing the perturbation
to grow. Although the early stages of this phenomenon
have been investigated extensively (Brouillette 2002),
few studies have explored the late-time behavior, i.e.,
when the mixing layer becomes turbulent. Though the
theory of incompressible turbulence in a single fluid is
mature, little is known about the turbulence generated in
accelerated multi-material flows, such as that produced
in the Richtmyer-Meshkov instability. Compressibility
and mass diffusion are important effects, so that the den-
sity and fluid properties vary based on the fluctuations of
the flow field. It is therefore far from obvious whether
the fundamental theories of turbulence and modeling are
valid in these situations.
Because of the challenges in carrying accurate and
controlled experiments of such phenomena, numerical
simulations are ideal tools to analyze these flows. How-
ever, a significant difficulty in the simulation of turbu-
lent flows with shocks and interfaces exists, in that the
requirements to represent shock waves, interfaces and
turbulence numerically are contradictory: while numeri-
cal dissipation is added to stabilize the solution near dis-
continuities (shocks, interfaces), this dissipation tends to
overwhelm the physical dissipation. To overcome this
drawback, a hybrid methodology is followed, in which
a shock- and interface-capturing scheme is used near
shock waves and interfaces, while a non-dissipative cen-
tral difference scheme is employed elsewhere.
In the present paper, a hybrid WENO/central differ-
ence method capable of solving the compressible mul-
ticomponent Navier-Stokes equations with mass and
heat transfer is introduced. A non-dissipative central
scheme is used in smooth regions while a weighted es-
sentially non-oscillatory (WENO) scheme is employed
near shocks (Larsson et al. 2007); a shock sensor spec-
ifies regions in which each scheme is applied. In the
present multifluid extension, an additional sensor based
on the WENO weights of the mass fraction detects re-


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


gions in which the fluid composition is discontinuous.
In these regions, a fully conservative interface-capturing
WENO scheme is employed, in which spurious pres-
sure oscillations and temperature errors at interfaces are
prevented. The resulting high-order accurate hybrid
method can therefore simulate turbulence, shock waves
and interfaces in a stable and accurate fashion. Three-
dimensional simulations multi-material mixing, with ap-
plication to the shock-gas curtain interaction of Balaku-
mar et al. (2008), are presented.

Governing equations

Thermodynamic relations.
The present work focuses on accurately simulate com-
pressible flows with multiple non-reacting gases. To
simplify the analysis, two gases, a and 3, are consid-
ered. The mass fraction of the first gas is z", and the
mass fraction of the second gas is z3 = 1 z" for bi-
nary systems. Thus, a single mass fraction is required,

Z ", 1 Z Z= (1)

The mixture molecular mass follows from the definition
of the mass fraction:
1 2 1
M+ -f (2)
M M M (2
In the present work, the gases are assumed to behave
ideally such that they obey the ideal gas law,


p pRT,


where


R,
R Rz + R (1 z).
M


Here, R, is the universal gas constant. It follows from
the definition of R that the specific heats at constant vol-
ume and pressure are given by


The ratio of specific heat is defined as


= Cp/Cv (6)

and thus depends on z. The following expression can
then be derived:
1 z M 1-z M
7 1 7a- 1 M + M (7)

Finally, defining the internal energy as

e = c,T, (8)


the ideal gas law can be re-written as
p
S- 1


C, = z + c(l z),


Cp Cpz +Cp3(1 -Z). (5)











Equations of motion.
The equations of motion for the multicomponent com-
pressible Navier-Stokes equations with heat and mass
transfer for N gases are:


0 + (PU )
ap a





' a +a- ('
at axj
OE a
8 xj


OTij
Oxj

S(UiTi
Oxj


(10a)


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


Thus, in the energy equation, the heat flux has two com-
ponents; the first is due to Fourier's law of heat conduc-
tion, and the second is the enthalpy diffusion term:


OT d
qj = -k + qj,

where k is the thermal conductivity.


(14)


(10b)
Binary systems.
(10) Equations (10) simplify significantly for binary sys-
tems:


qj),
(10d)


where a = 1,...., N 1. It is noted that the mass frac-
tions obey advection-diffusion equations, and that the
transport equations can be written in advection form. To
close the system, the equation of state (9) and appropri-
ate constitutive relations are required.
Constitutive relations.
The constitutive relations connect the diffusive fluxes
involved in transfer processes with the respective driving
potentials. When heat and mass transfer occur simul-
taneously in a fluid with non-zero velocity, the consti-
tutive relations are more complicated: diffusional mass
fluxes are not only created by composition gradients, but
also by pressure gradients (and body forces), as long
as the relevant forces act on the different fluid com-
ponents (Bird et al. 2002). As described by Eckert &
Drake (1972), in any originally isothermal mixture, the
diffusion-thermo effect creates an energy flux, which
will then create temperature differences. In general,
the thermal-diffusion and diffusion-thermo effects are
of a smaller order of magnitude than the effects due to
Fourier's and Fick's laws. Dufour and Soret effects are
neglected at this time in the present work.
For momentum diffusion, Newton's viscosity law is
used:


Tj = a u + ) + a
xj 9xij


2 a 9uk
23 6ij (11)
3 O xk


ap a
S+ (pa ) =0,

a(pz) a a
at + a- (pZUj) a
at Xj Oa-j

a(pu) 9
at xj .. .)


aE 0
-+-+
at axj

8x3


(15a)


(pDx


(15b)


(15c)


a
[u(E + p)] -= (uiT) (15d)
8xj
9T + 9hza + z
[k pD(ha +h ).
a T \ 5x1 a 5xje)
(15e)


The last term in equations (15e) can be written in several
enlightening forms:

1. Based on T and z:
a aT az
0 [kaOT+ pDT(ca- 09Z


2. Based on T and h:

a [(k cppD)T + pD Oh
3. Based onh andz:
3. Based on h and z:


0j [(
X C_


D) + pDT(c<
IU j


)t 9Zl


For mass diffusion, Fick's diffusion law is used:


a a za
ax


S D3 z3
3


For binary mixtures, the second term in equation (12) is
zero.
Mass transfer induces enthalpy diffusion in the energy
equation:
qj h j (13)


Transport coefficients.
In binary systems, the transport coefficients p, k and
D depend on the temperature; p and k additionally de-
pend on the fluid composition. The mixing rule for vis-
cosity used in the present work is that of Wilke, and for
thermal conductivity that of Wassiljewa with the Ma-
son and Saxena modification (Poling et al. 2000), which
have an analogous form (for ( = p, k):


(1 z)s'
.. + -z


(16)


z+(1-z "-.











where


,a3


[ 1( ,

+ 1+ )
(1 + T,3


In each fluid, the viscosity obeys a
Chapman-Enskog:


4 2

4]. (17)


given law, e.g.,


26.69 MT
QV, a2


where a is the collision diameter and Q, is the collision
integral:

Q, A(T*) B + Cexp(-DT*) + Eexp(-FT*),
(19)
where T* T/T, A 1.16145, B 0.14874, C
0.52487, D 0.77320, E = 2.16178, F = 2.43787
and Te e/k is the effective temperature characteristic
of the force potential function (k is Boltzmann's con-
stant).
In each fluid, the thermal conductivity obeys a given
law, e.g., Eucken:


c (1/ 9R)
S = 1 4 cc


c, (1 94
Af 1 4 (a


1)) (20)


For a binary mixture, Da = D". Thus, the mass
diffusivity coefficient is given by:

0.00266 T3/2
D D pV a)2 (21)

The collision integral for diffusion is given by

D = A(T)B + Cexp(-DT*) + Eexp(-FT*)
G exp(-HT*),
(22)

where A = 1.06036, B = 0.15610, C = 0.19300, D =
0.47635, E 1.03587, F 1.52996, G =1.76474,
H = 3.89411 and


2
1/ + ^


"a + Ta
2


Numerical method

The single-fluid Hybrid code.
The single-fluid Hybrid code (Larsson et al. 2007)
is based on the principle that turbulence and shock
waves are fundamentally different phenomena and
should thus be treated differently. Hence, to distin-
guish shock waves from smooth turbulent regions, the


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


method relies on a shock sensor based on vorticity
and dilatation that is similar to that of Ducros et al.
(1999). In smooth regions, a sixth-order accurate
explicit central differencing scheme is applied in split
form for improved nonlinear stability. In discontinuous
regions, a fifth-order accurate finite difference Weighted
Essentially Non-Oscillatory (WENO) scheme is used
(Jiang & Shu 1996) with a modified Harten, Lax and
van Leer (HLL) solver (Harten et al. 1983). The hybrid
nature of the code creates internal interfaces between
the central and WENO regions, the stability of which
was analyzed in Larsson & Gustafsson (2008). The
time-marching is handled using a standard fourth-order
accurate explicit Runge-Kutta scheme. The code has
been tested extensively (Johnsen et al. 2010a) and used
to study various problems in compressible turbulence
using direct numerical simulation (DNS), including the
interaction of isotropic turbulence and a normal shock
wave (Larsson & Lele 2009).

Numerical errors generated by shock-capturing
schemes in multi-material flows.
When applying high-order accurate shock-capturing
schemes to multiple fluids of different ratios of spe-
cific heats based on the Euler equations, numerical er-
rors are typically generated, e.g., spurious pressure os-
cillations, temperature errors and species conservation
errors. Though extensive work has been carried out
to prevent spurious pressure oscillations (Abgrall 1996;
Shyue 1998; Johnsen & Colonius 2006), temperature
and species conservation errors have mostly gone unre-
ported in the literature, primarily because they are either
passive effects (temperature errors) or are due to a lack
of resolution (species conservation errors). However, er-
rors in the temperature and mass fraction are detrimen-
tal to Navier-Stokes calculations in which diffusion is
present, because the diffusive terms depend on gradients
of the temperature and mass fraction: numerical errors
in these quantities generated in the convective terms then
lead to additional numerical errors in the diffusive terms.
Though it is true that diffusion will tend to damp these
Effects, they will nevertheless be present in the initial
3 transients and thus affect the solution.
k In the present effort, spurious pressure oscilla-
tions, temperature errors and species conservation er-
rors are prevented by following the procedure outlined
in Johnsen et al. (2010b). In order to prevent spurious
pressure oscillations, a transport equation for 1/ (7-1) is
solved in advection form, with a WENO reconstruction
of the average primitive variables (Johnsen & Colonius
2006). An additional transport equation for the mass
fraction is solved in conservative form, with the WENO
weights for pz consistent with those for p, thus prevent-
ing errors in the temperature and mass fraction (Johnsen








7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


0.06


0.04


0.02


0.0 1.0 2.0 3.0 4.0
tCHe/L
(a) Normalized Lo error in pressure.


U.60U +W
++

0.59 -

-1.0 -0.5 0.0 0.5 1.(
x/L
(b) Pressure.

1.50- ------------


1.50






SIII00


0.0
x/L


0.5 1.0


(c) Temperature.


Figure 1: Profiles at tcHe/L = 4 for the advection of a
material interface between helium and nitrogen. Solid
line: exact solution; pluses: fully conservative; open
triangles: Abgrall (1996); open squares: Shyue (1998);
open circles: present scheme.


et al. 2010b).
In order to illustrate these errors, the one-dimensional
advection of a material interface is considered. A top-
hat distribution of nitrogen in helium moves at the speed
of sound on a periodic domain. In the absence of physi-
cal diffusion, the distribution is expected to advect, and
there should be no changes in the pressure, velocity or
temperature. The density, pressure and temperature pro-
files are plotted after the distribution has travelled two
periods in Figure 1. Figure 2 shows time histories of the


10-
10-3
10 "
10'.
10-
10-


10-1




0.30


0I


o

S0


.20


1.0 2.0 3.0 4.0
tCHe/L
(b) Normalized Lo error in temperature.

* -- -~ ~ '


.UE11


0.0 1.0 2.0 3.0
tCHe/L
(c) Normalized error in mass of helium.


Figure 2: Time histories for the advection of a ma-
terial interface between helium and nitrogen. Pluses:
fully conservative; open triangles: Abgrall (1996); open
squares: Shyue (1998); open circles: present scheme.


errors in pressure, temperature and mass of helium. A
fully conservative scheme, that of Abgrall (1996) trian-
gles, that of Shyue (1998) and the present scheme are
included.
It is clear that the fully conservative scheme leads
to spurious pressure oscillations, which then propogate
to all other variables. The other three schemes are
specifically designed to prevent this drawback and thus
do not generate any errors in the pressure. The scheme
of Shyue (1998) generates errors in the temperature
up to 50% in this simple problem; when considering


+
. ++.

+ + o
+
+. ~ + +
++ r;.
+ --; -


x/L
(a) Density.


O
u
4


+ + + .+ + + +++ 4I +
S+ + + ++ + + +
S+ ++


. .........." H WH N ...... ............ ...".". .


r


^


""-1.0







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


the Euler equations only, this drawback is minor, since
temperature is (passi\ eh ) advected. However, for the
Navier-Stokes equations, any numerical errors in the
temperature are going to affect the diffusive terms and
thus propagate to the other variables. In the scheme of
Abgrall (1996), temperature errors are small; however
a significant drawback is that no explicit relationship
between temperature and mass fraction is provided.
The scheme of Shyue (1998) does not preserve the
mass of helium, as expected, nor does the scheme of
Abgrall (1996), but to a lesser degree. The present
multifluid algorithm generates no errors in the pressure,
temperature or mass of helium. Hence, no errors due to
the convective terms are generated.


The multifluid Hybrid code.

The Hybrid code can readily be extended to the
previously mentioned multifluid algorithm, though two
main items must be addressed. First, a new interface
sensor is developed. Following the idea of Hill & Pullin
(2'i1 4) for a shock sensor, the WENO weights of mass
fraction are used to determine regions in which the
WENO scheme must be used (i.e., regions with sharp
changes in mass fraction). Second, a finite volume
formulation of the WENO scheme must be used to
preserve the pressure equilibrium (Johnsen & Colonius
2006). However, the base hybrid code is written in
finite-difference form. Thus, a second-order error is
introduced when going from the point values to the
cell-average values in the WENO reconstruction. This
error is not expected to affect the solution significantly,
because it only affects regions in which WENO is ap-
plied, i.e., when the solution is discontinuous. However,
it is well known that any high-order accurate scheme
will reduce to first-order accurate at discontinuities, so
that the second-order accuracy achieved at the interface
is not expected to be detrimental to the overall solution.
Thus, DNS of turbulent multi-material flows can be
carried out using the present scheme.


Diffusion effects in the advection of a material inter-
face.

As shown previously, a material interface is simply
advected if there is no physical diffusion; there are no
changes in the pressure, temperature or velocity. How-
ever, if mass diffusivity, viscosity and thermal conduc-
tivity are included, the pressure, temperature and veloc-
ity no longer remain constant. This behavior is shown in
this section.

From the energy equation (15e) and using the mo-
mentum equation (15c), the following equation can be


derived for the internal energy:


at + a (peuj)
dt d).r


-p + 9(UiTij)
ax1 atry


(24)


The internal energy changes due to convection, pressure
work, viscous stresses, thermal conduction and enthalpy
diffusion (i.e., molecular diffusion). Thus, evolution
equations for temperature and pressure can be derived
using equations (8) and (9).

Temperature
The evolution equation for temperature can be de-
rived by substituting the definitions in equation (8)
into equation (24) and use the fact that c, is related
to z via equation (5). With the help of equation
(15b), an evolution equation for the temperature
can be written in the following form:


(IT iT \
P O-T + ujT
at Oxjr
+a k (T
x a ax


auk, a
axU1+ a(UTj)
iDT iX z
S- C pD) pD
T Oxj 9Oxj


+ (Ra R) T (pD az


(25)


It is therefore clear that, even for an initially
isothermal system with constant velocity, temper-
ature fluctuations will occur because of the compo-
sition gradients in the last term and will affect all
the other variables, since R" K R3 in general.

* Pressure
Before considering the evolution equation for pres-
sure, the evolution equation for 1/(7 1) is de-
rived. Using the relationship between z and 7, the
transport equation (15b) can be written in terms of
1/(7 1):


-Y 1 O} rY 1
1_ 1 ) D 1

1Ia pD a 1
p axj [r aj -Y 1


+ 2D


1
7 1


1 M 3 1

(26)


Thus, using the definition in equation (9), and equa-
tions (5), an evolution equation for pressure can be


+ a [kax+ ( c)pDT-0j
9x~j I x~j 9Oxij







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


written:


( 1 at a a jY 1 Xk
I-l)O[ +u3xJO 7-l P
at(UtTij)+ a ki (OT)
az ar

(c, c ) pD Oz OT
(a TPa T

+ (R- R3) Y T a pD .
a1 (T DQ.Ox
(27)

Similarly to the evolution equation for tempera-
ture, it is clear that, even for an initially isothermal
system with constant velocity, pressure fluctuations
will occur because of the composition gradients and
will affect all the other variables.

Extension to multiphase flows.
For flows of gases and liquids in which the only phase
interactions are described by the diffusive processes pre-
viously discussed, the present formulation is easily ex-
tended. Modifications must be made in the diffusivity
coefficients (mass diffusivity, viscosity, thermal conduc-
tivity) and in the equation of state (e.g., using a stiff-
ened or Tait-type equation). However, two phenom-
ena present significant numerical challenges, due to the
fact that they occur at the interface: surface tension and
phase change.

Results

The present interest lies in the turbulent multi-material
mixing that occurs in shock-accelerated flows. In or-
der to validate the computational results, high-quality
experimental data are required. Though many experi-
mental investigations of the Richtmyer-Meshkov insta-
bility have been carried out throughout the years, only
recently have detailed data of the turbulence statistics
been published. In particular, the experiments of Bal-
akumar et al. (2008), in which simultaneous Particle Im-
age Velocimetry (PIV) and Planar Laser Induced Fluo-
rescence (PLIF) measurements are carried out to provide
velocity and density fields, are considered. A Mach 1.22
shock traveling in air impacts a layer of SF6. The layer is
created by injecting 21 cylinders of SF6 and maintaining
them in equilibrium until the shock arrives. The problem
is nominally two-dimensional, until the shock reflects
off the end of the shock-tube and interacts with the now
multi-modal perturbation again (re-shock). Small per-
turbations in the other transverse direction exist and the
mixing layer becomes turbulent thereafter.
Preliminary two-dimensional simulations were car-
ried out to explore the parameter space and set up the


80


4 60


40


200 400 600 800 1000 1200
t (Ps)

Figure 3: Location of the center of mass of the SF6 layer
in the two-dimensional shock-curtain interaction. Solid
line: simulations; upward blue triangles: experiments
(no re-shock); downward green triangles: experiments
(with re-shock).


three-dimensional calculations, and meaningful compar-
isons were made to the experiments. For instance, the
center of mass of the SF6 layer was tracked with time
in the re-shock case and exhibited good agreement with
the experiments, as shown in Figure 3. As expected, the
SF6 layer is accelerated to the velocity downstream of
the transmitted shock; after re-shock, the center of mass
of the mixing layer is not expected to move significantly.
Balakumar et al. (2008) remark that it is not clear why
the center of mass deviates from rest at late times in the
experiments. In these simulations, diffusive effects are
included, but the results are not converged. The grid is
uniform and has 512 points across the width of the shock
tube, which has width 0.0762 m. Thus, there were ap-
proximately 24 points across a given initial perturbation.
Three-dimensional simulations are presented as well.
A small long-wavelength perturbation (5% of the initial
perturbation) is introduced in the second transverse di-
rection to render the problem three-dimensional, as ob-
served in the experiments (after re-shock). Iso-contours
of the mass fraction of SF6 are shown in Figure 4.
The flow is nominally two-dimensional until re-shock,
which occurs slightly before the second frame. There-
after, both the large and small initial perturbations be-
come much more prone to instability and their growth
evolves dramatically. However, the flow is not suffi-
ciently resolved to extract meaningful turbulence statis-
tics; in other words, the numerical dissipation may hin-
der the evolution of small scale structures. In these
three-dimensional calculations, diffusive effects are not
included. The resolution in the transverse direction is
256 points across the width of the shock tube. These
calculations required the code to be run in parallel, as
33 million points were used; thus, additional resolution




























3.2 "

( 3 -1A
(a) After the passage of the incoming shock (but before re-shock).


-0.5 ,1
2.5


2


-D S


(b) Just after re-shock


- "
0,5


2 -

1.5 1
(c) At a late time after re-shock

Figure 4: Isocontours of the mass fraction of SF6 for
the three-dimensional shock-curtain interaction.



requirements lead to the need for massive computational


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


resources and subsequent difficulties. Nevertheless, the
present code has been shown to scale well with the num-
ber of processors.


Conclusions

A hybrid WENO/central difference method has been
developed to solve the compressible multicomponent
Navier-Stokes equations with mass and heat transfer. A
non-dissipative central scheme is used in smooth regions
while a weighted essentially non-oscillatory (WENO)
scheme is employed near shocks; a shock sensor spec-
ifies regions in which each scheme is applied. In the
present multifluid extension, an additional sensor based
on the WENO weights of the mass fraction detects re-
gions in which the fluid composition is discontinuous.
In these regions, a fully conservative interface-capturing
WENO scheme is employed, in which spurious pres-
sure oscillations and temperature errors at interfaces are
prevented. The resulting high-order accurate hybrid
method can therefore simulate turbulence, shock waves
and interfaces in a stable and accurate fashion and is
used to study turbulent multi-material mixing in shock-
accelerated flows. In particular, the set up of the ex-
periments by Balakumar et al. (2008) is considered, in
which a shock in air interacts with a layer of SF6. Two-
and three-dimensional calculations are presented.
The future directions for this work will include fur-
ther comparisons with the experiments and more de-
tailed studies of the physics of turbulent multi-material
mixing. In particular, the role of diffusion and the effect
of the different terms in the equations will be investi-
gated. Further analysis will also be carried out to gen-
erate initial conditions more representative of the exper-
iments and the resolution requirements for the present
calculations. Finally, the present methodology will be
extended to flows of liquids and gases.


Acknowledgements

The authors wish to thank Sanjiva Lele and Frank Ham
for helpful conversations. This work was supported by
DOE-SciDAC2 (Grant DE-FC02-06-ER25787).


References

Abgrall R., How to prevent pressure oscillations in mul-
ticomponent flow calculations: A quasi-conservative ap-
proach, J. Comput. Phys., Vol. 125, pp. 150-160, 1996.

Balakumar B. J. Orlicz G. C., Tomkins C. D., and Pre-
stridge K. P., Simultaneous particle-image velocimetry-
planar laser-induced fluorescence measurements of
Richtmyer-Meshkov instability growth in a gas curtain







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


with and without reshock, Phys. Fluids, Vol. 20, 124103, Poling B. E., Prausnitz J. M., and O'Connell J. P., The
2008. Properties of Gases and Liquids, McGraw-Hill, New
York, 2000.
Bird R. B., Stewart W. E., Lightfoot E. N., Transport
Phenomena, John Wiley & Sons, New York, 2002. Shyue K. M., An efficient shock-capturing algorithm
for compressible multicomponent problems, J. Comput.
Brouillette M., The Richtmyer-Meshkov instability, Phys., Vol. 142, pp. 208-242, 1998.
Annu. Rev. Fluid Mech. Vol. 34, pp. 445-468, 2002.

Ducros F, Ferrand V, Nicoud F, Weber C., Darracq
D., Gacherieu C., Poinsot T., Large-eddy simulation of
the shock/turbulence interaction, J. Comput. Phys., Vol.
152,pp. 517-549,1999.

Eckert E. R. G. and Drake Jr. R. M., Analysis of Heat
andMass Transfer McGraw-Hill, New York, 1972.

Harten A., Lax P D., and van Leer B., On upstream
differencing and Godunov-type schemes for hyperbolic
conservation laws, SIAM Review, Vol. 25, pp 35-61,
1983.

Hill D. J. and Pullin D. I., Hybrid tuned center-
difference-WENO method for large eddy simulations in
the presence of strong shocks, J. Comput. Phys., Vol.
194, pp. 435-450, 2004.

Jiang G. S. and Shu C. W., Efficient implementation of
WENO schemes, J. Comput. Phys., Vol. 126, 202-228,
1996.

Johnsen E. and Colonius T., Implementation of WENO
schemes in compressible multicomponent flow prob-
lems, J. Comput. Phys. Vol. 219, pp. 715-732, 2006.

Johnsen E., Larsson J., Bhagatwala A. V, Cabot W. H.,
Moin P, Olson B. J., Rawat P. S., Shankar S. K., Sjo-
green B., Yee H. C., Zhong X., and Lele S. K., Assess-
ment of high resolution methods for numerical simula-
tions of compressible turbulence, J. Comput. Phys., Vol.
228, pp. 1213-1237,2010.

Johnsen E., Larsson J., and Ham F., Numerical errors
generated by shock-capturing schemes in compressible
multicomponent flows, J. Comput. Phys., in preparation.

Larsson J., Lele S. K., and Moin P., Effect of numeri-
cal dissipation on the predicted spectra for compressible
turbulence, CTR Ann. Res. Briefs, pp. 47-57, 2007.

Larsson J. and Gustafsson B., Stability criteria for hy-
brid difference methods, J. Comput. Phys., Vol. 227, pp.
2886-2898,2008.

Larsson J. and Lele S. K., Direct numerical simulation
of canonical shock/turbulence interaction, Phys. Fluids,
Vol. 21, 126101, 2009.




University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - Version 2.9.7 - mvs