Astro2020 Science White Paper
Understanding Galactic Star Formation with Next Generation X-ray Spectroscopy
and Imaging

Thematic Areas: $\square$ ✓ Planetary Systems $\square$ ✓ Star and Planet
Formation $\square$ Formation and Evolution of Compact Objects $\square$
Cosmology and Fundamental Physics $\square$ ✓ Stars and Stellar Evolution
$\square$ ✓ Resolved Stellar Populations and their Environments $\square$
Galaxy Evolution $\square$ Multi-Messenger Astronomy and Astrophysics

Principal Author:

Name: Scott J. Wolk Institution: Harvard–Smithsonian Center for Astrophysics
Email<EMAIL_ADDRESS>Phone: 617-496-7766

Co-authors: Rachel Osten (STSCI, JHU), Nancy Brickhouse (Harvard-Smithsonian
Center for Astrophysics), Moritz Günther (Kavli Institute, MIT), Laura A.
Lopez (The Ohio State University), Jeremy Drake (Harvard-Smithsonian Center
for Astrophysics), Benjamin F. Williams (University of Washington), Elaine
Winston (Harvard-Smithsonian Center for Astrophysics), Denis Leahy (University
of Calgary), Panayiotis Tzanavaris (NASA/GSFC-CRESST), David A. Principe
(Kavli Institute, MIT)

Abstract: This white paper is motivated by open questions in star formation,
which can be uniquely addressed by high resolution X-ray imaging and require
an X-ray observatory with large collecting area along good spectral
resolution. A complete census of star-forming regions in X-rays combined with
well matched infrared (IR) data will advance our understanding of disk
survival times and dissipation mechanisms. In addition, we will be able to
directly observe the effects of X-ray irradiation on circumstellar grain
growth to compare with grain evolution models in both high- and low-UV
environments. X-rays are native to stars at all phases of star formation and
affect planet-forming disks especially through flares. Moreover, X-rays trace
magnetic fields which weave through the flares, providing a unique, non-
gravitational feedback mechanism between disk and star. Finally, the bright
X-ray emission emanating from hot plasma associated with massive stars can
have large scale impacts on the topology of star-forming regions and their
interface with the interstellar medium (ISM).

## Recent Progress in Star Formation

Understanding the fundamental aspects of star formation and evolution has been
a major success of the “Great Observatories” era. Observations of star-forming
regions are among the most challenging in modern astronomy. These regions are
very crowded, with two-dimensional source densities exceeding 100 pc-2, while
their dynamic range is large, with pre-main sequence (PMS) brown dwarfs
residing within an arcsecond of massive O stars. The overall range in stellar
luminosity can easily cover 20 MV. These regions are often filled with both
hot gas – which creates spectroscopic difficulties – and cold dust – which
complicates the interpretation of the observed data and obscures some objects
completely.
Infrared (IR) missions such as Herschel, Spitzer and 2MASS have been critical
in identifying young stellar objects (YSOs) with excess emission from
envelopes and disks. Modern sub-millimeter observatories, such as ALMA, can
image some nearby disks in exquisite detail. Similarly JWST will reveal
remarkable detail on specific systems. High resolution optical imaging and
spectroscopy have already brought us new vistas of incredible systems such as
HH30. However, the optical data are limited to cases with fairly low dust
obscuration and are mostly tracking scattered light. Meanwhile, the IR
instruments only identify young stars while they are in possession of their
natal disks, or in some cases, second generation debris disks. In contrast,
X-ray indicators of stellar youth last for hundreds of million to billions of
years, covering a wide variety of masses and circumstances. X-ray data serve
not only as fingerprints of stellar youth, but the X-rays themselves are
produced from the stellar dynamo and give important constraints on stellar
magnetic field structure.

Figure 1: Graphic comparing the distance out to which Chandra ($\sim$ 600 cm2
effective area) and a modern, 10,000 cm2 effective area X-ray telescope are
sensitive enough to take a complete census (from brown dwarfs to massive
stars) of the young star clusters. While Chandra can only conduct a this
inventory out to the Orion cluster ($D\approx$400 pc), an improvement in
effective area with comparable spatial resolution would expand this reach to
several kpc.

## The Role of X-Ray Studies

Young stars are bright, ubiquitous X-ray sources. For high-mass stars, most of
the X-ray energy is coming from wind shocks and wind-wind collisions. The
winds then interact with the ISM, depositing thermal energy and creating
turbulence (e.g. MacLow & Klessen 2004). This turbulence inhibits other
ongoing star formation. For low-mass stars in the accretion phase, specific
X-ray line ratios act as temperature and density diagnostics which allow us to
directly infer the infall rate and constraining mechanisms (Brickhouse et al.
2010). This phase is also associated with enhanced flaring in which 100 MK
long-lived flares are likely responsible for rapid heating of protoplanetary
disks and can have deleterious effects on young planets. Moreover, X-ray
emission in general impacts circumstellar disk chemistry (Cleeves et al. 2015)
and eventually disk dispersal (Owen et al. 2011).

For lower-mass stars, X-rays are generated in the corona which has its origin
at the dynamo. Modeling based on recent data indicates that the dynamo takes
different forms depending upon the mass of the star (Yadav et al. 2016). For
the most part, X-ray luminosity is seen to decrease with age following an
exponential law with $L_{X}\propto t^{-0.5}$ (c.f. Favata & Micela 2003). This
is driven by the convective turnover rate and indicates that late M stars have
elevated X-ray fluxes for up to 5 billion years (Wright et al. 2013). This
long-term X-ray luminosity allows the identification of young stars long after
they lose their dust disks. Thus, stellar X-ray emission probes the evolution
of star formation within existing clusters, and and allows us to study the
star formation history of the whole Galaxy (Pillitteri et al. 2013, Getman et
al. 2014).

To date, Orion (at a distance of $\sim$400 pc) has been the main astrophysical
laboratory where all the available observatories have been brought to bear.
Studies of this cluster, and the relationship between its star-forming sub-
regions and various triggering mechanisms – from the youngest, deeply-embedded
protostars discovered by Herschel, to the $\zeta$ Ori sub cluster, newly
revealed in X-rays (Alves & Buoy 2012, Megeath et al. 2012, 2016, Getman et
al. 2005 etc.) – have demonstrated the power of using well-matched
observatories in concert. Multiwavelength facilities with comparable spatial
resolution offer a complete view of the star-formation process and its
relation to the ISM. In particular, X-rays are crucial to probe highly
embedded sources, which have unique X-ray signatures. This is demonstrated by
direct X-ray detections of protostars (c.f. Pravdo et al. 2009) as well as the
detection of 6.4 KeV iron fluorescence, which is from circumstellar disk
material being irradiated by stellar flares and exciting iron fluorescence
(Tsujimoto et al. 2005). Consequently, it is possible to apply the
reverberation mapping technique, used in studies of active galactic nuclei, to
resolve AU-scale details on stellar disks that are several kpc from the sun.

A complete multi-wavelength and multi-cluster census will be much more
impactful and informative than a study of a single cluster alone, no matter
how well-studied. For instance, data from Orion and the nearby clusters have
been used to constrain circumstellar disk survival times, a crucial factor to
halt planet formation. Such a study has been carried out for a relatively
small number of clusters within 1 kpc of the Sun (Mamajek 2009), but has not
been carried out for the more distant and more massive clusters. These kinds
of studies naturally lead into the study of transition disk timescales. The
term transition disk derives from the fact that such YSOs either lack an IR
excess at near- and mid- IR wavelengths (e.g. Wolk & Walter 1996) or have
spectral energy distribution (SED) slopes in-between those of typical star-
forming disks and pure photospheres (e.g. Allen et al. 2004). X-ray and high-
energy flux appears to have a negative effect on grain growth in clusters
(c.f. Winston et al. 2010). Meanwhile many of the structures observed in
broadband IR and optical observations seem to have their origin in an X-ray
bright hot plasma expanding out into the ISM and being constrained by existing
dust (c.f. Townsley et al. 2003 Wolk et al. 2002; Figure 2).

## X-ray Astronomy Beyond the Next Decade

While this white paper is directed at the 2020 Decadal committee, it is
already likely that there will be no new major X-ray mission launched before
$\sim$2030\. Therefore the scientific problems to be solved will be already 10
years beyond today. By 2030, the current facilities XMM-Netwon and Chandra
will be nearing end of life. While e-ROSITA will provide a deeper all-sky
survey of X-ray sources than currently exists, it is not designed as a
detailed followup instrument. XRISM will have limited sensitivity and spatial
resolution, despite the groundbreaking spectral resolution expected from its
microcalorimeter. Athena should launch in the early 2030’s but will lack the
high spatial resolution need to understand compact regions such as the
Kleinmann-Low Nebula even within 1 kpc (Figure 3). Arcsecond resolution is
also required to maximize sensitivity.

By 2030, many new missions will have made significant advances related to
stars; however, each will have its own biases or limitations. TESS and PLATO
will have surveyed millions of stars to find planets, and JWST and WFIRST will
have provided deep IR and optical imaging and spectroscopy in crowded, star-
forming fields down to well below the hydrogen burning limit. These
observations will provide constraints based primarily on disk-bearing stars,
leading to significant biases. In the radio regime, ALMA, SKA and the next
generation VLA will be making advances in the chemistry and structure of
planet-forming disks. Large aperture ground based optical telescopes will be
able to use Zeeman splitting to measure the surface magnetic fields in a
variety of stars. Gaia may appear to be the most appropriate mission to
identify YSOs in star-forming regions since, in theory, it will create a
three-dimensional model of all the stars in every star-forming region within
10 kpc. However, it should be noted that interstellar dust will greatly limit
that range in the plane of the galaxy, where star formation is usually found.
And moreover, the local dust and ionized emission, which surround young stars,
limits Gaia’s ability to map stars even in the closest clusters.

With the current generation of X-ray telescopes and those approved to be
operating before 2030, Orion is the only massive cluster close enough to
enable detailed study. However Orion is only one cluster. It is relatively
young, with an age of $\lesssim$20 Myr. As a single supercluster, it
represents a single set of initial conditions. Triggering mechanisms and total
mass should affect the supernova rate and in turn the final elemental
abundances and the energy budget returned to the ISM. It only represents a
single sample of an IMF. Even this is unclear as all the O stars seem to lie
either in the Orion nebular cluster or in the belt. On the other hand, there
are thousands of young stars below the belt in Lynds 1641 without a proximate
O or B star.

Figure 2: (Left) Composite image of 30 Doradus. Multi-million-degree gas
detected in X-rays by Chandra originating in shock fronts – similar to sonic
booms –formed by stellar winds and supernova explosions. In this rendering,
0.3-1.5 keV photons are red, 1.5-2.5 keV are green and 2.5-7 keV are blue.
This hot gas carves out gigantic bubbles in the surrounding cooler gas and
dust, shown here is IR dust emission from the Spitzer Space Telescope (green).

In one simple example, our current understanding of disk survival times is
limited by current X-ray sensitivity. Kuhn et al. (2014, 2015) studied
archival Spitzer and Chandra data of about 15 regions of massive star
formation. They find that IR-excess sources typically only represent about 50%
of the cluster. However they also find that the X-ray detection rates plummet
with distance. Most of the star-forming regions within the sensitivity limits
of Chandra and XMM-Newton are within 2 kpc, the furthest being NGC 1893 at 3.6
kpc. This represents a selection bias since further clusters are simply
ignored. Worse still, the incompleteness begins to become a problem for X-ray
sources below about 1 $M_{\odot}$ at 2 kpc. A greater effective area will
increase the grasp in the coming decade, but it will be important to also
ensure the sub-arcsecond imaging necessary to alleviate source confusion. As
the star-forming regions of interest move further away, the relative
separations of stars shrink. Moreover, as absolute stellar density increases
with star cluster mass, spatial resolution will be crucial with increasing
X-ray telescope sensitivity (cf. Figure 3).

Figure 3:  (Left) VLT JKL band image of the Becklin–Neugebauer Object in the
Kleinmann–Low Nebula. Several deeply embedded objects can be seen in red.
(Right) The same region observed by Chandra with 0.5′′ resolution. The blue
color is due to sources only detected in hard X-rays. This occurs partially
because young stars have harder emission and softer X-rays are absorbed by gas
near the star. Note that hard X-ray sources track the red stars well but not
perfectly. Note the 5′′ scale bar. This kind of analysis is only possible due
to the good match between Chandra’s subarcsecond imaging and the resolution of
VLT and most planned IR and optical missions in the 2020’s.

