ON CARBON BURNING IN SUPER ASYMPTOTIC GIANT BRANCH STARS

When a single star on the main sequence (MS) exhausts the supply of hydrogen in its core, the core contracts and its temperature increases, while the outer layers of the star expand and cool. The star becomes a red giant (e.g., Iben 1991; Stancliffe et al. 2009; Karakas & Lattanzio 2014). The subsequent onset of helium burning in the core causes the star to populate the horizontal branch for more metal-poor stars or the red clump for more metal-rich stars (Cannon 1970; Faulkner & Cannon 1973; Seidel et al. 1987; Castellani et al. 1992; Girardi 1999). After the star exhausts the supply of helium in its core, the carbon–oxygen (henceforth CO) core contracts while the envelope once again expands and cools along a path that is aligned with its previous red giant track. The star becomes an asymptotic giant branch (AGB) star (e.g., Hansen et al. 2004; Herwig 2005; Kippenhahn et al. 2012; Fishlock et al. 2014; Salaris et al. 2014).

The minimum mass for carbon ignition is usually referred to as M_up ≈ 7 ${M}_{\odot }$ and the minimum mass for neon ignition in the core is traditionally referred to as M_mas ≈ 10 ${M}_{\odot }$ (Becker & Iben 1979, 1980; García-Berro et al. 1997). Stars with zero-age MS (ZAMS) masses between ≈7 M_⊙ and ≈10 ${M}_{\odot }$ are designated as super-AGB stars (henceforth SAGB, Ritossa et al. 1996, 1999; Gil-Pons et al. 2005; Siess 2006, 2007, 2010; Poelarends et al. 2008; Doherty et al. 2010). Due to the slope of the stellar initial mass function inferred from observations (e.g., Jennings et al. 2012), single stars in this ZAMS mass range represent the population of stars that can produce the most massive white dwarfs, the most numerous supernovae, and possibly the least massive neutron stars (e.g., Doherty et al. 2015). SAGB stars may also make significant contributions to the Galactic inventory of isotopes such as ⁷Li, ¹⁴N, ²³Na, ^25–26Mg, ^26–27Al, and ⁶⁰Fe (Siess 2010; Ventura et al. 2013; Doherty et al. 2014a, 2014b).

After helium is exhausted in the core, stars ascending the SAGB develop partially electron degenerate CO cores ranging from ≈0.9 ${M}_{\odot }$ to ≈2.0 ${M}_{\odot }$ (pioneering studies of CO cores include Rakavy et al. 1967; Beaudet & Salpeter 1969; Boozer et al. 1973). Depending primarily on the ZAMS mass but also sensitively on the composition mixing model (Poelarends et al. 2008; Siess 2009), the ignition of carbon may not occur at all (for stars ≲7 ${M}_{\odot }$), may occur at the center of the star (for stars ≳10 ${M}_{\odot }$), or may occur somewhere off-center. In the off-center case, ignition is followed by the inward propagation of a subsonic burning front (Nomoto & Iben 1985; Timmes et al. 1994; García-Berro et al. 1997; Saio & Nomoto 1998). Trailing behind the burning front is a convective region that may extend outward ≈0.6 ${M}_{\odot }$; see Figure 1 for an illustration. When a steady-state, convectively bounded, subsonic, carbon-burning front (henceforth a "flame") propagates toward the center of the CO core, the flame leaves behind oxygen–neon–sodium (ONeNa) ashes. The inward-propagating carbon flame may or may not reach the center of the star, depending on the parameters adopted for composition mixing beyond the convective boundary set by mixing-length theory (MLT) (e.g., thermohaline, overshoot, semiconvection; Siess 2009; Stancliffe et al. 2009; Denissenkov et al. 2013). If the flame makes it to the center, then the original CO core is converted into an ONeNa core. Such SAGB stars can explode as electron capture supernovae (ECSNe) if their ONeNa core masses reach central densities in excess of the threshold density for the ²⁰Ne(e⁻, ν)²⁰F electron capture reaction (Miyaji et al. 1980; Nomoto 1984; Gutierrez et al. 1996; Jones et al. 2013), as may be the case for the Crab Nebula (Davidson et al. 1982; Nomoto et al. 1982; Wanajo et al. 2009) or for potentially explaining observations of sub-luminous SNe II-P (Smartt 2009). If the flame does not make it to the center, then the star is left with inner parts of the original CO core surrounded by a layer of ONeNa. Such hybrid white dwarfs may provide unusual SNe Ia progenitors (Siess 2009; Denissenkov et al. 2013; Chen et al. 2014; Wang et al. 2014).

Zoom In Zoom Out Reset image size

This paper explores the ignition and subsequent evolution of carbon burning in SAGB stellar models as a function of the ZAMS mass, initial rotation rate, and the magnitude of various mixing parameters such as convective overshoot, semiconvection, thermohaline, and angular momentum transport. We sample this multidimensional space with 2755 MESA stellar evolution models (Paxton et al. 2011; Paxton et al. 2013) that are evolved from the pre-MS to the end of carbon burning. All models have Z = 0.02 and a solar composition from Grevesse & Sauval (1998). Along the way we provide quasi-analytical models for interpreting the results. In Section 2 we discuss the input physics for our calculations, the composition mixing processes considered, and the MESA implementation of rotation and magnetic fields. In Section 3 we define the baseline parameters used for our calculations and discuss the grids used for our exploration into this three-dimensional (3D) parameter space. In Section 4 we present the results of our non-rotating models, an analytical approximation of the location of first carbon ignition, and the evolution of carbon-burning flames and flashes. In Section 5 we present the results of the effect of rotation and overshoot on the ignition, evolution, and death of carbon burning. In Section 6 we study the impact of the semiconvection, thermohaline, and angular momentum transport coefficients on the location of first carbon ignition on our results. In Section 7 we present the results of spatial and temporal convergence studies on our results, and in Section 8 we discuss our results and their implications.

Our numerical instrument is MESA version 6794. We use the included sagb_NeNa_MgAl.net reaction network, which follows 22 isotopes from ¹H to ²⁷Al to track hydrogen (pp chains, CNO, NeNa, and MgAl cycles), helium and carbon burning. The 51 thermonuclear reaction rates coupling these isotopes are from JINA reaclib version V2.0 2013-04-02 (Cyburt et al. 2010), energy-loss rates and their derivatives from thermal neutrinos are from the fitting formulae of Itoh et al. (1996), and electron screening factors for thermonuclear reactions in both the weak and strong regimes are from Dewitt et al. (1973), Graboske et al. (1973), and Alastuey & Jancovici (1978), with plasma parameters from Itoh et al. (1979). Poelarends et al. (2008) showed that increasing the mass loss rate could decrease the mass range for systems that will become ECSNe. We thus use a Reimers mass loss prescription (Reimers 1975) with η = 0.5 on the RGB and a Blöcker mass loss prescription (Bloecker 1995) with η = 0.05 on the AGB. The MESA inlists are publicly available.³

Analysis of a carbon-burning event requires knowledge of when and where carbon ignites. We identify the cell location of carbon burning, hence the ignition mass coordinate (M_f,s), by three criteria. First, we require ${\epsilon }_{\mathrm{nuc}}\gg {\epsilon }_{\nu }$ in a CO core. Second, we require that ⁴He is depleted in the ignition region as some of the lowest mass stars investigated would have a small amount of ¹²C burning near the boundary between the CO core and ⁴He shell. Finally, we require X(²⁰Ne) > X(²³Na) > X(²⁴Mg), which indicates that we have vigorous ¹²C +¹²C burning.

The end of a carbon-burning event is defined when no cell within 10% of the mass location of the flame, during the next time step, has _nuc $\gg$ _ν. We then define the final flame location (M_f,e), where carbon burning is quenched, as the minimum value of the mass location taken at the end of all the carbon-burning events. This is independent of any subsequent carbon flashes, as we record only the closest approach the carbon burning makes to the core.

2.1. Mixing

Treatment of convective processes within stellar interiors is essential for a physically accurate stellar model. We briefly discuss the composition mixing processes used in our calculations, how MESA models the mixing processes, and previous studies that guide our choices for our baseline mixing parameters. Values for our baseline parameters are given in Section 3.

We use the Schwarzschild criterion for convection along with the Cox implementation of MLT (Cox & Giuli 1968). The Schwarzschild criterion describes that a region is stable to convection if the temperature gradient of a piece of adiabatic matter is less than the temperature gradient of the stellar atmosphere:

$$\begin{eqnarray}&&{|\displaystyle \frac{d\ \mathrm{ln}\ T}{d\ \mathrm{ln}\ P}|}_{\mathrm{ad}}\lt {|\displaystyle \frac{d\ \mathrm{ln}\ T}{d\ \mathrm{ln}\ P}|}_{\mathrm{rad}}.\end{eqnarray} \tag{ 1 }$$

MLT has a free parameter, α_MLT, as described by Böhm-Vitense (1958). Values within $1.6\lesssim {\alpha }_{\mathrm{MLT}}\lesssim 2.2$ have been inferred by Noels et al. (1991), who compared observations of the α Centauri binary star system with stellar evolution models. Trampedach et al. (2014) also suggest $1.6\lesssim {\alpha }_{\mathrm{MLT}}\lesssim 2.05$ from calibrating one-dimensional (1D) stellar models to 3D radiation-coupled hydrodynamics simulations of convection in stellar surface layers.

Turbulent velocity fields have been suggested to decay exponentially beyond the Schwarzschild convective boundary defined by Equation (1) (Herwig et al. 1997; Ventura et al. 1998; Mazzitelli et al. 1999; Blöcker et al. 2000; Herwig 2000), leading to a diffusive treatment of mixing as

$$\begin{eqnarray}&&{D}_{\mathrm{OV}}={D}_{\mathrm{conv},0}\mathrm{exp}\left(-\displaystyle \frac{2z}{{f}_{\mathrm{ov}}{\lambda }_{P,0}}\right).\end{eqnarray} \tag{ 2 }$$

Here D_conv,0 is the convective diffusion coefficient at the convective boundary, z is the radial distance from the convective boundary, λ_P,0 is the local pressure scale height, and f_ov is an adjustable parameter. MESA offers the flexibility of allowing f_ov to be different for different convective regions (H burning, He burning, metal burning, and non-burning). However, we set f_ov to the same value in all convective regions (see Section 3 for the values chosen).

Semiconvection occurs when regions of the stellar interior are stable to the Ledoux criterion and unstable to the Schwarzschild criterion (Kippenhahn et al. 2012). This occurs when ∇_ad < ∇_T < ∇_L, where ∇_L = ∇_ad + B and B is the Brunt composition gradient. MESA treats semiconvection as a diffusive process (Langer et al. 1983, 1985; Heger et al. 2000; Zaussinger & Spruit 2013) with a diffusion coefficient

