**Rayleigh-Taylor and Richtmyer-Meshkov instabilities for fluids with
finite density differences**

S.I. Abarzhi ^+\ ,\*^, J. Glimm^\*\*^, K. Nishihara^\*^

^+^Center for Turbulence Research, Stanford University, CA 94305, USA;
^\*^Institute of Laser Engineering, Osaka University, 2-6 Yamada-oka,
Suita, Osaka 565-0871, Japan; ^\*\*^Department of Applied Mathematics,
State University of New York at Stony Brook and Brookhaven National
Laboratory, Upton NY 11973, USA

**Abstract.** We report nonlinear solutions describing the large-scale
coherent motion of bubbles and spikes in the Rayleigh-Taylor and
Richtmyer-Meshkov instability for fluids with finite density
differences. The non-local character of the interface dynamics is taken
into account with a multiple harmonic analysis. The theory yields a
non-trivial dependence of the bubble velocity and curvature on the
density ratio and reveals an important qualitative distinction between
the dynamics of Rayleigh-Taylor and Richtmyer-Meshkov bubbles.

## PASC N 47.20.-k, 47.20.-i, 47.54.+r, 52.35.Py, 52.57.-z

When a light fluid accelerates a heavy fluid, the misalignment of the
pressure and density gradients gives rise to an instability of the
interface, and produces eventually the turbulent mixing of the fluids,
\[1-3\]. This phenomenon is called the Rayleigh-Taylor instability (RTI)
if the acceleration is sustained \[1\] and the Richtmyer-Meshkov
instability (RMI) if the acceleration is driven by a shock or if it is
impulsive \[2\]. The RT/RM turbulent mixing is of extreme importance in
astrophysics, inertial confinement fusion, and many other applications
\[3\]. To obtain a reliable description of the mixing process, the
evolution of a large-scale coherent structure, the dynamics of
small-scale structures, and the cascades of energy should be understood.

The large-scale coherent structure is an array of bubbles and spikes
periodic in the plane normal to the direction of gravity or the initial
shock \[1,2,4\]. It appears in the nonlinear regime of RTI and RMI and
has a spatial period determined by the mode of fastest growth \[5\]. The
light (heavy) fluid penetrates the heavy (light) fluid in bubbles
(spikes). The density ratio is a determining factor of the instability
dynamics \[3,4\]. Singular aspects of the interface evolution (such as
the generation of vorticity and secondary instabilities, resulting in
the direct and inverse cascades of energy \[6-12\], a finite contrast of
the fluid densities and the non-linearity of the dynamics) cause
theoretical difficulties \[13-18\] and preclude elementary methods of
solution.

For fluids with highly contrasting densities (fluid-vacuum), the effect
of singularities on the interplay of harmonics and on the nonlinear
motion in RTI/RMI have been studied intensively over the decades
\[13-16\]; a new approach based on group theory has been developed
recently in \[15,16\]. The asymptotic theories \[13-16\] agreed with
experiments and simulations. For fluids with a finite density contrast,
the influence of singularities on the cascades of energy and the
large-scale coherent dynamics in RTI/RMI has yet to be elucidated
\[17,18\]. A few empiric models proposed in Refs.\[3,20\] cannot explain
the observations completely, \[21\], and their simplicity masks many
outstanding issues \[17,21\]. In a recent attempt \[22\] to reproduce
the results of the drag model \[20\] in a single-mode approximation, a
complete set of the boundary conditions were not satisfied (see below),
and the conservation laws were thus violated.

In this Letter we report multiple harmonic theoretical solutions for a
system of conservation laws, which describe the large-scale coherent
dynamics in RTI and RMI for fluids with a finite density differences in
general three-dimensional case. The analysis yields new properties of
the bubble front dynamics. In either RTI or RMI, the obtained
dependencies of the bubble velocity and curvature on the density ratio
differ qualitatively and quantitatively from those suggested by the
models \[3,20,22\]. For the first time, our theory reveals an important
qualitative distinction between the dynamics of the RT and RM bubbles.
Asymptotically, the RT bubble is curved, and its curvature has a strong
dependence on the density ratio, while the RM bubble flattens
independently on the density ratio. The difference in the curvatures of
the RT and RM bubbles inter-depends with that in their velocities. Our
results explain existing data \[7-12,23\].

Based on the experimental observations \[4,7,9\], we divide the fluid
interface into active and passive regions, similarly to Ref. \[24\]. In
the active regions (small scales), the vorticity is intensive, while the
passive regions (large scales) are simply advected. If the fluids have a
finite density contrast and the initial perturbation has small
amplitude, the cascades of energy are not extensive in the nonlinear
regime \[6-12\]. A considerable part of the fluid energy concentrates in
the large-scale coherent motion. The dynamics of the coherent structure
could be described by scalar fields. One can apply therefore the
spectral approach and group theory \[15,16,25\].

Let *t* be time, ![](media/image1.wmf) be the Cartesian coordinates, and
![](media/image2.wmf) be the fluid interface. In the region
![](media/image3.wmf) the fluid density is $\rho = \rho_{h}$ and the
velocity is $v = v_{h}$, while in the region ![](media/image4.wmf) the
fluid density is $\rho = \rho_{l}$ and the velocity is $v = v_{l}$.
Assuming $\nabla \cdot v = 0$, we derive from the conservation laws

