A new radiative cooling curve based on an up-to-date plasma emission code

K. M. Schure; D. Kosenko; J. S. Kaastra; R. Keppens; J. Vink

doi:10.1051/0004-6361/200912495

Free Access

Issue		A&A Volume 508, Number 2, December III 2009


Page(s)		751 - 757
Section		Interstellar and circumstellar matter
DOI		https://doi.org/10.1051/0004-6361/200912495
Published online		27 October 2009

A&A 508, 751-757 (2009)

A new radiative cooling curve based on an up-to-date plasma emission code

K. M. Schure¹ - D. Kosenko¹ - J. S. Kaastra^1,2 - R. Keppens^1,3,4 - J. Vink¹

1 - Astronomical Institute, University of Utrecht, Postbus 80000, 3508 TA Utrecht, The Netherlands
2 - SRON, Utrecht, The Netherlands
3 - Centre for Plasma Astrophysics, K.U. Leuven, Belgium
4 - FOM-Institute for Plasma Physics ``Rijnhuizen'', Nieuwegein, The Netherlands

Received 14 May 2009 / Accepted 24 September 2009

Abstract
This work presents a new plasma cooling curve that is calculated using the SPEX package. We compare our cooling rates to those in previous works, and implement the new cooling function in the grid-adaptive framework ``AMRVAC''. Contributions to the cooling rate by the individual elements are given, to allow for the creation of cooling curves tailored to specific abundance requirements. In some situations, it is important to be able to include radiative losses in the hydrodynamics. The enhanced compression ratio can trigger instabilities (such as the Vishniac thin-shell instability) that would otherwise be absent. For gas with temperatures below 10⁴ K, the cooling time becomes very long and does not affect the gas on the timescales that are generally of interest for hydrodynamical simulations of circumstellar plasmas. However, above this temperature, a significant fraction of the elements is ionised, and the cooling rate increases by a factor 1000 relative to lower temperature plasmas.

Key words: hydrodynamics - ISM: evolution - radiation mechanisms: thermal

1 Introduction

Optically thin plasmas can cool very efficiently, especially on timescales relevant to circumstellar evolution. A correct description of cooling is therefore important for simulations of e.g., stellar winds. The dominant cooling mode for plasmas within the temperature range of 10 ⁴-10⁷ K is by metal line transitions. In many cases, the plasma can be assumed to be in collisional ionisation equilibrium (CIE). A predefined cooling curve can then be used to account for radiative losses, where cooling rates depend only on the temperature of the plasma. In the calculation of the cooling curve, we can adjust the abundances and try to incorporate as many line transitions as are reliably known, to obtain a straightforward and efficient way of treating radiative cooling. Well-documented cooling curves used in HD and MHD computations are those for example developed by Dalgarno & McCray (1972), Raymond & Smith (1977), MacDonald & Bailey (1981), Sutherland & Dopita (1993), and Smith et al. (2008). Other approaches, such as that of Mellema & Lundqvist (2002), include time- and ion dependent cooling, with CSM abundances that depend on the evolutionary stage of the star, but treat the hydrodynamics only in one dimension.

We use the SPEX package (Kaastra et al. 1996) to calculate a new cooling curve, using up-to-date transition lines to calculate radiative losses for a plasma. This code is used as a spectral analysis code tailored to EUV and X-ray observations, in which energy band it is one of the most complete packages currently available. Because it includes a very complete prescription of line emission, cooling rates predicted by SPEX are higher than those of cooling curves available until now.

We implement the new cooling curve, as deduced from the SPEX CIE code, in the AMRVAC framework (Keppens et al. 2003), of which we use the hydrodynamics to demonstrate the influence of radiative losses on the gas dynamics. As an example, we use this code with the new cooling prescription to simulate the circumstellar medium of a massive star. Plasmas typical of this circumstellar medium are at fairly high temperatures and are mostly optically thin. Therefore, they can cool very efficiently.

The cooling at a certain point in the simulated plasma can be calculated from the local temperature and density in combination with the prescribed cooling curve. Different cooling curves can be chosen, which depend on metallicity and/or ionisation fraction (for low temperatures). We provide cooling rates for the individual elements, so that the cooling curve can be tailored to fit the required abundance ratios.

In Sect. 2, the SPEX code, the abundances, and the new cooling curve will be presented and compared to existing cooling curves. We show that the SPEX curve causes more cooling for temperatures from 10⁴ to 10⁶ K, a temperature range that is very important in the evolution of the circumstellar as well as the interstellar medium.

In Sect. 3, we discuss the range of validity of the CIE assumption and discuss the implications of photoionisation and shock heating on the cooling rates of the plasma.

