数理解析研究所講究録
第 1606 巻 2008 年 41-57

41

Vortex dynamics in incompressible
Richtmyer-Meshkov instability with surface
tensionChihiro Matsuoka*
Department of Physics,
Graduate School of Science and Technology,
Ehime University, Matsuyama, Ehime 790-8577, Japan

Abstract

Motion of a planar interface in incompressible Richtmyer-Meshkov (RM) instability withsurface tension is investigated numerically by use of the boundary integral method. It isshown that an interface rolls up without regularization of the interfacial velocity when the
Atwood number is relatively small. The comparison between the growth rate of an interface
with and without surface tension is made, and we present that the growth rate of the interface
with surface tension is equal to the one without surface tension at the as ymptotic stage whenthe Atwood number is identical to each other. A phenomenon known as ‘pinching’ in thephysics of drops is found in the final stage of calculations for relatively small Atwood numbers
and it is shown that this phenomenon is caused by a vortex pair induced on the interface.
We also present that when the surface tension coefficient is relatively large, stable oscillatory
motion appears for RM instability. When the gravity is taken into account, linearly stable but
asymptotically unstable motion can appear under a certain condition of three parameters;
the Atwood number, the gravity and the surface tension.

1

INTRODUCTION

When an inhomogeneous vorticity initially distributes on an interface between two fluids withdifferent densities and it is driven by some external force such as a shock wave, a corrugatedinterface eventually rolls up to a mushroom-like structure. This phenomenon, known as theRichtmyer-Meshkov (RM) instability 1, 2), is important in various areas such as astrophysicalsupemova 4), supersonic combustion and the inertial confinement fusion. Due to the existence ofdensity inhomogeneity in the system, the RM instability can be also a model of an inhomogeneousand nonuniform turbulence, together with the Rayleigh-Taylor (RT) instability 3).
When the surface tension effect exists on the interface, the RM and RT instabilities can describe
the pure capillary wave and capillary-gravity wave 5, 6), respectively in the limit of Atwood number$A=1$ (see Sec. 2 for the definition of ). These describe relatively low velocity flow without vortex
motion, i.e., shear velocity motion, however, there exists high velocity motion with surface tension
such as breakup of drops in two-fluid $(A\neq 0)$ inviscid systems 7). One of such systems is the impact
of meteorites on the primitive earth 8), in which meteorites in liquid state collide on the surface of$A$“matsuoka@phys.sci.ehime-u.ac.jp

42
the earth at high speed. This is an example of motion of an unstable fluid interface with density
stratification and surface tension. It seems that the surface tension effect is considerably small
in some real experiments for the RM instability, for instance, such that the interfacial instability
is generated by impulsively accelerating a tank containing two liquids reported by Jacobs et al.
, however, that can be important when the liquids used in experiments are replaced anotherones in which the surface tension effect cannot disregard.
As is well-known, the surface tension term is expressed by the derivative of curvature of theinterface, therefore, it contains the highest spatial derivative in the governing equations (see Sec. 2)
and this term governs interfacial motion in the high wave numbers together with the convective
term which gives rise to vortex motion. When the surface tension term is sufficiently large, theeffect may suppress the vortex motion such as the roll-up of the interface 12, 13). Hou et al. 14, 15)
numerically examined the long-time motion of vortex sheet with surface tension for the density
matched $(A=0)$ case and observed a phenomenon called ‘pinching singularity’ in the physics ofdrops 7).
Pinching or pinch-off is a phenomenon observed in capillary-driven motion of free surface
flows when the convective nonlinearity, i.e., vortex motion is taking into account in the system
7, 16).
The pinching phenomenon is important for the process which a single drop breaks up into
several drops and the pinching singularity or pinch-off is observed as a phenomenon at which theinterfacial curvature and the interfacial velocity diverge to infinity. This singularity is peculiar tosurface tension driven motion and it is considered 14, 15, 16, 20) as different from Moore’s curvature
singularity 17) which is observed in the inviscid vortex motion without surface tension. Hou et al.
14, 15, 20)
calculated this pinching singularity with high accurate numerical scheme and obtained
a result that the pinching region does not perfectly stick and there exists a narrow gap. This
pinching singularity is hard to appear for larger Atwood numbers (see Secs. 3 and 4).
In this article we examine the effect that the surface tension gives to the interfacial motionin the RM instability and show the results which do not appear in the zero surface tension case12, 13, 18, 19)
or in the density matched case 14, 15). The well-known result that the linear growthrate in the RM instability is proportional to time does not hold for finite surface tension case.
This result, which can be easily expected from the linear stability analysis, is also confirmed
in our numerical calculations. As the surface tension parameter becomes large, the roll-up of
the interface disappears and a stable oscillatory motion appears instead. The frequency of thisoscillation depends on the surface tension parameter for a fixed Atwood number and the motion
is not so regular as the standing wave motion for the RT instability 20).
When the gravity is taking into account, i.e., for the RT instability, linearly stable but asymptotically unstable motion appears when three parameters, the Atwood number, the gravity and thesurface tension coefficient satisfy a certain condition. This corresponds to a critical or marginallystable state at which the system changes from stable to unstable. The existence of such criticalmotion can be also expected from the linear stability analysis, however, it has never presented
how the motion is. We show that a kind of pinching phenomenon appears at the final stage ofthis critical motion (see Sec. 4.3). In Sec. 2 we present governing equations which are used forour numerical calculations. In Sec. 3 we present interfacial motion for the RM instability withvarious Atwood numbers and surface tension parameters using the scheme presented by Hou et al.
14). Stable oscillatory motion in the RM instability and motion taking into account the gravity,
i.e., motion in the RT instability are presented in Sec. 4. Section 5 is devoted to discussions andconclusion.
$9)_{-}11)$