$\ \left( {\overset{˙}{v}}_{h} + \nabla\left( v_{h}^{2}/2 \right) + g \right)\ \rho_{h}|_{\theta = 0} = \left( {\overset{˙}{v}}_{l} + \nabla\left( v_{l}^{2}/2 \right) + g \right)\ \rho_{l}|_{\theta = 0}$,
(1a)

$\ \left( \overset{˙}{\theta} + v_{h} \cdot \nabla\theta \right)\ \rho_{h}|_{\theta = 0} = \left( \overset{˙}{\theta} + v_{l} \cdot \nabla\theta \right)\ \rho_{l}|_{\theta = 0}$,
(1b)

where ![](media/image5.wmf), the dot indicates a partial
time-derivative, and ![](media/image6.wmf) is the gravity directed from
the heavy fluid to the light fluid with ![](media/image7.wmf). In
accordance with the condition Eq. (1a) the pressure is continuous at the
interface. Equation (1b) is the equation of continuity. If there is no
mass flux across the interface, the normal component of velocity is
continuous at the interface and

$\ \left( \overset{˙}{\theta} + v_{h} \cdot \nabla\theta \right)\ \rho_{h}|_{\theta = 0} = \left( \overset{˙}{\theta} + v_{l} \cdot \nabla\theta \right)\ \rho_{l}|_{\theta = 0} = 0$.
(1c)

The boundary conditions at the infinity close the set of the governing
equations

$v_{h}|_{z\text{=+}\infty} = v_{l}|_{z = - \infty} = 0$. (1d)

The length scale (i.e. the spatial period), the time-scale, and the
symmetry of the motion in Eqs.(1) are determined from the initial
conditions. We choose the value of the spatial period in a vicinity of
the wavelength of the mode of fastest growth set by surface tension and
viscosity in real fluids, \[5\]. For RTI, the time-scale is
$\tau_{\text{RT}}\text{\textasciitilde}\sqrt{\lambda/\text{Ag}}$ where
![](media/image8.wmf) is the Atwood number and ![](media/image9.wmf)
\[1\]. For RMI, ![](media/image10.wmf) in Eqs.(1) and
$\tau_{\text{RM}}\text{\textasciitilde}\lambda/v_{0}$, where $v_{0}$ is
the initial velocity value \[16,18,19\]. In order to be stable under
large-scale modulations, the coherent motion in Eqs.(1) must be
invariant under one of spatial symmorphic groups with translations in
the plane ![](media/image11.wmf) and inversion ![](media/image12.wmf)
and $y \rightarrow - y$ \[15,16,25\]. For a 3D flow, these are hexagonal
p6mm, square p4mm, rectangular p2mm, and few other groups, \[25\].

To describe the dynamics of the nonlinear bubble, we reduce Eqs.(1) to a
local dynamical system, as in \[15,16\]. All calculations are performed
in the frame of reference moving with velocity $v(t)$ in the
*z*-direction, where $v(t)$ is the velocity at the bubble tip in the
laboratory frame of references. For potential flow
$v_{h(l)} = \nabla\Phi_{h(l)}$, and we expand the potential
$\Phi_{h(l)}$ as a Fourier series. In the case of a 3D flow with
hexagonal symmetry p6mm,
$\Phi_{h} =$$\sum_{m = 1}^{\infty}{\Phi_{m}(t)}$$\left\lbrack \text{exp}\left( - \text{mk}\ z \right)\sum_{i = 1}^{3}{\text{cos}\left( {mk}_{i}r \right)}/3\text{mk} + z \right\rbrack$+cross
terms+$f_{h}(t)$ and
$\Phi_{l} =$$\sum_{m = 1}^{\infty}{{\overset{\sim}{\Phi}}_{m}(t)}$$\left\lbrack \text{exp}\left( \text{mk}\ z \right)\sum_{i = 1}^{3}{\text{cos}\left( {mk}_{i}r \right)}/3\text{mk} - z \right\rbrack$+cross
terms +$f_{l}(t)$, where $r = (x,y)$, $k_{i}$ are the vectors of the
reciprocal lattice with $k_{1} = k\left( \sqrt{3}/2,1/2 \right)$,
$k_{2} = k\left( \sqrt{3}/2, - 1/2 \right)$, $k_{3} = k_{1} - k_{2}$,
and $k = \left| k_{i} \right| = 4\pi/\sqrt{3}\lambda$, and $f_{h(l)}$
are time-dependent functions. The bubble tip is a stagnation point with
coordinates (0,0,0) in the unit cell. For $x\text{\textasciitilde}0$ and
$y\text{\textasciitilde}0$ the interface can be expanded as a power
series,
$z^{} = \sum_{N = 1}^{\infty}{\zeta_{N}(t)\left( x^{2} + y^{2} \right)^{N}}$+cross
terms, where $\zeta_{1}(t)$ is the time-dependent principal curvature at
the bubble tip. Using the foregoing expressions, we expand Eqs. (1) for
$x\text{\textasciitilde}0$ and $y\text{\textasciitilde}0$ and derive a
system of ordinary differential equations for the variables
$\zeta_{N}(t)$ and the moments
$M_{n}(t) =$$\sum_{m = 1}^{\infty}{\Phi_{m}\left( \text{km} \right)^{n}}$+cross
terms and
${\overset{\sim}{M}}_{n}(t) =$$\sum_{m = 1}^{\infty}{{\overset{\sim}{\Phi}}_{m}\left( \text{km} \right)^{n}}$+cross
terms, where *n* is an integer. For *N*=1, the conditions in Eqs.(1a-d)
take respectively the form

