Compressible Rayleigh–Taylor instabilities in supernova remnants

PHYSICS OF FLUIDS

VOLUME 16, NUMBER 12

DECEMBER 2004

Compressible Rayleigh–Taylor instabilities in supernova remnants X. Ribeyrea) and V. T. Tikhonchuk Centre Lasers Intenses et Applications, UMR 5107 CNRS-Université Bordeaux 1-CEA, Université Bordeaux 1, 351, Cours de la Libération, 33405 Talence Cedex, France

S. Bouquet Commissariat à l’Energie Atomique, DIF/Département de Physique Théorique et Appliquée, 91680 Bruyères le Chaˇtel, France

(Received 18 March 2004; accepted 3 September 2004; published online 11 November 2004) Partial modeling of the hydrodynamic evolution of the supernovae is one of the prominent applications of laboratory astrophysics. In particular the role of Rayleigh–Taylor instability (RTI) in supernova evolution needs to be explained. In this paper we analyze the compressibility effects on the RTI in the linear regime. We compare the compressible isothermal and stratified incompressible RTI growth rates and analyze the vorticity generation at the interface. We show that for several configurations the effect of compressibility can be significant in supernovae remnants. © 2004 American Institute of Physics. [DOI: 10.1063/1.1810182] I. INTRODUCTION

a hydrostatic quasiequilibrium in a constant deceleration field g, which is equivalent to a gravitational field in the rest frame. The characteristic height h of this atmosphere is ⬃K / g, where 冑K is the isothermal sound speed of the material. To take into account this scale length h in the RTI evolution one must consider the equations of compressible hydrodynamics.3,11,12 One needs to choose an appropriate value of the adiabat parameter for hydrodynamic modeling of astronomical objects. Numerical simulations of the supernovae evolution show that the compressibility is one of dominant parameters in the RTI growth rate.5 The recent work on the x-ray and energetic particles (cosmic-rays) emission in young supernova remnants allows compression ratios at shock13 fronts much higher than the standard Rankine–Hugoniot14 value, 4 for the adiabat index 5 / 3. Hence the model of an adiabatic ideal gas is not appropriate for such objects. One has to take care of these properties in the laboratory experiments and to use materials with relevant equations of state. The quasiisothermal shocks, with the polytropic coefficient close to one, allow higher compression ratios and they might provided more realistic models for young supernova remnants.15 Former studies of the compressibility effects on the RTI growth rate for isothermal fluids are inconclusive. According to Refs. 16 and 17, the compressibility has a destabilization effect and it accelerates the RTI growth. This is, however, in a contradiction with other works18,19 where the authors claim that the compressibility has a stabilizing influence on the RTI. Moreover, other papers20–22 show that the compressibility effect depends on the ratio between the sound velocities of two fluids, which, nevertheless, cannot be arbitrarily chosen but depends on the equations of state used. The origin of this controversy is partially due the fact that the comparison has not been achieved appropriately. For example, in Ref. 18 the compressible case is compared with the classical RTI (incompressible uniform fluid), while in Ref. 17 the compressible behavior is compared with the incompressible stratified case.

Present day high power/energy lasers open a possibility to create conditions similar to those found in several astronomical objects (supernovae, supernova remnants, gammaray bursts, giant planets, etc.). Laser experiments1,2 allow to improve the models describing specific astronomical objects and to validate the codes simulating processes deduced from observations. In this paper, we focus our interest on hydrodynamic evolution of supernovae (and their remnants) and, in particular, on the Rayleigh–Taylor instabilities (RTIs) that are known to play an important role in these phenomena. In 1975, Gull3 studied the formation of the optical filaments in young supernova remnants and attributed them to the RTI. The numerical simulations related to the supernova SN 1987A show that RTI has a major influence on the evolution of the exploding progenitor. Due to an inhomogeneity in the progenitor composition, the RTI arises at both the hydrogen– helium and helium–heavy metal transition zones.4 It might explain a large amount of mixed material (core and envelope mixture) which is necessary to explain the shape of the light curve5–7 of this supernova. The RTIs are also of a great concern in the laser fusion experiments.8 These instabilities are triggered by nonuniformities in laser irradiation or by the imperfections of target geometry. Similarity of RTI development suggests that laser experiments could be used to validate the astrophysical simulation codes. One can also reproduce certain specific features of astronomical objects with laser targets by using the scaling laws derived from the equations of ideal hydrodynamics in Refs. 9 and 10. In particular, these scalings account for compressibility effects, which play an important role in astrophysics. For supernova remnants evolution, the dense shell behind the contact discontinuity can be considered as being in a)

Electronic mail: [email protected]

1070-6631/2004/16(12)/4661/10/$22.00

4661

© 2004 American Institute of Physics

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

4662

Phys. Fluids, Vol. 16, No. 12, December 2004

Ribeyre, Tikhonchuk, and Bouquet

The goal of this paper is to study the RTI in several configurations and to compare them in terms of dimensionless parameters. We consider a compressible fluid, an incompressible stratified fluid, and incompressible uniform finite mass fluids with the same characteristic vertical scale length, which is given by the sound speed in each medium in the compressible case. All these cases can be characterized by a single parameter, i.e., the Atwood number. A relevant comparison of the growth rates of RTI and the corresponding velocities can be performed. II. LINEAR THEORY OF THE RAYLEIGH–TAYLOR INSTABILITY A. Euler equations

It has been shown9,10 that the transport effects, viscosity, and thermal conductivity, are unimportant in the linear stage of the RTI arising in both supernovae explosions and in laser targets evolutions. The starting point of our analysis is, therefore, the set of equations of ideal hydrodynamics that are given below