$t$

43

GOVERNING EQUATIONS

2

We consider a two-dimensional system, therefore, the interface is assumed to be a curve between
two fluids with different densities. The governing equations are equations which describes theinterfacial position $(x, y)$ and the pressure boundary condition, i.e., the Laplace-Young condition.
Due to the strong numerical instability caused by the surface tension term, the long-time computation is impossible with the explicit method which we have used for the vortex sheet motion
with zero surface tension parameter ). Instead, we adopt here the formulation which waspresented by Hou et al. 14) in order to calculate the interfacial motion with surface tension. We
briefly review the method here for convenience and completeness.
$13$

2.1

$18)_{-}19)$

The formulation for numerical calculations

The temporal evolution of the interface $X=(x(\beta, t), y(\beta, t))$ is given by(1)$X_{t}=Un+Tt$ ,whereis a Lagrangian parameter which parameterize the interface, is time, $U$ and $T$ are
the normal and tangential velocities of the interface, respectively, and the subscript denotes thedifferentiation with respect to the variable. The unit normal $n$ and the unit tangent are describedby the tangent angle to the interface as$\beta$$t$$t$$\theta$$n=(-\sin\theta, \cos\theta)$,$t=(\cos\theta, \sin\theta)$.

By Frenet’s formula for plane curves, we have$t_{s}=Kn$ ,$n_{s}=-Kt$ ,where $K$ is curvature of the interface and is arc length given byinterfacial position is specified withand as follows:$s$$s= \int\sqrt{x_{\beta}^{2}+y_{\beta}^{2}}d\beta$.

The$\theta$$s_{\beta}$$S_{\beta,t}$$\theta_{t}$

$=$

$=$

$T_{\beta}- \frac{\theta_{\beta}^{2}}{s_{\beta}}$

(2a)

,

$\frac{1}{s_{\beta}}(\frac{\theta_{\beta}}{s_{\beta}})_{\beta}+\frac{T}{s_{\beta}}\theta_{\beta}$

.

(2b)

We assume that the flow is incompressible and inviscid. Then the pressure boundary condition,i.e., the Laplace-Young condition which holds at the interface,(3)$p_{1}-p_{2}=\sigma K$

is rewritten into the form$\rho_{1}[\frac{\partial\phi_{1}}{\partial t}+\frac{1}{2}(\nabla\phi_{1})^{2}+g]-\rho_{2}[\frac{\partial\phi_{2}}{\partial t}+\frac{1}{2}(\nabla\phi_{2})^{2}+g]=\sigma K$

(4)

,

where
is the density of fluid $i(i=1,2),$is the pressure in fluidis the gravity, is thesurface tension coefficient andis the velocity potential in fluid which is related to each fluid
velocityas. Note that
holds in each fluid region by the incompressibility
condition.
$\rho_{i}$$i,$$p_{i}$$i$$\phi_{i}$$u_{i}$$u_{i}=\nabla\phi_{i}$$\Delta\phi_{i}=0$$g$$\sigma$44
Differentiating Eq. (4) with respect toof the second kind:$\beta$, we obtain the following Fredholm integral equation$\gamma_{t}=\sigma K_{\beta}+(\frac{T^{A}\gamma}{s_{\beta}}I_{\beta}+2A[s_{\beta}W_{t}\cdot t+\frac{1}{8}(\frac{\gamma}{s_{\beta}})_{\beta}^{2}-T^{A}W_{\beta}\cdot t+gy_{\beta}]$(5),, the Atwood number
is given byasis related with the circulation$A=(\rho_{2}-\rho_{1})/(\rho_{1}+\rho_{2})$ and we normalize the surface tension coefficient. The
velocity $W=(W_{x}, W_{y})$ induced by the vortex sheet is given by the Birkhoff-Rott equation 21)

where$\Gamma\equiv\phi_{1}-\phi_{2}$$\gamma$$A$$\gamma=\Gamma_{\beta}$$\sigmaarrow 2\sigma/(\rho_{1}+\rho_{2})$$W_{x}$$=$${\rm Re}[ \frac{1}{2\pi i}P.V.\int_{-\infty}^{\infty}\frac{\gamma(\beta’,t)d\beta’}{Z(\beta,t)-Z(\beta’,t)}]$$W_{y}$$=$$-{\rm Im}[ \frac{1}{2\pi i}$

P.V.

,

$/- \infty\infty\frac{\gamma(\beta’,t)d\beta’}{Z(\beta,t)-Z(\beta’,t)}]$(6),

where P.V. denotes the principal value and $Z=x+iy$ . The tangential velocity $T$ is decomposed
into the form
$T=T^{A}+W\cdot t$ ,

(7)