Cooling at lower temperatures cannot be accurately calculated with SPEX, because of incomplete data for $\lambda>2000$ Å. In Sect. 4, we briefly summarise our treatment of lower temperature plasmas, which follows Dalgarno & McCray (1972). The cooling rates of Smith et al. (2008), provide a detailed representation of especially the low-temperature (T <10⁴ K) cooling of plasmas. We select our low-temperature ionisation fraction to ensure the optimal correspondence with their cooling rates.

Table 1: Solar abundances from Anders & Grevesse (1989).

In Sect. 5, we show the effect that cooling has on the dynamics of the plasma, by using the SPEX cooling curve as implemented in the AMRVAC code. Because the cooling adds to the resolution requirements of the hydro-simulations, the grid-adaptive method used by AMRVAC maximises the performance while reaching sufficient resolution to enhance the higher density contrast and allow for instabilities to develop in the cooling shell. In Sect. 6, we summarise the important points and outline the major differences from earlier cooling curves.

2 SPEX cooling curve

To generate the cooling curve we use the SPEX package, version 2.00.11 (Kaastra & Mewe 2000). This code is tailored to fit spectra in the EUV and X-ray regime and can be used to calculate the luminosity of a given emission measure for plasmas at different temperatures and different choices of electron and ion temperatures. For our purpose, we use SPEX to calculate the emissivity of a plasma in collisional ionisation equilibrium (CIE). The code is based on the mekal line-emission model (Mewe et al. 1995), but many updates have since been made. The 15 elements presently included are H, He, C, N, O, Ne, Na, Mg, Al, Si, S, Ar, Ca, Fe, and Ni.

We calculate spectra for a grid of temperatures ranging from 10^3.8 to 10^8.16 K, with steps of $\log T=0.04$ K. Solar abundances are assumed, using the abundance table presented in Anders & Grevesse (1989), which is summarised here in Table 1 for elements that are taken into account. The emissivity $\Lambda_N$ (erg s^-1 cm³, normalised to $n_{\rm e} n_{\rm H}V=1$ cm^-3) is calculated by integrating the spectrum over an energy range from 10^-1 eV to 1 MeV. The SPEX package is currently one of the most complete spectral packages, and it is especially complete for temperatures higher than $\sim$ 10^4.86 K (see Sect. 4).

The cooling curve is tabulated in Table 2, and in Fig. 1 it is plotted against the CIE cooling curve that was calculated by Sutherland & Dopita (1993). Compared to their curve, and other CIE cooling curves currently available, we find that the cooling efficiency calculated by SPEX is higher. We attribute this difference to SPEX taking into account more line transitions in both X-ray and EUV than previous studies (5466 lines in SPEX compared to 1683 in Sutherland & Dopita 1993). Additionally, though Sutherland & Dopita (1993) use the same values of Solar abundance ratios, they use a slightly different ionisation balance, which they calculate with their MAPPINGS II code. SPEX uses the ionisation balance from Arnaud & Rothenflug (1985) for most elements, apart from the one from Arnaud & Raymond (1992) for iron (Mewe et al. 1995). A more complete list of Fe L complex lines is also included, based on the HULLAC atomic physics package (Klapisch et al. 1977), causing higher cooling efficiency around 10⁷ K. The differences in the treatment of Fe are shown in Mewe et al. (1995).

Additionally, we plot the cooling curve by MacDonald & Bailey (1981) (hereafter also referred to as MB81) in Fig. 1. We used this for comparison in our hydro-simulations in Sect. 5. This cooling curve is a composite of cooling rates by Raymond et al. (1976) for $10^6<T\le 10^8$ K and by Shapiro & Moore (1976) for 10⁴<T<10⁶ K. Below 10⁴ K, extrapolation of the cooling rates is used, towards zero at 100 K. Their treatment of cooling between 10⁴-10⁶ K is a non-equilibrium approach, applicable to a low density plasma that is overionised, as a result of rapid cooling ( $t_{\rm cool} < t_{\rm recombination}$ ). Time-dependent cooling, starting from ionisation equilibrium at 10⁶ K, gives rise to lower cooling rates, as discussed in Sect. 3. MB81 base their cooling rates on the inclusion of 12 elements (H, He, C, N, O, Ne, Mg, Si, S, Ca, Fe, Ni), and use slightly different abundances, which are not however responsible for the significant differences from our cooling rates.

$\begin{figure} \par\resizebox{8.8cm}{!}{\includegraphics{12495fg1}} \end{figure}$	Figure 1: Cooling curves compared: the higher cooling rates calculated with SPEX are mainly due to a more complete coverage of the line transitions, including Fe L and EUV lines.
Open with DEXTER

To apply the cooling curve in the AMRVAC code, we renormalise the cooling rates calculated with SPEX by multiplying by the electron to proton ratio ( $n_{\rm e}/n_{\rm H}$ ). We thus obtain a value for $\Lambda_{\rm hd}$ that we can use in our hydrodynamic application:

$\displaystyle \Lambda_{\rm hd}=\frac{n_{\rm e}}{n_{\rm H}}\Lambda_N.$

