EVOLUTION OF MASSIVE STARS WITH PULSATION-DRIVEN SUPERWINDS DURING THE RED SUPERGIANT PHASE

Mass loss due to stellar winds is one of the governing factors that determine the evolution of massive stars. Great improvements in our understanding of line-driven winds from massive hot stars (i.e., O-Type stars or Wolf-Rayet stars) have been achieved over the last three decades, and they are well reflected in recent stellar evolution models (see Puls et al. 2008, for a recent review). However, the driving mechanism for cool-star winds remains elusive. Note that most of the authors still rely on the empirical mass-loss rate given by de Jager et al. (1988, hereafter, JNH88), which is more than 20 years old, for modeling cool giant stars (see, however, Vanbeveren et al. 1998 and Salasnich et al. 1999 for stellar models that adopt alternative prescriptions). Late stages of massive star evolution during and after core helium burning are therefore relatively difficult to study, which is a major obstacle to theoretical investigation on supernova (SN) progenitors and their environments.

One of the key factors that should be taken into account for the study of the red supergiant (RSG) evolution is pulsation. It is well known that RSGs are unstable to radial pulsations driven by partial ionization of hydrogen in the outermost layers of the envelope (e.g., Li & Gong 1994), and pulsation is one of the most commonly invoked mechanisms for driving winds from cool stars (cf. Bowen 1988; Höfner et al. 2003; Neilson & Lester 2008). Although the detailed mechanisms for driving mass loss by pulsations are not well understood, simple arguments based on the gain of pulsation energy suggest that RSGs might experience a stronger mass loss for a higher growth rate of the pulsational instability (e.g., Appenzeller 1970). In particular, Heger et al. (1997) found that RSG pulsations can achieve large amplitudes that may be comparable to those in asymptotic giant branch stars when the pulsation periods become close to the thermal timescale of the envelope. They therefore suggested that such strong pulsations might lead to a “superwind” if the luminosity-to-mass ratio in an RSG were sufficiently high (see van Loon et al. 2005, 2008 for possibly related observations).

In the present study, we explore possible consequences of such a superwind in the evolution of SN progenitors. Our approach to this topic is as follows. First, we investigate the properties of RSG pulsations using our new stellar evolution models that follow the nonlinear evolution of pulsation (Section 2). This allows us to derive a relation between the pulsation growth rate and the structure of RSGs. Second, we assume that the mass-loss rate due to RSG winds should be significantly enhanced compared to that of JNH88, when the growth rate of pulsation becomes sufficiently large. The consequent evolution of massive stars can be very different from the standard picture given in previous work (Section 3). Then, we conclude the Letter by discussing implications for SN progenitors (Section 4).

RSG pulsations are not usually found in stellar evolution models, even if one uses a hydrodynamic code, for two reasons. First, in most stellar evolution calculations, the evolutionary time steps are much larger than the pulsation periods in RSGs, which are comparable to the dynamical timescale of the envelope. Second, pulsations tend to be numerically damped out with the implicit method adopted in stellar evolution codes. However, the study by Heger et al. (1997) shows that the predictions of linear stability analysis for basic pulsational properties (i.e., pulsation period and growth rate) can be accurately reproduced in nonlinear evolutionary calculations of RSGs by a hydrodynamic stellar evolution code, if the growth rate is sufficiently large and if the adopted time steps are sufficiently small. On the other hand, nonlinear evolutionary calculations do not require the condition of complete thermal and hydrostatic equilibrium of the star, in contrast to the case of linear stability analysis. We therefore decide to use an existing hydrodynamic stellar evolution code for our analysis of RSG pulsations. Details of our implicit hydrodynamic code are described in Yoon et al. (2006) and references therein.

To start, we construct non-rotating hydrostatic models at solar metallicity for 11 different initial masses (15, 16, 17, 18, 19, 20, 21, 23, 25, 30, and 40 M_☉). Since the hydrodynamic term in the momentum equation is switched off, pulsations do not appear in these model sequences. The H-R diagrams for some model sequences are shown in Figure 1.

Zoom In Zoom Out Reset image size

Then, we turn on the hydrodynamic term and restart the calculation from the reference models chosen at various evolutionary epochs during the supergiant phase for each model sequence. In total, 71 reference models are selected as starting points. In Figure 1, some of the starting points are shown, marked by filled circles and open diamonds. To guarantee enough time resolution for describing pulsations in the hydrodynamic calculations, we limit the time steps such that Δt ⩽ 0.01[R/(1000 R_☉)]²[(15 M_☉)/M] yr. This choice is based on the fact that pulsations of an RSG with M_init>15 M_☉ have a period on the order of 1000 days (e.g., Heger et al. 1997), and obey the relation P ∝ R²/M, where P is the pulsation period, R is the radius, and M is the total mass of an RSG (Gough et al. 1965).