A high-spatial resolution, high collecting area X-ray telescope with a wide-
field integral field spectrograph (IFS) will fundamentally change how we
approach star formation in the X-rays and keep pace with expected improvements
in technology already planned at other wavelengths. Given an effective area of
about 10,000 cm2, an observation of Orion could reach the same depth as the
10-day (850ks) COUP (Getman et al. 2005) program in under one day of
observing. Assuming the data were taken with an IFS, each of the nearly 1000
stars detected in the central cluster region (perhaps 5′ on a side) would have
enough signal to measure not just one temperature but to detect all critical
temperature- (and density-) sensitive lines such as Si xiii, Si xiv, Si xv, Mg
x, Mg xi, Mg xii, O vi, O vii, O viii. These data and the changes observed due
to flares would be irreplaceable for regrading the high energy exposure of the
inner disk and young (proto) exoplanets. These conditions are fundamental to
the local space environment and the formation and habitability of exoplanets.
Further, given the greater sensitivity provided by both 10,000 cm2 and arcsec
imaging, Orion would no longer be the only star-forming region we could
examine in such detail. Full surveys complete to one-tenth of a solar mass and
including high quality multi-temperature and density sensitive spectra could
be carried out on clusters, such as Carina (2.4 kpc), NGC 281 (2.1 kpc), NGC
3603 (6.9 kpc), Arches (8 kpc), and even 30 Doradus (50 kpc) in the Large
Magellanic Cloud. The duration of such a program would still be less than the
amount of time spent by Chandra on Orion alone.

## Summary

The next major X-ray facilities will not arrive until the 2030s. These
facilities will join an era of star and planet formation science that will
have been transformed by a fleet of missions and ground based facilities that
will come of age in the 2020s. As a community, we cannot allow this high
energy window into our Universe to close at a time when star and planet
formation will be a major focus of the next generation of multiwavelength
telescopes. X-rays will still play a fundamental role in our understanding of
cluster morphology, stellar evolution, and the survival timescales of planet-
forming disks and their dust evolution. Even the observed IMF will be out of
reach without X-ray facilities working together with their multiwavelength,
lower energy counterparts. To pursue these science cases, the community needs
major X-ray facilities that require effective area of $\approx$ 10,000 cm2, an
integral field spectrograph, and imaging quality on the order of an arcsecond.

References

Allen, L. E., et al. 2004, ApJS, 154, 363

Alves, J., & Bouy, H. 2012, A&A, 547, A97

Brickhouse, N. S., Cranmer, S. R., Dupree, A. K., Luna, G. J. M., & Wolk, S.
2010, ApJ, 710, 1835

Cleeves, L. I., Bergin, E. A., Qi, C., Adams, F. C., & Öberg, K. I. 2015, ApJ,
799, 204

Favata, F., & Micela, G. 2003, Space Science Reviews, 108, 577

Getman, K. V., et al. 2014, ApJ, 787, 108

Getman, K. V., Feigelson, E. D., Grosso, N., McCaughrean, M. J., Micela, G.,
Broos, P., Garmire, G., & Townsley, L. 2005, ApJS, 160, 353

Kuhn, M. A., Getman, K. V., & Feigelson, E. D. 2015, ApJ, 802, 60

Kuhn, M. A., et al. 2014, ApJ, 787, 107

Mac Low, M.-M., & Klessen, R. S. 2004, Reviews of Modern Physics, 76, 125

Mamajek, E. E. 2009, American Institute of Physics Conference Series, 1158, 3

Megeath, S. T., et al. 2016, AJ, 151, 5

Megeath, S. T., et al. 2012, AJ, 144, 192

Owen, J. E., Ercolano, B., & Clarke, C. J. 2011, MNRAS, 412, 13

Pillitteri, I., et al. 2013, ApJ, 768, 99

Pravdo, S. H., Tsuboi, Y., Uzawa, A., & Ezoe, Y. 2009, ApJ, 704, 1495

Townsley, L. K., Feigelson, E. D., Montmerle, T., Broos, P. S., Chu, Y.-H., &
Garmire, G. P. 2003, ApJ, 593, 874

Tsujimoto, M., Feigelson, E. D., Grosso, N., Micela, G., Tsuboi, Y., Favata,
F., Shang, H., & Kastner, J. H. 2005, ApJS, 160, 503

Winston, E., et al. 2010, AJ, 140, 266

Wolk, S. J., Bourke, T. L., Smith, R. K., Spitzbart, B., & Alves, J. 2002,
ApJL, 580, L161

Wolk, S. J., & Walter, F. M. 1996, AJ, 111, 2066

Wright, N. J., Drake, J. J., Mamajek, E. E., & Henry, G. W. 2013,
Astronomische Nachrichten, 334, 151

Yadav, R. K., Christensen, U. R., Wolk, S. J., & Poppenhaeger, K. 2016, ApJL,
833, L28


# Instability and Mixing of Gas Interfaces Driven by Cylindrically Converging
Shock Wave

Wei-Gang Zeng Department of Engineering Mechanics, School of Aerospace,
Tsinghua University, Beijing 100084, China.  Yu-Xin Ren  Department of
Engineering Mechanics, School of Aerospace, Tsinghua University, Beijing
100084, China.  Jianhua Pan<EMAIL_ADDRESS>Department of Engineering
Mechanics, School of Aerospace, Tsinghua University, Beijing 100084, China.
William. G. Lowrie Department of Chemical and Biomolecular Engineering,
Koffolt Labs, The Ohio State University, US.

###### Abstract

In the present paper, an efficient method to generate "pure" cylindrically
converging shock wave without a following contact surface is proposed firstly.
Then, the Richtmyer-Meshkov instabilities of two interfaces driven by the
generated cylindrically converging shock wave and the associated fluids’
mixing behaviors are numerically studied. The results show that the
instability of the interface is characterized by the growth of perturbation
amplitude before re-shock. However, the mixing of fluids is enhanced
dramatically after re-shock, which is manifested not only by the evolutions of
flow structures but also by the temporal behaviors of mixing parameters.
Further investigation shows that, although these two cases are of different
initial perturbations, their evolutions of mixing width and other mixing
parameters such as molecular mixing fraction, local anisotropy and density-
specific volume correlation could achieve the same laws of temporal behavior,
especially during the later stage after re-shock. These results to some extent
demonstrate that there also exist scaling law and temporal asymptotic
behaviors in the mixing zone for cylindrically converging shock wave driven
interface. Moreover, the analyses of turbulent kinetic energy spectrums in the
azimuthal direction at late stage also witness the $k^{-5/3}$ decaying law of
turbulent kinetic energy for the present inhomogeneity flows driven by
cylindrically converging shock wave, which further manifests that the fluids’
mixing is indeed enhanced at later time after re-shock.

††preprint: AIP/123-QED

##  I  INTRODUCTION

The interaction between shock and gases interface with initial perturbations
is complicated, which involves interfacial instabilities of various kinds and
the possible turbulent mixing at late stage. At the beginning of interaction,
the Richtmyer-Meshkov instability (RMI) dominates the flow, resulting in the
deposition of baroclinic vorticities and the growth of the perturbation. Then,
the accumulated vortices cause the formation of the primary Kelvin–Helmholtz
(KH) billows, which would lead to the dramatic mixing of fluids when the more
complex instabilities occur in these billows.

In the last three decades, the instabilities of shock driven gases interface
and the associated fluids’ mixing behaviors have been widely studied due to
their importance in Inertial Confinement Fusion (ICF) Lindl (1995), supernova
explosions Arnett (2000), and supersonic combustion Yang, Kubota, and Zukoski
(1993). Many theoretical models are proposed for the vortex generation Picone
and Boris (1988); Yang, Kubota, and Zukoski (1994); Samtaney and Zabusky
(1994); Niederhaus  _et al._ (2008) and growth rate of perturbation amplitude
Zhang and Sohn (1997); Vandenboomgaerde, Gauthier, and Mügler (2002);
Matsuoka, Nishihara, and Fukuda (2003); Sohn (2003) of RMI flows, most of
which are manifested not only by numerical simulations Samtaney and Zabusky
(1994); Niederhaus  _et al._ (2008); Li, Wang, and Guan (2019); Wang  _et al._
(2020a); Sun  _et al._ (2020); Liang  _et al._ (2020) but also by experiments
Wang  _et al._ (2015); Liu  _et al._ (2018); Zhai, Ou, and Ding (2019).
Recently, the essential development of fluids’ mixing at late stage attracts
great attention Zhou (2017). Therefore, the morphological behaviors of
turbulent mixing and the criterion of turbulent mixing transition are
extensively studied during past ten years Hill, Pantano, and Pullin (2006);
Thornber  _et al._ (2010); Tritschler  _et al._ (2014); Lombardini, Pullin,
and Meiron (2012); Thornber  _et al._ (2017); Olson and Greenough (2014);
Thornber  _et al._ (2019); Zhou  _et al._ (2019); Li  _et al._ (2019).
However, most of these investigations are mainly focused on interface
instabilities and the associated fluids’ mixing induced by planar shock wave.
By contrast, the studies on interface instabilities and the associated fluids’
mixing induced by converging shock wave, which would be more relevant to
engineering applications such as ICF, are much less abundant. Obviously, these
research topics are gaining increasing attention recently. For example, recent
studies have manifested that, due to the Bell-Plesset effect and Rayleigh-
Taylor effect, the perturbation amplitude of interface driven by converging
shock wave would grow in a different way both at early and later stages
Mikaelian (2005); Ding  _et al._ (2017); Liu  _et al._ (2019); Luo  _et al._
(2019); Suponitsky, Froese, and Barsky (2014); Rafei  _et al._ (2019).
However, the fluids’ mixing behaviors of interface driven by converging shock
wave at late stage remain to be open issues.

The motivation of the present paper is to investigate the instability and the
associated fluids’ mixing of gases interfaces driven by cylindrically
converging shock wave (CCSW) using high resolution finite volume (FV) method.
To this end, we adopt an efficient way to generate “pure” CCSW which can be
used as the incident shock wave for further numerical study of interaction
between CCSW and gases interface. Then, the RMI and the associated fluids’
mixing behaviors of two interfaces driven by such CCSW are numerically
investigated in detail. The evolutions of flow structure and fluids’ mixing
behaviors highlight that the instability of the interface is characterized by
the growth of perturbation amplitude before re-shock, while the fluids’ mixing
is enhanced dramatically after the interface being re-shocked. The enhanced
fluids’ mixing is also manifested by the exponential scaling laws of mixing
width as well as the temporal asymptotic behaviors of mixing parameters such
as molecular mixing fraction, local anisotropy and density-specific volume
correlation at later stage after re-shock. Moreover, the turbulent kinetic
(TKE) in the azimuthal direction also decays with a slop of $k^{-5/3}$ in a
relatively broader range of low wave numbers at later stage after re-shock,
which further confirms the enhanced fluids’ mixing since the inertial range is
extended during the developing process.

The remainder of this paper is organized as follows. The numerical framework
based on high resolution FV method for compressible two fluids is presented in
Section II. The method to generate “pure” CCSW and the verification for its
usage in CCSW/interface interaction are introduced in Section III. The
numerical simulations for the CCSW induced RMI flows and the corresponding
results which include the wave patterns, the fluids’ mixing behaviors and the
decaying law of TKE spectrums are well discussed in Section IV. And finally,
the conclusions remarks are given in Section V.

##  II NUMERICAL FRAMEWORK

###  II.1 Governing equations

Following our previous workZeng  _et al._ (2018a, b); Wang  _et al._ (2020b),
the integral form of governing equations for compressible two fluids with
consistent treatment of the convective terms at a material interface can be
written as

$\frac{\partial}{\partial
t}\int{\bm{Q}}d{\Omega}=\oint_{\partial\Omega}(\bm{F}_{c}-\bm{F}_{v})\cdot\bm{n}dS=\int{\bm{W}}d{\Omega}\
.$ (1)

