# Wavelength limits on isobaricity of perturbations in a thermally unstable radiatively cooling medium

## I INTRODUCTION

Thermal instability is thought to be a mechanism which causes a radiatively cooling medium (RCM) to break onto two phases: relatively dense cold filaments and clumps surrounded by rare hot plasma. This phenomenon is common for both cosmic and laboratory optically thin plasmas which lose their energy radiatively (see for more discussion ). A comprehensive linear theory of thermal instability has been developed in the pioneering paper by Field. The instability results from the fact that the radiative energy loss rate grows with density. A local increase of the density brought in a uniform equilibrium plasma, enhances cooling rate and stimulate temperature to fall. This leads in turn to subsequent progressive growth of the density in order to keep pressure constant.

Last decade, nonlinear aspects of thermal instability have been extensively investigated both analytically and numerically, and considerable progress in understanding of the dynamics of phase transitions and structure formation in bistable media was achieved.

The equations governing dynamics of RCM, besides the nonlinearities intrinsically connected with gas dynamics processes, contain also nonlinearities related to radiative energy loss rate, $`\dot{E}(\rho ,T)`$, which is square function of density and, as a rule, strongly nonlinear function of temperature. In general, nonlinear dynamics of such media is complex for analytical study, and several approximations were used to reduce the equations to more simple and useful form. This procedure is based on the presence of different characteristic time scales in RCM: characteristic radiative cooling time $`\tau _R`$, conductive time $`\tau _\chi `$, and acoustic time $`\tau _A`$. When these times differ considerably from each other one can find the fastest process which relaxes to its equilibrium state, and thus eliminate one of the governing equations. In this manner Doroshkevich and Zeldovich, assuming the medium to follow isobaric dynamical regime, have described dynamics of RCM by a nonlinear heat-diffusion equation

$$\frac{T}{t}=\frac{}{q}\chi (T)\frac{T}{q}\epsilon (T),$$

which belongs to a class of diffusion-reaction equations, here the reduced cooling function $`\epsilon (T)`$ represents the “reaction” term. In subsequent papers Meerson and Sasorov, and then Meerson et al.  have developed this approach to reduce the dynamics of RCM to more simple set of equations depending on interrelation between characteristic times. In the short-wavelength limit when acoustic time is the shortest, the dynamical regime is isobaric and the equations of motions are reduced to a diffusion-reaction equation in Lagrangian variables, and contrary, when acoustic time is the longest (the long-wavelength limit), the dynamics of RCM is described by equations of gas dynamics with equation of state determined by “thermal equilibrium”, i.e., by balance between radiative cooling and external heating.

The above approximations suggest that the wavelength of perturbations of hydrodynamical variables satisfies one of the two strong inequalities: $`\tau _A\tau _R`$ for the short-wavelength and $`\tau _A\tau _R`$ for the long-wavelength limits, respectively. Perturbations with wavelengths close to the characteristic length, i.e., $`\tau _A\stackrel{<}{}\tau _R`$ or $`\tau _A\stackrel{>}{}\tau _R`$, show complex behavior and require in general numerical study treated with the full set of equations. We address in this paper the question of what are the wavelength limits in which isobarical approximation can be used for description of thermal instability. We consider perturbations which initially are small (i.e., linear), and show that these limits are rather broad: at $`t>315\tau _A`$, perturbations with wavelengths much larger than wavelengths determined by Meerson et al. , i.e., $`\tau _A\tau _R`$, evolve isobarically. Our conclusion is therefore that at comparatively late stages long-wavelength perturbations can be described by a diffusion-reaction type equation for temperature and density, equivalent to the equation for isobaric short-wavelength perturbations. Initial stages of long-wavelength perturbations can be understood in terms of wave dynamics in a pressure-free medium, and are reasonably described by nonlinear advective-diffusion equation for velocity.

The paper is organized as follows: in Section II we present the basic equations and discuss limiting cases which allow analytical description, in Section III we show numerical results for long-wavelength perturbations. Section IV contains a summary.

## II  GOVERNING EQUATIONS, CHARACTERISTIC TIMES, LIMITS