⳵t␳ + div共␳vជ 兲 = 0,

共1兲

⳵tvជ + 共vជ · grad兲vជ = −

1 grad p + gជ , ␳

p = K␳ ,

共2兲 共3兲

where ␳, vជ , and p are, respectively, the density, the velocity, and the pressure. We restrict in this paper to an isothermal model where the fluid sound velocity C0 = 冑K is constant. This approximation is appropriate for many astronomical objects, for example, for the Kepler’s supernova remnants13 where the temperature near the contact discontinuity is a slowly varying function of the coordinate and time. B. Geometry of the problem

The compressible RTI is studied in a plane geometry with a gravitation field (or more generally, acceleration field) gជ directed along the z axis. A small amplitude perturbation of the interface z = 0 is considered in the plane perpendicular to the y axis: zint = ␩1共x , t兲 (Fig. 1). In the upper region, z ⬎ zint, all physical quantities are labeled with the superscript ⫹, i.e., vជ +, ␳+, p+, and for the lower region, z ⬍ zint, we use the superscript ⫺, i.e., vជ −, ␳−, p−. The initial state is the equilibrium which is given by vជ±0 = 0,

␳±0 共z兲,

p±0 共z兲,

共4兲

where p+0 共0兲 = p−0 共0兲. The perturbed quantities are designed by the subscript 1: vជ±1 共x,z,t兲,

␳±1 共x,z,t兲,

p±1 共x,z,t兲.

C. Compressibility and finite mass effects

In order to examine the compressibility effect on the RTI, it is important to compare the characteristics that are comparable for two different systems. For instance, we con-

FIG. 1. Geometry of the problem in two dimensions: a small perturbation ␩1共x , t兲 is applied at the interface z = 0 at t = 0.

sider the same density profiles and, consequently, the same masses for the two superimposed fluids from one system to another. We choose four specific initial states which are unstable with the respect to RTI. Three of them correspond to incompressible cases with various density profiles and the fourth one is compressible (Fig. 2). In the case of the classical incompressible medium, panel (A) in Fig. 2, the masses of the two fluids are infinite and the density profiles are uniform. In the compressible isothermal case [panel (B) in Fig. 2] the density profiles are exponential and the mass of the upper fluid is finite. In panel (C) of Fig. 2, a incompressible stratified configuration is chosen. The density profiles are identical to the ones in panel (B). Finally, in panel (D), a finite upper mass incompressible case is presented. The mass of the upper medium is chosen equal to the then upper mass for both cases (B) and (C). We will compare the development of the instability for each state. First goal, it should be pointed out that cases (A) and (B) are difficult to compare, because they have neither the same fluid mass nor the same density profiles. Our first step is therefore, to compare the compressible case (B) with the case (D) which is the incompressible finitemass case. Although, these systems have the same masses, they do not have the same density profiles. The second step corresponds to the comparison between case (B) and the incompressible stratified case (C). Since we have the same upper masses and the same density profiles, this comparison is actually the most relevant.

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

Phys. Fluids, Vol. 16, No. 12, December 2004

Compressible Rayleigh–Taylor instabilities

⳵tv1x = − K⳵x␦1 − g⳵x␩1,

⳵tv1z = − K⳵˜z␦1 .

4663

共7兲

Note that due to the jump conditions at the interface, p+共0兲 ± = p−共0兲 and ⳵t␩1 = v1z 共0兲 (A10), the constant K takes different values up and down from the interface. These equations have to be solved with the evanescent boundary conditions: ␦1 → 0 for ˜z → ± ⬁. The assumption of small perturbations allows one to use the superposition principle and to describe the perturbed fluid motion as an ensemble of linearly independent modes. We assume that any perturbed quantity f evolves in time and along the x axis proportionally to e␥t+ikx. Then for the ampli˜ , ␥ , k兲 one has tude of the density perturbation ¯␦1共z g − K⳵˜z2¯␦1 + 共␥2 + Kk2兲¯␦1 + g⳵˜z¯␦1 = − a共␥2 + Kk2兲, K

FIG. 2. Density profiles for four cases: (A)—incompressible uniform fluids (infinite mass above and below), (B)—compressible fluids, (C)— incompressible stratified fluids [finite mass, same characteristic as case (B)], (D)—incompressible finite mass layer above [same mass as in case (B) and (C)] and infinite mass below.

where a is the amplitude of the interface perturbation. Since the right-hand side of this equation does not depend on ˜, z one has to seek the solution in the following form: ¯␦ = −a / h + Semc˜z. Then we have an equation for m , 1 c Km2c − gmc = ␥2 + Kk2 , which has two solutions mc1,2 =

D. Linearized hydrodynamic equations

We have to derive the linear equations for the first order quantities. We cast the system of equations (1)–(3) in the new coordinate frame associated with the moving interface, zint = ␩1共x , t兲 :˜z = z − zint = z − ␩1共x , t兲, as it was proposed in Ref. 23. Since the function ␩1 accounts for the deformation of the interface, the condition ˜z = 0 is always valid on the interface. This procedure strongly facilitates the utilization of the jump conditions at the interface. Appendix A provides the details this transformation of coordinates.

The compressible barotropic isothermal fluid corresponds to case (B) in Fig. 2. The system (1)–(3) according to (4) at the zero order (equilibrium) reduces to relation between the pressure and the density: d˜z p0 = − ␳0g

共5兲

and also p0 = K␳0. This provides the exponential density profiles at the equilibrium:

冑

1 ␥2 2 +k . + 4h2 K

One can eliminate one constant from the evanescent boundary conditions. To satisfy this condition mc must be real and nonzero. Moreover, for the upper half space, ˜z ⬎ 0, one has m+c ⬍ 0, and for the lower half space m−c ⬎ 0: g mc±˜z ¯␦±共z − ± a, 1 ˜兲 = S±e K