where
is an arbitrary tangential velocity and its choice determines the frame of reference. The
normal velocity is given by $U=W\cdot n$ in our formulations.
.
The arbitrary tangential velocity
cannot take the same form between $\sigma=0$ and$\sigma=0$,
the arbitrary tangential velocity
is chosen as$Whe,nthe13,l8l9,22)$ surface tension coefficient$T^{A}$$U$$T^{A}$$\sigma\neq 0$$T^{A}$$T^{A}= \frac{\alpha}{2s_{\beta}}\gamma$(8),

in our calculations, whereis an artificial parameter such that, the arbitrary tangential velocity
is chosen as$\alpha\neq 0$$\alpha=\alpha(\mathcal{A})$for $A\neq 0$ . When

$T^{A}$$\sigma\neq 0$(9),$T^{A}=W \cdot t+/0\beta(\theta_{\beta}U-\frac{1}{2\pi}\int_{0}^{2\pi}\theta_{\beta}’Ud\beta’)d\beta$

in which
is a constantis set to be everywhere equal to its mean 14, 15), i.e.,with respect to for a fixed time .
By solving Eqs. (2a), (2b) and (5), we can recover $(x, y, \kappa)$ fromusing the relation$s_{\beta}= \frac{1}{2\pi}\int_{0}^{2\pi}s_{\beta}’d\beta’$$s_{\beta}$$t$$\beta$$(s_{\beta}, \theta, \gamma)$

where

2.2

$\kappa=\gamma/s_{\beta}$

$x$

$=$

$/0^{\beta_{s_{\beta}’\cos\theta(\beta’)d\beta’+\beta}}$

$y$

$=$

$/0^{\beta_{s_{\beta}’\sin\theta(\beta’)d\beta’}}$

’

(10)

’

is the true vortex sheet strength which is frame independent.

Spatial discretization and numerical scheme

As we stated in Sec. 1,
in Eq. (5) causes strong numerical instability called ‘stiffness’. In order
to avoid this instability, we decompose the Cauchy integral in Eq. (6) as$K_{\beta}$$P.V$ . $/- \infty\infty\frac{\gamma(\beta’,t)d\beta’}{Z(\beta,t)-Z(\beta^{l},t)}=P.V.\int_{-\infty}^{\infty}[\frac{1}{Z_{\beta}(\beta-\beta’)}$$+$$( \frac{1}{Z(\beta,t)-Z(\beta’,t)}-\frac{1}{Z_{\beta}(\beta-\beta’)})]\gamma(\beta’, t)d\beta’$

.

(11)

45
Then the normal velocity can be expressed as$U( \beta, t)=\frac{1}{2s_{\beta}}H[\gamma|(\beta, t)+E[\gamma](\beta, t)$

(12)

,

where
is the Hilbert transform of and the Fourier component offor high mode 23).
Expandingwith Fourier modes$H[\gamma]$$\gamma$$E[\gamma]$

is exponential decay$(\theta, \gamma)$$\theta=\sum_{n=-\frac{N}{2}+1}^{\frac{N}{2}}\hat{\theta}_{n}e^{in\beta}$,$\gamma=\sum_{n=-\frac{N}{2}+1}^{\frac{N}{2}}\hat{\gamma}_{n}e^{in\beta}$,

and substituting these into Eqs. (2b) and (5), we have$\hat{\theta}_{t}$$\hat{\gamma}_{t}$,$=$$\frac{|n|}{2s_{\beta}^{2}}\hat{\gamma}(n)+\hat{P}(n)$$=$$- \frac{\sigma n^{2}}{s_{\beta}}\hat{\theta}(n)+\hat{Q}(n)$

(13)

,

where $N$ is the grid number, i.e., the number of point vortices, andandare Fouriertransformed remainder terms in the original equations Eqs. (2b) and (5), respectively. The equation forin Eq. (2a) is solved without applying the Fourier transform.
Following Hou et al. 14), we use the second order Adams-Bashforth method and the CrankNicholson scheme in order to perform the temporal integration forin Eqs. (2a) and (13),respectively. Applying these explicit and implicit method simultaneously, we can recoverand , therefore,
and using the inverse Fourier transform after each temporal integration.
The spatial integrations of the first and the second terms on the right hand side of Eq. (11) arecalculated with the alternate point quadrature method 24) and the trapezoidal rule, respectively.
There also exists explicit stable numerical method with spectral accuracy 25), however, calculations
with that method break down before appearing the pinching singularity.
$\hat{P}(n)$$\hat{Q}(n)$$s_{\beta}$$s_{\beta}$$\theta$$s_{\beta},$$\gamma$

3

$x,$

$y$

$\kappa$

Pinching phenomenon in Richtmyer-Meshkov instability