### A Equations

The basic equations governing dynamics of a radiatively cooling gas were proposed originally by Field and include in the energy equation source and sink terms describing heating by external sources and radiative energy losses. In planar geometry the equations are written as

$$\frac{\rho }{t}+\frac{}{x}(\rho u)=0,$$
$`(1)`$   
$$\rho \left(\frac{u}{t}+u\frac{u}{x}\right)+\frac{P}{x}=\eta \frac{^2u}{x^2},$$
$`(2)`$   
$$\frac{1}{\gamma 1}\left(\frac{P}{t}+u\frac{P}{x}\right)+\frac{\gamma }{\gamma 1}P\frac{u}{x}+\rho (\rho ,T)\frac{}{x}\left(\chi \frac{T}{x}\right)\eta \left(\frac{u}{x}\right)^2=0,$$
$`(3)`$   
$$P\frac{}{\mu }\rho T=0,$$
$`(4)`$   
where $``$ is the net cooling function defined as $`(\rho ,T)=\rho \mathrm{\Lambda }(T)\mathrm{\Gamma }(\rho ,T)`$, with $`\mathrm{\Gamma }`$, as a rule, weekly dependent on $`\rho `$ and $`T`$, other notations have a common meaning.

### B Characteristic times

Three characteristic time scales are relevant for description of perturbations in RCM: acoustic time $`\tau _A=\lambda /c_0`$, conductive time $`\tau _\chi =\rho \lambda ^2/\mu \chi `$, and radiative cooling time $`\tau _R=T/\mathrm{\Lambda }\mu \rho `$, here $`c_0=\sqrt{\gamma P_0/\rho _0}`$ is the adiabatic sound speed. It is obvious that

$$\tau _\chi 3\lambda ^2/\mathrm{}c_0,$$

where $`\mathrm{}`$ is a free-path of gas particles, and thus in framework of hydrodynamical description, i.e., $`\lambda \mathrm{}`$,

$$\tau _\chi =\lambda \tau _A/\mathrm{}\tau _A.$$

For cooling rate per unit mass $`\rho \mathrm{\Lambda }(T)\sigma _ic_0\mathrm{\Delta }En/\mu `$, where $`\sigma _i`$ is the cross-section of inelastic collisions, $`\mathrm{\Delta }E`$, the energy lost by particles in a single inelastic collision, $`n`$ is the gas number density, characteristic cooling time is

$$\tau _R=\frac{T}{ϵ\mathrm{\Delta }E}\frac{\mathrm{}}{c_0},$$

where $`ϵ=\sigma _i/\sigma _e<1`$, the ratio of inelastic-to-elastic cross-sections. Since energy transfered in a single collision, as a rule, is small compared to thermal energy $`\mathrm{\Delta }E<T`$, we get from here $`\tau _R\mathrm{}/c_0`$.

### C Limits

Therefore, three different intervals of wavelengths are relevant, see Fig. 1: short wavelengths $`\lambda _H<\lambda <\lambda _F`$, intermediate wavelengths $`\lambda _F<\lambda <\lambda _A`$, and long wavelengths $`\lambda >\lambda _A`$, here $`\lambda _F=\sqrt{\chi _0T_0/\rho _0\mathrm{\Lambda }_0}`$, the Field length determined as a length where conductive heat transfer and radiative losses equate, i.e., $`\tau _R=\tau _\chi `$, $`\lambda _A=\tau _Rc_0`$ is the acoustic length where $`\tau _A=\tau _R`$ , and $`\lambda _H\mathrm{}`$ is the wavelength at which hydrodynamical description violates. It is seen from Fig. 1 that always the interrelation between characteristic lengths $`\lambda _H<\lambda _F<\lambda _A`$ is valid. From this point of view, the “direct crossover” regime corresponding to the case of large heat conduction, $`\lambda _A\lambda _F`$, is out of the hydrodynamical approximation.