共8兲

where m±c =

1. Compressible barotropic isothermal fluid

1 ± 2h

g ⫿ 2K±

冑冉 冊 g 2K±

2

+

␥2 + k2 . K±

The relation between the remaining two constants follows from the jumping condition at the interface. From Eq. ± 共0兲 one has (7) and the condition ⳵t␩1 = v1z S± = −

␥2 a. K±m±c

共10兲

From these relations one can reconstruct the velocity and the pressure profiles:

±

␳±0 共z˜兲 = ␳±0 共0兲e−z˜/h ,

±

where h± = K± / g is the characteristic height of system. At the first order, the hydrodynamic equations can be reduced to one equation for the relative density perturbation, ␦1 = ␳1 / ␳0, in the x, ˜z plane:

⳵2t ␦1

−

K共⳵2x

+

⳵˜z2兲␦1

g = − g⳵˜z␦1 − ⳵2t ␩1 + g⳵2x ␩1 , K

共6兲

共9兲

± ¯v1z ˜兲 = ␥aemc˜z, 共z

± ¯v1x ˜兲 = ia 共z

k␥ m±˜z e c , m±c

¯p±1 共z ˜兲 = K±␳±0 共z ˜兲¯␦±1 共z ˜兲. To find dispersion relation we use the continuity condition for the pressure perturbation: p+1 共0兲 = p−1 共0兲. By using Eqs. (8) and (10) we find the dispersion equation:17

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

4664

Phys. Fluids, Vol. 16, No. 12, December 2004

␥2 = g

␳±0 共0兲 − ␳−0 共0兲 . − ␳+0 共0兲/m+ + ␳−0 共0兲/m−

Ribeyre, Tikhonchuk, and Bouquet

共11兲

The unstable solutions correspond to ␥ positive. The perturbations of the pressure and velocity can be expressed in terms of the amplitude a of the interface perturbation ␩1共x , t兲 = ae␥t cos kx. The velocity perturbations are ± ˜,t兲 共x,z v1x

k ␥ ±z ␥t = − a ± emc˜+ sin kx, mc

± ˜,t兲 = ± a␥e⫿kz˜+␥t sin kx, 共x,z v1x ± ˜,t兲 = a␥e⫿kz˜+␥t cos kx, 共x,z v1z

冉

˜,t兲 = − ␳±0 a ⫿ p±1 共x,z

共17兲

冊

␥2 ⫿kz˜ e + g e␥t cos kx. k

The velocity field is incompressible and irrotational.

共12兲 3. Incompressible stratified fluid

± ˜,t兲 共x,z v1z

= a␥e

mc±˜+ z ␥t

共13兲

cos kx,

and the density and pressure perturbations are ± ˜,t兲 = − ␳±0 共z ˜兲 ␳1x 共x,z

冉

冊

a ␥2 m±˜z e c + g e␥t cos kx, K± m±c

˜,t兲 = − ␳±0 共z ˜兲a p±1 共x,z

冉

冊

␥2 m±˜z e c + g e␥t cos kx. m±c

Notice that this case corresponds to the irrotational flows, curl vជ 1 = 0. This can be verified directly from the formulas for the velocity perturbations (12) and (13).

2. Incompressible case

The limit of incompressible fluid can be obtained directly from the hydrodynamic equations by setting formally K → ⬁. This is the case (A) in Fig. 2. The system (1)–(3) under the conditions [Eq. (4)] and ␳ = const gives at zero order Eq. (5). At first order the system reduces to the equations for the pressure and velocity perturbations;

⳵xv1x + ⳵˜zv1z = 0, ␳0⳵tv1x = − ⳵x p1 − ␳0g⳵x␩1,

共14兲

␳0⳵tv1z = − ⳵˜z p1 .

共15兲

These equations also follow from Eqs. (6) and (7) in the limit K → ⬁, if we account for the fact that ␳1 → 0 and the quantity K␦1 = p1 / ␳0 remains finite. Note that Eq. (6) in the incompressible limit reduces to Bernoulli’s equation 共⳵2x + ⳵˜z2兲p1 + ␳0g⳵2x ␩1 = 0. The dispersion relation can be derived from these equations. It can be also obtained directly from Eq. (11). By setting K± → + ⬁ in Eq. (9) one finds m±c = ⫿ k and correspondingly

␥2 = kg

␳+0 共0兲 − ␳−0 共0兲 . ␳+0 共0兲 + ␳−0 共0兲

共16兲

This is the classical Rayleigh–Taylor growth rate for an incompressible uniform fluid. The spatial and temporal behavior of physical parameters follows from the compressible case Eqs. (12) and (13) in the limit, K± → + ⬁. Then ␳±1 共x ,˜z , t兲 = 0 and for other parameters one has

In this case we consider a fluid which is incompressible, but with the density that changes with the coordinate ␳0共z兲, at the equilibrium Eq. (4). This is the case (C) in Fig. 2. It cannot be reduced to the previous cases, because it is not compatible with the isothermal equation of state (3). The equation of state in this case is substituted by the condition of incompressibility, div vជ = 0. Then Eqs. (1) and (2) in the first order in the coordinate system associated with the perturbed interface take the following form:

⳵t␳1 + v1zd˜z␳0 = ⳵t␩1d˜z␳0 ,

共18兲

␳0⳵tv1x + ␳0g⳵x␩1 = − ⳵x p1 ,

共19兲

␳0⳵tv1z = − g␳1 − ⳵˜z p1 ,

共20兲

⳵xv1x + ⳵˜zv1z = 0,