In the above formula, $\bm{Q}=\begin{bmatrix}\rho&\rho Y_{i}&\rho\bm{u}&\rho
E&\theta\end{bmatrix}^{T}$ is the vector of quasi-conservative variables in
the control volume $\Omega$, where $\rho$ is the density of mixture, $Y_{i}$
is the mass fraction of specie $i$, $\bm{u}=[u\ v\ w]^{T}$ is the vector flow
velocity, $E$ is the total energy of mixture and $\theta=\dfrac{1}{\gamma-1}$
is the function of the specific heat ratio of mixture, $\gamma$. It should be
noted that, due to the quasi-conservative form of $\theta$, a source term
$\bm{W}=\begin{bmatrix}0&0&0&0&\theta\nabla{\cdot{\bm{u}}}\end{bmatrix}^{T}$
should be added to the right-hand-side of Eq. (1). In fact, the equation for
$\theta$ is introduced to achieve a consistent treatment of the material
interfaces and to remove non-physical oscillations in the vicinity of the
material interfaces, which was initially proposed by Abgrall Abgrall (1996)
and lately improved by Johnsen Johnsen and Ham (2012) for FV method using high
order reconstructions. Additionally, $\bm{F}_{c}$ and $\bm{F}_{v}$ are,
respectively, the inviscid flux and viscous flux on the control surface $S$
with the unit outward normal vector $\bm{n}$. Their definitions are given as
follows

$\bm{F}_{c}=\begin{bmatrix}\rho\bm{u}\\\ \rho Y_{i}\bm{u}\\\
\rho\bm{u}\bm{u}+p\textbf{I}\\\ \rho H\bm{u}\\\ \theta\bm{u}\end{bmatrix},$
(2)

$\bm{F}_{v}=\begin{bmatrix}0\\\ \rho D_{i}\nabla{Y_{i}}\\\
\underline{\bm{\tau}}\\\
\underline{\bm{\tau}}\cdot\bm{u}-\bm{q}_{c}-\bm{q}_{d}\\\ 0\end{bmatrix}.$ (3)

In the above two equations, I is the unit tensor, $H=E+p/\rho$ is the total
enthalpy of the mixture,
$\underline{\bm{\tau}}=2\mu\underline{\bm{S}}-\dfrac{2}{3}\mu(\nabla{\cdot{\bm{u}}})\textbf{I}$
is the viscous stress tensor, and
$\underline{\bm{S}}=\dfrac{1}{2}(\nabla{\bm{u}}+(\nabla{\bm{u}})^{T})$ is the
strain rate tensor. The heat conduction and diffusion flux are given as
follows

$\bm{q}_{c}=-\kappa\nabla{T},$ (4) $\bm{q}_{d}=-\sum_{l=1}^{2}\rho
h_{l}(D_{l}\nabla{Y_{l}}-Y_{l}\sum_{m=1}^{2}D_{m}\nabla{Y_{m}}),$ (5)

where, $T$ is the static temperature of the fluids’ mixture. For species $i=l$
or $m$, $h_{i}$ and $D_{i}$ are the individual enthalpy and the effective
binary diffusion coefficient Ramshaw (1990), respectively. Additionally, the
models for dynamic viscosity coefficient $\mu$, thermal conductivity $\kappa$
and the effective binary diffusion coefficient for the mixture are well
documented by Tritschler Tritschler  _et al._ (2014) and Shanka Shankar and
Lele (2014), and one can also refer to our previous paper Zeng  _et al._
(2018b) for details.

To close the governing equations, the equation of state (EOS) for the mixture
of ideal gases is adopted in the present paper. Its formula is given by

$p=\rho\overline{R}T,$ (6)

where $\overline{R}=\dfrac{R_{u}}{\overline{M}}$ is the gas constant of the
mixture with $R_{u}$ being the universal gas constant and
$\overline{M}=1/\underset{i=1}{\overset{2}{\sum}}\dfrac{Y_{i}}{M_{i}}$ being
the mean molecular mass of the mixture. As proposed by Johnsen Johnsen and Ham
(2012), the internal energy is related to pressure in the following form:

$E=\dfrac{p}{(\gamma-1)\rho}+\dfrac{1}{2}\left|\bm{u}\right|^{2},$ (7)

where the specific heat ratio of mixture, $\gamma$, is given by

$\dfrac{1}{\gamma-1}=\sum_{i=1}^{2}\dfrac{Y_{i}\overline{M}}{(\gamma_{i}-1)M_{i}}.$
(8)

In the above formulas, $M_{i}$ and $\gamma_{i}$ are, respectively, the
molecular mass and specific heat ratio of species $i$ .

###  II.2 Numerical Method

In the framework of FV method, the semi-discretized form of Eq. (1) is used to
update the cell-averaged physical states at the cell center, which is in the
following form

$\dfrac{d\overline{Q}}{dt}=-\dfrac{1}{\Omega}\sum_{nf=1}^{6}(\bm{F}_{c,nf}^{*}-\bm{F}_{v,nf}^{*})\triangle{S}_{nf}+\overline{W}.$
(9)

In the above formula, $\overline{Q}$ is the average state of $Q$ at each cell
center, $\bm{F}_{c,nf}^{*}=\bm{F}_{c}\cdot\bm{n}_{nf}$ and
$\bm{F}_{v,nf}^{*}=\bm{F}_{v}\cdot\bm{n}_{nf}$ are, respectively, the
numerical convective flux and viscous flux at the cell inteface $nf$ with area
of $\triangle{S}_{nf}$ and unit outward normal vector of $\bm{n}_{nf}$, and
$\overline{W}=\begin{bmatrix}0&0&0&0&\underset{nf=1}{\overset{6}{\sum}}\overline{\theta}(\bm{u}_{nf}\cdot{\bm{n}_{nf}}\triangle{S}_{nf})\end{bmatrix}^{T}$
is the average source term of the cell with $\overline{\theta}$ being the
corresponding average of the function of specific heat ratio for the mixture
and $\bm{u}_{nf}$ being the fluids’ velocity at cell interface.

Combined with the fourth-order MDCD reconstruction proposed by Wang Wang  _et
al._ (2013), the Harten-Lax-van Leer-Contac (HLLC) Riemann solver for quasi-
conservative form of governing equations of multi-fluids, which is initially
proposed by Abgrall Abgrall (1996), is used to calculate the numerical flux of
convection, $\bm{F}_{c,nf}^{*}$. Additionally, Green’s theorem is used to
integrate the numerical viscous flux, $\bm{F}_{v,nf}^{*}$. Moreover, to
evaluate the average source term, the formula proposed by Johnse Johnsen and
Ham (2012) for the fluids’ velocity at cell interface, $\bm{u}_{nf}$, is used,
which is highly consistent with the HLLC Riemann solver for the numerical flux
of convection. Once all the terms in Eq.(9) are evaluated, we update the flow
states temporally using the third-order total variation diminishing (TVD)
Runge-Kutta method proposed by Shu Jiang and Shu (1996). In terms of the
accuracy of the present simulation code, it has already been well demonstrated
in our published workZeng  _et al._ (2018a, b); Wang  _et al._ (2020b).

##  III Generating CCSW and Verification of CCSW/interface interaction

###  III.1 Generating CCSW

Generating cylindrically converging shock wave or spherically converging shock
wave is somewhat complicated, no matter in experiment Zhai  _et al._ (2010,
2012); Eliasson  _et al._ (2006); Dimotakis and Samtaney (2006);
Vandenboomgaerde and Aymard (2011); Kjellander, Tillmark, and Apazidis (2012);
Liverts and Apazidis (2016) or in numerical simulation Lombardini, Pullin, and
Meiron (2014a); Bhagatwala and Lele (2012). Based on Guderley’s theory of
converging shock wave Guderley (1942), Lombardini and Pullin successfully set
up the initial conditions for numerical study on the turbulent mixing driven
by spherical implosions Lombardini, Pullin, and Meiron (2014a, b). However,
there are some defects for Guderley’s method in generating converging shock
wave since the condition for the validity of Guderley’s theory is that the
shock wave must be strong enough. Another efficient method of generating
converging shock wave is based on the theory of shock tube, which is initially
used by Bhagatwala and Lele Bhagatwala and Lele (2012) to generate the
spherically converging shock wave. In this method, one just needs to increase
the pressure and density ratios at a specific radius to generate the "pure"
converging shock wave without a following contact surface. In this paper, we
follow the approach of Bhagatwala and Lele Bhagatwala and Lele (2012) to
generate CCSW.

According to the theory of shock tube, the region of low pressure gas
(initially in static state with thermal states of pressure $p_{1}$, density
$\rho_{1}$, and the ratio of specific heats $\gamma$) and the region of high
pressure gas (initially in static state with thermal states of pressure
$p_{4}$, density $\rho_{4}$, and the same ratio of specific heats $\gamma$)
are initially separated by the diaphragm. As soon as the diaphragm is broken,
a CCSW will propagate inward into the region of low pressure gas, while an
expansion fan will propagate outward into the region of high pressure gas. To
generate a "pure" CCSW with desired Mach number $M_{s0}$ and no following
contact surface, the pressure $p_{2}$ and density $\rho_{2}$ of the gas in the
intermediate region (the region between the CCSW and the expansion fan) at
initial time should meet the following set of implicit formulations Bhagatwala
and Lele (2012)

$\dfrac{\rho_{2}}{\rho_{1}}=\dfrac{(\gamma+1)M_{s0}^{2}}{(\gamma-1)M_{s0}^{2}+2},$
(10) $\dfrac{p_{2}}{p_{1}}=1+\dfrac{2\gamma}{\gamma+1}(M_{s0}^{2}-1),$ (11)
$\dfrac{\rho_{4}}{\rho_{1}}=\dfrac{\rho_{2}}{\rho_{1}}\left\\{1-\dfrac{a_{1}}{a_{4}}(\gamma-1)(\dfrac{\rho_{2}}{\rho_{1}}-1)\sqrt{\dfrac{1}{2\dfrac{\rho_{2}}{\rho_{1}}\left[\gamma+1-(\gamma-1)\dfrac{\rho_{2}}{\rho_{1}}\right]}}\right\\}^{-\dfrac{2}{\gamma-1}},$
(12)
$\dfrac{p_{4}}{p_{1}}=\dfrac{p_{2}}{p_{1}}\left\\{1-\dfrac{a_{1}}{a_{4}}(\gamma-1)(\dfrac{p_{2}}{p_{1}}-1)\sqrt{\dfrac{1}{2\gamma\left[2\gamma+(\gamma+1)(\dfrac{p_{2}}{p_{1}}-1)\right]}}\right\\}^{-\dfrac{2\gamma}{\gamma-1}}.$
(13)

In the above equations, $\phi_{4}/\phi_{1}$ and $\phi_{2}/\phi_{1}$, where
$\phi\in\\{\rho,p\\}$, are jump ratios across the initial diaphragm and the
desired CCSW, respectively, and $a_{i}=\sqrt{\gamma p_{i}/\rho_{i}}$ (with
$i=1,2$ or $4$) is the sound speed in the corresponding region.

In order to demonstrate the efficiency of the above method in generating CCSW,
we set up a two-dimensional (2D) numerical simulation and verify its results
using Guderley’s theory Guderley (1942). The initial radius of circular shock
(corresponding to the diaphragm position) for this simulation is $R_{0}=$0.5
$m$. To generate a converging shock wave with initial Mach number of
$M_{s0}=1.5$, the thermal states of initially static gas (air with specific
heat ratio $\gamma=1.4$) inside the circular shock wave is set to be:
$\rho_{1}=1.225\ kg/m^{3}$, $p_{1}=101325\ pa$, $a_{1}=340.3\ m/s$. The
corresponding thermal states of initially static air outside the circular
shock wave, $\rho_{4},\ p_{4}$ and $a_{4}$, can be given by iteratively
solving Eq.(10)-Eq.(13). Once these initial conditions are given, a converging
shock wave and a expansion fan will propagate inward and outward respectively
from the initial position. The thermal states of the air in the intermediate
region at the very beginning time, $\rho_{2},\ p_{2}$ and $a_{2}$, can also be
given by Eq.(10)-Eq.(13).

Obviously, the strength of the converging shock wave will be enhanced during
its propagation since the area of shock surface is decreased. The Mach number
($M_{s}$) of the converging shock wave at any radius $r$ during its inward
propagation in our numerical simulation can be derived from the following
formulation Anderson (2003)