In this section we show some numerical results for motion of the interface in the RM instability$xarrow kx,$ $yarrow ky$ and $tarrow kv\iota int$
with surface tension. For those calculations, we normalize$12,13)$
with the asymptotic linear growth rate
and the wave number in the system. When
, we adopt the implicit scheme stated in Sec. 2.2, while we use the explicit scheme (the fourthorder Runge-Kutta scheme) which was adopted in Ref. 13) when $\sigma=0$ . For the former case, thearbitrary tangential velocity
is chosen in the form of Eq. (9), while
is chosen in the form ofEq. (8) for the latter case. All calculations throughout this paper except the calculation for $\sigma=0$
in Fig. 4 (see Sec. 3.1) are.performed without regularization 27) of the Cauchy integral in Eq. (6).
Initial conditions for all calculations in the RM instability are taken that$\betaarrow k\beta,$$v_{lin^{26)}}$$k$$\sigma\neq 0$$T^{A}$$x(\beta, 0)=\beta$$T^{A}$,$y(\beta, 0)=0$ ,$\gamma(\beta, 0)=2\sin\beta$

,

(14)

for both of $\sigma=0$ and, and we set $g=0$ in Eq. (5) for the RM instability. The periodic
boundary condition is imposed and the interval $[-\pi, \pi]$ is basically divided with grid number$N=1024$ . For calculations such that the pinching singularity appears, we refine the grid number$N$ from $N=1024$ to $N=2048$ .
We set the upper fluid $(\rho_{2}, y>0)$ is heavier than the lower fluid$(\rho_{1}, y<0)$ for all calculations including the RT instability (see Sec. 4).
$\sigma\neq 0$

46

$x$

Fig. 1: Interfacial profiles, circulation and sheet strength for $A=0.2$ and $\sigma=0$ at $t=(a)$ and(d) 2.6, (b) and (e) 5.0, and (c) and (f) 6.2, where dashed and solid lines in (d) $-(f)$ denote thecirculation and the sheet strength , respectively. The regularized parameter is set to $\delta=0.1$
here.
$\Gamma$$\kappa$$\delta$$\Gamma$$\kappa$

47

Fig. 2: Interfacial profiles, circulation and sheet strength for $A=0.2$ and $\sigma=0.05$ at $t=(a)$and (d) 1.5, (b) and (e) 5.0, and (c) and (f) 7.15, where dashed and solid lines in (d) $-(f)$ denotethe circulation and the sheet strength , respectively. The circled region in (c) is magnified inFig. 3.
$\Gamma$$\kappa$$\Gamma$$\kappa$48
3.1

Pinching singularity with small surface tension coefficient49
comes large. When the Atwood number $A=0$ , the pinching singularity appear at four points [seeFig. 6 $(a)]$ .
Pinching is caused due to a vortex pair which exists at the tip of the closest two points. We
show strong vorticity with plus (b) and minus (a) sign with black and white circle, respectively
in the figure, where the corresponding peaks in are shown in Fig. 2 (f). This vortex pair has
slightly minus net vorticity and the other vortex pair which is in the symmetrical position withrespect to the y-axis has slightly plus net vorticity. Pinching singularity is a phenomenon whichis caused due to the competition between the inertial force and the surface tension, and it is
not observed for a vortex sheet with $\sigma=0$ , i.e., this singularity is different from the well-known
Moore’s curvature singularity 17) [see Fig.5 (a) and $(d)$ ]. When the surface tension coefficient islarge, pinching in the RM instability does not appear even though the Atwood number is small.
We show the growth rate, i.e., the velocity of bubble and spike for $\sigma=0,$ $\sigma=0.05$ and $\sigma=0.1$
with the same Atwood number $A=0.2$ in Fig. 4 (a). When $\sigma=0$ , there exists deviation fromwhich is the exact value of the growth rate of bubble and spike at $t=0$ . This deviation is), where we set to $\delta=0.1$ for thecaused by the existence of the finite regularized parameter$\sigma=0$$\delta=0$. When we set to
calculation of, this deviation disappears (see Fig. 7), however,the breakdown of computation arises at earlier time [see Fig. 5 (a) and $(d)$ ]. We see that theasymptotic growth rate of bubble and spike coincides for three when we fix an Atwood number.
This tendency is unchanged even though we choose another value of for $\sigma=0$ or larger for
as long as the Atwood number is not close to 1. Figure 4 (b) shows thecalculations ofabsolute value of maximum sheet strength (peak value) for the corresponding parameters with(a). The maximum sheet strength for $\sigma=0$ has the largest value in all and the value becomes
smaller as the surface tension coefficient becomes large. This tendency is also true for another
Atwood numbers.
$\kappa$$\pm 1$$\delta 18$$\sigma$$\delta$$\sigma$$\sigma\neq 0$$\kappa$$\sigma$$\sigma$

3.2

Interfacial profiles for various Atwood numbers and surface tension coefficients