共21兲

˜兲 is considered as a given function. Assuming the where ␳0共z harmonic perturbations ⬀e␥t+ikx, we have

␥¯␳1 + ¯v1zd˜z␳0 = ␥ad˜z␳0 ,

共22兲

¯1, − i␳0␥¯v1x + k␳0ga = − kp

共23兲

␳0␥¯v1z = − ¯␳1g − ⳵˜z¯p1 ,

共24兲

ik¯v1x + ⳵˜z¯v1z = 0.

共25兲

With Eqs. (22), (23), and (25) one finds the x component of the velocity, the pressure, and the density perturbations: i ¯v1x = ⳵˜z¯v1z, k

冉

¯p1 = −

␥ ␳0⳵˜z¯v1z − ␳0ga, k2 共26兲

冊

¯v ¯␳1 = ␤␳0 1z − a , ␥ where ␤ = −d˜z ln ␳0 is the stratification gradient. Inserting these relations in Eq. (24) one finds an equation for the z component of velocity:

冉

⳵˜z2¯v1z − ␤⳵˜z¯v1z − k2 1 +

冊

␤g ¯v = 0. ␥2 1z

共27兲

This equation agrees with the one derived before in Refs. 12 and 24. We consider here the case of a constant density gradient, ␤ = const. Then one can search for a solution in the exponential form ems˜z, and can find two roots. The appropriate solution can be chosen by applying the evanescent

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

Phys. Fluids, Vol. 16, No. 12, December 2004

Compressible Rayleigh–Taylor instabilities

± boundary conditions at infinity. Then the solution reads, ¯v1z ± = V±ems ˜z, where

ms±

␤± ⫿ = 2

冑

冉

冊

␤±2 2 ␤ ±g +k 1+ 2 . 4 ␥

共28兲

The condition on the interface is the same as for the com+ − 共0兲 = ¯v 1z 共0兲 = ⳵t␩1 Eq. (A10). That gives pressible case: ¯v 1z + − V = V = ␥a, hence ¯v1z = ␥ae

ms±˜z

共29兲

.

Calculating now the pressure perturbation from Eq. (26) and applying the continuity condition at the interface, ¯p +1 共0兲 = ¯p −1 共0兲, we find the dispersion relation

␥ 2 = k 2g

−

␳+0 共0兲 − ␳−0 共0兲 . ␳+0 共0兲ms+ + ␳−0 共0兲ms−

共30兲

In the limit of a homogeneous density ␤ → 0, from Eq. (28) one finds ms± = ⫿ k and the dispersion equation (30) reduces to the classical case Eq. (16). In Appendix B we demonstrate that ms+ = −ms− = −ms whatever ␥ and k be. Then we can rewrite the dispersion relation17 as

␥2 =

k2g ␳+0 共0兲 − ␳−0 共0兲 . ms ␳+0 共0兲 + ␳−0 共0兲

共31兲

In order to compare the incompressible stratified case with the compressible case we have to assume that the density profiles are the same, that is, ␤± = 1 / h± = g / K±. The results of such a comparison are presented in the following section. All physical parameters can be reconstructed from Eqs. (26) and (29) (␥ is assumed to be real and positive):

␥

+

+

+ ¯v1z ˜兲 = V+共e−k共z˜−h 兲 − ek共z˜−h 兲兲, 共z − ¯v1z ˜兲 = V−ekz˜ . 共z

A relation between the constants follows from the boundary condition at the interface: 2V+ sinh kh+ = V− = ␥a. Then from Eq. (26) one finds the pressure perturbation and assuming its continuity at the interface, ¯p +1 共0兲 = ¯p −1 共0兲, one obtains the dispersion relation for the incompressible finitemass case:

␥2 = kg

␳+0 共0兲 − ␳−共0兲 . ␳+0 共0兲tanh−1 kh+ + ␳−0 共0兲

共32兲

Evidently that in the limit of the thick layer, h+ → + ⬁, this equation reduces to the classical Rayleigh–Taylor form Eq. (16). From Eq. (26) we can find now the spatial and temporal profiles of all physical quantities. In the lower fluid the solutions are the same as given in Eq. (17). For the upper fluid we have + ˜,t兲 = ␥a 共x,z v1x

z ␥t cosh k共h+ − ˜兲 e sin kx, + sinh kh

+ ˜,t兲 = ␥a 共x,z v1z

z ␥t sinh k共h+ − ˜兲 e cos kx, + sinh kh

˜,t兲 = a␳+0 p+1 共x,z

冉

冊

z ␥2 cosh k共h+ − ˜兲 − g e␥t cos kx. sinh kh+ k

± z ␥t ˜,t兲 = ⫿ amse±ms˜+ 共x,z sin kx, v1x k

One sees that this is an irrotational flow.

± z ␥t ˜,t兲 = a␥e±ms˜+ 共x,z cos kx, v1z

III. ANALYSIS OF THE DISPERSION RELATIONS

␳±1 共x,z˜,t兲 = a␳±0 ␤±共e±ms˜z − 1兲e␥t cos kx,

冉

˜,t兲 = − a␳±0 ± p±1 共x,y,z

冊

␥2 mse±ms˜z + g e␥t cos kx. k2

Notice that this flow is no more irrotational. The instability results in the vorticity generation as long as there is a density gradient, ␤ ⫽ 0. 4. Incompressible and finite mass fluid

