7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
A lowdissipation 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, RichtmyerMeshkov instability
Abstract
Highspeed 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 NavierStokes equations with mass and heat transfer. A nondissipative central difference scheme
is used in smooth regions while a weighted essentially nonoscillatory (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 interfacecapturing WENO scheme is employed, in which spurious
pressure oscillations and temperature errors at interfaces are prevented. The resulting highorder accurate hybrid
method can therefore simulate turbulence, shock waves and interfaces in a stable and accurate fashion and is used
to study turbulent multimaterial mixing in shockaccelerated 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
threedimensional 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
Highspeed 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 RichtmyerMeshkov 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 latetime 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 multimaterial flows, such as that produced
in the RichtmyerMeshkov 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 interfacecapturing scheme is used near
shock waves and interfaces, while a nondissipative 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 NavierStokes equations with mass and
heat transfer is introduced. A nondissipative central
scheme is used in smooth regions while a weighted es
sentially nonoscillatory (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 interfacecapturing
WENO scheme is employed, in which spurious pres
sure oscillations and temperature errors at interfaces are
prevented. The resulting highorder accurate hybrid
method can therefore simulate turbulence, shock waves
and interfaces in a stable and accurate fashion. Three
dimensional simulations multimaterial mixing, with ap
plication to the shockgas 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 nonreacting 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 1z M
7 1 7a 1 M + M (7)
Finally, defining the internal energy as
e = c,T, (8)
the ideal gas law can be rewritten 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 NavierStokes 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 advectiondiffusion 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 nonzero 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
diffusionthermo effect creates an energy flux, which
will then create temperature differences. In general,
the thermaldiffusion and diffusionthermo 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 Oaj
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+(1z ".
where
,a3
[ 1( ,
+ 1+ )
(1 + T,3
In each fluid, the viscosity obeys a
ChapmanEnskog:
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 singlefluid Hybrid code.
The singlefluid 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 sixthorder accurate
explicit central differencing scheme is applied in split
form for improved nonlinear stability. In discontinuous
regions, a fifthorder accurate finite difference Weighted
Essentially NonOscillatory (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
timemarching is handled using a standard fourthorder
accurate explicit RungeKutta 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 shockcapturing
schemes in multimaterial flows.
When applying highorder accurate shockcapturing
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 NavierStokes 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/ (71) 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 onedimensional
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
103
10 "
10'.
10
10
101
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
NavierStokes 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
finitedifference form. Thus, a secondorder error is
introduced when going from the point values to the
cellaverage 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 highorder accurate scheme
will reduce to firstorder accurate at discontinuities, so
that the secondorder accuracy achieved at the interface
is not expected to be detrimental to the overall solution.
Thus, DNS of turbulent multimaterial 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 OT + 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)pDT0j
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
Il)O[ +u3xJO 7l 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 Taittype 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 multimaterial
mixing that occurs in shockaccelerated flows. In or
der to validate the computational results, highquality
experimental data are required. Though many experi
mental investigations of the RichtmyerMeshkov 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 twodimensional, until the shock reflects
off the end of the shocktube and interacts with the now
multimodal perturbation again (reshock). Small per
turbations in the other transverse direction exist and the
mixing layer becomes turbulent thereafter.
Preliminary twodimensional 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 twodimensional shockcurtain interaction. Solid
line: simulations; upward blue triangles: experiments
(no reshock); downward green triangles: experiments
(with reshock).
threedimensional 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 reshock 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 reshock, 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.
Threedimensional simulations are presented as well.
A small longwavelength perturbation (5% of the initial
perturbation) is introduced in the second transverse di
rection to render the problem threedimensional, as ob
served in the experiments (after reshock). Isocontours
of the mass fraction of SF6 are shown in Figure 4.
The flow is nominally twodimensional until reshock,
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
threedimensional 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 reshock).
0.5 ,1
2.5
2
D S
(b) Just after reshock
 "
0,5
2 
1.5 1
(c) At a late time after reshock
Figure 4: Isocontours of the mass fraction of SF6 for
the threedimensional shockcurtain 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
NavierStokes equations with mass and heat transfer. A
nondissipative central scheme is used in smooth regions
while a weighted essentially nonoscillatory (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 interfacecapturing
WENO scheme is employed, in which spurious pres
sure oscillations and temperature errors at interfaces are
prevented. The resulting highorder accurate hybrid
method can therefore simulate turbulence, shock waves
and interfaces in a stable and accurate fashion and is
used to study turbulent multimaterial 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 threedimensional 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 multimaterial
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
DOESciDAC2 (Grant DEFC0206ER25787).
References
Abgrall R., How to prevent pressure oscillations in mul
ticomponent flow calculations: A quasiconservative ap
proach, J. Comput. Phys., Vol. 125, pp. 150160, 1996.
Balakumar B. J. Orlicz G. C., Tomkins C. D., and Pre
stridge K. P., Simultaneous particleimage velocimetry
planar laserinduced fluorescence measurements of
RichtmyerMeshkov 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, McGrawHill, 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 shockcapturing algorithm
for compressible multicomponent problems, J. Comput.
Brouillette M., The RichtmyerMeshkov instability, Phys., Vol. 142, pp. 208242, 1998.
Annu. Rev. Fluid Mech. Vol. 34, pp. 445468, 2002.
Ducros F, Ferrand V, Nicoud F, Weber C., Darracq
D., Gacherieu C., Poinsot T., Largeeddy simulation of
the shock/turbulence interaction, J. Comput. Phys., Vol.
152,pp. 517549,1999.
Eckert E. R. G. and Drake Jr. R. M., Analysis of Heat
andMass Transfer McGrawHill, New York, 1972.
Harten A., Lax P D., and van Leer B., On upstream
differencing and Godunovtype schemes for hyperbolic
conservation laws, SIAM Review, Vol. 25, pp 3561,
1983.
Hill D. J. and Pullin D. I., Hybrid tuned center
differenceWENO method for large eddy simulations in
the presence of strong shocks, J. Comput. Phys., Vol.
194, pp. 435450, 2004.
Jiang G. S. and Shu C. W., Efficient implementation of
WENO schemes, J. Comput. Phys., Vol. 126, 202228,
1996.
Johnsen E. and Colonius T., Implementation of WENO
schemes in compressible multicomponent flow prob
lems, J. Comput. Phys. Vol. 219, pp. 715732, 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. 12131237,2010.
Johnsen E., Larsson J., and Ham F., Numerical errors
generated by shockcapturing 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. 4757, 2007.
Larsson J. and Gustafsson B., Stability criteria for hy
brid difference methods, J. Comput. Phys., Vol. 227, pp.
28862898,2008.
Larsson J. and Lele S. K., Direct numerical simulation
of canonical shock/turbulence interaction, Phys. Fluids,
Vol. 21, 126101, 2009.