$\ \left( {\overset{˙}{M}}_{1}/4 + \zeta_{1}{\overset{˙}{M}}_{0} - M_{1}^{2}/8 + \zeta_{1}g \right)\ \rho_{h} = \left( {\overset{\sim}{\overset{˙}{M}}}_{1}/4 - \zeta_{1}{\overset{\sim}{\overset{˙}{M}}}_{0} - {\overset{\sim}{M}}_{1}^{2}/8 + \zeta_{1}g \right)\ \rho_{l}$,
(2a)

$\ \left( {\overset{˙}{\zeta}}_{1} - 2\zeta_{1}M_{1} - M_{2}/4 \right)\ \rho_{h} = \left( {\overset{˙}{\zeta}}_{1} - 2\zeta_{1}{\overset{\sim}{M}}_{1} + {\overset{\sim}{M}}_{2}/4 \right)\ \rho_{l}$,
(2b)

$\ \left( {\overset{˙}{\zeta}}_{1} - 2\zeta_{1}M_{1} - M_{2}/4 \right)\ \rho_{h} = \left( {\overset{˙}{\zeta}}_{1} - 2\zeta_{1}{\overset{\sim}{M}}_{1} + {\overset{\sim}{M}}_{2}/4 \right)\ \rho_{l} = 0$,
(2c)

$M_{0} = - {\overset{\sim}{M}}_{0} = - v$. (2d)

In a single-mode approximation ($N = 1$, $m = 1$), only the first-order
amplitudes of the velocity potentials are retained in the expressions
for the moments, $M_{n} = \Phi_{1}k^{n}$ and
${\overset{\sim}{M}}_{n} = {\overset{\sim}{\Phi}}_{1}k^{n}$ with
*n*=0,1,2, and the system in Eqs.(1,2) is reduced to the expansion of
the Layzer-type \[14,15,16\]. Then, the Layzer-type expansion has a
solution, which conserves momentum, Eqs.(1a,2a), and mass, Eqs.(1b,2b)
and satisfies Eqs.(1d,2d). For this solution, in the nonlinear regime of
RTI, $t\text{>>}\tau_{\text{RT}}$, the bubble curvature is
$\zeta_{1} = \zeta_{L} = - \text{Ak}/8$ and its velocity
$v =$$v_{L\text{, RT}} = \sqrt{\text{Ag}/k}$. The velocity
![](media/image13.wmf) has the re-scaling ![](media/image14.wmf) and
agrees therefore with the empiric approach \[3\]. In the nonlinear
regime of RMI, $t\text{>>}\tau_{\text{RM}}$, the curvature of the
Layzer-type bubble is $\zeta_{1} = \zeta_{L} = - \text{Ak}/8$ and its
velocity $v = v_{L\text{, RM}} =$$\left( 2 - A^{2} \right)/\text{Akt}$.
However, in either RTI or RMI, the single-mode solution violates the
conditions of continuity of the normal component of velocity,
Eqs.(1c,2c), and requires mass flux across the interface. To evade this
issue, Ref. \[22\] has added an artificial term in the potential of the
light fluid, which violated the boundary conditions in Eqs.(1d,2d) and
introduced a non-uniform time-dependent mass flow at the infinity. The
solution \[22\] has agreed then with the assumptions of the drag model
\[20\]. In Refs.\[20,22\] the bubble curvature is
$\zeta_{1} = \zeta_{D} = \zeta_{L}|_{A = 1}$ and its velocity
$v =$$v_{D\text{, RT}} =$ $\sqrt{2A/(1 + A)}\sqrt{g/k}$ in RTI, and
$\zeta_{1} =$$\zeta_{D} = \zeta_{L}|_{A = 1}$ and
$v =$$v_{D\text{, RM}} =$$2/(1 + A)\text{kt}$ in RMI. We conclude that a
single-mode approximation does not satisfy the complete set of the
governing equations, Eqs.(1). The Layzer-type nonlinear solution (either
in our version or in the model \[22\]) violates the conservation laws,
Eqs.(1,2), and it is not therefore a physical solution. The reason of
this difficulty lies in a non-local character of the non-linearity in
Eqs.(1,2). A multiple harmonic analysis will resolve the issue.