Another possibility to compare with the compressible case is to consider a homogeneous but finite mass fluid. This is the case (D) in Fig. 2. We choose the thickness h of the upper fluid in such a way that it has the same mass as the fluid in the compressible case (B). Also they have the same density at the interface. It is easy to see that h+ = K+ / g. To find the dispersion equation for this case, we use the previous analysis for the incompressible stratified case, but we set the boundary condition for ¯v1z共h+兲 = 0 at the finite altitude ˜z = h+. We also suppress the stratification assuming ␤ = 0, then Eq. (28) can be simplified: ms± = ⫿ k. The solution to Eq. (27) with the correspondent boundary condition reads

4665

We have obtained the dispersion relations for four cases presented in Fig. 2: the incompressible uniform flow (the classical RTI) is described by Eq. (16). The general compressible case is described by Eq. (11). The case of an incompressible and stratified fluid is described by Eq. (30). Finally the case of an incompressible and finite mass fluid is presented in Eq. (32). For all cases the instability criterion is ␳+0 共0兲 ⬎ ␳−0 共0兲 and the equilibrium condition corresponds to the continuity of the pressure: p+0 共0兲 = K+␳+0 共0兲 = K−␳−0 共0兲 = p−0 共0兲. It is convenient to introduce the characteristic sound velocity Kˆ = 共K+K−兲1/4 and normalize the wave number by kˆ

冑

冑

= g / Kˆ and the growth rate by ␥ˆ = g / Kˆ. Then all considered cases depend only on one dimensionless parameter, which is the Atwood number, At =

␳+0 共0兲 − ␳−0 共0兲 . ␳+0 共0兲 + ␳−0 共0兲

Figure 3 shows numerical solutions to the dispersion equations for four chosen cases and for two Atwood numbers. (To solve these dispersion relations we use the numerical package “MAPLE.”) One can see that the compressibility

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

4666

Phys. Fluids, Vol. 16, No. 12, December 2004

Ribeyre, Tikhonchuk, and Bouquet TABLE I. Comparison of the asymptotic RT growth rates for the four fluids models.

Incompressible

␥RT = 冑Atkg

Compressible

␥C = k冑Kˆ冑2At/共1 − At2兲1/4

Compressible stratified

␥IS = k冑KˆAt/共1 − At2兲1/4

Incompressible finite mass

␥IFM = k冑Kˆ冑2At共1 − At2兲1/4/共1 + At兲

fluid. It is evident that the deviation from the ideal RT instability occurs at the wavelengths longer than the characteristic height hc or for the wave numbers smaller than the characteristic wave number kc = 2␲ / hc. Comparing the asymptotic behaviors for the compressible and incompressible stratified cases we find kc =

FIG. 3. Dependence of the normalized growth rate ␥ / ␥ˆ on the normalized wavenumber k / kˆ for four cases and for two Atwood numbers: At = 0.1 (a) and 0.9 (b). Solid line corresponds to the compressible case, dotted line— the incompressible case, dashed line—the incompressible stratified case, and dash-double-dotted line corresponds to the incompressible finite mass case. Dash-dotted line shows the ratio of the growth rates of the compressible and stratified incompressible cases.

effect is important when the density of two fluids are close, At Ⰶ 1. A. Asymptotic expressions for the growth rate

We consider the asymptotic solutions to the dispersion relations in two limits, k → ⬁ and k → 0. In the case of incompressible fluid we rewrite the dispersion relation as

␥RT = 冑Atkg.

共33兲

This is the result which all other dispersion equations produce in the short wavelength limit k → ⬁. On the contrary, in the long wavelength limit, k → 0 and for ␥ / k finite, the results depend strongly on the model. We compare all the models in the Table I. The comparison between these three fluid systems is justified because they have the same density jump (characterized by the Atwood number) and the same mass for the upper

g 冑 1 − A t2, 2Kˆ

␥C = ␥IS

冑

2 . At

共34兲

Thus the compressibility strongly enhances the RTI growth rate for small At Ⰶ 1 compared to the incompressible case if has the characteristic heights are the same. In the contrary, both the compressibility and finite mass effects suppress the development of the RTI compared to the classical RTI. For large Atwood numbers, where the density jump is large, the critical wave number becomes very small and the difference between the compressible and incompressible cases becomes less important. Notice that the growth rate in the finite mass case is independent of the gravity field g for k ⬍ kc. IV. VELOCITY FIELDS GENERATED BY THE RTI

To understand the development of the RT instability it is important to analyze the velocity fields and in particular the vorticity production. The flow velocity fields are shown for the four cases considered above in Fig. 4. These charts are constructed by using the velocity fields derived in the preceding section. It is well known that the classical RT instability in incompressible fluid produces the velocity field without divergence and rotation: div vជ = 0,

curl vជ = 0.

The same holds for the incompressible finite-mass case. The compressible case corresponds to divergent irrotational flow: curl vជ = 0 and div vជ ± =

a␥ m±˜+ e c z ␥t共k2 − m±2 c 兲cos kx, m±c

while the velocity field in the incompressible stratified case has no divergence, div vជ = 0, but has a nonzero vorticity: curl vជ ± =

a␥ ±m ˜+ e sz ␥t共k2 − ms2兲sin kx. k

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

Phys. Fluids, Vol. 16, No. 12, December 2004

Compressible Rayleigh–Taylor instabilities

4667

FIG. 6. Velocity shear at the interface for the compressible case (solid line) and for the incompressible stratified case (dashed line).

FIG. 4. Fluid velocity fields for four cases considered above: (A)—the incompressible uniform media, (B)—the compressible fluids, (C)—the incompressible stratified case, (D)—the incompressible finite mass media. The parameters are following: ␥ˆ t = 1, k = kˆ, kˆa = 0.05, and At = 1 / 3.

In order to better understand the effects of compressibility and stratification it is instructive to represent the continuity and Euler equations (1) and (2) in the following form: 1 div vជ = − dt␳ , ␳