(1)

The ratio $n_{\rm e}/n_{\rm H}$ is a temperature-dependent parameter that is also calculated by SPEX and is given in the last column of Table 2. The cooling rate of the plasma is given by

$\displaystyle L = \int n_{\rm H} n_{\rm H} \Lambda_{\rm hd}(T) {\rm d}V \quad({\rm erg~s}^{-1}).$

(2)

To allow the cooling curve to be adapted for different abundances, we list the contributions of the different elements to the cooling rates, as calculated with the SPEX code. This is tabulated in Table 4 (online only). The cooling rates are given again for solar abundances (Anders & Grevesse 1989). The contributions of the individual elements are plotted in Fig. 2 to illustrate the relative importance of the different elements to certain temperature ranges. The cooling rates are given for normalised values of $n_{\rm H} n_{\rm e}=1$ , and should be multiplied by $n_{\rm e}/n_{\rm H}$ when used in a hydrocode.

For different metallicities, we can easily construct new cooling curves from Table 4 (online only), in which the cooling contributions and ionisation fractions of the elements are given separately. Cooling rates could then be calculated by summing the contributions of the different elements, multiplied by the appropriate abundance ratio and ionisation fraction

$\displaystyle \Lambda_N(T)=\sum_i \frac{n_i}{n_i ({\rm solar})} \Lambda_N(X_i,T),$

(3)

where n_i is the abundance ratio (relative to H) of element X_i. In this way, new cooling curves can be created ad hoc when needed for simulations of plasmas with different metallicities.

$\begin{figure} \par\includegraphics[width=8.8cm]{12495fg2a} \includegraphics[width=8.8cm]{12495fg2b} \end{figure}$	Figure 2: Contributions of different elements to the cooling curve are given. Each of the plots shows a different set of elements. Important peaks are labelled with the name of the element. The total cooling curve (black solid line) is an addition of the individual elemental contributions.
Open with DEXTER

Table 2: Cooling curve for solar metallicity.

3 Assumptions and validity of CIE cooling rates

Collisional ionisation equilibrium is a valid assumption if an optically thin plasma is dominated by collisional processes and the cooling timescale is longer than either the ionisation or recombination timescales. Deviations from CIE can occur if there is additional photoionisation or there is non-equilibrium ionisation (NEI), i.e., the recombination or ionisation timescales are longer than the cooling/heating timescale. These circumstances will produce an overionised plasma if photoionisation is present or if the cooling timescale is shorter than the recombination timescale. In young objects, such as supernova remnants, the heating timescale is shorter than the ionisation timescale, giving rise to an underionised plasma.

NEI requires a time-dependent treatment of the cooling and is not easily implemented in standard hydrodynamical codes, but some specialised one-dimensional hydrodynamical codes for SNRs do exist (e.g., Kosenko et al. 2004; Badenes et al. 2003). On the other hand, photoionisation requires additional ionisation dependence on the photon-field, resulting in a higher ionisation than CIE for the same plasma temperature. The efficiency and reach of a photoionising source depends on both its spectrum and the ionisation parameter $\xi = \frac{L}{nR^2}$ , where R is the distance to the source. Depending on whether the elements are more or less effective coolants in their additionally ionised state, the cooling rates may be higher or lower, respectively.

Studies of the cooling rates of overionised plasmas, relevant to photoionised plasmas or the late radiative stages of supernova remnants, have been performed by various authors (e.g., Kafatos 1973; Shapiro & Moore 1976; Sutherland & Dopita 1993). The cooling rates in all of these cases are lower than that of a CIE cooling curve. For comparison, we plot the cooling curve by MB81 in Fig. 1. This curve uses the time-dependent cooling curve by Shapiro & Moore (1976) for the temperature range $4 < \log T < 6$ . The gas is initially assumed to be in ionisation equilibrium at T=10⁶ K, and cools time-dependently from there. Since the cooling timescale is shorter than the recombination timescale, the gas becomes overionised, leading to suppressed cooling rates in this regime. Similarly, in intergalactic plasmas where photoionisation is important, cooling rates tend to be suppressed (Wiersma et al. 2009). Alternatively, cooling rates may be much higher when a neutral gas is shock-heated, due to the initial ionisation of H and He (Cox & Raymond 1985). Unless full radiative transfer is taken into account, we therefore find CIE to be a reasonable assumption to use in circumstellar environments and other low density plasmas. Only when photoionisation is significant, which depends on the temperature of the star and the distance to the source, or when strong shocks with velocities higher than about 150 km s^-1 are present, deviation from CIE can become important. At higher shock speeds, the gas is heated to temperatures higher than $\log T \approx 5.7$ , at which cooling by H and He lines of C, N, and O, becomes suppressed. Even in our relatively strong shock between two subsequent wind-phases as presented in Sect. 5 however, the typical speed of the forward shock is less than 150 km s^-1 and the temperature is typically around 10⁵ K.