Figure 5 show interfacial and curvature profiles for a fixed Atwood number $A=0.2$ and various
with the regularized parameter $\delta=0$ . The time step is taken that $\Delta t=1.25x10^{-5}$ for $\sigma=0$
calculation, while $\Delta t=1.25\cross 10^{-4}$ for $\sigma=0.1$ and 0.3 calculations. The profiles for $\sigma=0$ and$\sigma=0.1$ are the ones at just before theirbreakdown times, where the refinement from $N=1024$$(0\leq t\leq 8.0)$ to $N=2048(8.0\leq t\leq 9.0)$ is performed for the calculation of $\sigma=0.1$ . We see thatthe interfacial profile (a) is smooth and the amplitude of it is small, however, the curvature (d) isvery large. This is well-known Moore’s curvature singularity 17) for $\delta=0$ and $\sigma=0$ . Note thatthe interface for $\sigma=0$ does not roll up unlike the case for finite . The circled region in (b) does
not stick perfectly as found in Fig. 3, where sharp peaks in the curvature profile (e) correspondto closest two points in that region. When the surface tension coefficientis large $(\sigma\geq 0.2)$ ,both of the curvature and the peak values in the sheet strength are small, and the circulation
almost does not change its shape from the initial state. The pinching phenomenon was not foundforwith all Atwood numbers. The calculation for $\sigma=0.3$ is stable and the breakdown ofcomputation does not occur even at past $t=15$ (see Sec. 4.2).
Figure 6 shows interfacial profiles for various Atwood numbers with a fixed surface tension
coefficient $\sigma=0.05$ . Time denoted in each figure caption corresponds to just before time atwhich the calculation breaks down, where the time step is taken that $\Delta t=1.25x10^{-4}$ for all
calculations. A refinement from $N=1024$ to $N=2048$ is performed for $A=0$ and $A=0.5$
calculations, where the former grid number is taken over $0\leq t\leq 4.0$ for $A=0$ and $0\leq t\leq 10.0$
for $A=0.5$ , while the latter grid number is taken over $4.0\leq t\leq 4.9$ for $A=0$ and $10.0\leq t\leq 11.5$
for $A=0.5$ , respectively. When $A=0$ , the pinching singularity appear at four points (circled$\sigma$$\delta$$\sigma$$\Gamma$$\kappa$$\sigma\geq 0.2$

50

(e)
$- 3$

-2-lOi
Lagrangianpara$\mathfrak{n}\epsilon$

$2$

ter

3

$\beta$Fig. 5: Interfacial profiles for $A=0.2$ with various $\sigma;(a)\sigma=0$ and $t=0.93,$ $(b)\sigma=0.1$ and$t=9.0$ , and (c) $\sigma=0.3$ and $t=12.0$ , where (d) $-(f)$ are curvature profiles corresponding to (a)$-(c)$ , respectively. The regularized parameter is set to $\delta=0$ for all calculations.
$\delta$

(b)
$-3$
$x$

$-2$

$-1$

$0$

$x$

1

2

3
$x$

Fig. 6: Interfacial profiles for various Atwood numbers; (a) $A=0$ and $t=4.9,$ $(b)A=0.5$ and$t=11.5$ , and (c) $A=1.O$ and $t=1.1$ , where $\sigma=0.05$ for all calculations. The circled regions in(a) and (b) do not stick perfectly.

51
regions in $(a))$ . The closest two points of the interface do not stick perfectly at all pinching points
as found in Fig. 3, however, the gap between two points is narrower than that for $A=0.2$ . The
profile (b) for $A=0.5$ is similar to the one in Fig. 5 (b), however, the gap between the closest two
points of the circled region is wider than the gap for $A=0.2$ .

4Linear stability, stable oscillatory motion in RichtmyerMeshkov instability and unstable motion in RaighleyTaylor instability

4.1

Linear stability analysis

When the surface tension coefficient becomes still larger, motion of the interface becomes stable.
Linear stability analysis suggests this fact. The linear growth rate of this system is calculated byequation (4) and the kinematic boundary condition,$\frac{\partial\eta}{\partial t}+\frac{\partial\phi_{i}}{\partial x}\frac{\partial\eta}{\partial x}=\frac{\partial\phi_{i}}{\partial y}$$(i=1,2)$(15)

where $y=\eta(x, t)$ is the deviation of the interface. Within the linear approximation, the deviation
is given as$\eta$$\eta\propto e^{i(\omega t-x)}$

,

(16)

in our normalization, in which$\omega=\sqrt{\frac{\sigma}{2}-Ag}$

,

(17)

[if we do not normalize the length with wave number , Eq. (16) and Eq. (17) are given as
and $\omega=(\sigma k^{3}/2-Agk)^{1/2}$ , respectively].
Equation (16) indicates that when the surface tension exists, the well-known linear growthrate in the RM instability $(g=0)$ which is proportional to time $t1,12,26$ ) does not hold any
longer and motion is stable within the linear stability analysis. When the surface tension is finite
but sufficiently small, the inertial force which produces the vorticity and governs the nonlinear
stage is more effective in the system, therefore, the linear solution (16) does not affect in theasymptotic stage. For that case, the interface becomes unstable and the roll-up or pinching
phenomenon appears at the final stage as we found in the previous section. On the other hand,stable oscillatory motion is possible for large surface tension coefficients. We show this stablemotion in Sec. 4.2. As we see from Eq. (17), the system can be unstable when $Ag>\sigma/2$ for $g\neq 0$ ,where we always assume that $Ag>0$ in our calculations. Motion of the interface for this case isshown in Sec. 4.3. When the marginally stable case that $Ag=\sigma/2$ , there appears strange motionthat is at rest for a long time and occurs a weak pinching phenomenon at the final stage. This
critical motion for the RT instability is also described in Sec. 4.3.
$k$$\eta\propto e^{i(\omega t-kx)}$

4.2

Asymptotic oscillatory motion in Richtmyer-Meshkov instability