ជ =−␻ ជ div vជ + 共␻ ជ · grad兲vជ + d t␻

1 grad ␳ ⫻ grad p, ␳2

ជ = curl vជ is the vorticity and dt is the full (Lagrangwhere ␻ ian) derivative. One concludes from the first equation that the velocity divergence is generated only due to the fluid compressibility. The source of the vorticity generation is due to the nonparallel gradients of the pressure and the density.

FIG. 5. The flow divergence, div vជ , (A) for the compressible case and the y component of the rotation, curl vជ , (B), for the incompressible stratified case. The parameters are the same as in the previous figure.

Thus in the barotropic fluid, where the pressure is a given function of the density, p = p共␳兲, no vorticity is produced. On the contrary, for the baroclinic fluid the pressure and density are not necessary parallel, the pressure inclination at the bended interface drives a vertical fluid motion.11 This effect is important for the wavelengths comparable to the stratification height. It corresponds to very different behavior at the nonlinear stage. The compressible irrotational flows tend to create long “finger structures” of a heavy fluid, while incompressible flows tend to develop “mushroom structures” with more efficient mixing.25 This difference in the flow velocity fields is illustrated in Fig. 5(a), where we represent the vorticity in the incompressible stratified fluid. We see two vortices rotate in opposite directions localized near the contact discontinuity. Another difference between the compressible and incompressible case concerns the velocity shear at the inter+ − 共0兲 − v1x 共0兲. For the compressible case face: ⌬v1x = v1x ⌬v1x = − ak␥

冉

1 m±c

−

1 m−c

冊

e␥t sin kx,

共35兲

while for the incompressible stratified case

␥ ⌬v1x = − 2msa e␥t sin kx, k

共36兲

the parallel velocities at the interface have the same magnitudes but the opposite signs. Figure 6 shows that the velocity shear is larger for the compressible fluid and achieves its maximum at the inflexion point of the interface. Hence the compressible case can generate more rapidly the Kelvin–Helmholtz instability than the incompressible stratified case. We show that in the case of the stratified incompressible fluid the vorticity is generated. For the other case, the development of the RT instability results in shear flow generation along the contact discontinuity. At the nonlinear stage of the RT instability this shear flow could become a subject of the Kelvin–Helmholtz instability. This would be a secondary source of vorticity.

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

4668

Phys. Fluids, Vol. 16, No. 12, December 2004

V. DISCUSSION AND CONCLUSION

We see that chosen four geometries allow a consistent comparison of the RTI development for the compressible and incompressible cases. A transformation of coordinates to the system associated with the perturbed surface allows to simplify the instability analysis and the dispersion relation can be casted in the same mathematical form for all cases that we studied. The compressibility decreases the RTI growth rate for ˆ / g, comwavelengths larger than the critical length hc ⬃ K pared to the classical RTI. But the compressibility increases ˆ / g if we compare this case the RTI growth rate for ␭ ⬎ hc ⬃ K with a stratified incompressible fluid, which is more appropriate. The compressibility effect depends inversely on the Atwood number and it is more important for the cases where the density jump is moderate. The growth rate is independent of the gravity field strength. Even under the conditions where the growth rates are similar for the compressional and stratified cases, the velocity fields are completely different. We have shown that the vorticity is generated at the linear stage of the instability in the incompressible stratified case but it does not in the compressional case. Also the velocity shear in the incompressible stratified case is smaller in comparison with the compressible case. It would be interesting to understand the manifestations of the divergence and the vorticity production in the nonlinear development of the RT instability. It is possible to estimate the characteristic RTI growth rate for a typical supernova remnant (SNR). We choose M 0 = 15 M 䉺, V0 = 6 ⫻ 108 cm/ s, T0 = 30 eV, and ␳m −24 3 = 10 g / cm and Tm = 1 eV, respectively, for the mass, the velocity, and the temperature of the ejecta, the density and the temperature of the interstellar medium. With these physical quantities it is possible to estimate the deceleration g of the ejecta and the At (Atwood parameter) in function of the radius of the supernova ejecta (Appendix C). For a radius of SNR r = Rs / 2 we estimate the time expansion t = 0.56ts (Appendix C) and we find At = 0.78 and g = 10−2 cm/ s2, where Rs = 2 ⫻ 1019 cm= 6.7 pc (1 pc= 3 ⫻ 1018 cm and M 䉺 = 2 ⫻ 1033 g) and ts = 950 years. Then we can calculate the critical length (34): hc = 1016 cm= Rs / 2000. Hence only perturbations with a very short wavelength, ␭ ⬍ hc, can be considered as incompressible. For perturbations ␭ ⬎ hc the compressibility effect on the RTI growth rate is important. For example, for the perturbation wavelength ␭ = r / 50= 2 ⫻ 1017 cm, we have ␥ct = 1.8 and ␥ist = 1.1. That is, the compressibility dominates the RTI excitation. We see that the RTI can be excited at the early stage of SNR expansion, in the pre-Sedov stage, where the swept mass of the interstellar medium is still smaller than the ejecta mass. The interface instability between the SNR ejecta and the interstellar medium is dominated by the compressibility effect though the stratification might also contribute if the densities of the interstellar medium and the ejecta are inhomogeneous. In this regime the growth time of the RTI is comparable with the characteristic time of the ejecta evolution. Therefore for a more accurate analysis one should consider the RTI development in the nonstationary regime.

Ribeyre, Tikhonchuk, and Bouquet

APPENDIX A: TRANSFORMATION OF THE COORDINATES