To find regular asymptotic solutions, describing the nonlinear evolution
of the bubble front in RTI or RMI, one should account for non-local
properties of the flow that has singularities, \[17,18\]. The
singularities determine the interplay of harmonics in the global flow as
well as in the local dynamics system, and affect therefore the shape of
the regular bubble. Assuming that the bubble shape, parameterized by the
principal curvature at its tip, is free, we find a continuous family of
regular asymptotic solutions for the local system. The family involves
all bubbles allowed by the symmetry of the global flow. For the regular
asymptotic solutions the interplay of harmonics is well captured. We
perform the stability analysis and choose the fastest stable solution in
the family as being physically significant.

For the Rayleigh-Taylor instability, we impose the conditions
${\overset{˙}{\zeta}}_{N} = {\overset{˙}{M}}_{n} =$${\overset{\sim}{\overset{˙}{M}}}_{n} = 0$
on Eqs.(2), keep the relation
$\zeta_{1} =$$- M_{2}/8M_{1} =$${\overset{\sim}{M}}_{2}/8{\overset{\sim}{M}}_{1}$
as a free parameter, and find in this way the velocity $v$, the
amplitudes $\Phi_{m},{\overset{\sim}{\Phi}}_{m}$, and the surface
variables $\zeta_{N}$ as functions on this parameter. For every *N*,
additional relations between the moments can be established \[15,16\].
At $N = 1$ and $m = 2$, with
$M_{0} =$$\Phi_{1} + \Phi_{2} =$$\left( 3M_{1}/2k \right)$
$\left( 1 - M_{2}/3M_{1}k \right)$ and
${\overset{\sim}{M}}_{0} =$${\overset{\sim}{\Phi}}_{1} + {\overset{\sim}{\Phi}}_{2} =$$\left( 3{\overset{\sim}{M}}_{1}/2 \right)$$\left( 1 - {\overset{\sim}{M}}_{2}/3{\overset{\sim}{M}}_{1}k \right)$,
one finds the bubble velocity as the function of the bubble curvature
and the Atwood number from Eqs. (2):

$v = \sqrt{\frac{g}{k}}\left( \frac{- 2A\left( \zeta_{1}/k \right)}{A\left( 3 - 8\left( \zeta_{1}/k \right) \right)^{2} - (1 - A)\left( \zeta_{1}/k \right)} \right)^{1/2}\left( 9 - \text{64}\left( \zeta_{1}/k \right)^{2} \right)$.
(3a)

For $N > 1$, the family of solutions could be found analogously. In the
interval $\zeta_{\text{cr}} < \zeta_{1} \leq 0$, the Fourier amplitudes
$\Phi_{m},{\overset{\sim}{\Phi}}_{m}$ exponentially decay with increase
in *m*, the lowest-order amplitudes are dominant, and higher order
corrections for solutions in (3a) are small. For
$\zeta_{1}\text{\textasciitilde}\zeta_{\text{cr}}$ the convergence is
broken, where $\zeta_{\text{cr}} \approx - k/6$ for $A \approx 1$ and
$\zeta_{\text{cr}} \approx 0$ for $A \approx 0$. For the fastest
solution in (3a) $\left( \partial v/\partial\zeta_{1} \right) = 0$ and
$\left( \partial^{2}v/\partial\zeta_{1}^{2} \right) < 0$. This is the
solution with

$\zeta_{1} = \zeta_{A\text{, RT}}$,
$v = v|_{\zeta_{1} = \zeta_{A,\text{ RT}}} = v_{A\text{, RT}}$. (3b)

Explicit analytical expressions for $\zeta_{A\text{, RT}}$ and
$v_{A\text{, RT}}$ are cumbersome and are not presented here. The
subscript *A* in (3b) emphasizes complex dependencies of the bubble
asymptotic velocity and curvature on the Atwood number. In the limiting
case of fluids with highly contrasting densities$A \approx 1$, the
solution (3b) takes the form

$\zeta_{A\text{, }\text{RT}} \approx - (k/8)\left( 1 - (1 - A)/8 \right)$,
$v_{A\text{, RT}} \approx \sqrt{g/k}\left( 1 - 3(1 - A)/\text{16} \right)$.
(3c)

In the other limiting case, $A \approx 0$, it can be expanded as

$\zeta_{A\text{, RT}} \approx - \left( 3k/\text{16} \right)A^{1/3}$,
$v_{A\text{, RT}} \approx (3/2)^{3/2}\sqrt{\text{Ag}/k}$. (3d)

To analyze stability of the obtained solutions, we perturb
$\zeta_{N},M_{n},{\overset{˙}{M}}_{n}$ slightly such that
$\zeta_{N} \rightarrow \zeta_{N} + \theta_{N}(t)$,
$M_{n} \rightarrow M_{n} + \Delta_{n}(t)$,
${\overset{\sim}{M}}_{n} \rightarrow {\overset{\sim}{M}}_{n} + {\overset{\sim}{\Delta}}_{n}(t)$.
We derive from Eqs.(2) $\theta_{N},$ $\Delta_{n},$
${\overset{\sim}{\Delta}}_{n}\text{\textasciitilde}\text{exp}\left( \omega\ (\ t/\tau) \right)$,
where the Liapunov's exponents
$\omega = \omega\left( \sqrt{\text{gk}},A,\left( \zeta_{1}/k \right) \right)$.
For stable solutions $\text{Re}\lbrack\omega\rbrack < 0$. Solutions with
$\zeta_{1}\text{\textasciitilde}\zeta_{\text{cr}}$ and
$\zeta_{1}\text{\textasciitilde}0$ appear unstable, while solutions with
$\zeta_{1}\text{\textasciitilde}\zeta_{A\text{, RT}}$ are stable.
Therefore, the physically significant solution in the RT family is the
solution (3b).