52
15
Figure 7 shows the growth rate of bubble and spikeat the linear stage up to $t=0.5)$ for various Atwood numbers. We do not use the regularized parameter in the calculation for $\sigma=0$ here, therefore, the deviation fromat $t=0$ observed in)Fig. 4 does not exist in Fig. 7. Each line forrepresents the growth rate for all from
to $\sigma=1.0$ at the Atwood number which is denoted$0\wedge^{--}t0.20.3t0$
in the figure except $A=1.0$ , i.e., the growth rate$-1.50$
does not depend on the surface tension coefficientwhen an Atwood number is fixed. When $A=1.0$ ,slight deviations begin to appear in the growth rate Fig. 7: Growth rate (velocity) of bubble andof spike from the neighborhood of $t=0.2$ for var- spike at the linear stage, where the solid,
ious . We choose $\sigma=0.1$ for $A=1.0$ in the dashed and dot-dashed lines denote $A=0.2$ ,figure, however, the growth rate of spike becomes 0.5 and 1.0 with finite surface tension paramsmaller in its absolute value if we take larger . The eters, respectively. The solid line with circle$\sigma=0$ .
growth rate of bubble is larger for smaller Atwood denotes $A=0.2$ withnumbers, while the growth rate of spike is larger for larger Atwood numbers (in its absolute value)
regardless of the value of surface tension coefficient. This tendency is also observed for $\sigma=013$ ).
When we compare at the same Atwood number $A=0.2$ , the growth rate of bubble for $\sigma=0$
is larger than the one for finite , while the growth rate of spike for $\sigma=0$ is smaller than thatfor finite in its absolute value. This deviation and Eq. (16) suggest that the linear growth rateof an interface with finite does not obey the well-known result in the RM instability that it is
proportional to time $t12,19$ ) at the linear stage. The calculation for $\sigma=0$ with $\delta=0$ breaks
down at $t=0.93$ (see Fig. 5) due to the appearance of curvature singularity, therefore, we cannot
compare two lines of $\sigma=0$ andfor $\delta=0$ up to the fully nonlinear stage such that theroll-up appears. However, the asymptotic growth rate for $\sigma=0$ coincides with that for finite in
the long-time computation when we use the regularized parameter as found in Fig. 4.
From equation (16), we see that the interface in the RM instability $(g=0)$ can oscillate withfrequency. Whenforis small, the inertial force exceeds the surface tension
effect and concentration of vorticity arises to the roll-up of the interface. When is large, thesurface tension stabilizes the system and as a result, a kind of periodic motion appears. Figure
8 (a) and (b) show this oscillation for $A=0.2$ with $\sigma=2.0$ . The spike which continues to growdownward turns over at $t=2.0$ , starts to grow upward and turns over again at the neighborhoodof $t=6.0$ . This oscillation is also confirmed by motion of bubble and spike in (b). The velocitiesof bubble and spike are equal to zero at $t=2.0,5.8$ and 9.5, those correspond to times at which
the amplitudes of bubble and spike become maximum in their absolutes values and the turningover occurs. The maximum sheet strength with various surface tension coefficients is depicted
in Fig. 8 (c), where we take the time step $\Delta t=1.25\cross 10^{-4}$ for calculations of $\sigma<1.0$ and$\Delta t=6.25x10^{-5}$ for calculations of, respectively. As we see from the figure, the periodof oscillation is longer for smaller surface tension coefficients. Time at which the maximum sheet
strength takes the minimum value corresponds to the turning over time of bubble and spike andthe time interval, i.e., the period which the minimum value appears, is almost constant when we
fix a surface tension coefficient.
Every motion for
was stable and sharp concentration in the sheet strength or thegrowth of high wave number components in the Fourier mode of the interfacial amplitude, which
are signs whether the calculation breaks down or not, were not found. We add that the frequencyof oscillation is little affected by the Atwood number, i.e., we obtain ‘similar’ motion with the same
frequency for larger Atwood numbers for a fixed surface tension coefficient when it is sufficiently$($$\overline{>-}$$\frac{g}{\Phi}0_{-}5$$\delta$$\pm 1$$\underline{\geq\omega}$$018Sbbi^{-A-}\overline{|}\sim^{-A^{-02(0^{-- 0_{0.4}}}}e^{--------}\ovalbox{\tt\small REJECT}_{---\sim\sim--}^{u_{\theta k_{6}^{--A^{-02}}}}A---05$$\underline{\varpi}$$\sigma\neq 0$$\Sigma^{-0.5}$$\sigma=0.05\ovalbox{\tt\small REJECT}$$\sigma$$\Leftrightarrow-1$$T_{\dot{t}}met$$\sigma$$\sigma$$\sigma$$\sigma$$\sigma$$\sigma\neq 0$$\sigma$$\delta$$\omega=\sqrt{\sigma}/2$$\sigma\neq 0$$\sigma$$\sigma$$\kappa$$\sigma\geq 1.0$$\kappa$$\sigma\geq 0.2$$\kappa$

53
large $(\sigma\geq 0.2)$ , where ‘similar’ indicates that time at which the turning over of an interface
occurs is almost the same, however, the profile of the interface or sheet strength is different foreach Atwood number. The oscillatory motion in the RM instability is not so regular as found in
the standing wave solution $(Ag<0)$ for the RT instability 20) and the frequency of oscillation is
not constant even though when $A=0$ . The period of oscillation depends on the surface tension
parameter and it becomes shorter as the surface tension parameter becomes larger.
$\kappa$

4.3

Unstable and critical motion in Raighley-Taylor instability

In this subsection we present some interfacial profiles for the RT instability for the comparison with
the RM instability. Initial conditions for the RT
instability $(g\neq 0)$ are taken that$x(\beta, 0)=0$ ,$y(\beta, 0)=-0.l\cos\beta$ ,$\gamma(\beta, 0)=0,(18)$for various , where the gravity is normalized aswith the wave number and the linear
growth rate in the system 13).
Figure 9 (a) shows the temporal evolution of
the interface for $A=0.2$ with $g=10$ and $\sigma=1.0$ ,where the time step is taken that $\Delta t=6.25\cross 10^{-5}$and the refinement from $N=1024(0\leq t\leq 4.0)$
to $N=2048(4.0\leq t\leq 5.0)$ is performed. We seethat the roll-up appears in spite of that the surface tension coefficient is sufficiently large, where
the circled regions which the pinching phenomenon
appears do not stick perfectly. Motion in the RM
instability with the corresponding parameters excluding is stable oscillation depicted in Fig. 8 (c)and the roll-up or the pinching phenomenon does
not appear. This suggests that the gravity is effective for concentration of vorticity. The roll-up
gradually disappears as the gravity becomes small
when the Atwood number and the surface tension
coefficient are flxed.
Figure 9 (b) shows interfacial profiles for $A=$
1.0 with $g=0.5$ and $\sigma=0.5$ , where the time
step
is taken that $\Delta t=1.25\cross 10^{-4}$ and therefinement from $N=$ 1024 $(0\leq t\leq 12.0)$ to$N=2048(12.0\leq t\leq 13.0,$ $t=13.0$ is the breakdown time here) is performed. The Atwood number$A=1.O$ with finite
indicates that the consideredsystem is a drop rather than a vortex sheet. The
single-valuedness of the interfacial shape breaks at$t=13.0$ , however, the roll-up does not appear.
Generally, the roll-up does not appear for $A=1.0$
regardless of the value of (including $g=0$ ) or thesurface tension coefficient. When we rise the value
$g$