Since recombination, collisional ionisation, and cooling have similar density dependences, the cooling rates that we calculate are independent of density, as long as the gas can be regarded as optically thin.

4 Low temperature cooling

The SPEX package is tailored to the X-ray regime, and includes continuum but no line emission for wavelengths longer than $\lambda>2000$ Å. Therefore, it may underestimate cooling for temperatures below $T=hc/\lambda k \sim10^{4.86}$ K. At temperatures below 10⁴ K, cooling is much less efficient, and in the case of the typical CSM that we are interested in, most cooling takes place above this threshold.

For low temperatures (T <10⁴ K), we use the cooling rates from a number of reactions as given by Dalgarno & McCray (1972). In this regime, fractional ionisation and molecular line transitions can play an important role. We can set the ionisation fraction $f_{i}=n_{\rm e}/n_{\rm H}$ , which is a free parameter, to accomodate for different circumstances in the circumstellar medium (CSM). This fraction strongly depends on specific conditions concerning shocks, radiation, and the history of the plasma. The cooling rates are included in Table 3 (online only), and below we repeat the equations that we use to calculate them. The values for an ionisation fraction of 10^-4 correspond most closely to the values that Smith et al. (2008) find for a Solar abundance plasma.

The elements that are taken into account for low-temperature cooling are O, C, N, Si, Fe, Ne, and S. The cooling is provided partly by excitation of singly charged ions with the thermal electrons. The electron-collision transitions that are takeninto account are given by the following equations (with the temperature T in Kelvin) (Dalgarno & McCray 1972):

$\displaystyle \Lambda_{\rm e}({\rm C^+})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}} f_i \right)T^{-1/2}\left[ 7.9\times10^{-20}{\rm e}^{-92/T}\right.$
		$\displaystyle + \left. 3.0 \times 10^{-17} {\rm e}^{-61900/T}\right]$	(4)
$\displaystyle \Lambda_{\rm e}({\rm Si^+})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}} f_i \right)T^{-1/2}\left[ 1.9\times10^{-18}{\rm e}^{-413/T}\right.$
		$\displaystyle + \left. 3.0\times 10^{-17}{\rm e}^{63600/T}\right]$	(5)
$\displaystyle \Lambda_{\rm e}({\rm Fe^+})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}} f_i \right)T^{-1/2}\left[4.8\times 10^{-18}{\rm e}^{2694/T}\right.$
		$\displaystyle +1.1\times10^{-18}\left({\rm e}^{-554/T}+1.3{\rm e}^{-961/T}\right)$
		$\displaystyle + \left. 7.8\times 10^{-18}{\rm e}^{3496/T}\right]$	(6)
$\displaystyle \Lambda_{\rm e}({\rm O})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}} f_i \right)1.74\times10^{-24}T^{1/2}$
		$\displaystyle \times\left[(1-7.6T^{-1/2}){\rm e}^{-228/T}\right.$
		$\displaystyle \left.+0.38(1-7.7T^{-1/2}){\rm e}^{-326/T}\right]$
		$\displaystyle +9.4\times10^{-23}T^{1/2}{\rm e}^{-22700/T}$	(7)
$\displaystyle \Lambda_{\rm e}({\rm O^+})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}} f_i \right)1.5\times10^{-17}T^{-1/2}{\rm e}^{-38600/T}$	(8)
$\displaystyle \Lambda_{\rm e}({\rm N})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}} f_i \right)8.2\times10^{-22}T^{1/2}$
		$\displaystyle \times(1{-}2.7\times10^{-9}T^2){\rm e}^{-27700/T}$	(9)
$\displaystyle \Lambda_{\rm e}({\rm S^+})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}} f_i \right)8.4\times10^{-18}T^{-1/2}{\rm e}^{-21400/T}.$	(10)

For low fractional ionisation, collisions with neutral hydrogen can add substantially to the cooling. The following transitions are taken into account:

$\displaystyle \Lambda_{\rm H}({\rm Si^+})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}}\right)7.4\times10^{-23}T^{-1/2}{\rm e}^{-413/T}$	(11)
$\displaystyle \Lambda_{\rm H}({\rm Fe^+})$	=	$\displaystyle \left(\frac{n_i}{n_{\rm H}}\right)1.1\times10^{-22}$
		$\displaystyle \times\left[{\rm e}^{-554/T}+1.4{\rm e}^{-961/T}\right].$	(12)

The total cooling that we take into account for T < 10⁴ K is a summation of contributions from Eqs. (4) to (12), plus contributions from collisions between electrons and hydrogen atoms as tabulated in Table 2 in Dalgarno & McCray (1972).

5 Application: test problem