$\dfrac{p_{b}}{p_{1}}=1+\dfrac{2\gamma}{\gamma+1}(M_{s}^{2}-1),$ (14)

where, $p_{b}$ is the pressure behind the inward propagating CCSW at the
corresponding radius.

Additionally, based on the theory of Guderley Guderley (1942), during the
inward propagation for a pure converging shock wave initially placed at the
radius of $R_{0}$, its radius at given time $t$, $r(t)$, can be addressed as

$r(t)=R_{0}(1-\dfrac{t}{t_{0}})^{\alpha}.$ (15)

In the above equation, $t_{0}$ is the total propagation time from the initial
radius to the center for converging shock wave, which is about 0.865 $ms$ in
the present simulation. Additionally, for pure CCSW propagating in gas with
$\gamma=1.4$, the Guderley exponent $\alpha\approx 0.835$. Based on Eq.(15),
one can directly derive the theoretical formula for Mach number of converging
shock wave during its inward propagation Si  _et al._ (2015)

$M_{s}=\dfrac{1}{a_{1}}\dfrac{\alpha
R_{0}}{t_{0}}\left(\dfrac{r}{R_{0}}\right)^{\left(\alpha-1\right)/\alpha},$
(16)

where, $a_{1}$ is the sound speed of gas inside the CCSW.

Based on the above formulations, both $r(t)$ and $M_{s}$ would be sensitive to
the value of $\alpha$. Particularly, when the shock is propagating near to the
center, a marginal change of $\alpha$ would result in a dramatical variation
of $r(t)$ and $M_{s}$. In order to take the effect of $\alpha$ into account,
we additionally chose another two values of Guderley exponent
$\alpha^{\pm}=\alpha(1\pm 5\%)$, with $\alpha=0.835$, for our comparation
study. Fig.1 and Fig.2 compare the evolution of nondimensionalized radius
$ln\left[{r(t)}/{R_{0}}\right]$ versus nondimensionalized time
$ln\left(1-{t}/{t_{0}}\right)$ and the evolution of Mach number $M_{s}$ versus
$r$ of CCSW in our numerical simulation with the corresponding results
obtained from the theory of Guderley Guderley (1942), respectively. As shown
in Fig.1 and Fig.2, the evolutions of nondimensionalized radius and Mach
number in our simulation agree well with those of Guderley’s theory when the
Guderley exponent is in marginal range of $[\alpha^{-},\alpha^{+}]$, which
could demonstrate that the method used by Bhagatwala and Lele Bhagatwala and
Lele (2012) to generate the spherically converging shock wave is also
efficient in generating CCSW.

Figure 1:  Evolution of nondimensionalized radius
$ln\left[{r(t)}/{R_{0}}\right]$ versus nondimensionalized time
$ln\left(1-{t}/{t_{0}}\right)$. Figure 2:  Evolution of Mach number $M_{s}(r)$
versus $r$ .

###  III.2 Verification of CCSW/interface interaction: Studies on the
amplitude growth of CCSW driven interfaces

One of most important features for RMI flow is the growth of perturbation
amplitude. Due to the impulsive acceleration of incident shock wave, the
initial perturbation will be stretched, resulting in the possible linear
growth of perturbation amplitude at early stage. Initially studied by
Richtmyer Richtmyer (1960) and later confirmed experimentally by Meshkov
Meshkov (1969), the growth rate of perturbation amplitude for interface driven
by planar shock wave is extensively studied during past three decades Zhang
and Sohn (1997); Vandenboomgaerde, Gauthier, and Mügler (2002); Matsuoka,
Nishihara, and Fukuda (2003); Sohn (2003). However, the perturbation growth of
interface driven by converging shock wave remains to be an attractive topic
with many open issues. Based on the CCSW generated by the method mentioned
above, using 2D numerical simulations, we try to explore some features of
perturbation growth for CCSW driven interface before re-shock in this
subsection. In addition to numerical simulations, the theoretical model
proposed by Mikaelian Mikaelian (2005) is also used to predict the amplitude
growth of CCSW driven interface for further mutual confirmation.

In the present numerical studies of amplitude growth, the initial perturbation
of gases interface is in the following form

$\eta_{0}=r_{0}-a_{0}\left[1.0-cos(n\varphi)\right],$ (17)

where, $r_{0}$ is the initial radius of the outer interface (radius of crest
at initial time), $a_{0}$ is the initial amplitude, $n$ is the azimuthal mode
number, and $\varphi$ is the azimuthal angle. To make our results be more
general, we take the effects of initial mode and initial amplitude into
account. Therefore, three cases with different initial azimuthal mode number
and/or initial amplitude are studied. The parameters of initial perturbations
for the three cases are well listed in Table I.

Table 1:  Parameters of initial perturbations for study of amplitude growth Case Index | $r_{0}\ (m)$ | $a_{0}\ (cm)$ | n
---|---|---|---
I | 0.38 | $1.0$ | $8$
II | 0.38 | $2.0$ | $8$
III | 0.38 | $1.0$ | $12$

Moreover, the CCSW with desired Mach number of $M_{s0}=1.5$ is initially
placed at radius of $R_{0}=0.4m$. The gas in the region between gases
interface and the CCSW is air with initial thermal states of $\rho_{1}=1.225\
kg/m^{3}$, $p_{1}=101325pa$, $a_{1}=340.3m/s$ and $\gamma_{1}=1.4$. Once more,
the initial thermal states of gas (air) outside the CCSW is iteratively solved
based on Eq.(10)-Eq.(13). Additionally, the gas inside the initial interface
is a mixture of air ($20\%$ in mass fraction) and sulphur hexafluoride (SF6,
with specific heat ratio $\gamma_{2}=1.1$) under the condition that both
temperature and pressure are in equilibrium states at the gases interfaces.
The velocity for the gases in the whole computational domain is zero at the
initial time. It would take $0.04\ ms$ for the CCSW to initially strike the
gases interface when its Mach number approximately reaches 1.51. For boundary
conditions, a viscous circular wall with radius of $1cm$ is placed around the
center. Additionally, the circular boundary at the outside of the
computational domain is large enough for the propagations of all possible
waves during the durations of simulation. Consequently, a zero-gradient
boundary condition is used for the outer boundary. A body-fitted mesh is used
for all three simulations, with 4096 cells in azimuthal direction and 2440
cells in radial direction. The amplitude of interface in numerical simulations
is then given by
$a(t)\\_\textsc{CFD}\approx\left[r_{crest}(t)-r_{trough}(t)\right]/2$, where
$r_{crest}(t)=max\\{r_{Y_{SF_{6}}=0.5}\\}$ and
$r_{trough}(t)=min\\{r_{Y_{SF_{6}}=0.5}\\}$ are, respectively, the radius of
the outer interface (crest) and inner interface (trough) of simulation results
as shown in Fig.3.

Figure 3:  Crest radius and trough radius in simulations.

As an alternative to numerical simulations, linear models are also widely used
for predicting the perturbation development of CCSW driven gases interface
Mikaelian (2005); Ding  _et al._ (2017); Liu  _et al._ (2019); Luo  _et al._
(2019); Suponitsky, Froese, and Barsky (2014). Based on the pioneered work of
Bell Bell (1951) and Plesset Plesset (1954), Mikaelian modelled the amplitude
growth rate of interface (with small ratio of initial amplitude to the initial
wave length) driven by CCSW as follows Mikaelian (2005)

$\frac{d^{2}a(t)}{dt^{2}}=-2\frac{\dot{r}}{r}\frac{da(t)}{dt}+(nA-1)\frac{\ddot{r}}{r}a(t).$
(18)

In the above formula, $n$ is the mode number in azimuthal direction,
$A=(\rho_{in}-\rho_{out})/(\rho_{in}+\rho_{out})$ is the Atwood number with
$\rho_{in}$ and $\rho_{out}$ respectively being the density inside and outside
the already shocked interface. Additionally, $a(t)$ and $r$ are, respectively,
the amplitude and average radius of the gases interface at time $t$. There is
no analytical formula for $r(t)$, consequently, it is approximately given by
$r(t)=\left[r_{crest}(t)+r_{trough}(t)\right]/2$ under the assumption that the
perturbation amplitude of interface is small. Obviously, for $t>t_{0}^{+}$
(where $t_{0}^{+}$ is the time when the incident CCSW passes through the
trough of interface), Eq.(18) can be integrated in the following form Ding
_et al._ (2017); Liu  _et al._ (2019); Luo  _et al._ (2019)

$a(t)-a_{0}^{+}=\dot{a}_{0}^{+}{\overline{r}_{0}}^{2}\int_{t_{0}^{+}}^{t}\dfrac{1}{r^{2}}dt+(nA-1)\int_{t_{0}^{+}}^{t}\dfrac{\int_{t_{0}^{+}}^{t}a(\tau)r\ddot{r}d\tau}{r^{2}}dt,$
(19)

where $a_{0}^{+}\approx a_{0}(1-V_{si}/V_{is})$ and
$\dot{a}_{0}^{+}=a_{0}^{+}(nA-1)V_{si}/\overline{r}_{0}$ are, respectively,
the initial amplitude and growth rate of amplitude at time $t_{0}^{+}$, and
$\overline{r}_{0}\equiv r_{0}-a_{0}\approx r_{0}$ for small $a_{0}$.
Additionally, in the above formulations, $V_{is}$ is the velocity of CCSW when
it initially strikes the gases interface and $V_{si}$ is the velocity jump of
the shocked gases interface. There are two terms on the right hand side of
Eq.(19), of which the first term denotes the Bell-Plesset (BP) effect on the
growth of amplitude and the second term denotes Rayleigh-Taylor (RT) effect on
the growth of amplitude Ding  _et al._ (2017); Liu  _et al._ (2019); Luo  _et
al._ (2019). These two terms are corresponding to the effects of the first and
second term on the right hand side of Eq.(18), respectively.

For the given $a_{0}^{+}$ and $\dot{a}_{0}$ based on initial conditions, the
perturbation amplitude of Mikaelian’s model can be obtained by directly
integrating Eq.(18) using standard 4th order Runge-Kutta method. In order to
identify the BP effect and RT effect on the growth of perturbation more
clearly, we integrate Eq.(18) in two ways. In what follows, $a(t)\\_$BP stands
for the theoretical amplitude calculated by only integrating the first term of
Eq.(18) on the right hand side, while $a(t)\\_$BPRT stands for the theoretical
amplitude calculated by integrating both terms of Eq.(18) on the right hand
side.

Obviously, the BP effect is fully resulted from the geometric convergence
effects on flows, since this term will be zero for the planar shock driven
interface Mikaelian (2005). Additionally, as we will discuss below, Rayleigh-
Taylor effect will not be enhanced by the geometric convergence at early
stage. However, it, indeed, plays an important role in the growth of
perturbation amplitude at later stage. Fig.4 shows the evolutions of
perturbation amplitude obtained from numerical simulations as well as
theoretical model of Mikaelian for all three cases. As shown in Fig.4, the
results of Mikaelian’s theoretical model agree well with our numerical results
at early stage for all three cases (see Region I for each case). Additionally,
there are slight differences between $a(t)\\_$BP and $a(t)\\_$BPRT during this
stage, which means that BP effect is the dominant factor for the growth of
perturbation amplitude at early stage. However, at later stage (see Region II
for each case), the evolutions of $a(t)\\_$BP will increasingly diverge from
the numerical results, while the evolutions of $a(t)\\_$BPRT can still mimic
the numerical results although there are some differences between them. These
results indicate that the RT effect plays an important role in in the
development of perturbation at later stage before re-shock (even overwhelms
the BP effect at the very later stage before re-shock since the amplitudes for
all three cases are decreasing at the end of Region II).

(a) Case I

(b) Case II

(c) Case III

Figure 4:  Comparisons of amplitude growth between numerical results and
theoretical results

Moreover, according to the theoretical model of Mikaelian [see Eq.(19)], at
the early stage when the BP effect dominates the amplitude growth, the
nondimensionalized amplitude $a(t)/a_{0}$ can be approximately given by

$\dfrac{a(t)}{a_{0}}\approx(1-\dfrac{V_{si}}{V_{is}})[1+C_{1}(nA-1)V_{si}\overline{r}_{0}].$
(20)