For the Richtmyer-Meshkov instability, the multiple harmonic regular
asymptotic solutions in Eqs. (2) with $g = 0$ and
$M_{n},{\overset{\sim}{M}}_{n}$,$v\text{\textasciitilde}1/t$ and
time-independent $\zeta_{N}$ could be found in a similar way. The RM and
RT families have a number of common properties such as the dominance of
the lowest order amplitudes, convergence, and critical solutions.
However the asymptotic dynamics of the bubble front in RMI appears quite
different from that in RTI. In the RM family the velocity has the form
for *N*=1

$v = \frac{1}{\text{kt}}\frac{\left( 3 + \text{10}A\left( \zeta_{1}/k \right) - \text{128}A\left( \zeta_{1}/k \right)^{3} \right)\left( 9 - \text{64}\left( \zeta_{1}/k \right)^{2} \right)}{\left( A\left( 3 - 8\left( \zeta_{1}/k \right) \right)^{2} - (1 - A)\left( \zeta_{1}/k \right) \right)}$.
(4a)

In contrast to the RT case, in the RM family
$\partial v/\partial\zeta_{1} < 0$ in the interval
$\zeta_{\text{cr}} < \zeta_{1} \leq 0$, and the fastest solution in Eq.
(4a) corresponds to a bubble with a flattened shape:

$\zeta_{1} = \zeta_{A\text{, RM}} = 0$,
$v = v|_{\zeta_{1} = \zeta_{A\text{, RM}}} = v_{A\text{, RM}} = 3/\text{Akt}$.
(4b)

Higher order corrections for the value of $v_{A\text{, RM}}$ are
reasonably small. The stability analysis shows that small deviations
from the solutions in the RM family have a power-law time-dependence
$\text{\textasciitilde}t^{\beta}$ with exponents
$\beta = \beta\left( A,\ \zeta_{1}/k \right)$. For stable solutions
$\text{Re}\lbrack\beta\rbrack < 0$, and solutions with a finite
curvature appear unstable, while flattened bubbles with
$\zeta_{1}\text{\textasciitilde}\zeta_{A\text{, RM}}$ are stable.
Therefore, the physically significant solution in the RM family is the
solution (4b).

We apply the foregoing analysis for 3D flows with other symmetries and
2D flows. In either RTI or RMI, the nonlinear dynamics of 3D highly
symmetric flows (p6mm, p4mm, \[25\]) coincide except for the difference
in the normalization factor *k*, \[26\]. A nearly isotropic shape of the
bubble is the reason of this universality. The 2D and 3D results are
similar qualitatively. In main order,
$\left\lbrack v_{A\text{, RT}} \right\rbrack_{\text{ 3D}}/\left\lbrack v_{A\text{, RT}} \right\rbrack_{\text{ 2D}} \approx \sqrt{3}$,
$\left\lbrack \zeta_{A\text{, RT}} \right\rbrack_{\text{ 3D}}/\left\lbrack \zeta_{A\text{, RT}} \right\rbrack_{\text{ 2D}} \approx 3/4$,
while
$\left\lbrack v_{A\text{, RM}} \right\rbrack_{\text{ 3D}}/\left\lbrack v_{A\text{, RM}} \right\rbrack_{\text{ 2D}} \approx 2$,
$\left\lbrack \zeta_{A\text{, RM}} \right\rbrack_{\text{ 3D}}/\left\lbrack \zeta_{A\text{, RM}} \right\rbrack_{\text{ 2D}} \approx 1$,
similarly to \[15,16\] for *A*=1. For 3D low-symmetric flows, the
asymptotic analysis shows a tendency of 3D bubbles to conserve isotropy
in the plane, and a discontinuity of the 3D-2D dimensional crossover.