The principal difference between the short and intermediate-wavelength intervals is in more efficient conductive heat transfer compared to radiative losses for short-wavelength perturbations, which leads to erasing of a monochromatic perturbation (a normal mode with fixed wavenumber $`k=2\pi /\lambda `$, e.g., $`\delta Te^{ikx}`$) with $`\lambda <\lambda _F`$ (see). A non-monochromatic localized isobaric perturbation (a wave packet of a set of normal modes, e.g., $`\delta T\delta T_ke^{ikx}𝑑k`$) with characteristic size $`L<\lambda _F`$ shows decrease in amplitude at initial time, which then changes to increase due to thermal instability of long-wavelength modes presenting in the perturbation. Since in these intervals $`\tau _A`$ is the shortest characteristic time, pressure balance sets in first over the perturbation equations reduce to the isobaric form

$$\frac{T}{t}+L(T,P)\frac{}{m}\left(\stackrel{~}{\chi }\frac{T}{m}\right)=0,$$
$`(5)`$   
where

$$L(T,P)=\frac{\gamma 1}{\gamma }\frac{\mu }{}\rho $$

is the reduced cooling function,

$$\stackrel{~}{\chi }=\frac{\gamma 1}{\gamma }\frac{\chi }{PT},$$

the reduced heat conduction coefficient, $`Tdm=dx`$, $`P=`$const.

In the long-wavelength limit, acoustic time $`\tau _A`$ is much larger than characteristic radiative cooling time $`\tau _R`$, and conductive heat transfer is negligible: $`\tau _\chi \lambda ^2/\chi \tau _A\tau _R`$, see Fig. 1. In this limit the equations can be reduced to more simple form describing gaseous medium with a given (anomalous) equation of state. Indeed, when $`\tau _A\tau _R`$ the energy balance sets in on the shortest time $`\tau _R`$, and thus temperature and pressure relax rapidly to the values determined from the equation of steady energy balance

$$(\rho ,T)=0.$$
$`(6)`$   
Note, that Eq. (6) with cooling function $`\mathrm{\Lambda }(T)`$ shown in Fig. 2 and $`\mathrm{\Gamma }=`$const maps $`(\mathrm{\Lambda },T)`$ plane on $`(P,\rho )`$ plane to give the equilibrium equation of state $`P=P(\rho )`$ of a van der Waals type. Equilibria with temperature in the interval where $`d\mathrm{ln}\mathrm{\Lambda }(T)/d\mathrm{ln}T<1`$ are thermally unstable, and correspond to anomalous part of the $`P(\rho )`$-curve with $`dP/d\rho <0`$. These arguments allowed Meerson and Sasorov to reduce the momentum equation of RCM in the long-wavelength limit to the one with a given (“equilibrium”) equation of state $`P=P(\rho )`$

$$\rho \left(\frac{u}{t}+u\frac{u}{x}\right)+\frac{P}{x}=\eta \frac{^2u}{x^2},$$
$`(7)`$   
where $`P=P(\rho )`$ is determined from Eq. (6). For cooling rate with $`d\mathrm{ln}\mathrm{\Lambda }(T)/d\mathrm{ln}T>1`$ below the temperature $`T=T_L`$ and above $`T=T_U`$, and $`d\mathrm{ln}\mathrm{\Lambda }(T)/d\mathrm{ln}T<1`$ in the intermediate interval, the equation of state has a form similar to the van der Waals equation with two rising (stable) branches of the curve $`P(\rho )`$ in the low density $`\rho <\rho _U`$ (high temperature $`T>T_U`$) and high density $`\rho >\rho _L`$ (low temperature $`T<T_L`$) ranges, and a falling (unstable) branch of the $`P(\rho )`$-curve in the intermediate range $`\rho _U<\rho <\rho _L`$ ($`T_L<T<T_U`$). $`N`$-shaped curves in Figures 4 and 6 show the “equilibrium” equation of state $`P(\rho )`$ corresponding to Eq. (6) with the cooling function $`\mathrm{\Lambda }(T)`$ depicted in Fig. 2 (see, also). This equation of state is commonly used for description of pattern formation in a bistable medium in the long-wavelength limit. (For the interstellar medium $`T_L10^2`$ K and $`T_U10^4`$ K, see).