We use the AMRVAC code with the new cooling prescription to simulate the circumstellar medium (CSM) of a massive star. The adaptive mesh refinement (AMR) ensures a computationally efficient approach to problems that have localised needs of high resolution, but where lower resolution suffices in substantial parts of the grid. The standard Euler equations are solved in a conservative way, where we apply an HLLC solver (Toro 1999) to capture more accurately the contact discontinuity.

We implement the stages of the evolution of a massive star in which a red supergiant (RSG) wind is followed by the much faster Wolf-Rayet (WR) wind. A RSG wind is typically slow, cool, and dense because of the low escape velocity from the surface of the giant star. The WR wind, however, is fast and hot, originating from the compact stellar core. The WR wind therefore sweeps up the RSG wind into a shell (Dwarkadas 2007; García-Segura et al. 1996; Garcia-Segura & Mac Low 1995), hereafter called the WR shell, whose thickness and stability is strongly affected by the inclusion of radiative cooling.

We set up the simulation in 2.5D on a spherical grid in a meridional plane over a full angle $\pi$ , using reflective boundary conditions and symmetry around the $\phi$ -axis. The radial extent is about 3 pc and we use a base resolution of $60\times48$ gridcells with 4-6 refinement levels. For each refinement level, the resolution is doubled. The values that we use for the wind parameters are the following: $\dot M_{\rm RSG}=1.54\times 10^{-5}~ M_\odot$ yr^-1, $v_{\rm RSG}=4.7$ km s^-1, $T_{\rm RSG}=10^3$ K, $\dot M_{\rm WR}=9.7\times 10^{-6}~ M_\odot$ yr^-1, $v_{\rm WR}=1.7\times10^3$ km s^-1, $T_{\rm RSG}=10^4$ K.

For the RSG phase, 4 refinement levels are used, and for the WR phase, we use 4-7 refinement levels and check for differences in the development of the instabilities and the resolution dependence of the radiative losses. The maximum resolution corresponds to $1.56\times 10^{15}$ cm $\times\ 0.0156^\circ$ .

In general, all plasmas are prone to radiative losses. Depending on the circumstances however, these can induce significant changes in the hydrodynamics, or have only minor quantitative effects. In the situation that we simulate, without cooling the maximum compression reached in the shell is a factor of four; the value expected for a strong shock in a mono-atomic gas with an adiabatic index of $\gamma=5/3$ . With cooling, the shell can become much more compressed and, as a result, much thinner. When the compression ratio exceeds a factor of about 21, the shell can become subject to the Vishniac thin shell instability (Vishniac 1983; Mac Low & Norman 1993; Ryu & Vishniac 1987).

The Vishniac thin shell instability, which actually is an overstability, occurs when a thin shell is perturbed while being driven by a high pressure region from the inside and encountering ram pressure on the outside. The high pressure from the WR wind acts in a perpendicular direction to the shell, whereas the ram pressure acts in a parallel direction to the velocity. Whereas the velocity is radial, the interior pressure in the case of a perturbed shell can create a flow along the shell surface, driving matter from advanced regions to regions lagging behind. The more tenuous regions then become more susceptible to deceleration by ram pressure, reversing the cycle.

Above a temperature of 10⁵ K, the plasma also becomes susceptible to the thermal instability. This occurs for shock velocities above 150 km s^-1, when $({\rm d}\ln \Lambda)/({\rm d}\ln T) < 2$ (e.g., Sutherland & Dopita 1993; Field 1965; Bertschinger 1986). For typical circumstellar wind parameters, the plasma temperature in the cooling shell remains below this regime.

The cooling of the plasma is calculated at each timestep and the thermal energy is correspondingly updated. The temperature of the plasma is first determined. The corresponding cooling rate is then found by interpolating between the appropriate values in the cooling curve. The energy loss rate per unit volume is calculated as follows:

$\displaystyle \frac{{\rm d}e}{{\rm d}t}=-\Lambda_{\rm hd}n_{\rm H}^2.$

(13)

The energy losses are subtracted from the thermal energy only, preventing numerical precision problems from producing negative energy.

Although the cooling rates are higher when we use the SPEX-calculated cooling curve, in both cases the thickness of the WR shell is not so much limited by the cooling curve used, as by the adopted numerical resolution. It is important however, that the shell becomes affected by instabilities in the presence of cooling, whereas it remains stable in the absence of cooling. We compare the differences for runs with different resolution and different cooling prescription. For comparison, we use the cooling curve as prescribed in MacDonald & Bailey (1981), since this is the cooling curve that is used in the ZEUS3D v3.4 code, as used by some of the authors before, for simulations of a WR wind in Schure et al. (2008).

$\begin{figure} \par\resizebox{8.8cm}{!}{ \includegraphics{12495fg3}} \end{figure}$

Figure 3:

WR wind colliding with RSG wind. The density of the CSM is shown about 1750 years after the initiation of the WR wind, showing the early development of the Vishniac thin shell instability. The left panel shows the CSM when no radiative losses are taken into account. The middle panel shows the case where the MB81 cooling prescription is applied, and the right panel shows the CSM with the cooling prescription as presented in this work. The black line indicates the contact discontinuity. The number of refinement levels used for the resolution is 5, i.e., $6.25 \times 10^{15}$ cm $\times \ 0.0625^\circ$ .

Open with DEXTER

Figure 3 shows the density of the circumstellar medium after the WR-phase has just set in for about 1750 years, and shows the early developments of the Vishniac thin shell instability. The fast wind from the Wolf-Rayet star has swept up the red supergiant wind of the previous evolutionary stage and four regions can be identified. From inside (R=0) out: the free streaming fast WR wind, the shocked WR wind. The black line indicates the contact discontinuity, behind which the shocked, and farther out the free streaming, RSG wind can be found. Only part of the grid is shown to enlarge the interesting region.

The left panel shows the CSM for the colliding winds in the absence of cooling. The shell is smooth and thick, with a maximum compression ratio that is similar to that of a strong shock, i.e., a factor 4. The forward shock is partly outside the plotted region of the grid and extends to a radius of about $1.2 \times 10^{18}$ cm. The middle and right panels illustrate results when cooling is present, for the cooling rate given by MacDonald & Bailey (1981) in the middle panel, and for the cooling rate calculated with the SPEX code in the right panel. The instability is induced by random density perturbations of the 1% level, which may induce little differences between runs. If we do not include the initial perturbations, the code retains the original spherical symmetry to high accuracy, and the instabilities do not arise until very late. The extrapolation of the MB81 cooling curve towards lower temperatures results in higher cooling rates below 10⁴ K. This gives rise to more cooling in the compressed RSG shell and a more patchy structure as can be seen in the middle panel of Fig. 3. For temperatures above 10⁴ K, cooling rates are higher for the SPEX curve. This applies to the region of shocked WR wind material and results in this area being more compressed in the right panel. Consequently, the WR shell in the middle panel extends a little farther out.

Total energy losses are higher with our new cooling prescription. The comparison is made in Fig. 4. The thermal energy per unit volume is plotted for the different cases in Fig. 4. The inclusion of cooling reduces the thermal energy substantially. We note again that the initial random perturbation will influence the later evolution, explaining some of the inter-resolution differences. However, the key point is the systematic reduction between the SPEX versus MB81 cooling prescriptions.

$\begin{figure} \par\resizebox{8.8cm}{!}{ \includegraphics{12495fg4}} \end{figure}$

Figure 4:

Thermal energy in the stellar winds. The black solid line shows the energy when no cooling is taken into account. The black dotted curve shows the energy when the MB81 cooling rates are used, and the black dashed curve shows the system for the cooling rates presented in this work. The gray lines show the resolution dependence of the thermal energy, where nl indicates the number of refinement levels used.

Open with DEXTER

$\begin{figure} \par\resizebox{8.8cm}{!}{ \includegraphics{12495fg5}} \end{figure}$

Figure 5:

The Vishniac thin shell instability depends on numerical resolution. At high resolution, the shell becomes thinner and thus the wavelength shorter. The different panels correspond to refinement levels of 4, 6, and 7, corresponding to maximum resolution of $1.25 \times 10^{16}$ cm $\times$ $0.125^\circ$ , $3.125\times 10^{15}$ cm $\times \ 0.03125^\circ$ , and $1.56\times 10^{15}$ cm $\times$ $0.0156^\circ$ . The simulation with 5 refinement levels is shown in Fig. 3.

Open with DEXTER

The wavelength of the Vishniac thin shell instability depends on the thickness of the shell. In Fig. 5, we show the effect of resolution on the development of the instabilities with the new cooling curve. For higher resolution, the shell becomes thinner, and as a result the wavelength of the instability shortens. Although the WR shell is not resolved in terms of the dynamics, we see that for thermal losses, convergence seems to be reached for the resolutions used.

A typical WR phase lasts for about 10⁵ years. In Fig. 6, we show the development of the WR shell later in the evolutionary stage, at which time the WR phase has lasted for about 22 700 yr. It is clear that the presence of cooling majorly alters the CSM evolution, by losing energy by radiation, and by the irregular shell that develops as a result of this. At later times, shells from earlier evolutionary phases may intervene with the CSM evolution and simulations need to take these phases into account (see e.g., Dwarkadas 2007; van Marle et al. 2007).

$\begin{figure} \par\includegraphics[width=8.8cm,clip]{12495fg6} \end{figure}$	Figure 6: Unstable WR shell at a time of about 22 700 yr after its transition to a WR star. The shell becomes highly irregular due to the Vishniac thin shell instability. The black line indicates the contact discontinuity. This simulation used 6 refinement levels for the resolution.
Open with DEXTER