$$\begin{eqnarray}&&{D}_{\mathrm{sc}}={\alpha }_{\mathrm{sc}}\left(\displaystyle \frac{K}{6{C}_{p}\rho }\right)\displaystyle \frac{{{\rm{\nabla }}}_{T}-{{\rm{\nabla }}}_{\mathrm{ad}}}{{{\rm{\nabla }}}_{L}-{{\rm{\nabla }}}_{T}},\end{eqnarray} \tag{ 3 }$$

where K is the radiative conductivity, C_p is the specific heat at constant pressure, and ${\alpha }_{\mathrm{sc}}$ is an adjustable dimensionless parameter describing the speed with which convective mixing may occur at the boundary defined by Equation (1). Ongoing efforts to calibrate such semiconvection models include multidimensional numerical simulations of double-diffusive convection (Spruit 2013; Zaussinger & Spruit 2013) and comparing massive star models with observations (Yoon et al. 2006).

Thermohaline mixing occurs when ∇_T − ∇_ad ≤ B ≤ 0. These are regions stable against convection, according to the Ledoux criterion, but have an inversion of the mean molecular weight (Kippenhahn et al. 1980). This type of mixing forms elongated fluid parcels, sometimes called "salt-fingers." MESA treats thermohaline mixing as a diffusion process (Ulrich 1972; Kippenhahn et al. 1980; Brown et al. 2013; Zemskova et al. 2014):

$$\begin{eqnarray}&&{D}_{\mathrm{th}}={\alpha }_{\mathrm{th}}\displaystyle \frac{3\;K}{2\rho {C}_{p}}\displaystyle \frac{B}{\left({{\rm{\nabla }}}_{T}-{{\rm{\nabla }}}_{\mathrm{ad}}\right)},\end{eqnarray} \tag{ 4 }$$

where ${\alpha }_{\mathrm{th}}$ is a dimensionless parameter, related to the aspect ratio of the salt-fingers. Estimates for this parameter lie in the range 1 ≲ ${\alpha }_{\mathrm{th}}$ ≲ 667 (Kippenhahn et al. 1980; Charbonnel & Zahn 2007; Cantiello & Langer 2010; Stancliffe 2010; Wachlin et al. 2011), with some multidimensional simulations suggesting this parameter is significantly overestimated (Denissenkov 2010; Traxler et al. 2011; Denissenkov & Merryfield 2011; Lattanzio et al. 2015).

2.2. Rotation and Magnetic Fields

MESA implements rotation by making the assumption that the angular velocity, ω, is constant over isobars; see Paxton et al. (2013) for the implementation of rotation into MESA. Such an assumption is often referred to as the shellular approximation (Zahn 1992; Meynet & Maeder 1997), and allows the stellar structure equations to be solved in 1D for a rotating star. For this study, rotation is initialized by imposing a solid body rotation law at ZAMS, where the total luminosity equals the nuclear burning luminosity.

The transport of angular momentum and material due to rotationally induced instabilities is followed using a diffusion approximation (e.g., Endal & Sofia 1978; Pinsonneault et al. 1989; Heger et al. 2000, 2005; Maeder & Meynet 2003, 2004; Suijs et al. 2008) for the dynamical shear instability (DSI), secular shear instability (SSI), Eddington–Sweet circulation (ES), Goldreich–Schubert–Fricke instability (GSF), and Spruit–Tayler dynamo (ST). See Heger et al. (2000) for a description of the physics of the different instabilities and the calculation of the respective diffusion coefficients.

Berger et al. (2005) investigated the Ca line profiles of a sample of DA white dwarfs, concluding that their rotational velocities are generally less than 10 km s⁻¹. These values are significantly less than the values determined by the rotating, nonmagnetic models of Langer et al. (1999). Internal magnetic torques as proposed by Spruit (1998) have been suggested as an effective mechanism to spin down the cores of these white dwarf progenitors during the giant phase. Suijs et al. (2008) showed that magnetic torques as calculated in Spruit (2002) produce rotational velocities in better agreement with the observed values of Berger et al. (2005). We therefore include internal magnetic fields and the Spruit–Tayler dynamo angular momentum mechanism for our rotating MESA models.

Magnetic fields are implemented in MESA using the formalism of Heger et al. (2005), where a magnetic torque due to a dynamo (Spruit 2002) allows angular momentum to be transported inside the star. The radial component, B_r, and the azimuthal component, B_ϕ, of the magnetic field are modeled as

$$\begin{eqnarray}&&{B}_{\phi }\sim r{(4\pi \rho )}^{1/2}{\omega }_{{\rm{A}}}\end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray}&&{B}_{r}\sim {B}_{\phi }{({rk})}^{-1},\end{eqnarray} \tag{ 6 }$$

where r is the radial coordinate, ρ the density, ω_A the Alfvén frequency, and k the wavenumber. These magnetic fields then provide a torque

$$\begin{eqnarray}&&S=\displaystyle \frac{{B}_{r}{B}_{\phi }}{4\pi }\end{eqnarray} \tag{ 7 }$$

which acts to slow down the star's rotation rate by decreasing the amount of differential rotation inside the star (Heger et al. 2005).

The initial rotation is normalized against the critical rotation rate for the star ${{\rm{\Omega }}}_{\mathrm{crit}}=\sqrt{(1-L/{L}_{\mathrm{edd}}){cM}/{R}^{3}}$, where c is the speed of light, M the mass of the star, R the stellar radius, L the luminosity and L_edd the Eddington luminosity. The initial magnetic field is B_r = B_ϕ = 0 for all our rotating models. Effects of rotationally induced mass loss are not included.

We define a set of baseline parameters and construct a number of grids surrounding that baseline set to investigate the sensitivity of carbon burning in the SAGB models to variation in the baseline parameters. Choices in the numerical values of the baseline parameters are based on the current understanding of the canonical values for SAGB and other stars when using MESA. Choices in the range of values explored were designed such that we could explore the various competing factors involved in stellar evolution models.

Table 1 lists the baseline mixing, spatial resolution, and temporal resolution parameters. The parameter f_c is the ratio of turbulent viscosity to the diffusion coefficient (Heger et al. 2000), f_μ is the sensitivity to μ-gradients (Heger et al. 2000), while δ_mesh controls the spatial resolution by determining the relative magnitude of changes between the adjacent cells (Paxton et al. 2011), and w_t controls the temporal resolution by modulating the magnitude of the allowed changes between time steps (Paxton et al. 2011). Baseline values for δ_mesh and w_t were based on computational requirements. However, see Section 7 for a discussion of their relative effects on our results. The other baseline parameters listed in Table 1 are discussed in Section 2.1.

Table 1.  Baseline Parameters

Parameter	Value
Mixing Length Theory (αMLT)	2.0000
Semiconvection (αsc)	0.0100
Thermohaline (αth)	1.0000
Overshoot (fov)	0.0160
Angular Momentum (${\eta }_{\mathrm{am}}$)	1.0000
Turbulent Viscosity (fc)	0.0333
μ-gradient Sensitivity (fμ)	0.0500
Mesh Delta Coefficient (δmesh)	0.5000
Variance Control Target (wt)	0.0001

Download table as:  ASCIITypeset image

3.1. Mass–rotation–overshoot Grid

We explore the (ZAMS mass, initial rotation, overshoot) parameter space with three slices through this 3D data cube. Table 2 lists the start, stop, and step values for two of the three quantities while holding the third quantity constant. The number of SAGB models is 1326 in the mass–rotation plane, 546 in the rotation–overshoot plane, and 546 in the mass–overshoot plane for a total of 2418 SAGB models. Figure 2 illustrates these three orthogonal slices though this 3D parameter space.

Zoom In Zoom Out Reset image size

Table 2.  Mass–Rotation–Overshoot Grid