Based on the obtained results, we expect the following dynamics of the
bubble front in the Rayleigh-Taylor and Richtmyer-Meshkov instabilities
for fluids with a finite density contrast in the case of a small initial
perturbation. In RTI, the bubble curvature $\zeta(t)$ and velocity
$v(t)$ grow as
$\text{\textasciitilde}\text{exp}\left( t/\tau_{\text{RT}} \right)$ in
the linear regime, $t\text{<<}\tau_{\text{RT}}$, and reach finite values
$\zeta_{1} \approx \zeta_{A\text{, RT}}$ and
$v \approx v_{A\text{, RT}}$ asymptotically,
$t\text{>>}\tau_{\text{RT}}$. The parameters $\zeta_{A\text{, RT}}$ and
$v_{A\text{, RT}}$ depend strongly on the Atwood number, Eq.(3b). For
$A < 1$ the bubbles are slower and less curved compared to the case of
$A = 1$, \[15\]; as $A \rightarrow 0$ the bubbles flatten and their
velocity approaches zero, Eq.(3d). In RMI, the bubble curvature and
velocity change as
$\zeta_{1}(t)\text{\textasciitilde} - k\ t/\tau_{\text{RM}}$ and
$v(t) - v_{0}\text{\textasciitilde}$$- v_{0}t/\tau_{\text{RM}}$ in the
linear regime, $t\text{<<}\tau_{\text{RM}}$; for
$t\text{\textasciitilde}\tau_{\text{RM}}$ the curvature $\zeta_{1}$
reaches an extreme value, dependent on the initial conditions; as
$t/\tau_{\text{RM}} \rightarrow \infty$ the bubble flattens,
$\zeta_{1} \rightarrow \zeta_{A\text{, RM}}$, and decelerates,
$v \rightarrow v_{A\text{, RM}}$, (4b). For $A < 1$ the bubbles move
faster than those for $A = 1$, \[16\], and for $A \leq 1$ the bubbles
are flat. The flattening of the bubble front is a distinct feature of
RMI universal for all $A$. It follows from the fact that RM bubbles
decelerate. As $A \rightarrow 0$ in Eq. (4b), the velocity
$v_{A\text{, RM}}\text{kt} \rightarrow \infty$ and this suggests that
for fluids with very close densities the bubble velocity has a much
faster time-dependence, such as $\text{\textasciitilde}t^{a}$ where
$- 1 < a < 0$, in a qualitative agreement with experiments in \[9\],
where for $A \approx 0$ the dependence $\text{\textasciitilde}1/t$ was
shown to underestimate the velocity data.

Figures 1, 2, and 3 compare our results with the models \[20,22\] and
with the foregoing Layzer-type solution, which agrees with Ref. \[3\].
For $0 < A \leq 1$ the asymptotic dynamics is different in RTI and RMI:
the bubble velocity reaches a constant value in RTI and decreases with
time in RMI. Therefore the distribution of pressure around the bubble is
distinct in RTI and RMI. This distinction should lead to a different
shape of the bubble front. Our theory adequately describes the effect:
in RTI the curvature $\zeta_{A\text{, RT}}$ has a strong dependence on
the Atwood number, while in RMI $\zeta_{A\text{, RM}} = 0$ for all $A$,
Fig. 1. The models \[3,20,22\] do not predict any difference in shape of
the RT and RM bubbles. In RTI, for $A = 1$ the values of
$\zeta_{A\text{, RT}} = \zeta_{L} = \zeta_{D}$ \[15\], however, for
$A \leq 0\text{.}6$ the quantitative difference among the values
$\zeta_{L}$, $\zeta_{D}$ and $\zeta_{A\text{, RT}}$ is significant,
Fig.1. In the limit $A \rightarrow 0$, the ratio between the curvatures
$\zeta_{D}k\text{\textasciitilde}$$\text{const}$ from the models
\[20,22\] and $\zeta_{A\text{, RT}}k\text{\textasciitilde}$
$- A^{\frac{1}{3}}$ from Eqs.(3b-d) approaches
$\left( \zeta_{D}/\zeta_{A\text{, RT}} \right) \rightarrow \infty$,
while for the Layzer-type solution
$\zeta_{L}/\zeta_{A\text{, RT}} \rightarrow 0$. Emphasize that for
$A = 0$ the RT instability does not develop and the bubble curvature
should remain zero for all $t$, in agreement with Eqs.(3b-d). In RMI,
the bubble in (4b) flattens asymptotically for all $A$, while the models
\[20,22\] suggest an $A$-independent finite value of the bubble
curvature $\zeta_{D} = \zeta_{L}|_{A = 1}$, Fig. 1.

Our results are in good agreement with existing data \[7-12,23\]. In
RTI, the 2D numerical simulations in Ref. \[23\] agree quantitatively
with the solution Eqs.(3b-d). In RMI, the bubble flattening has been
observed in many 2D and 3D observations \[9-12,23\], in particular, in
the experiments in Ref. \[9\], accurately reproduced by the simulations
in Ref.\[12\]. The measurement of the bubble curvature requires however
an improved diagnostics of the interface dynamics. Owing to the lack of
data, we cannot perform a detailed quantitative comparison with all
experiments and simulations here.

We conclude that in either RTI or RMI the bubble curvature is a
sensitive diagnostic parameter. The value of
$\lambda\left| \zeta_{1} \right|$ defines how flat or narrow the bubble
is for a given length scale. This dimensionless shape parameter is
related to the bubble velocity and determines the flow drag. The value
of the drag force in the chaotic RTI and RMI is still a subject to
controversy \[21\], and the models \[3,20\] do not predict the bubble
shape and width correctly, \[4\]. Our theory resolves the problem for
highly nonlinear RTI and RMI. The major issues to check by experiments
and simulations are the following. Does the curvature of the RT bubble
have a strong dependence on the Atwood number or reach an
*A*-independent value? Does the curvature of the RM bubble
asymptotically vanish or approach a finite value?