In the above formula, $\overline{r}_{0}$ and
$C_{1}=\int_{t_{0}^{+}}^{t}\dfrac{1}{r^{2}}dt$ would become approximately
constant since the initial amplitude is small. For all three cases in the
present study, $V_{si}$ and $V_{is}$ are the same since the Mach number of the
incident CCSW and the initial densities inside and outside the interface are
the same. Consequently, during this stage, $a(t)/a_{0}$ should be
approximately proportional to the wave number of initial interface. Fig.5
compares the evolutions of nondimensionalized amplitude for three numerical
simulations (versus nondimensionalized time $t/t_{wall}$, where
$t_{wall}=1.06\ ms$ is the approximate time when the CCSW strikes the inner
wall boundary). As shown in Fig.5, due to the same wave numbers of initial
perturbations, the nondimensionalized amplitudes of Case I and Case II grow
almost in the same way at the early stage. On the other hand, due to a larger
wave number of initial perturbation, the nondimensionalized amplitude of Case
III is larger than those of Case I and Case II during the corresponding stage.

Figure 5:  Evolutions of nondimensionalized amplitude for all three cases in
numerical simulations.

According to the above analyses, we can see that, at early stage, the BP
effect is the dominant factor of the perturbation growth for CCSW driven
interfaces, and the amplitude growth obtained from the theoretical prediction
of Mikaelian’s model agrees well with our numerical simulation results.
However, at later stage, the RT effect becomes important and would even
overwhelm the BP effect. Moreover, during the later stage, although there are
some differences between the results of theoretical prediction (including RT
effect) and our numerical simulations, the evolutions of $a(t)\\_$BPRT can
mimic the trend of our numerical results on the whole. In fact, such
differences of amplitude growth have also been observed between the results of
Mikaelian’s theoretical model and shock tube experiments Ding  _et al._
(2017); Luo  _et al._ (2019). There are several factors which can account for
the above differences between the evolutions of $a(t)\\_$BPRT and the results
of numerical simulations/experiments. The main one is that the amplitude of
perturbation becomes larger at later stage as the flows evolve. Consequently,
the assumption of Mikaelian’s model, which requires that the amplitude of
perturbation should be small, is not quite valid anymore. Another one is that,
in practice, we can not have the analytical formula for $r(t)$. Consequently,
the approximate formulation of $r(t)$ could introduce errors, especially for
the term of RT effect since it involves the second-order derivative of $r(t)$.
In summary, all the results mentioned above, to some extent, not only could
demonstrate the features of amplitude growth of perturbation for CCSW driven
interfaces, but also can manifest that the CCSW derived from the
aforementioned theory of shock tube indeed can be used for the studies of CCSW
induced RMI flows.

##  IV FLUIDS’ MIXING OF CCSW DRIVEN INTERFACES

The fluids’ mixing in the mixing zone of gases interface is another crucial
topic for shock-driven inhomogeneous flows. Better understandings of the
mixing behaviors of shock-driven flows can shed light on the mechanisms of
turbulent mixing Zeng  _et al._ (2018b) as well as the turbulence modeling for
such kind of flows Livescu  _et al._ (2009); Cabot and Zhou (2013). Recently,
the mixing behaviors of planar shock driven interface are widely studied Hill,
Pantano, and Pullin (2006); Thornber  _et al._ (2010); Tritschler  _et al._
(2014); Lombardini, Pullin, and Meiron (2012); Thornber  _et al._ (2017);
Olson and Greenough (2014). However, the corresponding behaviors of gases
interface driven by converging shock wave, which are more important for some
scientific disciplines such as ICF and supernova explosions, remain to be
further investigated. In this section, we follow the approach of implicit
large eddy simulations (ILES) to have a primary study on the mixing behaviors
of three-dimensional (3D) gases interfaces driven by CCSW. Moreover, to take
the effects of mode of initial perturbation on the mixing behaviors into
account and to make our results be more general, two cases with different
amplitudes and mode numbers are studied.

###  IV.1 Problems setup

As remarked by Mikaelian Mikaelian (2005), the fluids’ mixing of gases
interface driven by converging shock wave would not happen dramatically if the
wave length of initial perturbation is much larger than the amplitude of
initial perturbation. In order to enhance the fluids’ mixing at later stage
after the gases interfaces are re-shocked, the initial perturbation for both
cases in the present study are set to be a linear combination of "egg-carton"
Hill, Pantano, and Pullin (2006), which, in the cylindrical coordinate system,
can be addressed as

$\eta(\varphi,\
z)=\dfrac{r_{0}-a_{0}\times\left|cos(n_{\varphi}\varphi)cos(n_{z}z)\right|+\delta-r}{2\delta}.$
(21)

In the above formula, $r_{0}$ and $a_{0}$ are, respectively, the crest radius
and the amplitude of the initial interface; $n_{\varphi}$ and $n_{z}$ are,
respectively, the azimuthal mode number and axial mode number;
$r=\sqrt{x^{2}+y^{2}}$ is the radius away from the center and $\varphi$ is the
azimuthal angle. Combined with an initial diffusion layer (with thickness
$\delta=0.01cm$) of the form proposed by Latini Latini, Schilling, and Don
(2007a), the mass fraction of SF6 inside the CCSW can be addressed as