Let us consider the function ␩共x , y , t兲 which defines the position of the interface zint = ␩. We introduce a new independent variable ˜z given by ˜z = z − zint = z − ␩共x , y , t兲. The physical quantity f共x , y , z , t兲 is denoted ˜f 共x , y ,˜z , t兲 in the new coordinate system and we have, f共x , y , z , t兲 =˜f 共x , y ,˜z , t兲. Consequently, total derivatives of both functions are the same, df = df˜, and one has the following relations between the partial derivatives:

⳵i f = ⳵i˜f − ⳵i␩⳵˜z˜f ,

i = x,y,t;

⳵z f = ⳵˜z ˜f .

共A1兲 共A2兲

We rewrite the hydrodynamic equations by using these transformations. To simplify the notations we drop the tilde on all physical quantities:

⳵t␳ + div共␳vជ 兲 − ⳵t␩⳵z␳ − ⳵z共␳vជ 兲grad ␩ = 0,

共A3兲

␳⳵tvជ + ␳共vជ · grad兲vជ − ␳共⳵t␩ + vx⳵x␩ + vy⳵y␩兲 − ⳵zvជ − ⳵z p grad ␩ = gជ ␳ − grad p, p = K␳ .

共A4兲 共A5兲

For a static solution at zero order we find: vជ 0 = 0, grad p0 = ␳0gជ , ␩0 = 0, and p0 = K␳0. By linearizing the system around this state one has at the first order

⳵t␳1 + div共␳0vជ 1兲 − ⳵z␳0⳵t␩1 = 0,

共A6兲

␳0⳵tvជ1 − ⳵z p0 grad ␩ = gជ ␳1 − grad p1 ,

共A7兲

p1 = K␳1 .

共A8兲

This linearization is valid as long as the second order terms are small compared to the first order, that is,

⳵tvជ Ⰷ 共vជ · grad兲vជ .

共A9兲

Let a be the amplitude of wave, ␶ be the time period of oscillation of the fluid particles in the wave, and ␭ be the wavelength of oscillation. Then v ⬃ a / ␶, the time derivative of the velocity is ⬃v / ␶, and the space derivative is ⬃v / ␭. The condition (A9) yields a / ␶2 Ⰷ 1 / ␭共a / ␶兲2 or a Ⰶ ␭, i.e., the amplitude of the perturbation must be small compared to the wavelength. The equations in the new coordinate system provide also the jump conditions at the interface ˜z = 0. By integrating Eq. + (A6) across the interface we find: ␳+0 共0兲关⳵t␩1 − v1z 共0兲兴 − − + − = ␳0 共0兲关⳵t␩1 − v1z共0兲兴. Since ␳0 共0兲 ⫽ ␳0 共0兲 the following conditions must be satisfied at the fluid interface: + − 共0兲 = v1z 共0兲 = ⳵t␩1 . v1z

共A10兲

From the z component of the Euler equation (A7) we find the continuity condition for the pressure at the interface: p+0 共0兲 = p−0 共0兲 and p+1 共0兲 = p−1 共0兲 as well.

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

Phys. Fluids, Vol. 16, No. 12, December 2004

Finally, restricting ourselves to a two-dimensional perturbation, zint = ␩共x , t兲 we drop the y dependence write the linearized system of equations in the following form:

⳵t␳1 + ⳵x共␳0v1x兲 + ⳵z共␳0v1z兲 − ⳵z␳0⳵t␩1 = 0,

Compressible Rayleigh–Taylor instabilities

As expected, Rs is the radius of a sphere in the interstellar medium with density ␳m and with the mass equal to the mass of the ejecta. With these equations we find the acceleration g and the Atwood number in terms of the radius:

␳0⳵tv1x − ⳵x␩1⳵z p0 = − ⳵x p1 , ␳0⳵tv1z = − g␳1 − ⳵z p1 .

By inserting ␥ from the dispersion relation (30) in Eq. (28) we can find an expression for ms+:

␤+ − 2

冑

␳+0 ms+ − ␳−0 ms− ␤+2 2 , + k − ␤+ 4 ␳+0 − ␳−0

which can be represented in the following form: ms+2 + ␤+

␳+0 共m+ − ms−兲 = k2 . ␳+0 − ␳−0 s

Similar expression one finds also for ms− ms−2 + ␤−

␳−0 共m+ − ms−兲 = k2 . ␳+0 − ␳−0 s

By subtracting these two equations one has ms+2 − ms−2 = 共ms+ − ms−兲

␤+␳−0 − ␤−␳+ . ␳+0 − ␳−0

Having in mind that ms+ and ms− have opposite signs and that the boundary condition at the equilibrium is ␤−␳+0 = ␤+␳−0 , one has ms+ = −ms− = −ms. Then one can rewrite the dispersion relation (30) in the form (31).

APPENDIX C: SUPERNOVAE REMNANT (SNR) EVOLUTION

We consider a supernova ejecta with total mass M 0 and averaged velocity V0. The ejecta expands into an interstellar medium with the density ␳m. From the law of momentum conservation we have26

冉

M 0V 0 = M 0 +

冊

4␲␳m 3 r v, 3

where r and v are the radius and the averaged velocity of the ejecta at the time t. By using the definition of velocity, v = dr / dt, and integrating by quadratures, we find the SNR velocity and radius as function of time (we suppose that the velocity is uniform in the SNR):

冉 冊

t r 1 r = + ts Rs 4 Rs

4

,

g共r兲 =

3V20 共r/Rs兲2 , Rs 关1 + 共r/Rs兲2兴3

At共r兲 =

1 − 共r/Rs兲3 . 1 + 共r/Rs兲3

The instability of the SNR can develop for r ⬍ Rs where Atwood number is positive.

APPENDIX B: SIMPLIFICATION OF THE DISPERSION RELATION FOR A STRATIFIED FLUID