For the bubble velocity, the experiments and simulations in RTI \[7,8\]
have been in a reasonable agreement with the dependencies
$v_{L\text{, RT}}$ and $v_{D\text{, RT}}$ in the case of
$A\text{\textasciitilde}1$, \[3, 20\]. On the other hand, as shown in
the foregoing, the single-mode solutions $v_{L\text{, RT}}$ and
$v_{D\text{, RT}}$ violate the conservation laws. Figure 2 explains this
apparent paradox. For finite values of $A$ and for
$\zeta_{1}\text{\textasciitilde} - 1/\lambda$, the bubble velocity is
relatively insensitive to details of the interface dynamics. In the
interval $0\text{.}3 \leq A \leq 1$, the quantitative distinction among
the values of $v_{A\text{, RT}}$, $v_{D\text{, }\text{RT}}$ and
$v_{L\text{, RT}}$ is
$\text{\textasciitilde}\text{10} - \text{20}\text{\%}$. In the limit of
$A \rightarrow 0$, the distinction increases and
$v_{A\text{, RT}}:v_{D\text{, RT}}:$
$v_{L\text{, RT}} =$$(3/2)^{3/2}:\sqrt{2}:1$. In RMI, the bubble
velocities given by the drag models $v_{D\text{, RM}}$ in \[20,22\], and
the Layzer-type solution $v_{L\text{, RM}}$ differ significantly from
$v_{A\text{, RM}}$, (4b), Fig.3. The difference is however hard to
distinguish in experiments. Since the velocity decays as
$\text{\textasciitilde}C/k\ t$, the bubble displacement,
$\text{\textasciitilde}C\ \text{log}\left( t/\tau_{\text{RM}} \right)/k$,
is comparable with an experimental error, and the coefficient $C$ in
this dependence is impossible to evaluate. In contrast, the bubble
curvature is a reliable diagnostic parameter as discussed in the
foregoing.

Our theory describes the principal influence of the density ratio on the
large-scale nonlinear dynamics of the RT and RM bubbles. The analysis is
based on the assumptions that the flow dynamics is governed by a
dominant mode, the transfers of energy to smaller or larger scales are
not extensive, and the vorticity does not change the time-dependence of
the large-scale coherent motion. If these conditions are broken (for
example, for fluids with very close densities, $A \approx 0$), the
potential approximation may not give a correct time-dependence for the
asymptotic bubble motion. We address these important issues in the
future.

In conclusion, we have found a multiple harmonic solution for a complete
system of conservation laws, describing the large-scale coherent
dynamics in RTI/RMI in general three-dimensional case. The theory yields
new dependencies of the bubble curvature and velocity on the density
ratio in RTI and RMI, and predicts in a consistent manner an important
quantitative difference between the dynamics of the RT and RM bubbles.

**Acknowledgements**

S.I.A. thanks the Center for Turbulence Research, Stanford University;
S.I.A. and K.N. are grateful to the Japanese Society for the Promotion
of Science; J.G. thanks US Dept. of Energy, grants DE-FG02-90ER25084 and
DE-AC02-98CH10886, the Army Research Organization grant DAAD190110642,
the NSG grant DMS0102480 and the LANL contract C26730001014L

**References**

1.  Lord Rayleigh, Scientific Papers, II (Cambridge University Press,
    Cambridge, England, 1990); R.M. Davies, G.I. Taylor, Proc.Roy.Soc. A
    200, 375 (1950).

2.  R.D. Richtmyer, Commun. Pure & Appl. Math. 13, 297 (1960); E.E.
    Meshkov, Sov. Fluid Dyn.4,101(1969).

3.  D.H. Sharp, Physica D 12, 3 (1984).

4.  M. Schneider, G. Dimonte, B. Remington, Phys. Rev. Lett. 80, 3507
    (1998).

5.  S. Chandrasekhar, Hydrodynamic and hydro-magnetic stability, 3d ed.
    (Dover Publications Inc., New York, 1981), pp.428-477; R. Menikoff,
    R.C. Mjolness, D.H. Sharp, C. Zemach, Phys. Fluids 20, 2000, 1977;
    K.O. Mikaelian, Phys. Rev. E 47, 375 (1993).

6.  S.B. Dalziel, P.F. Linden, D.L. Youngs, Jour. Fluid Mech., 399, 1
    (1999).

7.  X.Y.He, R.Y.Zhang, S.Y.Chen, G.D.Doolen, Phys.Fluids 11, 1143
    (1999).

8.  G. Gardner et al, Phys. Fluids 31, 447 (1988); D.L. Youngs, Phys.
    Fluids A 3, 1312 (1991).

9.  J.W. Jacobs and J.M. Sheeley, Phys. Fluids 8, 405 (1996).

10. Yu.A. Kucherenko, et al., Izv. Akad. Nauk SSSR, Mekh. Zhidk. Gaza 6,
    157 (1978); S.G. Zaytsev, et. al., Sov. Phys. Doklady 35, 159
    (1990).

11. R. Holmes, J.W. Grove, D.H. Sharp, J. Fluid Mech.301, 51 (1995);
    R.L. Holmes et al., J. Fluid Mech. 389, 55 (1999); N.V. Volkov, A.E.
    Maier, A.P. Yalovets, Tech. Phys. Lett 27, 20 (2001).

12. A.D. Kotel'nikov, J. Ray, N.J. Zabusky, Phys. Fluids 12, 3245
    (2000).