$g$

$garrow g/(kv_{lin}^{2})$

$g$

$g$

$\Delta t$

$g$

$g$

Fig. 8: Stable oscillatory motion in the RM
instability for $A=0.2$ with large $\sigma;(a)$ theinterfacial profiles and (b) the growth rate
(velocity) of bubble and spike for $\sigma=2.0$ ,where the solid and dashed lines in (b) denote bubble and spike, respectively, and (c)for various
the maximum sheet strengthsurface tension coefficients $(\sigma=0.05-2.0)$ .
$\kappa$

54
of fixing the Atwood number and the surface tension parameter, the spread into the lateral direction of the ‘head’ part at the final stage graduallydisappears as found in the line at $t=10.0$ in the figure, and the breakdown of computation occurs
at earlier time.
The profiles presented in Fig. 9 are examples
of unstable motion which satisfies $Ag>\sigma/2$ in, interfacial motion is staEq. (17). Whenble oscillation and the oscillation can be regarded asa kind of standing wave which appears for $Ag<0$ .
When the critical situation that $Ag=\sigma/2$ , there
appears motion which does not almost move for along time and grows rapidly at the final stage. We
show this critical motion in Fig. 10, where the time
step is $\Delta t=1.25\cross 10^{-4}$ and we refine $N=512$$(0\leq t\leq 28.0)$ to $N=1024(28.0\leq t\leq 32.0)$ ,and further to $N=2048(32.0\leq t\leq 34.2)$ . Theinterface and the sheet strength (therefore, also
the circulation ) do not change their initial shapes
Eq. (18) up to the neighborhood of $t=28.0$ , and after this long linear stage, the interface grows rapidly
as found in Fig. 10 (a) and (b). The maximum sheet
strength or the growth rate of bubble and spikeat the linear stage grows with some oscillation little by little as found in (c). We take parameters$A=0.2,$ $g=5.0$ and $\sigma=2.0$ here, however, suchmotion stated above appears for any parameterswhich satisfy the relation $Ag=\sigma/2$ , although theinterval of its linear stage differs by the value of .
Generally, the interval of the linear stage is longer
for smaller .
$g$$\mathcal{A}g<\sigma/2$$\kappa$$\Gamma$$\kappa$

$g$

$g$

5

Discussions and conclusion

We have investigated motion of the density stratified interface with surface tension. Calculations
performed here have spectral, i.e., exponential ac-

curacy and we can regard obtained results as exact solutions of the Euler equation within machineaccuracy. We presented with the accuracy thatthe roll up of the interface can occur for relatively
small Atwood numbers without the regularized parameter in both of the RM and RT instabilities.
We have also shown that the pinching phenomenon
caused by a vortex pair induced on the interface
appears for both instabilities and the gap between
closest two points becomes wider as the Atwood
number becomes large.
For large surface tension coefficients, we have$\delta$