An example of such calculations is given in Figure 2. As the pulsation amplitude grows gradually, the surface velocity eventually reaches the sound speed. Since our code cannot describe shock waves, we investigate the pulsation properties only with the subsonic results. Figure 1 shows the pulsation period P and the growth rate η of the fundamental mode at different reference points on the H-R diagram. Here, the growth rate is defined as η = |v(t₀ + P)/v(t₀)|, where v(t₀) denotes a local maximum of the surface velocity located at a certain time t = t₀. According to this definition, pulsational instability means η>1.

Zoom In Zoom Out Reset image size

The growth rate obtained here should be considered only indicative, given the limitations of our adopted numerical method. For example, in our calculations, no pulsation is detected from some reference points. This does not mean that the considered RSGs are stable against the pulsational instability, but that the growth rate is too small for our code to follow. Also, in our calculations, instant adjustment of convective flux during pulsation is assumed, which is not physical given the similar timescales of convection and pulsation.³ Despite these uncertainties, our results are qualitatively in good agreement with those of previous studies: the growth rate is generally higher for a larger luminosity-to-mass ratio (L/M) and/or for a smaller thermal timescale of the envelope. Specifically, as shown in Figure 3, we find the following relation between the pulsation growth rate and the stellar structure from our results with the hydrodynamic model sequences:

$$\begin{equation} \log \eta = C_1\left(\frac{R}{{R_\odot }}\right)^2\left(\frac{M}{{M_\odot }}\right)^{-1} \left(\frac{\tau _\mathrm{KH,env}}{\mathrm{yr}}\right)^{-0.315} - C_2, \end{equation} \tag{ 1 }$$

where τ_KH,env (≔GM_tot M_env/LR) denotes the Kelvin–Helmoltz timescale of the hydrogen-rich envelope, C₁ = 9.219 × 10⁻⁶, and C₂ = 0.0392844.

Zoom In Zoom Out Reset image size

We confirm that RSGs are generally more susceptible to the pulsation instability for larger initial masses. Figure 1 indicates that for M_init>20 M_☉, pulsations can become strong enough to be detected with our code even before the star reaches the Hayashi line, while pulsations are only found after neon burning for M_init = 17 M_☉. We do not find any signature of pulsation for M_init = 15 M_☉.⁴ Therefore, if a superwind could be induced by strong pulsations, it would occur in an earlier evolutionary epoch for a higher mass star.

If the mass-loss rate were dramatically enhanced due to strong pulsations, how would it affect the evolution of supergiants and what would be the observational consequences?

To address this issue in more detail, we perform an experiment as the following. We assume that the mass-loss rate of a supergiant is significantly enhanced compared to the rate of JNH88, if pulsations are strong enough to be detected with our hydrodynamic stellar evolution code, i.e., if η>1 in Equation (1):

$$\begin{equation} \dot{M} = \eta ^{\alpha }\dot{M}_\mathrm{JNH88}, \end{equation} \tag{ 2 }$$

where $\dot{M}_\mathrm{JNH88}$ is the mass-loss rate of JNH88 and α is a free parameter. New evolutionary model sequences are constructed with this modified mass-loss rate. Here, the hydrodynamic term is switched off again, and normal evolutionary time steps are used.

It should be kept in mind that the empirical mass-loss rate of RSGs of JNH88 might already reflect the role of pulsations. However, as implied by the growth rates shown in Figure 1, only very unstable massive stars with M_init ≳ 19 M_☉ are significantly affected by this modification. Given the paucity of RSGs with M_init  >  20–25 M_☉ in our Galaxy and nearby galaxies (e.g., Massey et al. 2009), our assumption that the mass-loss history of such very unstable RSGs can be qualitatively different from that of most RSGs in the sample of JNH88 does not necessarily cause an inconsistency.

Figure 4 shows the evolutionary tracks and the mass-loss rates of model sequences with α = 1.75 (this choice for α is rather arbitrary), compared to the corresponding ones with α = 0, for two different masses (M_init= 20 and 25 M_☉). Two important points are worthy of note here.

Zoom In Zoom Out Reset image size

First, as the star gradually loses mass, the L/M ratio and the pulsation growth rate increase accordingly. If the value of α is large enough (as adopted here), this leads to a runaway increase in the mass-loss rate as shown in the lower panel of the figure. Second, the mass-loss rate decreases very rapidly when the star moves away from the pulsationally unstable regime on the H-R diagram. The peak of the mass-loss rate is reached when the helium mass fraction at the center decreases to 0.54 (0.031) and the total mass to 9.35 M_☉ (6.50 M_☉), for 25 M_☉ (20 M_☉) model sequence. As shown for the 20 M_☉ sequence, such a mass-loss history of runaway increase followed by sudden decrease can be repeated as the star moves in and out the unstable regime. The 20 M_☉ star finally ends its life as a yellow supergiant of M_tot = 6.1 M_☉, with only a small amount of hydrogen of about 0.5 M_☉. The most likely outcome of the death of such a star would be a Type IIb SN. The 25 M_☉ star has a hydrogen envelope of 0.22 M_☉ at core helium exhaustion and would produce either a Type Ib or a IIb SN depending on the subsequent history of mass loss.