13. P.R. Garabedian, Proc. Roy. Soc. A 241, 423 (1957); G. Birkhgoff
    & D. Carter, J. Math. Mech. 6, 769, (1957); Vanden Broeck, Phys.
    Fluids 27, 1090 (1984); J.A. Zuffiria, Phys. Fluids 31, 3199
    (1988); S. Tanveer, Proc. Roy. Soc. London A 441, 501 (1993); G.
    Hazak, Phys. Rev. Lett. 76, 4167 (1996).

14. D. Layzer, Astrophys. Jour. 122, 1 (1955).

15. S.I. Abarzhi, Phys.Rev.Lett.89, 337 (1998); Phys. Rev. E 59, 1729
    (1999).

16. S.I. Abarzhi, Physics Letters A 294, 95 (2002), Phys. Rev. E 66, art
    036301 (2002).

17. D.I. Meiron, G.R. Baker, S.A. Orszag, J. Fluid Mech. 114, 238
    (1982); D.W. Moore, Proc. R. Soc. London A 365, 105 (1979); M.
    Shelley, J. Fluid Mech. 119, 507 (1992); T.Y. Hou, J.S. Lowengrub,
    M.J. Shelley, Phys. Fluids 9,1933 (1997); S.J. Cowley, G.R.
    Baker, S. Tanveer, J. Fluid Mech. 378, 233 (1999).

18. C. Matsuoka, K. Nishihara, Y. Fukuda, Phys. Rev. E 67, art. 082303
    (2003).

19. J.G. Wouchuk, K. Nishihara, Phys. Plasmas 3, 376 (1996), Phys.
    Plasmas 4, 1028 (1997); S.W. Haan, Phys. Fluids B3, 2349 (1991);
    A.L. Velikovich, G. Dimonte, Phys. Rev. Lett. 76, 3112 (1996).

20. U. Alon, J. Hecht, D. Offer, D. Shvarts, Phys. Rev. Lett. 74, 534
    (1995); D. Oron, U. Alon, D. Offer, D. Shvarts, Phys. Plasmas 8,
    2883 (2001).

21. G. Dimonte, Phys. Plasmas 7, 2255 (2000).

22. V.N. Goncharov, Phys. Rev. Lett.88, 4502 (2002).

23. J. Glimm et al., SIAM J. Sci. Comput. 21, 2240 (2000); J. Glimm et
    al., SIAM J. Sci. Comput. 24, 208 (2002); S.I. Abarzhi, A. Lin, J.
    Glimm, Phys. Fluids, subm. 2002.

24. H. Aref, G. Tryggvason, Phys. Rev. Lett. 62, 749 (1989).

25. A.V. Shubnikov & V.A. Koptsik, Symmetry in science and art (Plenum
    Press, New-York-London, 1979). Basic elements for the group p6mm are
    periodicity in the (*x*, *y*) plane, 6-fold axis of rotation and
    mirror planes of reflection normal each other and parallel to the
    *z*-axis. Groups p4mm and p2mm have 4-fold and 2-fold axis of
    rotation. Symmorphic groups have no glide planes and screw axes.

26. For hexagonal symmetry $k = 4\pi/\lambda\sqrt{3}$, for square
    symmetry $k = 2\pi/\lambda$, for tubular flow
    $k = 2\beta_{0}/\lambda$, where $\beta_{0}$ is the first zero of the
    Bessel function $J_{0}$.

# Figure Captions

**Fig.1**: Dependence of the bubble curvature on the Atwood number *A*
for 3D highly symmetric flows in the Rayleigh-Taylor and
Richtmyer-Meshkov instabilities; *k* is the wavevector,
$\zeta_{L} = - \text{Ak}/8$ is the Layzer-type solution,
$\zeta_{D} = - k/8$ is given by the models \[20,22\],
$\zeta_{A\text{, RT}}$ and $\zeta_{A\text{, RM}} = 0$ are the solutions
(3b) and (4b) respectively.

**Fig.2**: Dependence of bubble velocity *v* on the Atwood number *A*
for 3D highly symmetric flows in the Rayleigh-Taylor instability; *k* is
the wavevector and *g* is gravity,
$v_{L\text{, RT}} = \sqrt{\text{Ag}/k}$ is the Layzer-type solution,
$v_{D\text{, RT}} = \sqrt{2\text{Ag}/(1 + A)k}$ is given by the models
\[20,22\], and $v_{A\text{, RT}}$ is the solution (3b).

**Fig.3**: Dependence of the bubble velocity *v* on the Atwood number
*A* for 3D highly symmetric flows in the Richtmyer-Meshkov instability;
*k* is the wavevector, *t* is time,
$v_{L\text{, RM}} = \left( 2 - A^{2} \right)/\text{Akt}$ is the
Layzer-type solution, $v_{D\text{, RM}} =$$2/(1 + A)\text{kt}$ is given
by the models \[20,22\], $v_{A\text{, RM}} = 3/\text{Akt}$ is the
solution (4b).

![](media/image15.bmp)

Figure 1

![](media/image16.bmp)

Figure 2

![](media/image17.bmp)

Figure 3