ms+ =

4669

v=

V0 , 1 + 共r/Rs兲3

where Rs = 共3M 0 / 4␲␳m兲1/3 and ts = Rs / V0.

1

A. Remington, D. Arnett, R. P. Drake, and H. Takabe, “Modeling astrophysical phenomena in the laboratory with intense lasers,” Science 284, 1488 (1999). 2 B. A. Remington, J. Kane, R. P. Drake, S. G. Glendinning, K. Estabrook, R. London, J. Castor, R. J. Wallace, D. Arnett, E. Liang, R. McCray, A. Rubenchik, and B. Fryxell, “Supernova hydrodynamics experiments on the Nova laser,” Phys. Plasmas 4, 1994 (1997). 3 S. F. Gull, “The x-ray, optical and radio properties of young supernova remnants,” Mon. Not. R. Astron. Soc. 171, 263 (1975). 4 T. Ebisuzaki, T. Shigeyama, and K. Nomoto, “Rayleigh–Taylor instability and mixing in SN 1987A,” Astrophys. J. Lett. 344, L65 (1989). 5 E. Müller, B. Fryxell, and D. Arnett, “Instability and clumping in SN 1987A,” Astron. Astrophys. 251, 505 (1991). 6 L. Wang, J. C. Wheeler, P. Höflich, A. Khokhlov, D. Baade, D. Branch, P. Challis, A. V. Filippenko, C. Fransson, P. Garnavich, R. P. Kishneer, P. Lundqvist, R. McCray, N. Panagia, C. S. J. Pun, M. M. Phillips, G. Sonnbo, and N. B. Suntzeff, “The asymmetric ejecta of supernova 1987A,” Astrophys. J. 579, 671 (2002). 7 D. Arnett, Supernovae and Nucleosynthesis, edited by J. P. Ostricker (Princeton University Press, Princeton, 1996). 8 S. Nakai and H. Takabe, “Principles of inertial confinement fusion-physics of implosion and concept of inertial fusion energy,” Rep. Prog. Phys. 59, 1071 (1996). 9 D. Ryutov, R. P. Drake, J. Kane, E. Liang, B. A. Remington, and W. M. Wood-Vasey, “Similarity criteria for the laboratory simulation of supernova hydrodynamic,” Astrophys. J. 518, 821 (1999). 10 D. D. Ryutov, R. P. Drake, and B. A. Remington, “Criteria for scaled laboratory simulations of astrophysical MHD phenomena,” Astrophys. J., Suppl. Ser. 127, 465 (2000). 11 N. A. Inogamov, “The role of Rayleigh–Taylor and Richmyer–Meshkov instabilities in astrophysics: An introduction,” in Astrophysics and Space Physics Reviews, Vol. 10, Part 2, edited by R. A. Sunyaev (Institute for Space Research, Russian Academy of Sciences, Moscow, 1999). 12 S. Chandrasekhar, Hydrodynamic and Hydromagnetic Stability (Oxford University Press, Oxford, 1968). 13 A. Decourchelle, D. C. Ellison, and J. Ballet, “Thermal x-ray emission and cosmic-ray production in young supernova remnant,” Astrophys. J. Lett. 543, L67 (2000). 14 L. D. Landau and E. M. Lifchitz, Fluid Mechanics (Pergamon, New York, 1959). 15 J. M. Blondin and D. C. Ellison, “Rayleigh–Taylor instabilities in young supernova remnants undergoing efficient particle acceleration,” Astrophys. J. 560, 244 (2001). 16 P. O. Vandervoort, “The character of the equilibrium of a compressible, inviscid fluid of varying density,” Astrophys. J. 134, 699 (1961). 17 I. B. Bernstein and D. L. Book, “Rayleigh–Taylor instability of selfsimilar spherical expansion,” Astrophys. J. 225, 633 (1978). 18 G. M. Blake, “Fluid dynamic stability of double radio sources,” Mon. Not. R. Astron. Soc. 156, 67 (1972). 19 D. Sharp, “An overview of Rayleigh–Taylor instability,” Physica D 26, 950 (1984). 20 M. S. Plesset and D. Hsieh, “General analysis of the stability of superposed fluids,” Phys. Fluids 7, 1099 (1964). 21 L. Baker, “Compressible Rayleigh–Taylor instability,” Phys. Fluids 26, 950 (1983). 22 Y. Elbaz, A. Rikanato, D. Oron, and D. Shvarts, “Compressibility effects on the Rayleigh–Taylor instability,” in the Proceedings of Eighth Workshop on Physics of Compressible Turbulent Mixing, 2001. 23 W. G. Mathews and G. R. Blumenthal, “Rayleigh–Taylor stability of compressible and incompressible radiation-supported surfaces and slabs: Ap-

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp

4670

24

Phys. Fluids, Vol. 16, No. 12, December 2004

plication to QSO clouds,” Astrophys. J. 214, 10 (1977). Lord Rayleigh, “Investigation of the character of the equilibrium of an incompressible heavy fluid of variable density,” Proc. London Math. Soc. 14, 170 (1883).

Ribeyre, Tikhonchuk, and Bouquet 25

C. L. Gardner, J. Glimm, O. McBryan, R. Menikoff, D. H. Sharp, and Q. Zhang, “The dynamics of bubble growth for Rayleigh–Taylor unstable interface,” Phys. Fluids 31, 447 (1988). 26 V. Sobolev, Cours d’Astrophysique Théorique (Mir, Moscow, 1975).

Downloaded 02 Dec 2004 to 147.210.28.65. Redistribution subject to AIP license or copyright, see http://pof.aip.org/pof/copyright.jsp