When the mass-loss rate reaches the maximum, the growth rate given by Equation (1) becomes as high as 11.3 and 10.6 for 25 and 20 M_☉ models, respectively. We checked if nonlinear evolutionary calculations also give such high growth rates, given that the relation of Equation (1) is only based on the result with η ≲ 8.0 (Figure 3). We find that, at the reference points marked by the filled circles in Figure 4, the hydrodynamic calculations give η= 9.9 and 9.0 for 25 and 20 M_☉ models, respectively, which are comparable to the values given by Equation (1).

If a pulsationally driven superwind (PDSW) phase could be induced by strong pulsations, this would have very important implications for SN progenitors. As implied by our model sequences presented above, such a PDSW can significantly reduce the critical zero-age main-sequence mass (M_crit) above which a huge fraction of the stellar envelope is removed before core collapse, thus producing no Type II-P SN. Based on our results, we suggest M_crit∼ 19–20 M_☉ with a PDSW (see Figures 1 and 4), while the models with the JNH88 mass-loss rate predict SN II-P progenitors up to ∼25 M_☉. Note that M_crit could be even further reduced by rotation, as implied by the result of Heger et al. (1997).

This might provide a plausible solution to the so-called RSG problem, i.e., the observed lack of Type II-P progenitors with M_init ≳ 16.5 ± 1.5 M_☉(Smartt et al. 2009). Such observation could also result from the presence of dusty circumstellar material around the stars, as an obscured progenitor would be estimated to have a lower initial mass in pre-SN images. However, Smartt et al. (2009) also noted that the number of SNe Type II-P in their sample is consistent with the expected number of stars in the range 8.5–16.5 M_☉ assuming a Salpeter initial mass function. This could be a possible indication of the fact that stars with higher initial mass do actually lose most of their hydrogen envelope, which cannot be easily understood with the canonical mass-loss rate of JNH88.

The PDSW is expected to significantly affect the circumstellar medium around an RSG (B. van Veelen 2011, in preparation). When the star dies, the shock produced by the collision between the SN ejecta and the circumstellar material transforms kinetic energy into thermal energy. This energy can be radiated away at different wavelengths, resulting in brightening the SN remnant for long times (Chevalier 1977). If the collision occurs directly after the SN explosion, the emission can even alter the spectrum and light curve of the SN. This scenario is usually invoked to explain the subclass of Type IIn SNe (Schlegel 1990; Filippenko 1997). To reproduce the light curve of the most luminous Type IIn SNe (SNe IIn) like SN 2006gy and SN 2006tf, very massive shells are required, which need to be ejected in eruptive events a few years before core collapse (e.g., van Marle et al. 2010). Pulsational pair instability (Woosley et al. 2007) and luminous blue variable like eruptions (e.g., Smith et al. 2007) have been discussed to explain the presence of shells of 10–25 M_☉ around the progenitors of luminous SNe IIn. A PDSW phase is an unlikely explanation for such extreme environments. However, it might be interesting for those Type IIn where the required mass in the stellar vicinity is on the order of a few solar masses or less.

In this context, it is interesting to consider the circumstellar medium around the RSG VY CMa. The stellar surrounding appears shaped by episodic mass ejections which occurred about 500–1000 year ago (Smith et al. 2009). The complex morphology of the circumstellar matter (CSM) is suggestive of a possible interaction between convection and pulsation (see Figure 13 in Smith et al. 2009), which for these stars is predicted to occur on similar timescales (e.g., Heger et al. 1997). The mass-loss rate of (1–2) × 10⁻³ M_☉ yr^-1 derived by Smith et al. (2009) is comparable to the one expected during the PDSW phase (see Figure 4). Even if the mass-loss prescription for PDSW we used is somewhat arbitrary, we want to stress that the energy available from the growth of pulsation is enough to drive mass-loss rates up to ∼10⁻² M_☉ yr^-1. Smith et al. (2009) state that an extreme RSG like VY CMa would produce a Type IIn event like SN 1988Z if it were to explode in its current state. Therefore, the occurrence of PDSWs in RSG might explain moderately luminous SN IIn if the enhanced mass loss takes place less than about ∼10³ years before core collapse. With the mass-loss prescription of Equation (2) with a sufficiently large α (∼2), this would occur in a narrow range around ∼19 M_☉ for non-rotating stars. An accurate determination of the expected rate of Type IIn due to PDSW requires a more realistic mass-loss prescription.

To conclude, if a PDSW phase could be induced by strong pulsations, we would expect a substantial change in the late evolution of single massive stars. The mass-loss rate during the RSG phase would increase dramatically for stars with M_init ≳ M_crit compared to the JNH88 rate, and the PDSW phase should start earlier in the evolution of more massive stars. The resulting pre-SN structure of these stars is affected as well as their CSM. For single stars, this suggests the following sequence in SN types as a function of increasing initial mass: II-P → II-L → IIn → IIb → Ib → Ic. Detailed estimate for the mass range of each SN requires a better prescription for the pulsationally enhanced mass loss, which deserves further study.

We are grateful to Stan Woosley, Bob van Veelen, and Norbert Langer for many useful discussions. S.C.Y. is supported by the DOE SciDAC Program (DOE DE-FC02-06ER41438).