6 Conclusion and discussion

SPEX is a very complete package that calculates the emissivity of a plasma at a certain temperature, taking into account the abundances. We calculate the radiative losses for a plasma in CIE, and compare these with the curve by Sutherland & Dopita (1993), which was calculated with the same basic assumptions. We find that cooling rates calculated by SPEX tend to be higher than those of Sutherland & Dopita (1993). This is not surprising, since more lines have since been added and the SPEX package is, to our knowledge, the most complete package for X-rays and EUV. For low temperatures, we reuse the cooling prescriptionby Dalgarno & McCray (1972). Using that, and the SPEX results, we provide a recommended cooling curve for Solar abundances, which can be adjusted to different ionisation fractions at low temperatures. Additionally, changes can be decided for each cooling curve according to specific requirements on the abundances, by using the separate cooling rates of the individual elements.

In our application, we calculate the CSM of a massive star that experiences a RSG phase followed by a WR phase. The colliding winds form a shell that, in the presence of cooling, becomes very thin and susceptible to the Vishniac thin shell instability. Taking into account the cooling of the plasma is very important when accurately simulating systems such as these.

Acknowledgements

We thank the referee, J.C. Raymond, for his helpful comments that helped us to improve the manuscript. This study has been financially supported by a Vidi grant from the Netherlands Organisation for Scientific Research (NWO). This work was sponsored by the Stichting Nationale Computerfaciliteiten (National Computing Facilities Foundation, NCF) for the use of supercomputer facilities, with financial support from the Nederlandse Organisatie voor Wetenschappelijk Onderzoek (Netherlands Organisation for Scientific Research, NWO).

References

Anders, E., & Grevesse, N. 1989, Geochim. Cosmochim. Acta, 53, 197 [NASA ADS] [CrossRef]
Arnaud, M., & Raymond, J. 1992, ApJ, 398, 394 [NASA ADS] [CrossRef]
Arnaud, M., & Rothenflug, R. 1985, A&AS, 60, 425 [NASA ADS]
Badenes, C., Bravo, E., Borkowski, K. J., & Domínguez, I. 2003, ApJ, 593, 358 [NASA ADS] [CrossRef]
Bertschinger, E. 1986, ApJ, 304, 154 [NASA ADS] [CrossRef]
Cox, D. P., & Raymond, J. C. 1985, ApJ, 298, 651 [NASA ADS] [CrossRef]
Dalgarno, A., & McCray, R. A. 1972, ARA&A, 10, 375 [NASA ADS] [CrossRef]
Dwarkadas, V. V. 2007, ApJ, 667, 226 [NASA ADS] [CrossRef]
Field, G. B. 1965, ApJ, 142, 531 [NASA ADS] [CrossRef]
Garcia-Segura, G., & Mac Low, M.-M. 1995, ApJ, 455, 160 [NASA ADS] [CrossRef]
García-Segura, G., Langer, N., & Mac Low, M.-M. 1996, A&A, 316, 133 [NASA ADS]
Kaastra, J. S., & Mewe, R. 2000, in Atomic Data Needs for X-ray Astronomy, ed. M. A. Bautista, T. R. Kallman, & A. K. Pradhan, 161
Kaastra, J. S., Mewe, R., & Nieuwenhuijzen, H. 1996, in UV and X-ray Spectroscopy of Astrophysical and Laboratory Plasmas: Proceedings of the Eleventh Colloquium on UV and X-ray (Tokyo: Universal Academy Press), ed. K. Yamashita, & T. Watanabe, 411
Kafatos, M. 1973, ApJ, 182, 433 [NASA ADS] [CrossRef]
Keppens, R., Nool, M., Tóth, G., & Goedbloed, J. P. 2003, CoPhC, 153, 317 [NASA ADS]
Klapisch, M., Schwob, J. L., Fraenkel, B. S., & Oreg, J. 1977, J. Opt. Soc. Amer. (1917-1983), 67, 148 [NASA ADS] [CrossRef]
Kosenko, D. I., Sorokina, E. I., Blinnikov, S. I., & Lundqvist, P. 2004, Adv. Space Res., 33, 392 [NASA ADS] [CrossRef]
Mac Low, M.-M., & Norman, M. L. 1993, ApJ, 407, 207 [NASA ADS] [CrossRef]
MacDonald, J., & Bailey, M. E. 1981, MNRAS, 197, 995 [NASA ADS]
Mellema, G., & Lundqvist, P. 2002, A&A, 394, 901 [NASA ADS] [EDP Sciences] [CrossRef]
Mewe, R., Kaastra, J. S., & Liedahl, D. A. 1995, Legacy, 6, 16
Raymond, J. C., & Smith, B. W. 1977, ApJS, 35, 419 [NASA ADS] [CrossRef]
Raymond, J. C., Cox, D. P., & Smith, B. W. 1976, ApJ, 204, 290 [NASA ADS] [CrossRef]
Ryu, D., & Vishniac, E. T. 1987, ApJ, 313, 820 [NASA ADS] [CrossRef]
Schure, K. M., Vink, J., García-Segura, G., & Achterberg, A. 2008, ApJ, 686, 399 [NASA ADS] [CrossRef]
Shapiro, P. R., & Moore, R. T. 1976, ApJ, 207, 460 [NASA ADS] [CrossRef]
Smith, B., Sigurdsson, S., & Abel, T. 2008, MNRAS, 385, 1443 [NASA ADS] [CrossRef]
Sutherland, R. S., & Dopita, M. A. 1993, ApJS, 88, 253 [NASA ADS] [CrossRef]
Toro, E. F. 1999, Riemann solvers and numerical methods for fluid dynamics - A practical introduction - 2nd edn. (Berlin: Springer)
van Marle, A. J., Langer, N., & García-Segura, G. 2007, A&A, 469, 941 [NASA ADS] [EDP Sciences] [CrossRef]
Vishniac, E. T. 1983, ApJ, 274, 152 [NASA ADS] [CrossRef]
Wiersma, R. P. C., Schaye, J., & Smith, B. D. 2009, MNRAS, 393, 99 [NASA ADS] [CrossRef]