Fig. 9: Temporal evolution of an interface in
the Rayleigh-Taylor instability for (a) $A=$
0.2, $g=10$ and $\sigma=1.0$ , and (b) $A=1.0$ ,$g=0.5$ and $\sigma=0.5.$ ’ where the solid, dashed,
dot-dashed and solid with circle lines in (a)denote $t=2.8,3.6,4.0$ and 5.0, respectively,
while those lines in (b) denote $t=6.0,8.0$ ,10.0 and 13.0, respectively.

55
shown that there exists stable oscillatory motion for the RM instability. Such oscillatory motion
is also possible for the RT instability if the relation that $Ag<\sigma/2$ is satisfied for the Atwood
number , the gravity and the surface tension , however, the motion is more regular than the
one for the RM instability, as found in the standing wave motion for $Ag<0$ . When, wehave found marginally stable motion which needs a very long time in order to reach the nonlinear
growth appears. This marginally stable motion is possible for any parameter values of the Atwood
number , the gravity and the surface tension , and the length of the linear stage depends onthe value of the gravity . We mention here that the oscillatory motion in the RM instability is
not linear motion and the stability of the system is determined due to the competition between
the surface tension and the inertial force which effect should appear in the higher order. Detailed
nonlinear analysis will be presented elsewhere.
$A$

$g$

$\sigma$$\mathcal{A}g=\sigma/2$

$A$

$g$

$\sigma$$g$

Acknowledgments
The author thanks Professor K. Nishihara for his fruitful discussions. The author also thanks
Professor R. Krasny and Professor T. Sakajyo for useful comments and advice for numericalcalculations. This work is partially supported by the Laser Engineering, Osaka University.

REFERENCES
1)

R.D. Richtmyer, Comm. Pure. Appl. Math. 13, 297 (1960), E. E. Meshkov, Fluid Dyn. 4, 101(1969).

2)

A.L. Velikovich, Phys. Fluids 8, 1666 (1996).

3$)$

A. Celani, A. Mazzino and L. Vozella, Phys.Rev.Lett. 96, 134504 (2006), W. Cabot and A.
W. Cook, Nature Phys. 2, 562 (2006), M. Chertkov, Phys.Rev.Lett. 91, 115001 (2003).

4

$)$

5

$)$

6

$)$

7

$)$

J. Eggers, Rev. Mod. Phys. 69865 (1997) [references therein].

8$)$

T. Kadono and M. Arakawa, Icarus 173295 (2005).

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

W. D. Arnett, J. N. Bahcall, R. P. Kirshner and S. E. Woosley, Ann. Rev. Astron. Astrophys.
27, 629 (1989).

G. B. Whitham, in linear and Nonlinear Waves, (John Wiley&Sons, 1999).
H. Okamoto and M. Shoji, in The mathematical theory of permanent progressive water-waves,(World Scientific, 2001).

10) M. A. Jones and J. W. Jacobs, Phys. Fluids 9, 3078 (1997).

11) P. R. Chapman and J. W. Jacobs, Phys. Fluids 18, 074101 (2006).
12) C. Matsuoka, K. Nishihara and Y. Fukuda, Phys. Rev. $E67$ , 036301 (2003), 68, 029902(E)(2003).

13) C. Matsuoka and K. Nishihara, Phys. Rev. $E73$ , 026304 (2006) [references $therein|,$ $74$ ,049902(E) (2006).

5614) T. Y. Hou, J. S. Lowengrub and M. J. Shelley, J. Comput. Phys. 114, 312 (1993).
15) T. Y. Hou, J. S. Lowengrub and M. J. Shelley, Phys. Fluids 9, 1933 (1997).
16) D. Leppinen and J. R. Lister, Phys. Fluids 15, 568 (2003).

17) D. W. Moore, Proc. Roy. Soc. London, Ser. A 365, 105 (1979).
18) C. Matsuoka and K. Nishihara, Phys. Rev. $E73,055304(R)$ (2006).
19) C. Matsuoka and K. Nishihara, Phys. Rev. $E74$ , 066303 (2006).

20) H. D. Ceniceros and T. Y. Hou, Math. Comput. 67, 137 (1998).

21) P. G. Saffman, in Vortex Dynamics, (Cambridge, 1992), pp. 141.
22) G.R. Baker, D. I. Meiron and S.A. Orszag, J. Fluid Mech. 123, 477 (1982).
23) J. T. Beale, T. Y. Hou and J. Lowengrub, SIAM J. Numer. Anal. 33, 1797 (1996).

24) A. Sidi and M. Israeli, J. Sci. Comput. 3, 201 (1988).
25) G. Baker and A. Nachbin, SIAM J. Sci. Comput. 19, 1737 (1998).

26) J. G. Wouchuk and K. Nishihara, Phys. Plasmas 3, 3761 (1996), 4, 1028 (1997).

27) R. Krasny, J. Fluid Mech. 184, 123 (1987).

57

Fig. 10: Temporal evolution of an interface(a), sheet strengthand the maximum
strength
sheet
at the initial stage (c) for $A=$
0.2, $g=5.0$ and $\sigma=2.0$ , where the solid,
dashed, dot-dashed and solid with circle lines
in (a) and (b) depict $t=22.0,28.0,32.0$ and34.2, respectively.
$\kappa(b)$