We will mostly concentrate in this paper on the dynamics of perturbations in an RCM which is in the upper marginally stable point: $`TT_U`$, $`d\mathrm{ln}\mathrm{\Lambda }/d\mathrm{ln}T1`$, thus even small (linear) perturbations break the initial state. In this state an RCM can be prepared, e.g., by adiabatically slow decrease of the intensity of external heating source in an RCM kept at fixed pressure $`P_U`$, or adiabatically slow increase of pressure at fixed intensity of heating sources. As at the equilibrium point $`T=T_U`$ the compressibility is infinite, i.e., $`dP/d\rho =0`$, small perturbations around this point evolve to form localized structures as was demonstrated by Meerson and Sasorov. In the next order positive perturbations of density decrease pressure, providing anomalous (negative) compressibility

$$P(\rho )=P_U\frac{a^2}{2}\delta \rho ^2,$$
$`(8)`$   
where $`\delta \rho =\rho \rho _U`$, which amplifies the nonlinearity connected with the advective term in momentum equation (3). It is clear that nonlinear term in Eq. (8) does not saturate the instability, and therefore both the advective term and anomalous negative compressibility can be stabilized only by a strong nonlinearity in the equation of state $`P(\rho )`$, which gives rise to the second stable branch. This saturated state corresponds therefore to a new state of a medium separated onto two phases. The difference between densities of two phases being in pressure equilibrium is considerable, as a rule about two orders of magnitude, and only numerical study can be applied to these strongly nonlinear saturation stages. In order to study dynamics of long-wavelength localised structures in a fluid with negative compressibility Meerson and Sasorov have approximated the falling branch of the equation of state as $`P(\rho )=P_U+\mathrm{const}/\rho `$, and described the so-called unlimited instability.

Initial stages, when $`|\delta \rho |\sqrt{2P_U/a^2}`$, can be described by the advective-diffusion equation

$$\rho \left(\frac{u}{t}+u\frac{u}{x}\right)=\eta \frac{^2u}{x^2}.$$
$`(9)`$   
Numerical simulations show that equation (9) gives qualitatively reasonable description of velocity profiles up to $`t\stackrel{<}{}(215)\tau _A`$ depending on wavelength. Within this time interval pressure varies from $`P_U`$ to minimal value $`P=P_U\mathrm{\Delta }P_m`$, and arised pressure gradients force then velocity to steepen sharply. It forms, in turn, a short-wavelength region, where thermal regime separates the “equilibrium” equation of state $`P=P(\rho )`$ and isobaric distribution sets in rapidly, and then velocity relaxes gradually to profiles adequate to isobaric dynamics. We argue therefore, that at larger times long-wavelength perturbations in RCM evolve isobarically, and can be described by equation (5).

## III NUMERICAL RESULTS

To confirm this conclusion we show here two numerical models which follow the evolution of perturbations from initial linear state to final asymptotic behaviour. In numerical simulations we have used a 6-parametric finite difference scheme of Eqs. (1)–(4) which is free of nonphysical effects of numerical viscosity and absolutely stable (the detailed description of the scheme will be given elsewhere).

To make presentation more clear we use dimensionless variables with length, time and velocity normalized, respectively, to the initial size $`\lambda `$ of a perturbation (which is assumed to be localized), to the dynamical time $`\tau _A=\lambda /c_0`$, and to the adiabatic sound speed $`c_0`$; hydrodynamical variables $`\rho ,P`$ and $`T`$ are normalized to their unperturbed values. In these units cooling time is equal to $`\stackrel{~}{\tau }_R=\tau _R/\tau _A`$. As we concern here dynamics of long-wavelength perturbations, the characteristic size is $`\lambda >\lambda _F`$, (in numerical examples given in the paper $`\lambda _F=0.001\lambda `$). Dimensionless cooling function $`\mathrm{\Lambda }(T)`$ has been taken in the form