Online Material

Table 3: Cooling curve for solar metallicity for T<10⁴ K, as computed from Dalgarno & McCray (1972), see Sect. 4.

Table 4: Cooling rates per element for solar metallicity.

Footnotes

... code: Tables 3 and 4 are only available in electronic form at http://www.aanda.org
... SPEX: http://www.sron.nl/spex
... abundances: http://www.phys.uu.nl/ schure/cooling

All Tables

Table 1: Solar abundances from Anders & Grevesse (1989).

Table 2: Cooling curve for solar metallicity.

Table 3: Cooling curve for solar metallicity for T<10⁴ K, as computed from Dalgarno & McCray (1972), see Sect. 4.

Table 4: Cooling rates per element for solar metallicity.

All Figures

$\begin{figure} \par\resizebox{8.8cm}{!}{\includegraphics{12495fg1}} \end{figure}$	Figure 1: Cooling curves compared: the higher cooling rates calculated with SPEX are mainly due to a more complete coverage of the line transitions, including Fe L and EUV lines.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.8cm]{12495fg2a} \includegraphics[width=8.8cm]{12495fg2b} \end{figure}$	Figure 2: Contributions of different elements to the cooling curve are given. Each of the plots shows a different set of elements. Important peaks are labelled with the name of the element. The total cooling curve (black solid line) is an addition of the individual elemental contributions.
Open with DEXTER
In the text

$\begin{figure} \par\resizebox{8.8cm}{!}{ \includegraphics{12495fg3}} \end{figure}$	Figure 3: WR wind colliding with RSG wind. The density of the CSM is shown about 1750 years after the initiation of the WR wind, showing the early development of the Vishniac thin shell instability. The left panel shows the CSM when no radiative losses are taken into account. The middle panel shows the case where the MB81 cooling prescription is applied, and the right panel shows the CSM with the cooling prescription as presented in this work. The black line indicates the contact discontinuity. The number of refinement levels used for the resolution is 5, i.e., $6.25 \times 10^{15}$ cm $\times \ 0.0625^\circ$ .
Open with DEXTER
In the text

$\begin{figure} \par\resizebox{8.8cm}{!}{ \includegraphics{12495fg4}} \end{figure}$	Figure 4: Thermal energy in the stellar winds. The black solid line shows the energy when no cooling is taken into account. The black dotted curve shows the energy when the MB81 cooling rates are used, and the black dashed curve shows the system for the cooling rates presented in this work. The gray lines show the resolution dependence of the thermal energy, where nl indicates the number of refinement levels used.
Open with DEXTER
In the text

$\begin{figure} \par\resizebox{8.8cm}{!}{ \includegraphics{12495fg5}} \end{figure}$	Figure 5: The Vishniac thin shell instability depends on numerical resolution. At high resolution, the shell becomes thinner and thus the wavelength shorter. The different panels correspond to refinement levels of 4, 6, and 7, corresponding to maximum resolution of $1.25 \times 10^{16}$ cm $\times$ $0.125^\circ$ , $3.125\times 10^{15}$ cm $\times \ 0.03125^\circ$ , and $1.56\times 10^{15}$ cm $\times$ $0.0156^\circ$ . The simulation with 5 refinement levels is shown in Fig. 3.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=8.8cm,clip]{12495fg6} \end{figure}$	Figure 6: Unstable WR shell at a time of about 22 700 yr after its transition to a WR star. The shell becomes highly irregular due to the Vishniac thin shell instability. The black line indicates the contact discontinuity. This simulation used 6 refinement levels for the resolution.
Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.