$Y_{\text{SF}_{6}}=\left\\{\begin{matrix}Y_{\text{SF}_{6},0}&if\ \eta\geqslant
1.0&\\\ Y_{\text{SF}_{6},0}\times(1.0-e^{\left|\eta\right|^{8}ln\beta})&if\
0\leqslant\eta<1.0&\\\ 0&if\eta<0&\end{matrix},\right.$ (22)

where $\beta$ is the machine zero. According to the above formulation, at
initial time, the mass fraction of $\text{SF}_{6}$ for the mixture of
$\text{SF}_{6}$ and air inside the interface (the diffusion layer) is
$Y_{\text{SF}_{6},0}$ ($Y_{\text{SF}_{6},0}=0.75$ for both cases), while the
gas is pure air outside the interface. The parameters of the initial
perturbations for both cases are listed in Table II. Furthermore, we initially
place the CCSW with desired Mach number of $M_{s0}=1.5$ at radius of
$R_{0}=0.4\ m$ for both cases. The thermal states of air in the region between
the initial gases interface (diffusion layer) and the CCSW are the same as
that we set in the second part of Section III, which will result in the same
thermal states of air outside the initial CCSW. Moreover, for the present 3D
simulations, we extend the axial width, $L_{z}$, to be 0.128 $m$. The
configurations of initial flow field are shown in Fig.6. As for boundary
conditions, a viscous cylindrical wall with radius of $1\ cm$ is placed around
the center, and the periodic boundary condition is used along axial direction.
Additionally, in order to avoid the effects of waves reflected by the outer
boundary, the cylindrical outer boundary is far away from the flow structure
evolving region (the region with fine grid). To reduce computational costs, a
hyperbolic mesh stretching is applied between the fine-grid domain and the
outer boundary along the radial direction. The whole computational domain is
discretized by a body-fitted mesh with total cells of $1024\times 640\times
128$ (azimuthal cells $\times$ radial cells $\times$ axial cells).

Table 2:  Parameters of initial perturbations for CCSW induced interfacial fluids’ mixing Case Index | $r_{0}\ (m)$ | $a_{0}\ (mm)$ | $n_{\varphi}$ | $n_{z}$
---|---|---|---|---
1 | 0.38 | $5.0$ | $32$ | $8$
2 | 0.38 | $10.0$ | $16$ | $4$

(a) Case 1

(b) Case 2

Figure 6:  Configurations of initial flow field for both cases.

Due to the moderate grid resolution, we hereby remark that the numerical
dissipation would overweigh the physical one for the present simulations.
Consequently, the present studies can be categorized as a class of ILES, in
which the equations are implicitly filtered by the discretization and the
numerical dissipation is treated as a surrogate for an explicit subgrid-scale
model Schilling and Latini (2010). As remarked by Grinstein Grinstein,
Gowardhan, and Ristorcelli (2013) and Attal Attal and Ramaprabhu (2015),
although ILES could only resolve the length scales of turbulent mixing driven
by advection and convective stirring, many studies show that they are indeed
suited to (moderate) high-Reynolds-number flows in which shocks and interfaces
are present Thornber  _et al._ (2010); Gowardhan and Grinstein (2011); Zhou
_et al._ (2014); Thornber  _et al._ (2011, 2012). Moreover, our previous study
on reshocked heavy gas curtain shows that the fluids’ mixing behaviors of
coarse grid (corresponding to the results of ILES) agree well with those of
fine grid (corresponding to the results of direct numerical simulation)
statistically Zeng  _et al._ (2018b), which, to some extent, manifests that
our numerical method is appropriate for ILES.

###  IV.2 Results and Discussions

####  IV.2.1 Wave patterns and flow structures evolutions

The wave patterns of RMI induced by converging shock are more complicated than
those of RMI induced by planar shock, which would affect the resulting
evolutions of flow field to some extent. Fig.7 shows the typical wave patterns
of CCSW induced RMI (results of Case 1 in the axial view). At early stage
($t=0.240\ ms$) shown in Fig.7(a), the reflected expansion fan (REF), which is
associated with the incident CCSW and generated at the initial time, would
propagate outward all the time. Additionally, as the results of interaction
between the incident CCSW and initial interface, the reflected shock wave
(RSW) will always propagate outward while the first transmitted shock wave
(FTSW) will propagate inward initially. Obviously, the REF and the RSW would
not impose much effects on the later evolutions of the flows since they
propagate outward all the time and would not be reflected by the outer
boundary which is large enough. However, the FTSW will be reflected by the
inner wall boundary and re-impact the gases interface, which will result in
the second transmitted shock wave (STSW) and the reflected rarefaction wave
(RRW) at later time ($t=1.740\ ms$) as shown in Fig.7(b).

To identify the propagations of shock waves and their effects on the
evolutions of flow structures more clearly, Fig.8 and Fig.9 show the details
of flow evolutions and the propagations of shock waves for Case 1 and Case 2,
respectively. As shown in Fig.8(a) and Fig.9(a), at early stage ($t=0.840ms$),
the gases interfaces will move inward with growth of perturbation amplitudes
since the FTSWs will induce inward radial velocities. However, the inward
propagating FTSWs will be reflected by the inner wall boundary. Then, the
reflected FTSWs will propagate outward and begin to re-impact the gases
interfaces, which are well shown in Fig.8(b) and Fig.9(b). After the reflected
FTSWs re-impact the gases interfaces, the gases interfaces will move outward
since the resulting STSWs will propagate outward and induce outward main
radial velocities [see Fig.8(c)-(d) and Fig.9(c)-(d)]. Moreover, based on the
evolutions shown in Fig.8 and Fig.9, we can see that the flow structures are
characterized by the growth of perturbation amplitudes and the fluids’ mixing
is not intensive at early stage (Therefore, we just show one quarter of the
flow fields). However, the fluids’ mixing is dramatically enhanced by the
second RMI after the gases interfaces are re-shocked [see the SF6 mass
fraction iso-surface shown in Fig.8(c)-(d) and Fig.9(c)-(d)].

The morphological patterns of waves and the motions of flow structures (the
evolutions of positions of inner and outer interfaces) mentioned above for
both cases are depicted quantitatively in Fig.10. As shown in Fig.10, there
are some unique features for the CCSW induced RMI flows. One is that the FTSW
will move faster as it propagates inward due to the deformation (decreasing
area) of shock surface. This is quite different from the RMI flows induced by
planar shock wave, since, for planar shock driven RMI flows, the transmitted
shock wave will propagate forward with a nearly constant velocity Latini,
Schilling, and Don (2007b); Orlicz  _et al._ (2009). Another feature is that
the movements of inner and outer gases interfaces are nonlinear versus time
before re-shock, while, for planar shock driven RMI flows, the shocked
interface will move forward with an approximately constant velocity as well
Latini, Schilling, and Don (2007b); Orlicz  _et al._ (2009); Balakumar  _et
al._ (2008). Obviously, the nonlinear movements of inner and outer gases
interfaces before re-shock will result in the nonlinear growth of perturbation
amplitudes or mixing zone width at the very beginning (see Fig.11 shown
below). Moreover, we hereby remark that the model of Mikaelian mentioned in
the second part of Sections III would not be applicable to the nonlinear
growth of perturbation amplitudes for these two cases, since, for both cases,
their ratios of initial amplitude to initial wave length are not small enough.

(a) $t=$0.240 $ms$

(b) $t=$1.740 $ms$

Figure 7:  Typical wave patterns of the CCSW induced RMI.

(a) $t=$0.840 $ms$

(b) $t=$1.440 $ms$

(c) $t=$1.740 $ms$

(d) $t=$2.040 $ms$

Figure 8:  Evolutions of flow structures and shock waves propagations for Case
1.

(a) $t=$0.840 $ms$

(b) $t=$1.440 $ms$

(c) $t=$1.740 $ms$

(d) $t=$2.040 $ms$

Figure 9:  Evolutions of flow structures and shock waves propagations for Case
2.

(a) Case 1

(b) Case 2

Figure 10:  Morphological wave patterns and quantitative evolutions of flow
structures.

####  IV.2.2 Scaling law of mixing width

The self-similar scaling laws of mixing width for the RMI flow induced by
planar shock wave are extensively investigated Thornber  _et al._ (2010);
Tritschler  _et al._ (2014); Lombardini, Pullin, and Meiron (2012); Thornber
_et al._ (2017, 2011, 2012). For the planar shock induced RMI flows without
re-shock, Dimonte Dimonte, Frerking, and Schneider (1995) initially
demonstrated that the mixing width is an exponential function versus time
$\delta(t)\sim t^{\sigma}$. For the planar shock induced RMI flow with re-
shock, later studies of Thornber and Young Thornber  _et al._ (2011, 2012)
show that the mixing width after re-shock scales as
$\delta(t)\sim(t-t^{{}^{\prime\prime}}_{0})^{\sigma_{r}}$, where
$t^{{}^{\prime\prime}}_{0}$ is a virtual time and always set to be the time of
reshocking instant. Moreover, the exponents ${\sigma}$ and $\sigma_{r}$ for
the above scaling laws would be varied depending on some factors such as Mach
number of incident shock wave Lombardini, Pullin, and Meiron (2012). Recently,
the evolutions of mixing width for the RMI flows induced by converging shock
wave are paid more attention Lombardini, Pullin, and Meiron (2014a); Rafei
_et al._ (2019); Boureima, Ramaprabhu, and Attal (2018), while the
corresponding scaling laws are not fully reported and remain to be open
issues. Consequently, in this part, we try to gain further insights into the
scaling laws of mixing width for RMI flows induced by CCSW.

In the cylindrical coordinate system, the mixing width for RMI flows induced
by CCSW can be defined as

$\delta(t)=4\int_{r_{min}^{IMZ}}^{r_{max}^{IMZ}}<Y_{\text{SF}_{6}}>_{\varphi
z}(1-<Y_{\text{SF}_{6}}>_{\varphi z})dr.$ (23)

In the above formulation, $r_{min}^{IMZ}$ and $r_{max}^{IMZ}$ are,
respectively, the minimum radius and maximum radius of the inner mixing zone
(IMZ) in which the average mass fraction of sulphur hexafluoride
$<Y_{\text{SF}_{6}}>_{\varphi z}\in[0.05,0.70]$. Additionally, for arbitrary
scalar $\phi$, $<\phi>_{\varphi z}$ is the ensemble average of $\phi$ on the
cylindrical shell (in $\varphi z$ plane), which is defined as

$<\phi>_{\varphi z}(r,t)=\frac{1}{2\pi
L_{z}}\iint\phi(r,\varphi,z,t)d{\varphi}dz.$ (24)

The evolutions of mixing width versus time for both cases and their
corresponding evolutions on Log-Log scale are shown in Fig.11(a) and
Fig.11(b), respectively. As depicted in Fig.(9), the evolutions of mixing
width for both RMI flows induced by CCSW almost follow the same scaling laws.
At early stage, the scaling law of mixing width is $\delta(t)\sim t^{0.65}$,
which is quite similar to the scaling law of mixing width for planar shock
driven RMI flow with 3D broadband perturbations on an initial interface
(Thornber  _et al._ , 2010). Additionally, at the earlier stage after re-
shock, the mixing width for both cases scales as
$\delta(t)\sim(t-t^{{}^{\prime\prime}}_{0})^{0.25}$ (with re-shocking instant
$t^{{}^{\prime\prime}}_{0}\approx 1.5ms$). Actually, this scaling law of
mixing width for CCSW induced RMI flows is also widely reported for RMI flows
driven by planar shock during the earlier stage after re-shock Zeng  _et al._
(2018b); Tritschler  _et al._ (2014); Thornber  _et al._ (2011). However, at
the later stage after re-shock for present two CCSW driven RMI flows, the
mixing width seems to scale as $\delta(t)\sim t^{\sigma}$ again, while the
scaling exponent becomes ${\sigma}\approx 1.75$. This recovering scaling law
of $\delta(t)\sim t^{\sigma}$ at later stage after re-shock seems to be unique
for CCSW driven RMI flows since no other similar results have been reported
for planar shock driven RMI flows. It should be noted that, although we
witness the same scaling laws of mixing width for the present two CCSW induced
RMI flows, the scaling exponents at corresponding stages would be varied for
other cases since they would largely depend on the mode of initial
perturbations, impulsive Mach numbers and other factors Lombardini, Pullin,
and Meiron (2012).

(a) Normal evolutions versus time

(b) Evolutions versus time on a Log-Log scale

Figure 11:  Temporal evolutions of mixing width.

####  IV.2.3 Temporal asymptotic behaviors of mixing parameters

As mentioned above, although the flow structures are characterized by the
growth of perturbation amplitude before re-shock, the fluids’ mixing is
dramatically enhanced after re-shocked. To figure out the level of fluids’
mixing as well as the isotropy/homogeneity properties of the mixing zone, we
quantitatively investigate the temporal asymptotic behaviors of molecular
mixing fraction $\Theta$, local anisotropy $a_{i,vol}$ and density-specific
volume correlation $b_{vol}$ in this subsection.

The molecular mixing fraction can characterize the relative amount of
molecularly mixed fluid within the mixing layer. It can be interpreted as the
ratio of molecular mixing to large-scale entrainment by convection motion.
Following the definition of Youngs Youngs (1991, 1994), in the cylindrical
coordinate system, the formulation of molecular mixing fraction can be
expressed as

$\Theta(t)=\frac{\int_{r_{min}^{IMZ}}^{r_{max}^{IMZ}}<Y_{\text{SF}_{6}}(1-Y_{\text{SF}_{6}})>_{\varphi
z}dr}{\int_{r_{min}^{IMZ}}^{r_{max}^{IMZ}}<Y_{\text{SF}_{6}}>_{\varphi
z}(1-<Y_{\text{SF}_{6}}>_{\varphi z})dr}.$ (25)

The temporal asymptotic behaviors of the molecular mixing fraction for present
two CCSW induced RMI flows are shown in Fig.12. As shown in Fig.12, the ratios
of molecular mixing to large-scale entrainment by convection motion are
relatively small for both cases at the early stage before re-shock, while they
increase as the instabilities evolve. Additionally, the molecular mixing
between fluids for both cases is sharply increased, to some extent, after the
gases interfaces are re-shocked by the reflected FTSW (after $t\approx
1.5\times 10^{-3}s$). Moreover, at the later stage after re-shock, the
evolutions of molecular mixing fraction for both cases become asymptotic, with
a final value being 0.93 approximately. Actually, the above asymptotic
behavior of molecular mixing fraction for the present CCSW induced RMI flows
highly resembles that for some planar shock driven RMI flows which are well
reported numerically and experimentally Zeng  _et al._ (2018b); Orlicz,
Balasubramanian, and Prestridge (2013)

Figure 12:  Temporal evolutions of the molecular mixing fraction.

For mixing flows, the anisotropy and inhomogeneity are of significance since
both of them are important to large-eddy and Reynolds-averaged Navier-Stokes
modeling Cabot and Zhou (2013). To figure out the properties of local
anisotropy and inhomogeneity of fluids’ mixing for the present two CCSW driven
RMI flows, we investigate the temporal asymptotic behaviors of the volume-
averaged anisotropy $a_{i,vol}$ and the volume-averaged density-specific
volume correlation $b_{vol}$, respectively. Their formulations are given as
follows

$a_{i,vol}=\frac{1}{\Delta
r}\int_{r_{min}^{IMZ}}^{r_{max}^{IMZ}}<\frac{\left|u_{i}^{{}^{\prime\prime}}\right|}{\left|u_{r}^{{}^{\prime\prime}}\right|+\left|u_{\varphi}^{{}^{\prime\prime}}\right|+\left|u_{z}^{{}^{\prime\prime}}\right|}-\frac{1}{3}>_{\varphi
z}dr,$ (26) $b_{vol}=\frac{1}{\Delta
r}\int_{r_{min}^{IMZ}}^{r_{max}^{IMZ}}(<\frac{1}{\rho}>_{\varphi
z}<\rho>_{\varphi z}-1)dr.$ (27)

In the above two equations, $\Delta r=r_{max}^{IMZ}-r_{min}^{IMZ}$ is the
length of the inner mixing zone in the radial direction. Additionally, for
arbitrary scalar $\phi$, its fluctuating part, $\phi^{{}^{\prime\prime}}$, is
given by

$\phi^{{}^{\prime\prime}}=\phi-\overline{\phi},$ (28)

where $\overline{\phi}=<\rho\phi>_{\varphi z}/<\rho>_{\varphi z}$ is the
ensemble Favre average of $\phi$ on the cylindrical shell. Moreover, for
$i=r,\phi$ and $z$, $u_{i}^{{}^{\prime\prime}}$ respectively indicates the
fluctuating part of radial, azimuthal and axial velocity, and, $a_{i,vol}$
respectively denotes the corresponding volume-averaged anisotropy in radial,
azimuthal and axial direction.

Based on Eq.(26), the volume-averaged anisotropy $a_{i,vol}$ would range from
$-\frac{1}{3}$ to $\frac{2}{3}$, manifesting the ratio of the TKE in a
specific direction to the total TKE. Larger value of $a_{i,vol}$ implies
larger fluctuation or TKE in the corresponding direction. Moreover, the
positive value of $a_{i,vol}$ indicates that the TKE in the corresponding
direction is dominant, while the negative value of $a_{i,vol}$ implies that
the TKE in the corresponding direction is less important. For isotropically
mixing flow, $a_{i,vol}$ should be almost zero in all directions. Fig.13 shows
the temporal evolutions of the volume-averaged anisotropy in all three
directions for both cases. As shown in Fig.13, the volume-averaged anisotropy
in each direction would achieve a final asymptotic value for both case.
However, the temporal asymptotic behavior of $a_{i,vol}$ for the present CCSW
induced RMI flows is different from that for planar shock driven RMI flows.
For the planar shock driven RMI flows, the magnitude of $a_{i,vol}$ in all
directions would achieve a very small asymptotic value at later stage Zeng
_et al._ (2018b); Orlicz, Balasubramanian, and Prestridge (2013), implying
that the flows would become much less anisotropic. However, for the present
CCSW driven RMI flows during the duration of simulations, the temporal
asymptotic values in three directions follow the law of
$a_{r,vol}>0>a_{\varphi,vol}>a_{z,vol}$ and the magnitudes of $a_{r,vol}$ and
$a_{z,vol}$ are much larger than zero to some extent. These results indicate
that, during the duration of simulations, the fluids mixing would always be
anisotropic for the present CCSW driven RMI flows. Actually, the above results
are consistent with the corresponding evolutions of TKE. For CCSW driven RMI
flows, the TKE in the radial direction and the TKE in the azimuthal direction
are always more important than the TKE in the axial direction, simply because
the converging and expanding effects will continuously perturb the flow field
in these two directions.

Figure 13:  Temporal evolutions of the anisotropy.

The volume-averaged density-specific volume correlation is critical in second-
moment turbulence modeling for variable density flows Livescu  _et al._
(2009). According to Eq.(27), $b_{vol}$ is a non-negative parameter. For
nearly homogeneous flow, $b_{vol}$ would become very small. However, if the
fluids’ mixing is spatially inhomogeneous, the value of $b_{vol}$ would be
large. Fig.14 shows the temporal evolutions of the volume-averaged density-
specific volume correlation for the present two CCSW induced RMI flows. As
shown in Fig.14, the volume-averaged density-specific volume correlation for
both cases is decreasing on the whole. Moreover, at the later stage after re-
shock, the volume-averaged density-specific volume correlation for both cases
would asymptotically achieve the same relatively small value of 0.04. These
results imply that the fluids’ mixing for the present CCSW driven flows would
become much less inhomogeneous at later stage after re-shock.

Figure 14:  Temporal evolutions of the density-specific volume correlation.

####  IV.2.4 Turbulent kinetic energy spectrums

As mentioned by Tritschler Tritschler  _et al._ (2014), a fully isotropic
mixing zone is never obtained for shock induced fluids’ mixing flows, although
the fluid’s mixing would become less anisotropy and less inhomogeneous at
later stage. However, the theory of TKE spectrum for homogeneous isotropic
turbulence is often used as the theoretical framework for the numerical
analyses of shock induced RMI flows Tritschler  _et al._ (2014). According to
this theory, there would be a broadened inertial range in TKE spectrums once
the fluids’ mixing is enhanced. Therefore, alternative to the aforementioned
analyses in physical space, the enhanced fluids’ mixing after re-shock is
analyzed in Fourier representation for both cases in this subsection.
Additionally, due to the moderately small cell number in the axial direction,
the 2D TKE spectrums analyses in $\varphi z$ plane would conceal the
characteristics of TKE spectrums at relatively high wave numbers in the
azimuthal direction (which would be more important for the present CCSW
induced RMI flows since, as mentioned above, the TKE in the azimuthal
direction would be more important to some extent). Consequently, only the TKE
spectrums in the azimuthal direction after re-shock are analyzed for both
cases in this subsection.

In the cylindrical coordinate system, the average TKE spectrum for the inner
mixing zone in the azimuthal direction is given by

$E^{a}(k_{\varphi},t)=\frac{1}{\Delta
r}\int_{r_{min}^{IMZ}}^{r_{max}^{IMZ}}<E(r,k_{\varphi},z,t)>_{\varphi z}dr,$
(29)

where $E(r,k_{\varphi},z,t)$ is the TKE spectrum of wave number $k_{\varphi}$
in azimuthal direction at time $t$ on specific radial $r$ and axial position
$z$. The formulation of TKE spectrum is given by

$E(r,k_{\varphi},z,t)=\widehat{u_{r}^{{}^{\prime\prime}}}\widehat{u_{r}^{{}^{\prime\prime}}}^{*}+\widehat{u_{\varphi}^{{}^{\prime\prime}}}\widehat{u_{\varphi}^{{}^{\prime\prime}}}^{*}+\widehat{u_{z}^{{}^{\prime\prime}}}\widehat{u_{z}^{{}^{\prime\prime}}}^{*}.$
(30)

In the above equation, for arbitrary scalar $\phi$, $\widehat{\phi}$ denotes
its Fourier transform in the azimuthal direction, and $\widehat{\phi}^{*}$
indicates the corresponding complex conjugate of $\widehat{\phi}$. Fig.15
shows the average TKE spectrums at three instants of the later stage after re-
shock for the present two CCSW driven RMI flows. As depicted in Fig.15, at the
relatively early stage after re-shock $(t=2.440ms)$, the decaying law of
$k^{-5/3}$ for the TKE spectrums of these two cases only locates in a
narrowband of wave numbers. However, as the flows evolve, the TKE spectrums
for both cases would decay with a slope of $k^{-5/3}$ in a broader range of
wave numbers. These results imply that the inertial range is extended during
the developing process and manifest, to some extent, that the fluids’ mixing
is enhanced at later time after re-shock.

(a) Case 1

(b) Case 2

Figure 15:  Average TKE spectrums in the azimuthal direction.

##  V CONCLUSIONS

Based on the theory of shock tube, we successfully generate "pure" CCSW
without a following contact surface. Additionally, studies on the amplitude
growth of gases interfaces manifest that the generated CCSW is efficient for
studying on CCSW induced RMI flows. Then, the instabilities and fluids’ mixing
behaviors of two gases interfaces driven by CCSW are numerically investigated
using high resolution FV method. The morphological wave patterns and the
evolutions of flow structures imply that the instabilities of the interfaces
are characterized by the growth of perturbation amplitude before re-shock,
while the fluids’ mixing is dramatically enhanced after re-shock. Detailed
analyses of the fluids’ mixing parameters show that the evolutions of mixing
width and other mixing parameters could achieve the same laws of temporal
behavior for the present two CCSW induced RMI flows, which indicates the
existences of scaling law and temporal asymptotic behaviors for fluids’ mixing
parameters in the mixing zone. Additionally, due to the converging/expanding
effects on the flow fields, the motions of shock wave and inner/outer gases
interfaces for the present CCSW driven RMI flows are nonlinear versus time,
which is quite different from the results of planar shock wave induced RMI
flows. These nonlinear developments of flow fields would lead to some unique
features for the fluids’ mixing of CCSW driven RMI flows, one of which is
that, at later stage after re-shock, the fluids’ mixing would be less
isotropic than that of planar shock wave induced RMI flows although both of
them would reach final temporal asymptotic behaviors in each direction.
Further analyses of TKE spectrums in the azimuthal direction at later stage
after re-shock also witness the $k^{-5/3}$ decaying law of TKE spectrums for
the present CCSW driven RMI flows. Both the temporal behaviors of mixing
parameters and the decaying law of TKE spectrums manifest that the fluids’
mixing is indeed enhanced at later time after re-shock.

## AVAILABILITY OF DATA

Raw data were generated at the TH-2 Supercomputer. Derived data supporting the
findings of this study are available from the corresponding author upon
reasonable request.

###### Acknowledgements.

The work of this paper is supported by 2016YFA0401200 of the National Key
Research and Development Program of China and U1430235 of the National Natural
Science Foundation of China.

## References

  * Lindl (1995) J. Lindl, “Development of the indirect-drive approach to inertial confinement fusion and the target physics basis for ignition and gain,” Physics of plasmas 2, 3933–4024 (1995).
  * Arnett (2000) D. Arnett, “The role of mixing in astrophysics,” The Astrophysical Journal Supplement Series 127, 213 (2000).
  * Yang, Kubota, and Zukoski (1993) J. Yang, T. Kubota, and E. E. Zukoski, “Applications of shock-induced mixing to supersonic combustion,” AIAA journal 31, 854–862 (1993).
  * Picone and Boris (1988) J. M. Picone and J. P. Boris, “Vorticity generation by shock propagation through bubbles in a gas,” Journal of Fluid Mechanics 189, 23–51 (1988).
  * Yang, Kubota, and Zukoski (1994) J. Yang, T. Kubota, and E. E. Zukoski, “A model for characterization of a vortex pair formed by shock passage over a light-gas inhomogeneity,” Journal of Fluid Mechanics 258, 217–244 (1994).
  * Samtaney and Zabusky (1994) R. Samtaney and N. J. Zabusky, “Circulation deposition on shock-accelerated planar and curved density-stratified interfaces: models and scaling laws,” Journal of Fluid Mechanics 269, 45–78 (1994).
  * Niederhaus  _et al._ (2008) J. H. J. Niederhaus, J. A. Greenough, J. G. Oakley, D. Ranjan, M. H. Anderson, and R. Bonazza, “A computational parameter study for the three-dimensional shock-bubble interaction,” Journal of Fluid Mechanics 594, 85–124 (2008).
  * Zhang and Sohn (1997) Q. Zhang and S. I. Sohn, “Nonlinear theory of unstable fluid mixing driven by shock wave,” Physics of Fluids 9, 1106–1124 (1997).
  * Vandenboomgaerde, Gauthier, and Mügler (2002) M. Vandenboomgaerde, S. Gauthier, and C. Mügler, “Nonlinear regime of a multimode Richtmyer-Meshkov instability: A simplified perturbation theory,” Physics of Fluids 14, 1111–1122 (2002).
  * Matsuoka, Nishihara, and Fukuda (2003) C. Matsuoka, K. Nishihara, and Y. Fukuda, “Nonlinear evolution of an interface in the Richtmyer-Meshkov instability,” Phys Rev E 67, 036301 (2003).
  * Sohn (2003) S. I. Sohn, “Simple potential-flow model of Rayleigh-Taylor and Richtmyer-Meshkov instabilities for all density ratios,” Physical Review E 67, 026301 (2003).
  * Li, Wang, and Guan (2019) D. Li, G. Wang, and B. Guan, “On the circulation prediction of shock-accelerated elliptical heavy gas cylinders,” Physics of Fluids 31, 056104 (2019).
  * Wang  _et al._ (2020a) G. Wang, Y. Wang, D. Li, and B. Guan, “Numerical study on shock-accelerated gas rings,” Physics of Fluids 32, 026102 (2020a).
  * Sun  _et al._ (2020) P. Sun, J. Ding, S. Huang, X. Luo, and W. Cheng, “Microscopic Richtmyer-Meshkov instability under strong shock,” Physics of Fluids 32, 024109 (2020).
  * Liang  _et al._ (2020) Y. Liang, Z. Zhai, X. Luo, and C.-y. Wen, “Interfacial instability at a heavy/light interface induced by rarefaction waves,” Journal of Fluid Mechanics 885 (2020).
  * Wang  _et al._ (2015) X. Wang, D. Yang, J. Wu, and X. Luo, “Interaction of a weak shock wave with a discontinuous heavy-gas cylinder,” Physics of Fluids 27, 064104 (2015).
  * Liu  _et al._ (2018) L. Liu, Y. Liang, J. Ding, N. Liu, and X. Luo, “An elaborate experiment on the single-mode Richtmyer-Meshkov instability,” Journal of Fluid Mechanics 853 (2018).
  * Zhai, Ou, and Ding (2019) Z. Zhai, J. Ou, and J. Ding, “Coupling effect on shocked double-gas cylinder evolution,” Physics of Fluids 31, 096104 (2019).
  * Zhou (2017) Y. Zhou, “Rayleigh-Taylor and Richtmyer-Meshkov instability induced flow, turbulence, and mixing. II,” Physics Reports 723, 1–160 (2017).
  * Hill, Pantano, and Pullin (2006) D. J. Hill, C. Pantano, and D. I. Pullin, “Large-eddy simulation and multiscale modelling of a Richtmyer-Meshkov instability with reshock,” Journal of Fluid Mechanics 557, 29–61 (2006).
  * Thornber  _et al._ (2010) B. Thornber, D. Drikakis, D. L. Youngs, and R. J. R. Williams, “The influence of initial conditions on turbulent mixing due to Richtmyer-Meshkov instability,” Journal of Fluid Mechanics 654, 99–139 (2010).
  * Tritschler  _et al._ (2014) Tritschler, V.K., Olson, B.J., Lele, S.K., Hickel, S., Hu, and X. and, “On the Richtmyer-Meshkov instability evolving from a deterministic multimode planar interface,” Journal of Fluid Mechanics 755, , 429–462 (2014).
  * Lombardini, Pullin, and Meiron (2012) M. Lombardini, D. I. Pullin, and D. I. Meiron, “Transition to turbulence in shock-driven mixing: a mach number study,” Journal of Fluid Mechanics 690, 203–226 (2012).
  * Thornber  _et al._ (2017) B. Thornber, J. Griffond, O. Poujade, N. Attal, H. Varshochi, P. Bigdelou, P. Ramaprabhu, B. Olson, J. Greenough, and Y. Zhou, “Late-time growth rate, mixing, and anisotropy in the multimode narrowband Richtmyer-Meshkov instability: The $\theta$-group collaboration,” Physics of Fluids 29, 105107 (2017).
  * Olson and Greenough (2014) B. J. Olson and J. A. Greenough, “Comparison of two- and three-dimensional simulations of miscible Richtmyer-Meshkov instability with multimode initial conditions,” Physics of Fluids 26, 101702 (2014).
  * Thornber  _et al._ (2019) B. Thornber, J. Griffond, P. Bigdelou, I. Boureima, P. Ramaprabhu, O. Schilling, and R. Williams, “Turbulent transport and mixing in the multimode narrowband Richtmyer-Meshkov instability,” Physics of Fluids 31, 096105 (2019).
  * Zhou  _et al._ (2019) Y. Zhou, T. T. Clark, D. S. Clark, S. Gail Glendinning, M. Aaron Skinner, C. M. Huntington, O. A. Hurricane, A. M. Dimits, and B. A. Remington, “Turbulent mixing and transition criteria of flows induced by hydrodynamic instabilities,” Physics of Plasmas 26, 080901 (2019).
  * Li  _et al._ (2019) H. Li, Z. He, Y. Zhang, and B. Tian, “On the role of rarefaction/compression waves in Richtmyer-Meshkov instability with reshock,” Physics of Fluids 31, 054102 (2019).
  * Mikaelian (2005) K. O. Mikaelian, “Rayleigh-Taylor and Richtmyer-Meshkov instabilities and mixing in stratified cylindrical shells,” Physics of Fluids 17, 094105 (2005).
  * Ding  _et al._ (2017) J. Ding, T. Si, J. Yang, X. Lu, Z. Zhai, and X. Luo, “Measurement of a Richtmyer-Meshkov instability at an air-SF6 interface in a semiannular shock tube,” Physical review letters 119, 014501 (2017).
  * Liu  _et al._ (2019) L. Liu, J. Ding, Z. Zhai, and X. Luo, “Richtmyer-Meshkov instability of a sinusoidal interface driven by a cylindrical shock,” Shock Waves 29, 263–271 (2019).
  * Luo  _et al._ (2019) X. Luo, M. Li, J. Ding, Z. Zhai, and T. Si, “Nonlinear behaviour of convergent Richtmyer-Meshkov instability,” Journal of Fluid Mechanics 877, 130–141 (2019).
  * Suponitsky, Froese, and Barsky (2014) V. Suponitsky, A. Froese, and S. Barsky, “Richtmyer-Meshkov instability of a liquid-gas interface driven by a cylindrical imploding pressure wave,” Computers & Fluids 89, 1–19 (2014).
  * Rafei  _et al._ (2019) M. E. Rafei, M. Flaig, D. Youngs, and B. Thornber, “Three-dimensional simulations of turbulent mixing in spherical implosions,” Physics of Fluids 31, 114101 (2019).
  * Zeng  _et al._ (2018a) W.-G. Zeng, J.-H. Pan, Y.-X. Ren, and Y.-T. Sun, “Numerical study on the turbulent mixing of planar shock-accelerated triangular heavy gases interface,” Acta Mechanica Sinica 34, 855–870 (2018a).
  * Zeng  _et al._ (2018b) W. G. Zeng, J. H. Pan, Y. T. Sun, and Y. X. Ren, “Turbulent mixing and energy transfer of reshocked heavy gas curtain,” Physics of Fluids 30, 064106 (2018b).
  * Wang  _et al._ (2020b) Q. Wang, R. Deiterding, J. Pan, and Y.-X. Ren, “Consistent high resolution interface-capturing finite volume method for compressible multi-material flows,” Computers & Fluids , 104518 (2020b).
  * Abgrall (1996) R. Abgrall, “How to prevent pressure oscillations in multicomponent flow calculations: a quasi conservative approach,” Journal of Computational Physics 125, 150–160 (1996).
  * Johnsen and Ham (2012) E. Johnsen and F. Ham, “Preventing numerical errors generated by interface-capturing schemes in compressible multi-material flows,” Journal of Computational Physics 231, 5705–5717 (2012).
  * Ramshaw (1990) J. D. Ramshaw, “Self-consistent effective binary diffusion in multicomponent gas mixtures,” Journal of Non-Equilibrium Thermodynamics 15, 295–300 (1990).
  * Shankar and Lele (2014) S. K. Shankar and S. K. Lele, “Numerical investigation of turbulence in reshocked Richtmyer-Meshkov unstable curtain of dense gas,” Shock Waves 24, 79–95 (2014).
  * Wang  _et al._ (2013) Q. Wang, Y. Ren, Z. Sun, and Y. Sun, “Low dispersion finite volume scheme based on reconstruction with minimized dispersion and controllable dissipation,” Science China Physics, Mechanics and Astronomy 56, 423–431 (2013).
  * Jiang and Shu (1996) G.-S. Jiang and C.-W. Shu, “Efficient implementation of weighted eno schemes,” Journal of Computational Physics 126, 202–228 (1996).
  * Zhai  _et al._ (2010) Z. Zhai, C. Liu, F. Qin, J. Yang, and X. Luo, “Generation of cylindrical converging shock waves based on shock dynamics theory,” Physics of Fluids 22, 041701 (2010).
  * Zhai  _et al._ (2012) Z. Zhai, T. Si, X. Luo, J. Yang, C. Liu, D. Tan, and L. Zou, “Parametric study of cylindrical converging shock waves generated based on shock dynamics theory,” Physics of Fluids 24, 026101 (2012).
  * Eliasson  _et al._ (2006) V. Eliasson, N. Apazidis, N. Tillmark, and M. Lesser, “Focusing of strong shocks in an annular shock tube,” Shock waves 15, 205–217 (2006).
  * Dimotakis and Samtaney (2006) P. Dimotakis and R. Samtaney, “Planar shock cylindrical focusing by a perfect-gas lens,” Physics of Fluids 18, 031705 (2006).
  * Vandenboomgaerde and Aymard (2011) M. Vandenboomgaerde and C. Aymard, “Analytical theory for planar shock focusing through perfect gas lens and shock tube experiment designs,” Physics of Fluids 23, 016101 (2011).
  * Kjellander, Tillmark, and Apazidis (2012) M. Kjellander, N. Tillmark, and N. Apazidis, “Energy concentration by spherical converging shocks generated in a shock tube,” Physics of fluids 24, 126103 (2012).
  * Liverts and Apazidis (2016) M. Liverts and N. Apazidis, “Limiting temperatures of spherical shock wave implosion,” Physical review letters 116, 014501 (2016).
  * Lombardini, Pullin, and Meiron (2014a) M. Lombardini, D. I. Pullin, and D. I. Meiron, “Turbulent mixing driven by spherical implosions. Part 1. Flow description and mixing-layer growth,” Journal of Fluid Mechanics 748, 85–112 (2014a).
  * Bhagatwala and Lele (2012) A. Bhagatwala and S. K. Lele, “Interaction of a converging spherical shock wave with isotropic turbulence,” Physics of Fluids 24, 035103 (2012).
  * Guderley (1942) K. Guderley, “Starke kugelige und zylindrische verdichtungsstosse in der nahe des kugelmitterpunktes bnw. der zylinderachse,” Luftfahrtforschung 19, 302 (1942).
  * Lombardini, Pullin, and Meiron (2014b) M. Lombardini, D. I. Pullin, and D. I. Meiron, “Turbulent mixing driven by spherical implosions. Part 2. Flow description and mixing-layer growth,” Journal of Fluid Mechanics 748, 113–142 (2014b).
  * Anderson (2003) J. D. Anderson,  _Modern compressible flow: with historical perspective_ (McGraw-Hill, 2003) p. 90.
  * Si  _et al._ (2015) T. Si, T. Long, Z. Zhai, and X. Luo, “Experimental investigation of cylindrical converging shock waves interacting with a polygonal heavy gas cylinder,” Journal of Fluid Mechanics 784, 225–251 (2015).
  * Richtmyer (1960) R. D. Richtmyer, “Taylor instability in shock acceleration of compressible fluids,” Communications on pure and applied mathematics 13, 297–319 (1960).
  * Meshkov (1969) E. Meshkov, “Instability of the interface of two gases accelerated by a shock wave,” Fluid Dynamics 4, 101–104 (1969).
  * Bell (1951) G. I. Bell, “Taylor instability on cylinders and spheres in the small amplitude approximation,” Los Alamos Scientific Laboratory, Los Alamos, NM, Technical Report No. LA-1321 (1951).
  * Plesset (1954) M. Plesset, “On the stability of fluid flows with spherical symmetry,” Journal of Applied Physics 25, 96–98 (1954).
  * Livescu  _et al._ (2009) D. Livescu, J. Ristorcelli, R. Gore, S. Dean, W. Cabot, and A. Cook, “High-Reynolds number Rayleigh-Taylor turbulence,” Journal of Turbulence , N13 (2009).
  * Cabot and Zhou (2013) W. Cabot and Y. Zhou, “Statistical measurements of scaling and anisotropy of turbulent flows induced by Rayleigh-Taylor instability,” Physics of Fluids 25, 015107 (2013).
  * Latini, Schilling, and Don (2007a) M. Latini, O. Schilling, and W. S. Don, “Effects of WENO flux reconstruction order and spatial resolution on reshocked two-dimensional Richtmyer-Meshkov instability,” Journal of Computational Physics 221, 805–836 (2007a).
  * Schilling and Latini (2010) O. Schilling and M. Latini, “High-order WENO simulations of three-dimensional reshocked Richtmyer-Meshkov instability to late times: dynamics, dependence on initial conditions, and comparisons to experimental data,” Acta Mathematica Scientia 30, 595–620 (2010).
  * Grinstein, Gowardhan, and Ristorcelli (2013) F. Grinstein, A. Gowardhan, and J. Ristorcelli, “Implicit large eddy simulation of shock-driven material mixing,” Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 371, 20120217 (2013).
  * Attal and Ramaprabhu (2015) N. Attal and P. Ramaprabhu, “Numerical investigation of a single-mode chemically reacting Richtmyer-Meshkov instability,” Shock Waves 25, 307–328 (2015).
  * Gowardhan and Grinstein (2011) A. Gowardhan and F. Grinstein, “Numerical simulation of Richtmyer-Meshkov instabilities in shocked gas curtains,” Journal of Turbulence , N43 (2011).
  * Zhou  _et al._ (2014) Y. Zhou, F. F. Grinstein, A. J. Wachtor, and B. M. Haines, “Estimating the effective reynolds number in implicit large-eddy simulation,” Physical Review E 89, 013303 (2014).
  * Thornber  _et al._ (2011) B. Thornber, D. Drikakis, D. Youngs, and R. Williams, “Growth of a Richtmyer-Meshkov turbulent layer after reshock,” Physics of Fluids 23, 095107 (2011).
  * Thornber  _et al._ (2012) B. Thornber, D. Drikakis, D. Youngs, and R. Williams, “Physics of the single-shocked and reshocked Richtmyer-Meshkov instability,” Journal of Turbulence , N10 (2012).
  * Latini, Schilling, and Don (2007b) M. Latini, O. Schilling, and W. S. Don, “High-resolution simulations and modeling of reshocked single-mode Richtmyer-Meshkov instability: Comparison to experimental data and to amplitude growth model predictions,” Physics of Fluids 19, 024104 (2007b).
  * Orlicz  _et al._ (2009) G. Orlicz, B. Balakumar, C. Tomkins, and K. Prestridge, “A Mach number study of the Richtmyer-Meshkov instability in a varicose, heavy-gas curtain,” Physics of Fluids 21, 064102 (2009).
  * Balakumar  _et al._ (2008) B. Balakumar, G. Orlicz, C. Tomkins, and K. Prestridge, “Simultaneous particle-image velocimetry-planar laser-induced fluorescence measurements of Richtmyer-Meshkov instability growth in a gas curtain with and without reshock,” Physics of Fluids 20, 124103 (2008).
  * Dimonte, Frerking, and Schneider (1995) G. Dimonte, C. E. Frerking, and M. Schneider, “Richtmyer-Meshkov instability in the turbulent regime,” Physical review letters 74, 4855 (1995).
  * Boureima, Ramaprabhu, and Attal (2018) I. Boureima, P. Ramaprabhu, and N. Attal, “Properties of the turbulent mixing layer in a spherical implosion,” Journal of Fluids Engineering 140 (2018).
  * Youngs (1991) D. L. Youngs, “Three-dimensional numerical simulation of turbulent mixing by Rayleigh-Taylor instability,” Physics of Fluids A: Fluid Dynamics 3, 1312–1320 (1991).
  * Youngs (1994) D. L. Youngs, “Numerical simulation of mixing by Rayleigh-Taylor and Richtmyer-Meshkov instabilities,” Laser and particle beams 12, 725–750 (1994).
  * Orlicz, Balasubramanian, and Prestridge (2013) G. Orlicz, S. Balasubramanian, and K. Prestridge, “Incident shock Mach number effects on Richtmyer-Meshkov mixing in a heavy gas layer,” Physics of Fluids 25, 114101 (2013).