$$\mathrm{\Lambda }(T)=\left\{\alpha _1+\alpha _2\mathrm{tan}[\frac{\pi }{2}(T^{}T)]\right\}^1,$$
$`(10)`$   
as shown in Fig. 2. This function represents typical cooling function for bistable media, and generates equation of state of the van der Waals type when “thermal equilibrium”, i.e., $`\rho \mathrm{\Lambda }(T)=\mathrm{\Gamma }(\rho ,T)`$ is reached; we assume through the paper $`\mathrm{\Gamma }(\rho ,T)=\mathrm{\Lambda }(1)`$=const. The unperturbed state was taken near the upper unstable point. We present here results with initial conditions taken in the form of small velocity perturbations corresponding to a localized symmetric inflow to the center of a grid zone

$$u(x,0)=\{\begin{array}{cc}u_0\mathrm{sin}[k(xx_0)],|xx_0|2\pi /k,\hfill & \\ 0,\mathrm{otherwise},\hfill & \end{array}$$
$`(11)`$   
while perturbations of $`\rho `$, $`T`$, and $`P`$ were zero. This choice of initial conditions seems to be most appropriate in study of long-wavelength perturbations. Due to radiative losses initial conditions corresponding to perturbations of density, temperature and pressure are inevitably more affected by the procedure of preparation such perturbations. We therefore, do not consider here alternative dynamics which would start from small perturbations of density, temperature and pressure.

For adopted type of initial conditions velocity perturbations disturb density, temperature and pressure, and radiative cooling drives them to relax in a short time to values determined by thermal equilibrium (6). Arised gradients of pressure redistribute then velocities, such that the resultant velocity field tends to set the flow in isobaric regime. At late stages the flow is approximately isobaric through the perturbed region. The boundary conditions used correspond to a “closed” box, i.e., zero velocity and gradients of density and pressure. In most models the grid zone size was 20 times of the size of a localized perturbation, however test runs with larger size (80 times of the perturbation) show results to be insensitive to the grid size. This is a direct consequence of the fact that the compressibility in the unperturbed state (which is taken to be close to the upper unstable point) is very high, and a localized perturbation propagates slowly: the acoustic velocity $`c_U=\sqrt{(dP/d\rho )_U}0.25`$. Moreover, the initial perturbation brings most fraction of perturbed gas into a state with negative compressibility where the perturbation does not propagate. The time required for perturbations to reach a boundary is therefore at least 80 dynamical times, which is larger than the total computational time in all cases we considered. This means, that the boundary conditions we used weakly affect dynamics of perturbations. (Note, that systems whose unperturbed states lie on the low-density stable branch well below the upper unstable point, can be brought into thermal instability only by perturbations of finite, i.e., not small, amplitude. We do not consider such cases in this paper.)

Fig. 3 shows dynamics of a perturbation with initial amplitude of velocity $`u_0=0.03`$ and with $`\stackrel{~}{\tau }_R=0.33`$ (i.e., $`\lambda =3\lambda _A`$), thus it corresponds to the long-wavelength limit. On the upper panel velocity profiles are plotted. At initial times, $`t10`$, as gas in central regions is involved in motion the velocity amplitude decreases. In gases with normal equation of state, i.e., $`dP/d\rho >0`$, perturbation decays onto two acoustic waves propagating in opposite directions, and a standing and damping entropy mode \[in case of perturbations given by eq. (11), decayed waves represent two pieces of sinusoidal waves with velocity amplitude of $`u_0/2`$\]. Radiative losses with $`\stackrel{~}{\tau }_R1`$ modify the inflow, transfering the internal (compressed) gas onto the anomalous part of equation of state ($`dP/d\rho <0`$), and thus make the perturbation locally growing, i.e., non-decaying. At times $`10<t<15`$ the velocity profile steepens and motion concentrates in a rather narrow region with increasing density, as shown in the bottom panel.

Relaxation of pressure to the equilibrium (anomalous) equation of state $`P(\rho )`$ occurs at times $`t\stackrel{~}{\tau }_R`$ (these stages are not shown). Increasing density in central regions amplifies pressure gradient through the perturbation, and thus forces the velocity amplitude to increase. It results in more sharp steepening of velocity profiles compared to the pressure-free case as seen in Fig. 3a for time interval $`10<t<15`$: the profile for $`t=15`$ clearly shows two inflection points outside the central region very close to peak values, which is not what can be expected for pressure-free motions, or motions with normal equation of state. In turn, the inflow with increasing velocity amplitude forms very narrow region (a “droplet” with sizes $`\mathrm{\Delta }x\lambda `$, as seen in Fig. 3b), through which pressure equilibrium sets in rapidly. As a result, at these stages structure of the perturbation can be understood as distinguished onto two regions: first, the “long-wavelength” part outside the “droplet” where the equilibrium equation of state dominates, and second, the “short-wavelength” part (the “droplet” itself) where motion tends to set in the pressure equilibrium. The two regimes manifest themselves clearly in Fig. 4, where loci $`P\rho `$ connecting pressure and density through the computational zone for different times, i.e., $`P(x,t)`$ and $`\rho (x,t)`$ profiles with eliminated $`x`$, are plotted: the short piece of line at the upper equilibrium point shows $`P\rho `$ locus for $`t=10`$. The locus corresponding to $`t=15`$ is readily seen to consist of two parts: a small piece of a $`P\rho `$ curve starting from the unperturbed state with decreasing pressure and density almost parallel to the growing branch of the equation of state (corresponding to external, long-wavelength regions of the perturbation where density and pressure decrease), which then turns around and goes with increasing density and almost constant pressure (this part corresponds to internal, short-wavelength regions). Thus, at these stages pressure separates the equilibrium curve $`P(\rho )`$, and subsequent evolution is characterized by settling the overall flow down to isobaric regime. At $`t\stackrel{>}{}20`$ most fraction of density perturbation concentrates in the isobaric “droplet” with characteristic size of $`|xx_0|=0.2`$ and in the conductive interface, which are in pressure equilibrium with gas outside (variation in pressure through the perturbed region is less than 0.001, and main contribution comes from surrounding regions where density and pressure lie on the equilibrium equation of state). (Strongly speaking, in models with “closed” box boundary conditions we used, pressure depends on time. However, as mentioned above, in all considered cases asymptotic behaviour of perturbations establishes on time interval shorter than the time required for a perturbation to reach boundaries. One can think, therefore, that isobaric regimes asymptotically established in numerical models with “closed” box boundaries give qualitatively correct description of the situation in systems with free boundary and fixed external pressure.)

After the short-wavelength droplet is formed, velocity decreases to $`u<0.01`$ with profiles corresponding to a weak inflow, so more gas continues to increase mass of the droplet. Subsequent evolution follows the isobaric regime with $`P\rho `$ loci practically coincident with that reached at $`t=20`$, as shown in Fig. 4. Note that in models with fixed equilibrium equation of state $`P=P(\rho )`$, see, the established inflow velocity is, in general, higher than in our calculations with energy equation (3) solved explicitly. This is due to the fact that in models with the equilibrium equation of state, $`P=P(\rho )`$, pressure reaches minimal value in the interface layer, and thus at equal conditions (e.g., viscosity and conductivity) arised pressure gradients generate higher acceleration.

In Fig. 5 and 6 we show results for a perturbation with larger spatial size: $`\lambda =30\lambda _A`$, i.e., $`\stackrel{~}{\tau }_R=0.033`$. Initial conditions are same as in the previous case. The perturbation is seen to reveal similar behavior as that of smaller size, however, the dynamics is more violent. The difference is mainly in time scales and the amplitude of velocity: perturbations of longer wavelength develop higher velocities and pass same stages quicker. The enhanced steepening of the velocity profile due to negative compressibility occurs at $`t\stackrel{<}{}2`$. Here the difference between enhanced and normal steepening is more obvious than in previous case: rather flat increase in velocity in outer regions, and very sharp growth near the symmetry plane. At $`t3`$ the inflow with amplitude $`u0.08`$ forms a narrow short-wavelength droplet with almost constant pressure inside \[as seen from $`P\rho `$ loci on $`(P,\rho )`$ plane, Fig. 6\]. However, as the velocity amplitude developed is large, pressure in the droplet is higher than the unperturbed value by $`\mathrm{\Delta }P\rho u^20.01`$, and in very central parts of the droplet overpressure reaches the value $`\mathrm{\Delta }P0.03`$ due to higher densities \[note, that $`P\rho `$ locus in the central overpressured region goes along the equilibrium equation of state $`P(\rho )`$ corresponding to the second stable brach\]. On next stages, $`3<t\stackrel{<}{}4`$, the overpressure leads to expansion of the short-wavelength region (with almost isobaric pressure distribution inside) to relax to the equilibrium with external regions. This relaxation occurs more violent than in the previous case, due to larger variations in pressure, $`|\mathrm{\Delta }P|0.010.03`$ around the unperturbed state. As a result, subsequent growth of gas mass in the short-wavelength (droplet) region proceeds in a damped oscillatory regime, as seen on $`(P,\rho )`$ plane in Fig. 6.

The dependence of dynamics of perturbations on wavelength (or equivalently on $`\tau _R`$) can be understood qualitatively from simple estimations of steepening of velocity profiles. The characteristic formation time of a short-wavelength region (a droplet) due to steepening can be estimated from $`at_d^22\lambda `$, where $`a`$ is the acceleration caused by the negative compressibility,

$$a\frac{|P_m|}{\rho _0}\frac{\mathrm{\Delta }P_m}{\rho _0\lambda },$$
$`(12)`$   
$`\mathrm{\Delta }P_m`$ is the maximal deviation of pressure from the unperturbed value due to density increase. This gives in dimensionless units

$$t_d\sqrt{\frac{2}{\mathrm{\Delta }P_m}}.$$
$`(13)`$   
The deviation of pressure from the unperturbed value increases along the equilibrium $`P(\rho )`$-curve until radiative cooling time $`\tau _R`$ is less than the acoustic time for the shortest wavelength $`\lambda _s`$ generated due to steepening. Thus, $`\mathrm{\Delta }P_m`$ can be determined from the condition $`\stackrel{~}{\tau }_R=\lambda _s`$. To estimate $`\lambda _s`$ we assume that perturbation of density caused by initial inflow (11) in dimensionless units is $`\delta \rho _0(x)u_0(x)`$. Then, due to thermal instability it grows as $`\delta \rho (x,t)u_0\mathrm{exp}(t/\tau _R)`$. On the other hand, one can expect that $`\lambda _s\lambda \mathrm{exp}(t/\tau _R)`$, thus giving $`\delta \rho u_0/\lambda _s`$. Substituting here $`\mathrm{\Delta }P_mc_A^2\delta \rho `$, where $`c_A^2=|dP/d\rho |`$ estimated on the anomalous branch of equation of state \[$`c_A^20.1`$ at $`\rho 1.5`$\]

$$\mathrm{\Delta }P_m=\frac{u_0c_A^2}{\lambda _s},$$
$`(14)`$   
We obtain therefore for $`t_d`$:

$$t_d\sqrt{\frac{2\stackrel{~}{\tau }_R}{u_0c_A^2}},$$
$`(15)`$   
which provides with reasonable accuracy an estimate for the time when the numerical solutions turn from the regime with the equilibrium equation of state, $`P(\rho )`$, to the isobaric regime.

It is clear from these estimates and from numerical results presented, that at small values of the ratio $`u_0c_A^2/\stackrel{~}{\tau }_R`$ (less than 0.1) long-wavelength pertubations develop quiescently as isobaric. In dimensional units the condition for long-wavelength perturbations to evolve in a quiescent isobaric regime is written as

$$\frac{\tau _A}{\tau _R}\stackrel{<}{}0.1\frac{c_0}{u_0}\left(\frac{c_0}{c_A}\right)^2.$$
$`(16)`$   
When $`u_0c_A^2/\stackrel{~}{\tau }_R>0.1`$, high pressure gradients developed at initial “equilibrium” stages \[i.e., along $`P(\rho )`$ equation of state\] give rise to very fast steepening of velocity profiles, and as a result to formation of the overdense and overpressured central short-wavelength droplet, which then relaxes to the equilibrium in an oscillatory regime with variations in pressure of $`|\mathrm{\Delta }P|\stackrel{<}{}0.03`$ around the unperturbed value. Regimes with the equilibrium equation of state, $`P(\rho )`$, can exist for long-wavelength perturbations only at initial stages, $`t\stackrel{<}{}t_d`$. Subsequent evolution implies generation of a short-wavelength core, which in turn results in separation of pressure of the equilibrium curve $`P(\rho )`$. Formation of such a short-wavelength region is an inevitable consequence of the “unlimited” instability when narrow localized structure develops in a medium on the anomalous branch of the equilibrium equation of state<sup>4</sup>. In Fig. 7 we show the dependence of pressure in the central droplet on time for several values of $`\stackrel{~}{\tau }_R`$ and fixed $`c_A^2=0.1`$ at $`\rho =1.5`$. The amplitude of pressure oscillations is clearly seen to increase for smaller $`\stackrel{~}{\tau }_R`$: approximately, at initial times, when the central droplet follows the equilibrium equation of state, $`\mathrm{\Delta }P_0\mathrm{\Delta }P_mu_oc_A^2/\stackrel{~}{\tau }_R`$. The period of oscillations is $`𝒯\lambda _s\sqrt{\rho _d/\mathrm{\Delta }P}`$, where $`\rho _d=5`$ is the density on the second stable branch; at initial times $`𝒯2.24\stackrel{~}{\tau }_R^{3/2}/\sqrt{u_0c_A^2}`$, however then it grows due to damping of $`\mathrm{\Delta }P`$.

## IV SUMMARY

We have considered nonlinear evolution of initially small (linear) perturbations in a radiatively cooling medium in the long-wavelength limit, i.e., when acoustic time is larger than radiative cooling time. We have not been assuming pressure and density to satisfy the equilibrium equation of state, but studied the problem numerically with the radiative energy loss term in energy equation treated explicitly. For initial conditions in the form of small localized velocity perturbation around a uniform equilibrium close to the upper unstable point, we have shown that early stages of evolution tend to proceed along the equlibrium equation of state, $`P(\rho )`$, with velocity profiles to steepen as in pressure-free case. At intermediate stages density increases enough and anomalous (negative) compressibility becomes important and results in more sharp steepening of velocity profiles, which in turn forms a narrow (short-wavelength) region where pressure separates the value determined by the equilibrium equation of state and sets in approximate isobaric distribution rapidly. One may say that thermal instability works as a regulator which prevents large pressure gradient to arise. The maximal amplitude of pressure developed in the perturbation when it separates the equilibrium curve $`P(\rho )`$, is proportional to the amplitude of initial velocity perturbation and inversely proportional to the radiative cooling time \[as described by Eq. (14)\]

$$\mathrm{\Delta }P_mP_0\frac{u_0}{c_0}\frac{\tau _A}{\tau _R}\left(\frac{c_A}{c_0}\right)^2.$$
$`(17)`$   
At late stages, gas in the short-wavelength region reaches the high-density (stable) equilibrium branch of the equation of state, and subsequent evolution corresponds to the isobaric inflow of low-density gas onto the high-density core (the droplet). We have shown thus, that nonlinear evolution of long-wavelength perturbations after the time interval

$$t_d\frac{c_0}{c_A}\sqrt{2\tau _A\tau _R\frac{c_0}{u_0}},$$
$`(18)`$   
reaches quiescent isobaric regime, if the wavelength is restricted from above as

$$\lambda <0.1\lambda _A\frac{c_0}{u_0}\left(\frac{c_0}{c_A}\right)^2.$$
$`(19)`$   
In these limits, long-wavelength perturbations in thermally unstable media can be described as isobaric by the reduced equation (5). Perturbations with wavelengths longer than the limit (19), after $`t=t_d`$ when pressure separates the equilibrium equation of state, develop localized oscillating structures with initial pressure amplitude determined by Eq. (14) and damped then due to viscosity.

###### Acknowledgements.

YS acknowledges financial support from Russian Foundation of Basics Research (grant 94-02-05016-a), and from Italian Consiglio Nazionale delle Ricerche within NATO Guest Fellowship programme 1996 (Ann. No. 219.29). We acknowledge anonimous referees for critical remarks.