Variable	Start	Stop	Step	Constant
MZAMS	6.0	11.0	0.1	fov = 0.016
${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$	0.0	0.5	0.02
${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$	0.0	0.5	0.02	MZAMS = 8.0
fov	0.0	0.02	0.001
MZAMS	6.0	11.0	0.2	${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25
fov	0.0	0.02	0.001

Download table as:  ASCIITypeset image

For each slice within the 3D data cube we chose one quantity to be held fixed. In the mass–rotation plane this is the overshoot value, and our choice of f_ov = 0.016 is based on the canonical value for this overshoot model (Herwig 2000). In the rotation–overshoot plane we held the initial mass fixed at M = 8 ${M}_{\odot }$, which was selected because, based on the non-rotating models, we expected we could induce a range of behaviors from non-ignition, to off-center ignition, to central ignition. In the mass–overshoot plane, we held the initial rotation fixed at ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}=0.25$, purely as a middle-ground value between non-rotating models and our upper bound of ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.5. Our choice for the range of values covered by each grid was based on a requirement to have a comprehensive sample over the canonical baseline values for SAGB stars.

3.2. Mixing Coefficients Grid

We investigate the (semiconvection, thermohaline, overshoot, angular momentum diffusion) parameter space on the location of the first carbon ignition in an 8 ${M}_{\odot }$ ZAMS, ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25 model with selected points in this 4D data cube. Table 3 lists these quantities and their selected values. We choose to limit the range of ${\alpha }_{\mathrm{th}}$ in this grid to span only the lower values discussed in the literature. ${\eta }_{\mathrm{am}}$ is a scale factor that alters the strength of angular momentum diffusion in each cell; see Paxton et al. (2013) for details. The total number of SAGB models in this grid is 4⁴ = 256. This grid permits assessment of the relative strengths of each mixing process.

Table 3.  Mixing Coefficients Grid

Variable	Values
${\alpha }_{\mathrm{sc}}$	0.0	10−3	10−2	10−1
${\alpha }_{\mathrm{th}}$	0.0	0.1	1.0	10.0
fov	0.000	0.001	0.016	0.020
${\eta }_{\mathrm{am}}$	0.0	0.5	1.0	1.5

Download table as:  ASCIITypeset image

3.3. Spatial and Temporal Convergence Grid

Finally, we examine the spatial and temporal convergence properties of a subset of our SAGB models. The MESA parameter ${\delta }_{\mathrm{mesh}}$ broadly controls the spatial resolution and ${w}_{t}$ broadly relates to the temporal resolution. For ${\delta }_{\mathrm{mesh}}$ = 0.5 there are ≈5000−10,000 spatial points from center to surface. Spatial resolutions necessary to capture carbon-burning flames and flashes are discussed in Section 7. For ${w}_{t}$ = 10⁻⁴ a temporal resolution of ≈10 years is achieved. Table 4 lists these quantities and their values. The total number of SAGB models in this resolution sensitivity grid is 3⁴ = 81.

Table 4.  Spatial and Temporal Convergence Grid

Variable	Values
MZAMS	7	8	9
${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$	0.0	0.25	0.5
${\delta }_{\mathrm{mesh}}$	0.1	0.5	1.0
${w}_{t}$	10−5	10−4	10−3

Download table as:  ASCIITypeset image

We begin by considering a series of non-rotating stellar models using our baseline mixing parameters as described in Section 2. Figure 3 shows the Kippenhahn plots of these non-rotating stars as representative samples for all the stellar models. In the 7 ${M}_{\odot }$ case, carbon ignites off-center at ≈0.6 ${M}_{\odot }$. However, it rapidly quenches and does not propagate toward the core. The 7.5 ${M}_{\odot }$ case undergoes a series of carbon flashes, with each flash igniting closer to the core, but it is unable to form a steady-state flame. In the 8 ${M}_{\odot }$ case, an off-center ignition occurs at ≈0.15 ${M}_{\odot }$. A flame propagates inwards and through a series of distinct flashes as it approaches the core, and almost reaches the center. For the 9 ${M}_{\odot }$ case, carbon ignites at the center. In both the 8 and 9 ${M}_{\odot }$ models, secondary flashes at ≈0.5 ${M}_{\odot }$ are due to off-center carbon burning. For all cases the core is undergoing significant cooling, primarily through photo-neutrino and plasma neutrino losses prior to the first ignition of carbon (Nomoto 1984, 1987; Ritossa et al. 1996, 1999). In the 7.5, 8.0, and 9.0 ${M}_{\odot }$ cases the stars undergo a series of subsequent carbon flashes that travel away from the core. For the 9.0 ${M}_{\odot }$ case with no overshoot (bottom right in Figure 3) the flame ignites off-center, contrary to the center ignition of the 9 ${M}_{\odot }$ model with f_ov. The star undergoes a flash and then a flame, which propagates all the way to the center. The model with closest morphology is the 8 ${M}_{\odot }$ case, which has a flash and flame, but carbon burning does not reach the center. The difference between models with overshoot and without will be discussed further in Section 6.

Zoom In Zoom Out Reset image size

The He and CO core masses for the non-rotating models are shown in Figure 4 as a function of the ZAMS mass. Stellar models with M_ZAMS < 7 ${M}_{\odot }$ do not ignite carbon (Becker & Iben 1979, 1980) and are not shown. The carbon core mass increases linearly with ZAMS mass over the range shown. For M_ZAMS ≲ 7.5 ${M}_{\odot }$, the He forms a radiative, geometrically thin, burning shell accreting material onto the CO core at a rate of ≈10⁻⁶$\dot{{M}_{\odot }}\ {\mathrm{yr}}^{-1}$. Between ≈7.5 ${M}_{\odot }$ and ≈8.0 ${M}_{\odot }$, Figure 4 shows that the He shell transitions from a geometrically thin shell to a geometrically thick shell. For M_ZAMS ≲ 8.0 ${M}_{\odot }$, the geometrically thick He shell mass grows linearly with ZAMS mass over the range shown, accreting material onto the CO core at ≈$3\times {10}^{-6}\dot{{M}_{\odot }}\ {\mathrm{yr}}^{-1}$ at 8.0 ${M}_{\odot }$, and increasing to $\approx 2.0\times {10}^{-5}\dot{{M}_{\odot }}\ {\mathrm{yr}}^{-1}$ for the 11 ${M}_{\odot }$ stellar models.

Zoom In Zoom Out Reset image size

As the star evolves, the He shell grows in mass, reaching ≈1.7 ${M}_{\odot }$ for the 7 ${M}_{\odot }$ stars and up to 3.2 for the 11 ${M}_{\odot }$ stars. Once the He shell reaches its peak size, the ⁴He depletes quickly, leaving a CO core. Shortly afterwards the ⁴He shell begins rapidly accreting material onto the CO core, eventually depleting itself as well. At the lowest ZAMS masses (≲7 ${M}_{\odot }$), ignition occurs after the ⁴He shell has accreted onto the CO core and the star has finished its second dredge-up (2DU) (Becker & Iben 1979). In the transition region where we have a series of flashes (≈7−7.8 ${M}_{\odot }$), the stars are undergoing their 2DU while igniting carbon. In the higher mass systems, ≈7.9−8.2 ${M}_{\odot }$, where we have steady-state flames or central carbon ignition, the star ignites carbon before the 2DU and before any significant accretion onto the CO core can occur. Above 8.3 ${M}_{\odot }$ we have dredge-out events (Ritossa et al. 1999), where the ⁴He shell grows an outwardly moving convection zone that merges with the inwardly moving convective envelope.

At the base of the flame, we can ask whether a packet of convective material can penetrate into the region ahead of the flame, transferring energy that will decrease the flame's lifetime. A simplified derivation (L. Bildsten 2015, private communication) assumes that the length scale l of the flame front (the distance over which the temperature decreases from a peak inside the flame to the background value) is much less than the local pressure scale height, l $\ll$ H. This implies that the pressure is constant across the subsonic flame front. A convective packet will move from a region of high temperature to one of lower temperature, at constant pressure. Assuming adiabatic motion, a fluid packet keeps its original temperature. Thus, the buoyancy felt by this convective packet is

$$\begin{eqnarray}&&a=g\displaystyle \frac{d\mathrm{ln}T}{{dr}}\delta r,\end{eqnarray} \tag{ 8 }$$

where a is the buoyancy acceleration, g is the local acceleration due to gravity, and δr is the distance the packet moves ahead of the flame. Simplifying,

$$\begin{eqnarray}&&a\approx -\displaystyle \frac{g}{l}\delta r,\end{eqnarray} \tag{ 9 }$$

and solving this harmonic oscillator equation we find

$$\begin{eqnarray}&&\delta r\approx \displaystyle \frac{{v}_{c}}{N},\end{eqnarray} \tag{ 10 }$$

where v_c is the convective velocity and N is the Brunt–Väisälä frequency. Using the local scale height $H={c}_{s}^{2}/g$, where c_s is the local sound speed,

$$\begin{eqnarray}&&\displaystyle \frac{\delta r}{l}=\displaystyle \frac{{v}_{c}}{{c}_{s}}{\left(\displaystyle \frac{H}{l}\right)}^{1/2}.\end{eqnarray} \tag{ 11 }$$

In a typical carbon flame we have ${v}_{c}/{c}_{s}\approx {10}^{-4}$ and H/l ≈ 10, thus δr $\ll$ l and a convective fluid packet cannot penetrate the flame front.

4.1. First Ignition of Carbon

Figure 5 shows the mass location of the first ignition as a function of the stellar model ZAMS mass. For our choice of baseline mixing parameters, the lowest mass star to ignite carbon is a 7 ${M}_{\odot }$ model and ignition occurs off-center at a mass coordinate of ≈0.65 ${M}_{\odot }$. As the ZAMS mass increases, the location of the off-center first ignition moves steadily inwards in mass (Siess 2007). For ZAMS masses larger than 8.4 ${M}_{\odot }$, carbon ignites at the center. Figure 5 also shows the local density at the location of first ignition. All stellar models that ignite carbon off-center, 7 ${M}_{\odot }$ ≤ M_ZAMS ≤ 8.4 ${M}_{\odot }$, do so at a nearly constant density of log ρ ≈ 6.2, or ρ ≈ 1.5 × 10⁶ g cm⁻³. For stellar models that undergo central carbon ignition, $8.4\;{M}_{\odot }$ < M_ZAMS ≤ 11 ${M}_{\odot }$, the density at ignition monotonically decreases. Models that do not ignite carbon will eventually form a CO white dwarf (WD). Those models that undergo off-center carbon ignition but where the burning does not reach the center will form a hybrid CO+ONe WD. Model stars that ignite carbon at the center will eventually form an ONeNa WD, which may explode as an ECSNe.

Zoom In Zoom Out Reset image size

To a first approximation, at ignition the nuclear burning timescale and thermal diffusion timescales are equal (Timmes & Woosley 1992):

$$\begin{eqnarray}&&{\tau }_{\mathrm{diff}}\sim \displaystyle \frac{\sigma }{\rho \ {C}_{p}}\qquad {\tau }_{\mathrm{burn}}\sim \displaystyle \frac{E}{{\epsilon }_{\mathrm{nuc}}},\end{eqnarray} \tag{ 12 }$$

where σ is the thermal conductivity, C_p is the specific heat capacity at constant pressure, E is the local thermal energy, and _nuc is the screened nuclear energy generation rate during carbon burning including energy losses due to neutrino cooling. For a given temperature, density, and composition both C_p and E may be calculated from an equation of state. For carbon burning, _nuc takes the form (e.g., Woosley et al. 2004)

$$\begin{eqnarray}&&{\epsilon }_{\mathrm{nuc}}=6.7\times {10}^{23}\ {{X}}^{2}{(}^{12}{\rm{C}})\ {\rho }_{6}\ {f}_{\mathrm{screen}}\ {\lambda }_{\mathrm{12,12}}\ -\ {\epsilon }_{\nu },\end{eqnarray} \tag{ 13 }$$

where X(¹²C) is the carbon mass fraction, ρ₆ is the density divided by 10⁶ g cm⁻³, λ_12,12 is the unscreened nuclear reaction rate for carbon burning, and f_screen is the electron screening factor. Using a MESA 501 isotope reaction network that includes neutrino losses with an initial composition of X(¹²C) = 0.3 and X(¹⁶O) = 0.7 to calculate ${\epsilon }_{\mathrm{nuc}}$ over the relevant range in the ρ–T plane, we find that positive values of _nuc may be approximated by the power law

$$\begin{eqnarray}&&{\epsilon }_{\mathrm{nuc}}\approx 1.6\times {10}^{7}{\left(\displaystyle \frac{T}{7\times {10}^{8}}\right)}^{23}{\left(\displaystyle \frac{\rho }{2\times {10}^{6}}\right)}^{1.2}.\end{eqnarray} \tag{ 14 }$$

Results for carbon ignition for any large reaction network, including the 501 isotope network used here, are expected to be similar to those for the smaller 22 isotope net used in the SAGB models because both networks have the same set of key isotopes and reaction rates crucial for carbon burning. (e.g., Timmes et al. 2000).

Fitting the other quantities in Equation (12) in a similar manner, we find the following expressions for the nuclear burning timescale and thermal diffusion timescales:

$$\begin{eqnarray}&&{\tau }_{\mathrm{burn}}=5.1\times {10}^{9}{\left(\displaystyle \frac{T}{7\times {10}^{8}}\right)}^{-32}{\left(\displaystyle \frac{\rho }{2\times {10}^{6}}\right)}^{-0.8}\end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray}&&{\tau }_{\mathrm{diff}}=4.0\times {10}^{9}{\left(\displaystyle \frac{T}{7\times {10}^{8}}\right)}^{-2.4}{\left(\displaystyle \frac{\rho }{2\times {10}^{6}}\right)}^{-1}.\end{eqnarray} \tag{ 16 }$$

These two timescales are shown in Figure 6. Equating the two timescales gives

$$\begin{eqnarray}&&{\left(\displaystyle \frac{T}{7\times {10}^{8}}\right)}^{29.6}{\left(\displaystyle \frac{\rho }{2\times {10}^{6}}\right)}^{-0.2}=1.3.\end{eqnarray} \tag{ 17 }$$

At the threshold of vigorous carbon burning, T ≈ 7 × 10⁸ K, this expression gives a unique and constant ignition density

$$\begin{eqnarray}&&{\rho }_{\mathrm{ign}}\approx 2.1\times {10}^{6}\;{\rm{g}}\;{\mathrm{cm}}^{-3},\end{eqnarray} \tag{ 18 }$$

which is consistent with the constant ignition density of ρ ≈ 1.5 × 10⁶ g cm⁻³ found in the MESA SAGB models that ignite carbon off-center, 7 ${M}_{\odot }$ ≤ M_ZAMS ≤ 8.4 ${M}_{\odot }$, of Figure 5. We find this result also holds for our rotating SAGB models.

Zoom In Zoom Out Reset image size

Figure 7(a) shows that the carbon core of a non-rotating 7 ${M}_{\odot }$ ZAMS star can be well approximated by an n ≈ 3.0 polytropic model, (b) an 8 ${M}_{\odot }$ ZAMS star by an n ≈ 3.7 polytrope, and (c) a 9 ${M}_{\odot }$ ZAMS star by an n ≈ 4.3 polytrope. The density structures for these polytropic models were calculated using an open-source tool.⁴ In addition, the polytropic relations offer an explanation for why the location of the first ignition moves steadily inwards in mass for the off-center ignition cases (see Figure 5). The density structures of the n = 3.0, 3.7, and 4.3 polytrope models are shown in Figure 8. The mass location for a fixed ignition density (dashed horizontal line) where ignition occurs (dashed vertical lines) moves monotonically inwards as the polytropic index increases, as the carbon core mass increases, as the ZAMS mass increases.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We now turn to the decrease in the central density for those SAGB models in Figure 4 that centrally ignite carbon. Homology relations between the central density ρ_c and the mass M (e.g., Hansen et al. 2004; Kippenhahn et al. 2012) for a chemically homogeneous star characterized by a mean molecular weight μ, constant opacity, ideal gas equation of state, and power-law energy generation rate ${\epsilon }_{\mathrm{nuc}}\propto {T}^{\nu }$ indicate

$$\begin{eqnarray}&&{\rho }_{c}\propto {\mu }^{\frac{3(4-\nu )}{\nu +3}}{M}^{\frac{2(3-\nu )}{\nu +3}}.\end{eqnarray} \tag{ 19 }$$

For carbon burning near ignition, Equation (14) shows ν ≈ 23, and Equation (19) then gives ρ_c ∝ M^−1.5. The negative exponent shows that the density at first ignition, for those 8.4 ${M}_{\odot }$ < M_ZAMS ≤ 11 ${M}_{\odot }$ models that undergo central carbon ignition, monotonically decreases as the mass of the carbon core increases with a slope that is consistent with the rate of decline shown in Figure 5.

4.2. Carbon Burning flashes and Transition to a Steady-State Flame

Local characteristics before, during, and after the first off-center carbon flash in a 7.5 ${M}_{\odot }$ model are shown Figure 9. Before the first ignition of carbon, the CO core is in its most compact, most electron degenerate configuration. The first off-center carbon burning flash is thus the most powerful; any subsequent flashes or steady-state flames take place under more extended, less degenerate configurations. In addition, the energy released during the first ignition decreases as the ZAMS mass increases because the CO core is not as compact and not as degenerate. Furthermore, the first carbon flash impacts the base of the convective envelope more strongly in lower mass SAGB models than in higher mass SAGB models (García-Berro et al. 1997; Siess 2006), partly because of their more compact configuration but also because the first flash occurs farther from the center in lower mass models (see Figure 5).

Zoom In Zoom Out Reset image size

At the start of carbon burning at the ignition density of ${\rho }_{\mathrm{ign}}\;\approx \;$ 1.5 × 10⁶ g cm⁻³, marked by the vertical dotted line in Figure 9, the energy generation rate rapidly rises. Figure 9 shows that ≈10% of the energy produced by nuclear reactions is lost to neutrinos, with the balance of the thermal energy transported by convection (Ritossa et al. 1996; Saio & Nomoto 1998). This is not the ≈50% expected for steady-state burning, as the 7.5 ${M}_{\odot }$ model is undergoing a time-dependent flash and not a steady-state flame. Were it not for carbon burning, the surface luminosity and radius of these SAGB models would continue to increase in step with an increasing helium-burning luminosity, similar to lower mass, M_ZAMS ≲ 7 ${M}_{\odot }$ AGB stars that do not experience carbon burning (García-Berro et al. 1997; Karakas & Lattanzio 2014). These dramatic events mainly impact the core and have only a modest effect on the structure of the outer envelope and surface luminosity.

Figure 9 shows that the rapid injection of energy causes the entire core to expand and cool (Ritossa et al. 1999). A convective region also develops above the ignition location in response to the injection of energy from carbon burning. The overall expansion of the core decreases the electron degeneracy parameter ($\eta ={\mu }_{e}/{kT}\to -\infty$ means non-degenerate where μ_e is the electron chemical potential, $\eta \to \infty$ means perfect degeneracy). The expansion extinguishes the flash partially because the nuclear energy generation rate drops below the critical luminosity (e.g., Siess 2006) and partially because the ignition density ρ_ign is pushed deeper into the stellar model. The first flash thus quenches and does not become a steady-state flame that propagates inwards toward the center.

After the first flash is quenched and the deposition of energy from nuclear reactions ceases, the core again contracts but to a less compact, less degenerate configuration. This contraction leads to the second ignition of carbon. The second flash (and any subsequent flashes) occurs at a deeper mass coordinate because the accretion of C from the He burning shell slightly increases the core mass, which moves the location of the critical density ${\rho }_{\mathrm{ign}}$ inwards where there is also fresh, previously unburned fuel. Similar evolution pathways can occur for other fuels and other masses: for example, neon burning in more massive models (e.g., Jones et al. 2013).

The evolution of the second and any subsequent flashes is notably different than the first flash as shown in Figure 3. When these later flashes develop, their convective region grows into regions previously occupied by previous flashes where carbon had been depleted. Thus the nuclear energy production rate, which from Equation (13) is proportional to X²(¹²C), is reduced. Furthermore, the subsequent core expansion and induced thermodynamic changes are significantly smaller. The first flash lasts the longest and is the most vigorous, while subsequent flashes have shorter lifespans, release less energy, and the time interval between flashes becomes shorter. Models that have a geometrically thin helium envelope show more flashes than models with a geometrically thick helium envelope (see Figure 4).

Each successive flash releases less energy, expands the core by a smaller amount, and moves the mass location of the critical density inward at a slower rate. This allows a flash to transition to a steady-state flame. Thus, after one or more flashes a steady state may be achieved. Combustion in steady-state flames is incomplete; only a small portion of the carbon burns. A condition of balanced power is set up where the rate of energy emitted as neutrinos from the base of the convective region equals the power available from the unburned fuel that crosses the flame front. The inward propagation of the flame by thermal conduction is limited by the temperature at the base of the convective shell, which cannot greatly exceed the adiabatic value. These two local conditions give a unique speed for the flame, with typical values of ≈0.1 cm s⁻¹. We verified the flame speed at several locations in the SAGB models and it is commensurate with previous local studies, with speeds of order 10⁻³–10⁻² cm s⁻¹ (e.g., Timmes et al. 1994; Ritossa et al. 1996; García-Berro et al. 1997; Siess 2006; Denissenkov et al. 2013). The flame lives ≈20,000 years on its journey toward the center.

The convective nature of the material behind the flame has two key consequences for its journey to the core. First, the temperature behind the flame front (i.e., toward the surface) is bounded. The ONeNa ashes of the burning are not allowed to assume an arbitrary value of temperature; rather convection fixes the temperature behind the flame front. The second feature of the convective material is the lack of abundance gradients behind the flame. That is, convection uniformly mixes the ashes of the partial burning.

For our standard mixing parameter settings, models with ZAMS masses in the range 7−7.6 ${M}_{\odot }$ do not achieve a steady-state flame. Instead, they undergo a series of flashes, where each flash occurs closer to the core (Siess 2009; Denissenkov et al. 2013, and see Figure 3). The number of flashes increases as the ZAMS mass increases, until the ZAMS mass exceeds 7.6 ${M}_{\odot }$ when the first flash transitions into steady-state flame. For stars between 7.7 and 8.4 ${M}_{\odot }$ the off-center steady-state flame begins closer to the center. The dependence of this mass range on the composition mixing parameters is discussed in Section 5.

4.3. Does the Burning Reach the Center?

Whether off-center carbon burning, either as a steady-state flame or as a series of time-dependent flashes, reaches the center in these models depends on the ZAMS mass and the adopted mixing parameters. When convective mixing operates within the Schwarzschild boundaries, a flame will propagate all the way toward the center (Nomoto & Iben 1985; García-Berro et al. 1997; Saio & Nomoto 1998; Siess 2009; Doherty et al. 2010). Convective overshoot is not strictly required for carbon flames to be quenched away from the center. Additional mixing processes such as thermohaline mixing have been shown to effectively quench the carbon flame away from the center while convection is operating within only the strict Schwarzschild boundaries (Siess 2009; Denissenkov et al. 2013). Mixing processes that, by design, extend beyond the MLT convective boundary take unburned carbon fuel ahead of the flame front and mix this fresh fuel with the ashes of the convective region behind the flame front. This starves the flame of fuel, with the nuclear energy production rate proportional to the square of the carbon abundance (see Equation (13)). For instance, instead of a fresh carbon mass fraction of 0.3 the carbon mass fraction near the ignition point may be depleted to 0.1 and polluted with enhanced abundances of ¹⁶O, ²⁰Ne, ²³Na, and ²⁴Mg. Examples of such convective boundary mixing processes include thermohaline and overshoot. For large enough mixing parameters, the flame either disintegrates and sputters in a series of fuel-starved flashes moving toward the center, or is extinguished before reaching the center (Siess 2009; Denissenkov et al. 2013). Where the flame is extinguished, if it is extinguished, can have repercussions for the composition of the subsequent white dwarf that is formed, and from there, possible consequences for supernova Type Ia models.

Figure 10 shows the location where carbon burning is extinguished (red curve, left y-axis) for the non-rotating models with baseline mixing parameters. The distance (in mass coordinates) that the flame traveled from birth to death is shown by the blue curve and the right y-axis. For 7.0 ≤ M_ZAMS ≤ 8.2, carbon burning does not reach the center. As the ZAMS mass increases, the flame or flashes get closer to the center, eventually reaching the center at 8.2 ${M}_{\odot }$. In terms of the mass traversed, the flame (or flashes) increases its travel distance from 7.0 to 7.5 ${M}_{\odot }$ and then decreases for higher mass models. This transition occurs as stars with masses between 7.0 and 7.5 ${M}_{\odot }$ undergo one or more flashes, where each flash does not travel, but each subsequent flash ignites closer to the core. For higher mass models a steady-state flame reaches the center. The dependence of these results on the chosen mixing parameters is discussed in Section 5.

Zoom In Zoom Out Reset image size

For the models between 7.8 ${M}_{\odot }$ and 8.2 ${M}_{\odot }$, which do not reach the center, they transition from a steady-state flame into a series of flashes, as seen in Figure 3 for the 8 ${M}_{\odot }$ case. These flashes ignite in regions where the X(¹²C) abundance has dropped due to the mixing. A subsequent flash thus requires a higher critical density (see Equation (18)) to ignite in the presence of a lower abundance of carbon. These flashes are able to drive a local expansion of the core and the critical density to move inwards into the core, analogous to the flashes seen in the <7.8 ${M}_{\odot }$ models. The flashes eventually stop forming once the core can no longer reach the critical density, leaving unburned carbon in the core.

In this section we investigate the impact of rotation and convective overshoot on carbon burning in the SAGB models using the range of values in Table 2 and baseline mixing parameters listed in Table 1. As an example of the rotation characteristics in the carbon core, Figure 11 shows the angular momentum diffusion coefficients, specific angular momentum, angular frequency, and magnetic field profiles at first carbon ignition for an 8.0 ${M}_{\odot }$ ZAMS model with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.2. Rotation is initialized by imposing a solid body rotation law at ZAMS. For ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.2 this corresponds to an initial angular frequency of ω = 1.522 × 10⁻⁵ rad s⁻¹ and a total angular momentum of L = 1.90 × 10⁵¹ erg s. At first carbon ignition, Figure 11 (middle panel) shows the carbon core has spun up by a factor of ≈30 to ω ≈ 6 × 10⁻⁴ rad s⁻¹ and rotates as a solid body. The largest angular momentum diffusion coefficient, by several orders of magnitude, at first carbon ignition is due to the Spruit–Tayler dynamo. The implied ≈1 MG radial component of the magnetic field is shown in the bottom panel. DSI and SSI are not shown due to their negligible contributions in this model.

Zoom In Zoom Out Reset image size

Figure 12 illustrates the evolution of the specific angular momentum from ZAMS to first carbon ignition for the 8.0 ${M}_{\odot }$ ZAMS model initialized at ZAMS with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25 and f_ov = 0.016. During the MS phase (model numbers 1000−1400) the specific angular momentum is uniformly distributed throughout the model. As the star ascends the RGB and then the AGB (model numbers 1400−3200) the specific angular momentum is extracted from the core and redistributed into the envelope, decreasing the specific angular momentum in the core by a factor of ≈100. Boundaries between convective and non-convective regions are distinguished by sharp transitions in the specific angular momentum, with non-convective regions having the least specific angular momentum. The first ignition of carbon, off-center in this case, occurs around model 3700 and is marked.

Zoom In Zoom Out Reset image size

While Poelarends et al. (2008) claim that mass loss will strongly affect the final outcome of SAGB, we note that in our models (with only one mass loss rate used) mass loss seems to have a minimal effect on the star up to carbon ignition. In Figure 12 the effect of mass loss is visible as the white region at the top of the figure and only significant near the end of the core helium-burning phase. While mass loss will extract angular momentum from the star, the core's angular momentum is unaffected.

Rotation during the MS supplies a prolonged source of hydrogen fuel that builds a slightly more massive helium core than non-rotating models (e.g., Heger et al. 2000; Maeder & Meynet 2000; Lagarde et al. 2012). The increase in core mass (≈0.05−0.2 ${M}_{\odot }$), shifts the effective temperature and luminosity. However, once core helium is depleted the two tracks converge on the AGB and there is little difference in the CO core masses between rotating and non-rotating stars. Figure 13 shows a portion of the Hertzsprung–Russel diagram of two 8 ${M}_{\odot }$ ZAMS models from H depletion to He depletion: one with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.0 (blue curve) and one with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.2 (red curve).

Zoom In Zoom Out Reset image size

5.1. Evolution of Core Rotation

Figure 14 (left panel) shows the evolution of the central rotation rate, ω_center, for different initial rotation rates for the 7.5 ${M}_{\odot }$ ZAMS mass model. During the MS phase there is little evolution in the central rotation rate, spinning down by ≈20%. Figure 14 (left panel) shows the transition region from core hydrogen burning to shell hydrogen burning and then to core helium burning. As a star leaves the core hydrogen-burning phase (log₁₀ Age ∼ 7.596 for the ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.1 model), the nuclear energy in the core decreases. This causes a contraction of the core and subsequently a spin-up of the core (Palacios et al. 2006). The convective core spins at a constant rate throughout the convective region, and spins faster than the outer non-convective layers. As the nuclear energy generated decreases, the convective region inside the core shrinks toward the center (Sills & Pinsonneault 2000). As it shrinks, magnetic fields in the outer radiative layers propagate inwards into regions that were previously convective. These fields act to slow down the core, removing the rotation differential between the core and the envelope.

Zoom In Zoom Out Reset image size

Eventually the convective region recedes entirely, allowing the magnetic fields to propagate through to the center, which causes the rapid spin-down of ω_center. At this point the core is still contracting, thus the core begins to spin up again.

There are now two possible outcomes, seen in Figure 14 (left panel): either the core has a second peak in core rotation (${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.1, 0.3, 0.5) or the core rotation plateaus (${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.2, 0.4), the outcome of which depends on the sign of dω/dM. Stars with dω/dM < 0 form a second peak in the core rotation profile. When convection restarts in the core, due to helium burning, the convective region expands outwards. As it does so it engulfs slower rotating material, which slows down the convective center of the star. Hence the second peak in core rotation occurs when convection restarts. For stars with dω/dM > 0, a plateau is reached in the core rotation profile. Here, when the convection restarts and expands outwards it engulfs faster rotating material and thus spins up. However, this forces the convective core to expand, thus as more material is engulfed into the convective region, it will begin to slow the core down. Hence the maximum ω_center occurs after convection has started. Typically, we find that the convective core will grow to ≈0.05 ${M}_{\odot }$ before the spin-up ceases. While we believe these qualitative aspects hold, we caution against inferring quantitative predictions from this, as we have found that the sign of dω/dM is dependent on model resolution.

Figure 14 (right panel) shows the evolution of ω_center at the end of the core helium-burning phase up to the start of carbon burning. Note the change in the scale of ω_center; the centers have spun down by ≈25% during the core helium-burning phase. First, we see a glitch in the rotation rate at the end of the core helium phase. This is due to the same process that occurs at the end of core hydrogen burning. The convective core shrinks, allowing the magnetic fields to propagate inwards, slowing the core down. This is countered by the core contracting and spinning up. We can also see that the initially faster rotating stars are evolving more slowly due to their ability to mix fresh fuel into the core (Heger et al. 2000; Maeder & Meynet 2000).

We see that as the star forms its CO core the core spins up, due to core contraction. Carbon ignition occurs at the glitches seen at log₁₀ ω_center ≈ −3.12 rad s⁻¹, when all stars have the same core rotation rate. These glitches occur because of the core expansion due to the carbon-burning events. Subsequent episodes of carbon burning can be seen as smaller rotation glitches in Figure 14 (right panel insert). We can also see that the faster the star initially rotates, the later the ignition occurs, due to their slower evolution.

5.2. Mass–rotation Plane

Figure 15 (top left) shows the location of the first carbon ignition in the mass–rotation plane for a fixed f_ov = 0.016. Only models with 7 ${M}_{\odot }$ ≲ M_ZAMS ≲ 8 ${M}_{\odot }$ feature off-center ignition, as models with M_ZAMS ≲ 7 ${M}_{\odot }$ do not ignite carbon, while models with ${M}_{\mathrm{ZAMS}}\gtrsim 8\ {M}_{\odot }$ feature central carbon ignition. Figure 15 (top right) shows the location where carbon flames and flashes are quenched. Only in models having off-center ignition within the relatively narrow range 8.0 ${M}_{\odot }$ ≲ M_ZAMS ≲ 8.2 ${M}_{\odot }$ do the flames or flashes reach the center.

Zoom In Zoom Out Reset image size

For a fixed ZAMS mass the ignition and quenching locations are mostly independent of ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ values between 0.0 and 0.5. This occurs because the transport of angular momentum from the core to the overlying layers is efficient during the giant branch phases of evolution. Thus, regardless of the initial rotation rate, by the time the carbon core forms, the central regions are rotating as a solid body with similar angular frequencies for a fixed ZAMS mass. Figure 15 (bottom left) shows the center angular frequency, ${\omega }_{\mathrm{center}}$, in the mass–rotation plane at the fixed baseline overshoot value. Note that ω_center only spans a factor of ≈2 over the entire plane; all models rotate with similar angular frequencies at a fixed ZAMS mass. We can quantitatively explain, to first order, the rate at which the carbon core spins up between formation of the carbon core and first ignition of carbon from angular momentum conservation and the mass–radius relationship of polytropes. When the rotating carbon core forms, its total angular momentum is

$$\begin{eqnarray}&&{L}_{i}={I}_{i}{\omega }_{i}\sim {c}_{i}{M}_{i}{R}_{i}^{2}{\omega }_{i},\end{eqnarray} \tag{ 20 }$$

where I is the moment of inertia, M_i is the mass of the carbon core, R_i is the radius of the carbon core, and c_i is a constant that depends on the density structure. At first ignition of carbon, the angular momentum of the more massive contracting CO core is

$$\begin{eqnarray}&&{L}_{f}={I}_{f}{\omega }_{f}\sim {c}_{f}{M}_{f}{R}_{f}^{2}{\omega }_{f}.\end{eqnarray} \tag{ 21 }$$

Conserving angular momentum over this phase of evolution gives

$$\begin{eqnarray}&&\displaystyle \frac{{\omega }_{f}}{{\omega }_{i}}\sim \displaystyle \frac{{M}_{i}{R}_{i}^{2}}{{M}_{f}{R}_{f}^{2}}.\end{eqnarray} \tag{ 22 }$$

Assuming the non-rotating polytropic mass–radius relation

$$\begin{eqnarray}&&R\sim {M}^{(1-n)/(3-n)}\end{eqnarray} \tag{ 23 }$$

applies at the first order, substitution into Equation (22) gives

$$\begin{eqnarray}&&\displaystyle \frac{{\omega }_{f}}{{\omega }_{i}}\sim \displaystyle \frac{{M}_{i}^{(5-3{n}_{i})/(3-{n}_{i})}}{{M}_{f}^{(5-3{n}_{f})/(3-{n}_{f})}},\end{eqnarray} \tag{ 24 }$$

where n_i is the polytropic index at the formation of the carbon core and n_f is the polytropic index at first carbon ignition. For the angular frequency of the core to increase, ω_f > ω_i, the polytropic index is restricted to be in the range 1 $\lt \;n\;\lt$ 3. For example, for an 8 ${M}_{\odot }$ ZAMS model with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.2 at formation of the carbon core, we find M_i ≈ 0.92 ${M}_{\odot }$ and n_i ≈ 1.5. Using a least-squares fitting program to generate the polytropic index n for a sequence of MESA profiles between formation of the carbon core and first carbon ignition (see Figure 7 for an example), we find the left-hand side and the right-hand side of Equation (24) agree to within a factor of ≈2 for the 7 ${M}_{\odot }$, 8 ${M}_{\odot }$, and 9 ${M}_{\odot }$ models shown in Figure 15. The center spins up, on average, by a factor of ≈40 between formation of the carbon core and first carbon ignition. Since the carbon core at first ignition rotates as a solid body with similar angular frequencies for a fixed ZAMS mass, the carbon core masses are nearly in the same state independently of the initial rotation rate, as shown in Figure 15 (bottom right).

We find the carbon cores are rotating with periods between 0.1 and 1.0 days at first carbon ignition on the AGB. During the red giant branch (RGB) phase of evolution we find the helium cores have periods of ≈2.5 days, again independently of the initial rotation rate. Mosser et al. (2012) measured rotational splittings in a sample of about 300 red giants observed during more than two years with Kepler. They found these splittings are dominated by core rotation. Periods range between 10 and 100 days with larger periods for red clump stars than RGB stars. They inferred a ZAMS mass range of 1.2−1.5 ${M}_{\odot }$, less massive than our rotating SAGB models.

Stars with masses <7 ${M}_{\odot }$ will eventually form CO WDs after removing their outer envelopes. Between 7 ${M}_{\odot }$ and ≲ 8 ${M}_{\odot }$ where the carbon burning does not reach the core, our models suggest these stars will form CO+ONe hybrid WDs. Stars with masses ≳8 ${M}_{\odot }$ form ONe WDs as the carbon flames will burn away the ¹²C. ECSNe are expected for stars with masses >9 ${M}_{\odot }$, due to the CO core mass being greater than the Chandrasekhar mass (Eldridge & Tout 2004).

5.3. Overshoot–rotation Plane

Figure 16 (top left) shows the location of the first carbon ignition in the overshoot–rotation plane for a fixed ZAMS mass of 8 ${M}_{\odot }$. For this case, overshoot is a dominant factor in setting the location of the first, off-center, carbon ignition. This first ignition of carbon can be made to occur at almost any mass coordinate within the carbon core of the 8 ${M}_{\odot }$ model by varying the overshoot parameter. The f_ov = 0.0 case, where convective mixing operates only within the Schwarzschild boundaries and does not extend beyond the MLT convective boundary, does not ignite carbon for any of the initial rotations rates. The smallest non-zero overshoot parameter in our grid, f_ov = 0.001, gives mass locations for the first ignition of carbon that are furthest from the center, closest to the outer boundary of the carbon core. Progressively larger values of the overshoot parameter generally move the location of the off-center ignition closer to the center.

Zoom In Zoom Out Reset image size

Figure 16 (top right) shows the quenching location, where the flame and flashes die, in the overshoot–rotation plane for a fixed ZAMS mass of 8 ${M}_{\odot }$. The flame and flashes approaches the center for nearly all the models; only in models with ${f}_{\mathrm{ov}}\leqslant 0.02$ does the burning become quenched relatively far from the center (Denissenkov et al. 2013). Similar to our analysis of the mass–rotation plane, Figure 16 (bottom left) shows that the carbon core rotation rate is approximately constant, to within a factor of ≈2, regardless of the initial rotation rate. There is evidence for a weak dependence of the overshoot parameter on the rotation rate. Figure 16 (bottom right) shows that the carbon core mass increases with increasing values of f_ov, again nearly independently of the initial ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$. A larger core mass, in turn, leads to the first carbon ignition occurring deeper in the star.

Comparing the results of the overshoot–rotation grid with the mass–rotation grid, we find that the carbon core mass is the quantity that most strongly determines the structure of the flame. For example, the boundary between cases that ignite off-center and those that do not ignite carbon (ZAMS masses ≈7 ${M}_{\odot }$) depends on whether the star can form a carbon core of ≈1.05 ${M}_{\odot }$, which is necessary to reach the critical density in Equation (18).

5.4. Mass–Overshoot Plane

Figure 17 (top left) shows the mass location of carbon ignition in the mass–overshoot plane at a fixed ZAMS rotation of ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25. For f_ov = 0.0 the minimum mass needed to ignite carbon is 9 ${M}_{\odot }$, ≈2 ${M}_{\odot }$ greater than the baseline case, and models up to 11 ${M}_{\odot }$ ignite carbon off-center. While no overshoot may be unphysical, even a small amount of overshoot moves the minimum ZAMS mass for ignition considerably, down to 7.8 ${M}_{\odot }$. Increasing f_ov decreases the ZAMS mass required for off-center carbon ignition and decreases the minimum ZAMS mass needed for central ignition of carbon.

Zoom In Zoom Out Reset image size

The width of the ZAMS mass range where the model stars ignite carbon off-center is approximately constant with respect to overshoot, ΔM_ZAMS/Δf_ov ≈ 1.6 ${M}_{\odot }$. That is, overshoot uniformly moves the ZAMS mass boundaries where off-center carbon ignition occurs. For example, the sloped contours in Figure 17 (top left) show that when f_ov = 0.0, models in the mass range ≈8.9 to 11 ${M}_{\odot }$ have off-center ignition (as found by, e.g., Siess 2006, 2007). When f_ov = 0.008 this mass range shifts to ≈7.4–9.4 ${M}_{\odot }$, and when f_ov = 0.016 the off-center carbon ignition range shifts to ≈7.2–8.8 ${M}_{\odot }$.

Figure 17 (top right) shows the final fate of the carbon burning flames and flashes, the quenching location in the mass–overshoot plane. With f_ov = 0.0, all flames and flashes reach the core (e.g., bottom right plot of Figure 3). As the overshoot parameter increases, the carbon burning is less likely to reach the center. For ${f}_{\mathrm{ov}}\lt$ 0.01, the flame either has a single flash (similar to the top left panel of Figure 3 for the 7 ${M}_{\odot }$ model) or undergoes a single flash then a steady-state flame (similar to the middle right panel of Figure 3 for the 8.2 ${M}_{\odot }$ model). Only when the overshoot parameter is large, ${f}_{\mathrm{ov}}\gt 0.01$, is there an intermediate evolution of a flash and then a steady-state flame that does not reach the core (similar to the middle left panel of Figure 3 for the 8.0 ${M}_{\odot }$ model).

Figure 17 (bottom left) shows the central angular frequency at the point of first ignition in the mass–overshoot plane with a fixed ZAMS rotation of ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25. We find that models that undergo central carbon ignition have lower angular frequencies than the off-center igniting stars. This is due to the carbon cores being spun up during the cooling phase, owing to the contraction of the carbon core (see Equation (24)). As the heavier stars ignite carbon earlier they have less time in which to be spun up than the lower mass stars have at first ignition. As before, only models with carbon core masses greater than ≈1.05 ${M}_{\odot }$ ignite carbon, as shown in Figure 17 (bottom right). As before, the maximum carbon core mass that ignites off-center carbon burning is ≈1.2 ${M}_{\odot }$, similar to that shown in Figures 15 and 16.

Figure 18 shows the final CO core mass for the non-igniting models (white region in Figure 17). The maximum mass for a CO core is 1.05 ${M}_{\odot }$; stars with heavier CO cores ignite carbon burning. We also find a trend for increasing overshoot to increase the final CO core mass, as noted previously. Doherty et al. (2015) reported a grid of models with the Monash stellar evolution code, over a range of metallicities, to investigate the fate of AGB and SAGB stars. Comparing our results with their solar metallicity results (their Figure 6), we find that for a given ZAMS mass our rotating CO core masses are 0.05−0.1 ${M}_{\odot }$ larger. Increasing the overshoot can mimic decreasing the metallicity in terms of the final CO WD mass.

Zoom In Zoom Out Reset image size

Table 5 shows overshoot has the most significant effect on the location of the initial flame, with no ignition for no overshoot for the 8 ${M}_{\odot }$, ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25 model, as long as the scale factor for the strength of angular momentum diffusion ${\eta }_{\mathrm{am}}$ > 0. As overshoot increases, the flame ignition occurs deeper in the star. This is due to changes in the ⁴He and CO core masses during the star's evolution. For instance, comparing Figure 3 bottom left and bottom right, we can see the effect of overshoot for the 9 ${M}_{\odot }$ model. Primarily, the model with overshoot ignites at the center, while the no-overshoot model ignites off-center. Without overshoot the CO core mass, the size of the helium shell, and the ignition location are comparable to the 7 ${M}_{\odot }$ model with overshoot.

Table 5.  Ignition Locations in Solar Masses for the Mixing Grid

		${\eta }_{\mathrm{am}}$
		0.0				0.5				1.0				1.5
${\alpha }_{\mathrm{sc}}$	${\alpha }_{\mathrm{th}}$	fov				fov				fov				fov
		0.000	0.001	0.016	0.020	0.000	0.001	0.016	0.020	0.000	0.001	0.016	0.020	0.000	0.001	0.016	0.020
0.000	0.00	...a	...	0.65	0.28	...	0.49	0.15	0.18	...	0.43	0.17	0.16	...	0.49	0.17	0.13
0.000	0.10	0.50	0.71	...	...	...	0.46	0.17	0.05	0.43	0.42	0.14	0.05	...	0.47	0.18	0.05
0.000	1.00	...	0.63	0.35	0.38	0.39	0.46	0.11	0.18	...	0.40	0.16	0.16	...	0.48	0.15	0.16
0.000	10.0	0.50	0.65	...	0.54	...	0.49	0.16	0.15	...	0.43	0.18	0.05	...	0.43	0.14	0.16
0.001	0.00	0.81	...	0.65	...	...	0.59	0.21	0.17	...	0.83	0.16	0.14	...	0.85	0.23	0.04
0.001	0.10	0.79	0.48	0.28	0.37	...	0.83	0.20	0.05	...	0.85	0.14	0.18	...	0.69	0.12	0.05
0.001	1.00	0.83	...	0.58	0.71	...	0.84	0.15	0.17	...	0.83	0.18	0.05	...	0.71	0.17	0.15
0.001	10.0	0.33	0.69	0.32	0.49	...	0.57	0.18	0.16	...	0.62	0.08	0.05	...	0.65	0.11	0.16
0.010	0.00	0.19	0.57	0.45	0.75	...	0.63	0.21	0.16	...	0.65	0.17	0.14	...	0.63	0.16	0.05
0.010	0.10	...	...	0.30	0.39	...	0.69	0.20	0.05	...	0.62	0.13	0.04	...	0.68	0.12	0.05
0.010	1.00	0.73	...	...	0.73	...	0.66	0.17	0.18	...	0.55	0.20b	0.05	...	0.85	0.16	0.14
0.010	10.0	0.49	0.57	0.66	0.36	...	0.56	0.18	0.04	...	0.63	0.08	0.17	...	0.64	0.10	0.05
0.100	0.00	0.28	...	0.39	0.76	...	0.66	0.12	0.17	...	0.84	0.15	0.14	...	0.83	0.16	0.06
0.100	0.10	0.54	...	0.28	0.37	...	0.52	0.20	0.08	...	0.58	0.14	0.17	...	0.83	0.11	0.05
0.100	1.00	0.51	...	0.57	0.73	...	0.68	0.15	0.18	...	0.83	0.20	0.06	...	0.55	0.17	0.13
0.100	10.0	0.69	...	0.52	0.55	...	0.49	0.18	0.04	...	0.58	0.08	0.15	...	0.60	0.12	0.16

Notes.

^aEllipses represents models with no ignition. ^bBaseline model.

Download table as:  ASCIITypeset image

There are two distinct populations in the thermohaline models: those with small values of ${\alpha }_{\mathrm{th}}$, which ignite a flame at M ≈ 0.5 ${M}_{\odot }$, and those with large ${\alpha }_{\mathrm{th}}$ values, which ignite a flash at M ≈ 0.8 ${M}_{\odot }$ (though only if f_ov = 0.001). In those models that ignite at M ≈ 0.5 ${M}_{\odot }$, thermohaline mixing has little impact on the location of flame ignition. However, there is some variation due to thermohaline mixing before the flame ignites. Before the flame ignites, when _nuc ≫ _ν, there is a region in the vicinity of the ignition point that undergoes weak ¹²C + ¹²C burning (with ${\epsilon }_{\mathrm{nuc}}\ll {\epsilon }_{\nu }$). This weak burning is able to drive a region of weak thermohaline mixing as a precursor to the vigorous carbon burning. Those models that ignite at M ≈ 0.8 ignite a flash, but under different conditions to those that we predict for the other flashes. Here, there is thermohaline mixing between the CO and helium shell, which allows us to form a small region (M ≈ 0.05 ${M}_{\odot }$) where X(¹²C) ∼ X(¹⁶O). This higher fraction of carbon allows the ignition to occur at a lower density. This flash then prevents an ignition occurring deeper in the CO core where X(¹²C) ≈ 0.3, which we assume in Equation (14).

Thermohaline mixing also affects the flame once burning has commenced; the mixing pulls ¹²C material from below the flame (Siess 2009). As the strength of thermohaline mixing increases, the sub-flashes, seen in Figure 3, middle left panel, merge into one continuous flame, due to the increased carbon abundance. However, as Siess (2009) showed, this mixing eventually starves the flame of fuel, preventing it from reaching the core.

The effect of semiconvection is almost negligible, over the range of values considered here. Semiconvection acts near regions of convection, but it only acts for short periods of time in our models, and thus has limited ability to change the composition of models before the formation of the CO core. It can, however, act during the carbon burning once the convective region has formed, mixing the burnt material with unburnt CO. Again this effect is small and plays a limited role in the evolution of the flame. In Figure 3, where we have secondary carbon flashes (top left, middle left, and middle right panels), those flashes that occur near the ⁴He shell can form brief semiconvective regions across the shell and into the convective envelope. This may provide a way to detect the product of the flashes in the surface abundances, though these flashes (and the semiconvective regions) are short-lived, which will limit the material transferred.

Changes in ${\eta }_{\mathrm{am}}$, which is a global scale factor on the strength of angular momentum mixing, primarily act by changing the strength of thermohaline mixing. As ${\eta }_{\mathrm{am}}$ increases, the amount of thermohaline mixing increases as well, which allows the mixing of material from the core into the flame (Siess 2009) to increase, though again this effect is small.

Zaussinger & Spruit (2013) found for a 15 ${M}_{\odot }$ model on the MS that the semiconvection mixing timescale is long (10¹⁰ years), which explains why the semiconvection has little effect on these systems. Siess (2009) showed that thermohaline mixing has a limited effect in the evolution up to the carbon ignition, due to the lack of the ³He(³He, 2p)⁴He reaction, which is necessary to set up the mean molecular weight inversion needed for thermohaline mixing. Brown et al. (2013) propose a model for mixing by fingering convection in the parameter regime relevant for stellar (and planetary) interiors that is supported by 3D direct numerical simulations.

The angular momentum diffusion term has limited impact due to its implementation as a global scale factor on the angular momentum mixing process (Paxton et al. 2013), thus the value itself is not physically motivated. However, we have varied it to test whether missing sources of angular momentum, like internal gravity waves (Kumar & Quataert 1997), would have an impact. Given that in Figure 11 we have shown that the individual diffusion coefficients can vary by 10 orders of magnitude, however, a change of ≈50% in ${\eta }_{\mathrm{am}}$ is insignificant. For additional sources of angular momentum mixing to have an impact, they must be able to affect the size of the CO core, like overshoot does, to have a detectable impact. Compositional changes (thermohaline and semiconvection) are too weak to have an appreciable impact on the ignition location due to their limited ability to change the CO core mass.

Accurately capturing the nuclear burning and thermal transport within a steady-state, convectively bounded, carbon burning front, or within the time-dependent carbon burning flashes, requires spatial resolutions ≲2 km (Timmes et al. 1994; Ritossa et al. 1996; García-Berro et al. 1997). Siess (2006) use as many as ≈50 grid points to describe the precursor flame between the convection region and the minimum in the luminosity profile below the convective region. Denissenkov et al. (2013) and Chen et al. (2014) use more than 100 mass zones. Figure 19 shows the profiles of temperature and nuclear energy generation rate of a carbon-burning flame in our 8 ${M}_{\odot }$ ZAMS model with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25 and f_ov = 0.016. Distances ≲3700 km lie ahead of the flame front, and distances between ≈3700 and ≈3900 km contain the region where thermal conduction dominates nuclear burning. Distances between ≈3900 and ≈4200 km contain the body of the flame front, which reaches a peak temperature of ≈7.5 × 10⁸ K and peak energy generation rates of ≈8.9 × 10⁶ erg g⁻¹ s⁻¹. Distances ≳4200 km contain the convectively bounded region of the flame. The critical temperature at which the heating due to nuclear reactions equals the energy diffused away by neutrino and conductive processes in the steady state is about T_crit ∼ 5.5 × 10⁸ K. The location of this critical temperature is marked in Figure 19 with a black cross. The profiles shown in Figure 19 also capture the flame structure with 1−2 km resolution with ≈400 mesh points. The flame structure propagates inward toward the center at speeds of ≈0.1 cm s⁻¹, consistent with the values reported in Timmes et al. (1994).

Zoom In Zoom Out Reset image size

Figure 20 shows the time step and Kelvin–Helmholtz thermal timescale, ${\tau }_{\mathrm{kh}}={{GM}}_{c}^{2}/{R}_{c}L$, for the CO core of the 8 ${M}_{\odot }$ ZAMS model with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25 and f_ov = 0.016. As helium is depleted in the core at model number ≈2800, the time step begins to decrease from ≈10,000 years due to the increase in nuclear burning. At first carbon ignition, model ≈3700, the time step is ≈10 years and it decreases to ≈1 year as the flame and flashes propagate toward the center. At model number ≈4500 the flashes have reached their closest approach to the center. The thermal timescale increases as the core increases in mass until the first ignition, where it then rapidly decreases due to the increased luminosity. The thermal timescale then peaks again shortly before the next ignition at model ≈4400. This time, however, the flame generates less energy and the thermal timescale is reduced by a smaller amount, compared to the first ignition. On average the flame lifetime is ≈10% that of the thermal timescale of the core.

Zoom In Zoom Out Reset image size

The location of first carbon ignition in the 7 ${M}_{\odot }$, 8 ${M}_{\odot }$, and 9 ${M}_{\odot }$ models as a function of spatial and temporal resolution is shown in Figure 21. Each model has the baseline f_ov = 0.016. Spatial resolution in MESA is generally controlled by ${\delta }_{\mathrm{mesh}}$, with smaller values providing an increase in the number of cells. Temporal resolution is loosely controlled by ${w}_{t}$, the allowed change in the size of variables during a time step, with smaller values decreasing the size of the time steps taken. See Paxton et al. (2011) for a detailed discussion of these two MESA control parameters.

Zoom In Zoom Out Reset image size

For zero rotation, Figure 21 (top left) shows that all values of ${\delta }_{\mathrm{mesh}}$ and ${w}_{t}$ give the same location of first carbon ignition, suggesting convergence has been attained. Increasing the spatial resolution has little impact on the location of the flame, while increasing the temporal resolution shows a slight decrease in the ignition location.

For ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.25 (Figure 21 top right), the location of first carbon ignition depends on the values of ${\delta }_{\mathrm{mesh}}$ and ${w}_{t}$. For the 7 ${M}_{\odot }$ case the highest resolution model (red-dashed), ${w}_{t}$ = 10⁻⁵ and ${\delta }_{\mathrm{mesh}}$ = 0.1, agrees with our baseline model (green-dotted–dashed), ${w}_{t}$ = 10⁻⁴ and ${\delta }_{\mathrm{mesh}}$ = 0.5, that there is no ignition. At ${w}_{t}={10}^{-3}$, as the spatial resolution decreases, the ignition point is pushed deeper into the star. At 8 ${M}_{\odot }$, most of the models have converged around an ignition point of 0.1–0.2 ${M}_{\odot }$, except for the lowest resolution model, which has a center ignition. Models with ${w}_{t}$ = 10⁻⁵ show little variation as ${\delta }_{\mathrm{mesh}}$ changes, while as ${w}_{t}$ increases in size the ${\delta }_{\mathrm{mesh}}$ term becomes more significant. These studies suggest that our baseline values for ${\delta }_{\mathrm{mesh}}$ and ${w}_{t}$ for off-center ignition are well within the convergence envelope. All models converge on a 9 ${M}_{\odot }$ star having center ignition.

At ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ = 0.5 (Figure 21 bottom), for the 7 ${M}_{\odot }$ case our baseline parameters agree with the highest resolution model about the lack of ignition. However, for all other values of ${w}_{t}$ and ${\delta }_{\mathrm{mesh}}$ there is a considerable spread in ignition points. For the 8 ${M}_{\odot }$ case the results have clustered around 0.05–2.0 ${M}_{\odot }$, except for the ${w}_{t}={10}^{-3}$ case, where the results have a spread of 0.5 ${M}_{\odot }$. All models agree again that the 9 ${M}_{\odot }$ case has a central ignition.

Overall, our baseline models agree within ≈0.1 ${M}_{\odot }$ with the highest resolution models we ran, for the ignition point. As the rotation rate of the star increases we can see that the choice of resolution terms becomes more significant and that there is a larger spread in possible values. Changing the temporal resolution has the most effect on the initial location of the flame. The choice of spatial resolution becomes more significant only as the temporal resolution decreases. Thus our choice of baseline parameters appears to be a good compromise in terms of precision of results and computational effort: decreasing ${w}_{t}$ increases the computational time by a similar amount, while decreasing ${\delta }_{\mathrm{mesh}}$ increases the memory requirements for the model. However, they also show a necessary requirement for carbon flame modelers to look critically at their choice of model resolution (Timmes et al. 1994; Ritossa et al. 1996; Siess 2006; Doherty et al. 2010; Denissenkov et al. 2013; Chen et al. 2014).

We have investigated the detailed and global properties of carbon burning in SAGB stars with 2755 stellar evolution models. These models consumed 200,000 core-hours (roughly three days per model) and yielded over 2 TB of decimated data (a limited number of MESA profiles were stored). To our knowledge this represents the largest block of computer resources used for a MESA survey to date. We note that every model ran from the pre-MS to the end of carbon burning (if carbon ignition was achieved) without failure and without intervention.

With these models, the location of first ignition, whether off-center or central, the quenching location of the carbon burning flames and flashes, the angular frequency of the carbon core, and the carbon core mass have been surveyed as a function of the ZAMS mass, initial rotation rate, and the magnitude of various mixing parameters such as convective overshoot, semiconvection, thermohaline, and angular momentum transport. We now compare our results to previous efforts and discuss methods for calibrating the f_ov parameter within a given overshoot implementation.

Georgy et al. (2013) found that rotation of a 9 M_⊙ model can increase the lifetime spent on the MS compared to that of a non-rotating model. This increase in MS lifetime is caused by rotational mixing supplying a sustained amount of fresh hydrogen into the convective core. They include modifications to the stellar structure equations due to centrifugal acceleration described by Kippenhahn et al. (1970) and Endal & Sofia (1976), assuming the angular velocity is constant on isobars (Zahn 1992). Georgy et al. (2013) also adopt an instantaneous method of overshoot with d_over/H_p = 0.10 applied to the H- and He-burning boundaries. For a non-rotating 9 ${M}_{\odot }$ model at the end of core He burning they find a ratio of convective core to total mass of ${M}_{\mathrm{cc}}/{M}_{\mathrm{tot}}\approx 0.10$ while our corresponding rotating model yields a larger value of M_cc/M_tot ≈ 0.15. We find a more modest difference between the non-rotating and rotating 8 ${M}_{\odot }$ models of ${M}_{\mathrm{cc}}/{M}_{\mathrm{tot}}\approx 0.14$ and M_cc/M_tot ≈ 0.15, respectively, as shown in Figure 13.

While our MESA models use a similar implementation for rotation, our calculations differ from Georgy et al. (2013) in that we include the effects magnetic torques, which aid in significantly inhibiting the spin-up of the convective core of the star during its evolution. For example, we find for an 8 ${M}_{\odot }$ ZAMS model with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ ≈ 0.5, an angular velocity at the center of the core of ${\mathrm{log}}_{10}\ {\omega }_{\mathrm{center}}\approx -3.4\ \mathrm{rad}\;{{\rm{s}}}^{-1}$ at the start of carbon ignition, compared to a more rapid value of ${\mathrm{log}}_{10}\ {\omega }_{\mathrm{center}}\approx -0.9\ \mathrm{rad}\;{{\rm{s}}}^{-1}$ when internal magnetic fields are neglected. Magnetic field torques that inhibit spin-up of the stellar model result in less massive convective cores due to the less efficient rotational mixing. We also find that magnetic torques can account for the less drastic shift in luminosity on the HR diagram for an 8.0 ${M}_{\odot }$ ZAMS model with ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ ≈ 0.2 (see Figure 13), contrary to the larger differences in the HR diagram tracks of the rotating model shown in Georgy et al. (2012).

For their stars past 2DU, Doherty et al. (2015) found the CO core mass for a 7 ${M}_{\odot }$ ZAMS model to be ≈0.8 ${M}_{\odot }$, increasing to 1.375 ${M}_{\odot }$ for their 9.5 ${M}_{\odot }$ ZAMS models. In mild contrast, our models predict the CO core mass for a 7 ${M}_{\odot }$ ZAMS model to be 1.05 ${M}_{\odot }$ and the highest mass star to produce a Chandrasekhar core to be a 9 ${M}_{\odot }$ ZAMS, assuming an overshoot of f_ov = 0.016. These differences are likely due to the treatment of the convective boundaries, with Doherty et al. (2015) using a search for convective neutrality rather than a convective-decay prescription, leading to differences in the size of the ⁴He and CO core masses.

Siess (2006) found for models without overshoot that the ZAMS mass range that ignites carbon off-center is 9−11.3 ${M}_{\odot }$. This is comparable to our f_ov = 0.0, rotating models (Figure 17), which yield a value of M_up ≈ 8.8 M_⊙ and M_mas > 11 M_⊙.

For a 9.5 ${M}_{\odot }$ ZAMS mass with no overshoot or thermohaline mixing, Denissenkov et al. (2013) found an off-center ignition mass of 0.665 ${M}_{\odot }$ with the flame proceeding to the center, consistent with our results. With f_ov = 0.007 they found that carbon ignites off-center but the flames and flashes do not reach the center. In contrast, we find in this case that the model star undergoes a central ignition. We speculate that this difference is due to Denissenkov et al. (2013) only including overshoot once the CO core has formed, where we include it from the pre-MS onwards. Thus, the CO core in the f_ov = 0.007 model of Denissenkov et al. (2013) will be smaller than in our models, and hence the ignition will occur off-center.

While our models are not completely comparable to those of Jones et al. (2013), who use f_ov = 0.014 except at the base of a burning region where they use f_ov = 0.007, they find an 8.2 ${M}_{\odot }$ model ignites off-center while models with M_zams > 8.8 ${M}_{\odot }$ ignite centrally. We find for models with 0.007 ≤ f_ov ≤ 0.014 that an 8.2 ${M}_{\odot }$ will ignite off-center while only models with M_zams > 9.4 ${M}_{\odot }$ will always ignite centrally.

Arguably the biggest uncertainty in stellar models is the treatment of convection. The overshooting parameter in particular, regardless of how it is implemented within a specific numerical instrument, critically influences all outputs of stellar evolution (e.g., Maeder 1975, 1976). Figure 17 in particular demonstrates that the properties of carbon burning in SAGB models are no exception, especially the range of ZAMS masses that experience off-center ignition. Testing on a small number of models suggests that the most significant location for overshoot is in regions of He burning, followed by carbon burning. Regions with H burning or no burning show little difference in ignition location with respect to changes in overshoot. The effects of convective overshoot on the stellar models considered in this work are in agreement with previous work by Siess (2007), who showed that for ${f}_{\mathrm{ov}}=0.016$ applied at the edge of the convective boundary, M_up can transition from 8.90 ± 0.10 ${M}_{\odot }$ to 7.25 ± 0.25 M_⊙ for Z = Z_⊙. We find a similar transition where ${M}_{\mathrm{up}}\approx 8.8\;{M}_{\odot }$ for our rotating, f_ov = 0.0 model, which shifts to a value of M_up ≈ 7.2 M_⊙ for f_ov = 0.016 (see Figure 17). Gil-Pons et al. (2007) found similar results upon investigating a grid of zero-metallicity stars with f_ov = 0.12 using an instantaneous overshooting formalism (Herwig et al. 1997), contrary to the diffusive approach used in this work. They find a value of ${M}_{\mathrm{up}}\;\approx$ 6.0 ${M}_{\odot }$ and M_mas ≈ 7.8 ${M}_{\odot }$. The adoption of instantaneous overshooting, as well as $Z\ll {Z}_{\odot }$, is likely to contribute to the modest discrepancy in values of M_up and M_mas.

Traditionally the value of the overshooting parameter for a given overshooting model is calibrated by fitting isochrones against the width of the terminal-age MS in color–magnitude diagrams, or the surface abundances, of young and intermediate-age clusters (e.g., Maeder 1976; Maeder & Mermilliod 1981; Mermilliod & Maeder 1986; Schaller et al. 1992; Herwig 2000; VandenBerg et al. 2006; Kamath et al. 2012). Photometry and spectroscopy of binary systems offer another avenue for calibration of overshooting because these measurements can provide the radii, effective temperatures, and masses. In addition both components of the binary need to lie on the same isochrone and fit their respective evolutionary tracks (Schroder et al. 1997; Pols et al. 1997; Ribas et al. 2000; Claret 2007; Meng & Zhang 2014; Stancliffe et al. 2015). High-precision high-cadence space photometry from the CoRoT and Kepler missions opens up a newer method for calibration of overshooting and other mixing processes in stellar interiors (Neiner et al. 2012; Montalbán et al. 2013; Tkachenko et al. 2014; Guenther et al. 2014; Aerts 2015).

MESA implements the time-dependent treatment of convective overshoot mixing of Herwig (2000) with the traditional calibration method leading to f_ov = 0.016. It is unknown if this value of f_ov in this specific overshoot model applies to masses other than the ones used for calibration, if it is consistent with values derived from binary systems or asteroseismology, or if it applies to advanced burning stages of stellar evolution. However, we have shown that for a dense grid of SAGB models taken to the end of carbon burning, utilizing our adopted baseline parameters, values of M_up ∼ 7.0 ${M}_{\odot }$ and ${M}_{\mathrm{mas}}\;\sim$ 8.4 ${M}_{\odot }$ are nearly independent of initial ZAMS rotational values of ${({\rm{\Omega }}/{{\rm{\Omega }}}_{\mathrm{crit}})}_{i}$ ∼ 0.0–0.5. While our SAGB models have been evolved from the pre-MS phase through the end of carbon burning, for models that do not ignite carbon, and those that ignite carbon off-center, the initial rotational rate may play a larger role in the final rotational rate of the WD that will eventually be born. For a given ZAMS mass and overshoot parameterization, we suggest that strong claims of quenching of carbon burning at an appreciable distance from the center to yield hybrid CO + ONeNa white dwarfs should be viewed with caution.

The authors thank Pavel Denissenkov for sharing his MESA inlists and Lars Bildsten for detailed discussions. We also thank the SPIDER collaboration for insightful discussions: K. Augustson, M. Browning, M. Cantiello, J. Fuller, R. Orvedahl, B. Paxton, J. Toomre, R. Townsend, and E. Zweibel. Finally, we thank the participants of the 2014 MESA Summer School for experimenting with some of the SAGB models: L. Arcavi, W. Ball, E. Bauer, P. Beck, J. Blumenkopf, J. Brown, T. Ceillier, D. Clausen, R. Connolly, J. Goldstein, A. Lauer, E. Leiner, J. McKeever, B. Mulligan, A. Nagy, A. Ordasi, J. Ostrowski, M. Renzo, V. Schmid, W. Strickland, T. Sukhbold, M. Sun, S. Triana, S. Valenti, T. Van Reeth, J. Vos, M. Vuckovic, D. Wilcox, and M. Windju. This project was supported by NASA under TCAN grant NNX14AB53G, by NSF under SI² grant 1339600, and by NSF under PHY 08-022648 for the Physics Frontier Center "Joint Institute for Nuclear Astrophysics—Center for the Evolution of the Elements" (JINA-CEE). C. E. F. acknowledges partial support from Arizona State University under the 2014 CLAS Undergraduate Summer Enrichment Award. The computing resources for the grid of MESA models were provided by the ASU Advanced Computing Center.