ON VARIATIONS OF PRE-SUPERNOVA MODEL PROPERTIES

MESA provides several controls to specify the mass resolution of a model. Sufficient mass resolution is required to accurately determine gradients of stellar structure quantities, but an excessive number of cells impacts performance. One parameter for changing the mass resolution in regions of rapid change is mesh_delta_coeff (${\delta }_{\mathrm{mesh}}$), which acts a global scale factor limiting the change in stellar structure quantities between two adjacent cells. Lower values of ${\delta }_{\mathrm{mesh}}$ increase the number of cells. Another parameter controlling mass resolution is the maximum fraction of a star's mass in a cell, max_dq. That is, the minimum number of cells in a stellar model is 1/max_dq. We use ${\delta }_{\mathrm{mesh}}\,=1.0$ and ${\mathtt{\max }}\_{\mathtt{dq}}={\rm{\Delta }}{M}_{\max }/{M}_{\star }(\tau )$, where ${M}_{\star }(\tau )$ is the mass of the stellar model in solar mass units at time τ and ${\rm{\Delta }}{M}_{\max }$ is a parameter we vary between 0.1 ${M}_{\odot }$ and 0.005 ${M}_{\odot }$. We choose to vary max_dq instead of mesh_delta_coeff to enable us to set a minimum level of mass resolution in the model.

MESA also offers a rich set of timestep controls. The parameter ${w}_{t}$ broadly controls the temporal resolution by modulating the magnitude of the allowed changes between individual timesteps. At a finer level of granularity, dX_nuc_drop_limit limits the maximum allowed change of the mass fractions between timesteps for mass fractions larger than dX_nuc_drop_min_X_limit. We use ${w}_{t}$ = 5 × 10⁻⁵, dX_nuc_drop_limit = 10⁻³, and dX_nuc_drop_min_X_limit = 10⁻³ for the evolution between the pre-MS to the onset of core Si-burning, where we loosen the criteria to allow larger timesteps; ${w}_{t}$ = 5 × 10⁻⁵, dX_nuc_drop_limit = 5 × 10⁻², and dX_nuc_drop_min_X_limit = 5 × 10⁻².

The timestep control delta_lg_XH_cntr_min regulates the time step as hydrogen is depleted in the core, which aids in resolving the transition from the ZAMS to the terminal-age main sequence (TAMS). Similarly, the timestep controls delta_lg_XHe_cntr_min, delta_lg_XC_cntr_min, delta_lg_XNe_cntr_min, delta_lg_XO_cntr_min, and delta_lg_XSi_cntr_min control the timestep as one of the major fuels is depleted in the core. These timestep controls are useful for obtaining, for example, convergence of mass shell locations, smoother transitions as a stellar model exits core H-burning in the HR diagram (i.e., the "Henyey Hook," Kippenhahn et al. 2012), and smoother trajectories in the central temperature ${T}_{{\rm{c}}}$-central density ${\rho }_{{\rm{c}}}$ plane. We use a mass fraction value of 10⁻⁶ for all these fuel depletion timestep controls. Additionally, we use MESA's default timestep controls for controlling changes in the hydrodynamics. At the point where hydrodynamics becomes important, during Si-burning, we find we are limited by the delta_lg_XSi_cntr_min control rather than by changes in the hydrodynamics.

MESA evolves models of massive stars from the pre-MS to the onset of core collapse with the nuclear burning fully coupled to the hydrodynamics using a single, possibly large, reaction network (see Paxton et al. 2015). This capability avoids the challenges of (a) operator splitting errors from evolving the hydrodynamics and nuclear burning independently; (b) stitching together different solution methods for different phases of evolution, for example, combining a reaction network with equilibrium methodologies such as quasi-static equilibrium (QSE) or nuclear statistical equilibrium (NSE), or with the use of an adaptive network that modifies itself based on the most populous isotopes present; or finally, (c) evolving a stellar model with a small reaction network while carrying along a larger reaction network that does not impact the energy generation rate or composition of the stellar structure, i.e., passive coprocessing. MESA's unified approach is not just a solution to issues of accuracy or self-consistency. It offers an improvement by providing a single-solution methodology—in situ reaction networks evolved simultaneously with the hydrodynamics.

We evolve each stellar model with one of the five nuclear reaction networks shown in Figure 1. We consider a small network, approx21_cr60_plus_co56.net (hereafter approx22.net), where each reaction pathway has been predetermined (i.e., a "hardwired" network). Such approximate networks, an α-chain backbone with aspects of H-burning, heavy ion reactions, and iron-group photodisintegration, are a traditional workhorse in massive star models (e.g., Weaver et al. 1978; Woosley & Weaver 1988; Heger et al. 2000; Heger & Woosley 2010; Chatzopoulos & Wheeler 2012). We also use four "softwired" networks, where after specifying the isotopes, all allowed reaction pathways between the isotopes are linked. Specifically, we consider mesa_79.net, which contains isotopes up to ⁶⁰Zn (black dots); mesa_127.net, which adds neutron-rich isotopes in the iron group (purple crosses); mesa_160.net, which adds more neutron-rich isotopes and a few proton-rich isotopes (green squares); and mesa_204.net, which adds isotopes lighter than ³⁶Cl (black squares), and includes the isotopes identified in Heger et al. (2001) as important for the electron fraction, ${Y}_{{\rm{e}}}$, in core-collapse models.

Zoom In Zoom Out Reset image size

The softwired reaction networks are designed to yield approximately the same final ${Y}_{{\rm{e}}}$ as a 3298 isotope reaction network in one-zone burn calculations of the thermodynamic history of 15 ${M}_{\odot }$ cores (Paxton et al. 2015). For example, in one-zone burn calculations the mesa_204.net reaction network gives a final ${Y}_{{\rm{e}}}=0.4032$, while the 3298 isotope reaction network gives ${Y}_{{\rm{e}}}=0.4039$.

All forward thermonuclear reaction rates are from the JINA reaclib version V2.0 2013-04-02 (Cyburt et al. 2010). Inverse rates are calculated directly from the forward rates (those with positive Q-value) using detailed balance, rather than using fitted rates. The nuclear partition functions used to calculate the inverse rates are from Rauscher & Thielemann (2000). Electron screening factors for both weak and strong thermonuclear reactions are from Alastuey & Jancovici (1978) with plasma parameters from Itoh et al. (1979). All the weak rates are based (in order of precedence) on the tabulations of Langanke & Martínez-Pinedo (2000), Oda et al. (1994), and Fuller et al. (1985). Thermal neutrino energy losses are from Itoh et al. (1996).

We treat convection using the mlt++ approximation (Paxton et al. 2013). This prescription takes the convective strength as computed by mixing length theory (MLT) and reduces the superadiabaticity in radiation dominated convection zones (e.g., near the iron opacity peak). This decreases the temperature gradient in these regions, and the artificial suppression can be viewed as treating additional, unmodeled, energy transport mechanisms in these regions (e.g., Jiang et al. 2015). Values of $1.6\lesssim {\alpha }_{\mathrm{MLT}}\lesssim 2.2$ have been inferred from comparing observations with stellar evolution models (Noels et al. 1991; Miglio & Montalbán 2005; Aerts et al. 2010), and 3D hydrodynamic simulations of the deep core and surface layers also suggest a similar range (Viallet et al. 2013; Trampedach et al. 2014; Jiang et al. 2015). These efforts inform our baseline choice of ${\alpha }_{\mathrm{MLT}}$ = 1.5, although this is by no means the only choice (Dessart et al. 2013).

Convective overshooting is assumed to be an exponential decay beyond the Schwarzschild boundary of convection (Herwig et al. 1997). The convective diffusion coefficient, ${D}_{\mathrm{conv},0}$, is measured at a distance ${f}_{\mathrm{ov},{\rm{D}}}$ inside the convection zone. Overshooting then extends beyond the edge of the convection zone a fraction ${f}_{\mathrm{ov}}$ of the local pressure scale height ${\lambda }_{{\rm{P}},0}$,

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

where ${D}_{\mathrm{OV}}$ is the overshooting diffusion coefficient at a radial distance z from the edge of the convection zone boundary. Our baseline choices for the overshoot mixing are ${f}_{\mathrm{ov}}=0.004$ and ${f}_{\mathrm{ov},{\rm{D}}}=0.001$. This choice is motivated by recent calibration of a 1D 25 ${M}_{\odot }$ pre-SN model to idealized 3D hydrodynamic simulations of turbulent O-burning shell convection (Jones et al. 2016). MESA offers the flexibility to have different values for convection driven by different types of nuclear burning (H, He, C, and others) and whether the burning occurs near the core or in a shell, but for simplicity, we take these values to be the same for all regions.

Thermohaline mixing is included in the models when ${{\rm{\nabla }}}_{{\rm{T}}}-{{\rm{\nabla }}}_{\mathrm{ad}}\leqslant B\leqslant 0$, where B is the Brünt composition gradient. The thermohaline mixing diffusion coefficient ${D}_{\mathrm{th}}$ is parameterized as

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

where ${\alpha }_{\mathrm{th}}$ is a dimensionless parameter, K is the radiative conductivity, and ${C}_{{\rm{p}}}$ is the specific heat at constant pressure (Ulrich 1972; Kippenhahn et al. 1980). Continuing attempts to calibrate ${D}_{\mathrm{th}}$ include multidimensional simulations of fingering convection under stellar conditions (Traxler et al. 2011; Brown et al. 2013; Garaud et al. 2015). We set ${\alpha }_{\mathrm{th}}$ = 2.0 except during core silicon burning, where we set ${\alpha }_{\mathrm{th}}$ = 0.0.

Finally, we include semiconvection when ${{\rm{\nabla }}}_{\mathrm{ad}}\lt {{\rm{\nabla }}}_{{\rm{T}}}\lt {{\rm{\nabla }}}_{{\rm{L}}}$, where ${{\rm{\nabla }}}_{{\rm{L}}}={{\rm{\nabla }}}_{\mathrm{ad}}+B$. We take the strength of semiconvective mixing ${D}_{\mathrm{semi}}$ as Langer et al. (1983) and Langer et al. (1985),

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

where ${\alpha }_{\mathrm{sc}}$ is a dimensionless parameter. Ongoing efforts to calibrate 1D semiconvection models include multidimensional simulations of double diffusive convection (Spruit 2013; Zaussinger & Spruit 2013) and comparing massive star models with observations (Yoon et al. 2006). While recent work suggests that to a first-order approximation, the effect of semiconvection can be neglected (Moll et al. 2016), we adopt a baseline choice of ${\alpha }_{\mathrm{sc}}\,=0.01$.

Figure 2 shows the MS and giant branch (GB) evolution of a 15 ${M}_{\odot }$ model at two resolutions ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ and ${\rm{\Delta }}{M}_{\max }=0.005$ ${M}_{\odot }$, evolved without mass-loss and with the approx22.net network. As stars evolve off of the pre-MS and onto the ZAMS, they begin burning hydrogen in a convective core, predominantly via the CNO cycle, but also via primordial ³He in the pp chain (Hansen et al. 2004; Iben 2013a, 2013b). This convective region expands outwards to consume approximately half the mass of the star. Once the primordial ³He is exhausted, at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{CHe}}-\tau )/\mathrm{yr})\approx 7.04$ in Figure 2, the convection region recedes and the core temperature and density increase. At this point in the evolution, most of the primordial ¹²C is already piled up in ¹⁴N as a result of the proton capture onto ¹⁴N reaction rate (Arnett 1996; Iliadis 2007). Owing to the strong temperature dependence of the CNO cycle, the energy production increases and the convective region expands slightly farther out in mass coordinate than the first peak at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{CHe}}-\tau )/\mathrm{yr})\approx 7.039$. The increased energy generation from the CNO cycle causes the star itself also to begin expanding at this point.

Zoom In Zoom Out Reset image size

The convective core shrinks during core H-burning. This recession leaves behind a composition gradient (μ-gradient, ${{\rm{\nabla }}}_{\mu }$) and lower levels of nuclear burning. Our models form a layered structure of semiconvective and convective regions, labeled an intermediate convection zone (ICZ) in Figure 2, at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{CHe}}-\tau )/\mathrm{yr})\approx 5.75$ at the approximate mass coordinate of the maximal extent of the core convection zone (≈6 ${M}_{\odot }$) (Langer et al. 1985; Heger et al. 2000; Hirschi et al. 2004). As the mass resolution or the initial mass of the stellar model increases, the structure and behavior the ICZ changes. For the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models, in Figure 2, when the mass resolution is increased from ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ to ${\rm{\Delta }}{M}_{\max }=0.005$ ${M}_{\odot }$, the fraction (in mass coordinates) of semiconvection in the ICZ decreases while the amount of convection increases.

Figure 2 shows that as the mass resolution increases the fine structure of the H-core convection zone, the maximal extent of the H-core extent and the edge of the overshoot region are better determined. All of which determine how much fresh H fuel is pulled into the core during the MS and how much fuel is available just outside the convective core to form the ICZ. At higher resolution less hydrogen is burned when the ICZ forms, which makes the star more compact as a consequence of the decreased burning. As less hydrogen has been burned, ${{\rm{\nabla }}}_{\mu }$ is smaller.

Semiconvection occurs when ${{\rm{\nabla }}}_{\mathrm{ad}}\lt {{\rm{\nabla }}}_{{\rm{T}}}\lt {{\rm{\nabla }}}_{{\rm{L}}}$ where ${{\rm{\nabla }}}_{{\rm{L}}}$ is a function of ${{\rm{\nabla }}}_{\mu }$. When ${{\rm{\nabla }}}_{\mu }$ is large enough, it can maintain a semiconvective region over the entire ICZ. If ${{\rm{\nabla }}}_{\mu }$ is small in a local region, then the semiconvective region may transition into a convective region (Langer et al. 1985). This occurs in Figure 2 in the higher resolution ${\rm{\Delta }}{M}_{\max }=0.005$ ${M}_{\odot }$ panel where ${{\rm{\nabla }}}_{\mu }$ is smaller, allowing more convection to develop inside the semiconvective region (Mowlavi & Forestini 1994). This is most clearly seen as the bifurcation of the ICZ in the lower resolution ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ (left panel), while in the right panel the convection zone is embedded inside the semiconvective region. Higher mass resolutions also capture the substructure in the ${{\rm{\nabla }}}_{{\rm{L}}}$, which become the seeds for the convective regions to form.

Langer et al. (1985) showed similar variations in the ICZ structure by varying the strength of the ${\alpha }_{\mathrm{sc}}$ parameter. As the efficiency of semiconvection increases, the fine structure increases. Our MESA models achieve a similar result by varying the mass resolution for a fixed ${\alpha }_{\mathrm{sc}}$.

A secondary convective zone (SCZ) then forms at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{CHe}}-\tau )/\mathrm{yr})\approx 4.5$ in both models, when H-shell burning begins. Figure 2 shows how the structure and behavior of the SCZ change with mass resolution. The low-resolution model (left panel) has layered semiconvection/convection that stays outside of the H-shell. In the higher resolution model (right panel), the SCZ is dominated by convection that propagates into the H-shell burning region. This merger of the SCZ and the H-shell leads to fresh fuel being injected into the H-shell and allows the star to live longer. If the star has mass-loss, then the longer H-burning lifetimes create larger He-cores that will last longer and thus increase the amount of mass lost (Section 3.2).

This difference in behavior arises because the SCZ forms within the mass coordinates previously occupied by the ICZ semiconvective/convective mixing region. The larger convectively mixed region during the ICZ in the higher resolution model implies that the μ-gradient left behind by the receding H-core is reduced or even zero. In turn, this reduction in the μ-gradient allows convection to form over a larger mass range in the SCZ. That is, a preceding semiconvective phase is not necessary to reduce the μ-gradient. This is shown in Figure 2 by the relative difference in mass covered by semiconvection and convection at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{CHe}}-\tau )/\mathrm{yr})\approx 4.5$.

The ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models without mass-loss are qualitatively similar, and Figure 3 shows a representative set of MESA models. When the ICZ forms, multiple layers of alternating semiconvection and convection form (Langer et al. 1985). At lower mass resolutions (left panel), this structure begins at $M\approx 17$ ${M}_{\odot }$ and propagates inwards to $M\approx 13$ ${M}_{\odot }$. As the mass resolution increases (right panel), the convective fingers start at a similar mass coordinate, but will propagate deeper into the star. Figure 3 shows an example of the fingers moving sufficiently inward to penetrate the H-burning core. There is a brief increase in the nuclear burning, as shown by the near step function in burning at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{CHe}}-\tau )/\mathrm{yr})\approx 5.5$ that arises because the ICZ injects unburnt hydrogen into the core. This leads to a longer MS lifetime for these models. This penetration of the core is a binary process, either the ICZ penetrates the core or it does not. Thus, one may expect step functions in the MS lifetime (Section 3.2).

Zoom In Zoom Out Reset image size

For the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models with mass-loss, Figure 4 shows a third possible state. As the mass resolution increases the SCZ, which formed at the TAMS, it decreases in extent. The left panel of Figure 4 shows that for a relatively coarse mass resolution of ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$, the SCZ forms multiple convection/semiconvection zones out to $M\approx 12$ ${M}_{\odot }$, while the right panel shows that for a finer mass resolution of ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$, the SCZ forms a single convection zone out to $M\approx 10$ ${M}_{\odot }$. The increased extent of the SCZ region allows more hydrogen fuel to be mixed inwards toward the core, and increases the hydrogen mass fraction by a few percent at the H-burning shell. This allows the H-shell to burn more energetically and for longer. Thus we would expect the models with lower spatial resolution to evolve slower to core collapse, but this effect is limited because of the small change in hydrogen mass fraction.

Zoom In Zoom Out Reset image size

Figure 5 shows the lifetime for models without mass-loss (top panel) and for models with mass-loss (bottom panel). The most striking feature is the sharp bifurcation of the lifetimes as a function of mass resolution, most prominently seen in the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ and 30 ${M}_{\odot }$ models without mass-loss at a spatial resolution of ${\rm{\Delta }}{M}_{\max }=0.05$ ${M}_{\odot }$. This is a direct consequence of the differences in behavior of the ICZ and SCZ during the MS and GB phases. Models where the ICZ penetrates the MS H-core, or where the SCZ penetrates the H-shell as in the 15 ${M}_{\odot }$ models, live up to ≈5% longer because fresh hydrogen fuel is injected. Hatching indicates models that do not reach core collapse, although all reach at least Si-burning. Thus we include them in our comparisons up to Si-burning.

Zoom In Zoom Out Reset image size

Stellar models with and without mass-loss have similar ages, with a trend for stars with mass-loss to have shorter lifetimes (e.g., El Eid et al. 2004). The spread in age is larger for models without mass-loss than for those with mass-loss. This is chiefly due to differences in the size of the timesteps, thus this is a numerical artifact. While all models have the same temporal control ${w}_{t}$, models with mass-loss require a smaller timestep to keep the variation per step within the limits specified by ${w}_{t}$. Smaller timesteps allow the mass lost to be treated more accurately, as the amount of mass-loss per step depends on the stellar parameters at the start of the step, thus smaller timesteps will better capture fast changes in the star's structure.

Compared to the pre-SN lifetimes from Paxton et al. (2011), Limongi et al. (2000), Woosley et al. (2002), and Hirschi et al. (2004), we find in general that the total stellar lifetimes considered agree to within $\approx \,\pm 5 \%$. The models of Woosley et al. (2002) tend to predict longer lifetimes. Specifically, for the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ model, they predict an $\approx 3 \%$ longer lifetime until core collapse. This value increases to $\approx 7 \%$ for the ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ model. Limongi et al. (2000) predict shorter lifetimes with $\approx 7 \%$ difference for their ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ model. The models of Paxton et al. (2011) and Hirschi et al. (2004) agree with our values to within $\lesssim 1 \%$ for ${M}_{\mathrm{ZAMS}}=15$ and 20 ${M}_{\odot }$ models, and spread to ±2% for ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$. Each set of pre-SN models considered here for comparison uses a slightly differing set of input physics. For example, our models use a more efficient value for convective overshoot than that in the Paxton et al. (2011) models and considers larger nuclear reaction networks. In addition, each set of pre-SN models considered uses different (and often unspecified) mass and temporal resolutions. We find considerable agreement among the total stellar lifetimes, with the largest discrepancies occurring at ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$.

Models with the same mass resolution show little variation in the final age with respect to changes in the number of isotopes in the reaction network. For example, the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models with and without mass-loss in Figure 5 show uniformity in the stellar age as a function of network size on either side of the bifurcation point. This is due to efficient mixing of the composition by convection; it does not make much difference where or when the SCZ penetrates the H-shell, for once it does, the H-shell is injected with similar amounts of fresh fuel (see Figure 2). As the initial mass of the stellar model increases, the magnitude of the spread in the stellar lifetime decreases. This is due to smaller variations in when and where the ICZ penetrates the H-core.

There are outlier models for the smallest network, approx22.net: the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ model without mass-loss and ${\rm{\Delta }}{M}_{\max }=0.02$ ${M}_{\odot }$, the ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ model without mass-loss and ${\rm{\Delta }}{M}_{\max }=0.005$ ${M}_{\odot }$, and the ${M}_{\mathrm{ZAMS}}\ =20$ ${M}_{\odot }$ model with mass-loss and ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$. These are caused by the SCZ penetrating deep into the star, from the edge of the H-burning shell to the radiative region at the center of the star. This injects enough unburnt hydrogen to prolong the lifetimes.

The ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models without mass-loss in Figure 5 show a lifetime bifurcation at a finer spatial resolution of ${\rm{\Delta }}{M}_{\max }=0.005$ ${M}_{\odot }$ than the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$. In addition, the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models with mass-loss show the bifurcation is inverted—increasing resolution decreases the stars lifetime—as a result of the relative extent of the SCZ region discussed for Figure 4.

For the ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ and 30 ${M}_{\odot }$ models, Figure 5 shows that the lifetime bifurcation point varies depending on whether or not there is mass-loss. The stars without mass-loss bifurcate at ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$ and ${\rm{\Delta }}{M}_{\max }=0.02$ ${M}_{\odot }$, respectively, while the mass-losing models bifurcate at ${\rm{\Delta }}{M}_{\max }=0.005$ ${M}_{\odot }$ and ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$, effectively one level of resolution higher. This variation may be due to our grid of mass resolutions being insufficient to resolve the bifurcation.

We encourage careful consideration of the mass resolution when modeling massive stars with MESA. Not only to the choice MESA's ${\delta }_{\mathrm{mesh}}$ parameter, which only increases the number of zones if MESA detects a large enough spatial gradient, but also the number of zones used throughout the model. Much of the variation seen in the MS and GB phases is set by the stellar structure at the start of the MS, which is early enough that ${\delta }_{\mathrm{mesh}}$ may not have sufficient time to act. From Figures 5 and for our choice of ${\alpha }_{\mathrm{sc}}$, we recommend a mass resolution of at least ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$, which equates to at least $\approx 2000$ mass cells during the pre-MS and MS evolution in order to reasonably resolve the ICZ and SCZ features. In addition, we recommend temporal controls that limit the timestep as major fuels are depleted in the core (see Section 2.1).

Figure 6 shows the He-core mass at the onset of convective core He-burning as a function of mass resolution ${\rm{\Delta }}{M}_{\max }$ and the number of isotopes in the nuclear reaction network. We define the He-core mass as the location where $X{(}^{4}\mathrm{He})\gt 0.1$ and $X{(}^{1}{\rm{H}})\lt 0.01$ and measured when the star reaches the base of the GB. They share many of the same features—bifurcations and outliers—discussed for the stellar lifetimes in Figure 5. For the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models without mass-loss the initial He-core mass ranges from ≈2.79 to 2.81 ${M}_{\odot }$, while models with mass-loss span ≈2.72–2.77 ${M}_{\odot }$. For the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models without mass-loss the initial He-core mass range between ≈4.6 and 4.7 ${M}_{\odot }$, while models with mass-loss span between ≈4.5 and 4.6 ${M}_{\odot }$. Our ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ models without mass-loss show initial He-core mass ranges of ≈6.8 to 7.2 ${M}_{\odot }$, while models with mass-loss span ≈6.54–6.66 ${M}_{\odot }$. The ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models without mass-loss show the initial He-core mass ranges between ≈9.2 and 9.8 ${M}_{\odot }$, while models with mass-loss span ≈8.63–8.78 ${M}_{\odot }$. In general, a longer core H-burning phase produces a higher He-core mass. The initial He-core masses in Figure 6 and their range as a function of mass resolution and number of isotopes in the reaction network are commensurate with the He-core masses found in numerous studies (e.g., Nomoto & Hashimoto 1988; Langer 1991; Wellstein & Langer 1999; Woosley et al. 2002; Petermann et al. 2015; Choi et al. 2016).

Zoom In Zoom Out Reset image size

Figure 7 shows the evolution of the central helium mass fraction, ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$, during the core He-burning phase of the 20 ${M}_{\odot }$ models without mass-loss. These central abundances, due to convective mixing, also represent the mass fraction values over the entire He-burning core. The core of these stars begin with ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})\approx 1$ and ${\rm{X}}{(}^{16}{{\rm{O}}}_{{\rm{c}}})\approx 0$ and, after ≈10⁶ yr, end with ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})\approx 0$ and ${\rm{X}}{(}^{16}{{\rm{O}}}_{{\rm{c}}})\approx 0.8$, with the remainder being mostly ¹²C.

Zoom In Zoom Out Reset image size

The left panel of Figure 7 shows the evolution of ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$ as a function of the number of isotopes in the nuclear reaction network, at a mass resolution of ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models without mass-loss. The slopes of the lines show that the models are evolving at similar rates. The positive offset in the mesa_204.net is caused by the core being slightly cooler at the start of core helium convection; for the approx22.net the temperature is 177 × 10⁶ K and for mesa_204.net it is 175 × 10⁶ K. Nevertheless, the mesa_204.net model reaches ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})\approx 0$ at the same time as the other models. The spikes in ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$ for the small approx22.net, expanded in the inset plot, are due to core breathing pulses (henceforth CBP, see Castellani et al. 1985). These occur when the ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$ value is low, such that if the core convection expands slightly outwards in mass, then a small entrainment of unburnt ⁴He leads to a large increase in the nuclear energy production (Straniero et al. 2003). For the approx22.net models, the edge of the convection core can move inwards and outwards ≈0.5 ${M}_{\odot }$ on timescales of ≈200 yr. The models with larger nuclear networks have significantly smaller CBPs as their convection cores move inwards and outwards ≈0.05 ${M}_{\odot }$, entraining a significantly smaller amount of additional ⁴He fuel.

The left panel of Figure 7 shows that the approx22.net almost doubles ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$ during a CBP. The nuclear energy generated from a CBP is sufficient to expand the stellar envelope, which can be seen as a blue loop in theoretical HR diagrams (e.g., Constantino et al. 2016). The impact of a CBP is also mirrored in the evolution of ${\rm{X}}{(}^{16}{{\rm{O}}}_{{\rm{c}}})$, and leads to a larger ${\rm{X}}{(}^{16}{{\rm{O}}}_{{\rm{c}}})$ at the end of core He-burning and thus a smaller ¹²C mass fraction. All these models at constant ${\rm{\Delta }}{M}_{\max }$ have between 1600 and 1700 zones and a similar timestep distribution at this evolutionary stage, thus the differences between CBPs are not due to changes in resolution. The presence of CBPs has been suspected as being a numerical artifact (Caputo et al. 1989; Boothroyd et al. 1993; Constantino et al. 2016), and we show that it can also be due to the choice of the nuclear network.

The right panel of Figure 7 shows the evolution of ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$ as a function of mass resolution for approx22.net. There is a wide variation in ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$. The ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ and ${\rm{\Delta }}{M}_{\max }\ =0.05$ ${M}_{\odot }$ models both have CBPs, with the ${\rm{\Delta }}{M}_{\max }=0.05$ ${M}_{\odot }$ models starting their CBPs at earlier times with a larger ${\rm{X}}{(}^{4}{\mathrm{He}}_{{\rm{c}}})$. This reinforces the notion that these approx22.net CBPs are purely numerical artifacts. Increasing the resolution above ${\rm{\Delta }}{M}_{\max }=0.05$ ${M}_{\odot }$ (which is also necessary for the ICZ to penetrate the H-core) is sufficient to prevent the CBP from forming. Models with CBPs have between 1500 and 2000 zones, while models without CBP have between 2000 and 7000 zones.

During core He-burning, the convective core grows, largely because the mass of the He core itself grows from the overlying H-burning shell. This growth of the He core has two effects. First, as the mass of the He core grows, so does the core luminosity, but the radius of the convective core stays nearly the same. This causes the density of the core to steadily increase as core He-burning proceeds. Second, as the mass of the He core rises, the ratio of the gas pressure to the total pressure decreases, which favors efficient convection in the core.

The core prior to C-ignition is composed of the ashes of core He-burning, mainly ¹²C and ¹⁶O in a ≈1:4 ratio for the choice of ¹²C($\alpha ,\gamma$)¹⁶O reaction rate from the JINA reaclib version V2.0 2013-04-02 (see Figure 7, and West et al. 2013). There are also other isotopes present in trace amounts such as the neutron-rich ${}^{\mathrm{21,22}}$Ne, ${}^{\mathrm{25,26}}$Mg from processing the ashes of the CNO elements during He-burning; about 1% of X(²⁰Ne) from ¹⁶O($\alpha ,\gamma$)²⁰Ne during He-burning; the light s-process (which we do not follow); and other heavy elements present since the pre-MS from the initial composition.

Figure 8 shows the ⁴He mass fraction profile during C-burning for a subset of the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models. The left panel plots the X(⁴He) profile as a function of the number of isotopes in the network for the models without mass-loss and a mass resolution of ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$. The inset plot is a zoom-in of the inner edge of the He shell. We find that the inner and outer mass location of the He shell agree within $\approx \,\pm 0.05$ ${M}_{\odot }$ across all network sizes. The right panel plots the X(⁴He) profile as a function of mass resolution for the models without mass-loss, which use the largest reaction network considered, mesa_204.net. The inset plot is a zoom of the inner edge of the He shell. We find that mass resolutions of ${\rm{\Delta }}{M}_{\max }=0.1$ and 0.05 ${M}_{\odot }$ do not sufficiently resolve the He-burning shell locations. For mass resolutions ${\rm{\Delta }}{M}_{\max }\ \lesssim 0.02$ ${M}_{\odot }$, we find convergence of the He shell inner and outer boundaries.

Zoom In Zoom Out Reset image size

Figure 9 shows the C-core mass at the onset of C-burning as a function of mass resolution and the number of isotopes in the nuclear reaction network. The core mass is defined as the mass location where $X{(}^{12}{\rm{C}})\gt 0.1$ and $X{(}^{4}\mathrm{He})\lt 0.01$, measured when the core reaches ${\mathrm{log}}_{10}({T}_{{\rm{c}}}/K)=8.95$. These maps share many of the same variations, chiefly the bifurcations and outliers, that are inherited from core H-burning (see Figure 5). For the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models without mass-loss, the initial C-core mass is in the range $\approx 2.49\mbox{--}2.57$ ${M}_{\odot }$, while models with mass-loss span $\approx 2.43\mbox{--}2.52$ ${M}_{\odot }$. The ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models without mass-loss have initial C-core masses spanning $\approx 4.0\mbox{--}4.3$ ${M}_{\odot }$, while models with mass-loss span $\approx 3.9\mbox{--}4.1$ ${M}_{\odot }$. Our ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ models without mass-loss show initial C-core masses of $\approx 5.7\mbox{--}6.3$ ${M}_{\odot }$, while models with mass-loss span $\approx 5.6\mbox{--}5.8$ ${M}_{\odot }$. The ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models without mass-loss have an initial C-core mass range of $\approx 7.8\mbox{--}8.7$ ${M}_{\odot }$, while models with mass-loss span $\approx 7.3\mbox{--}7.8$ ${M}_{\odot }$.

Zoom In Zoom Out Reset image size

When we compare this with Figure 6, the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models with initially more massive He core have less massive C-cores at carbon ignition. This is because less massive ZAMS models have a larger electron degeneracy in the C-core, a larger ${\rho }_{{\rm{c}}}$, and thus an enhanced screening factor for the ¹²C+¹²C reaction rate. For models more massive than the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$, those with larger He-cores have the largest C-core at carbon ignition.

All the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models undergo convective C-ignition at the center of the star followed by a series of three convective flashes, where each additional flash ignites at the approximate maximum mass location of the previous convective C-flash. All the ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ and 30 ${M}_{\odot }$ models show C-ignition under radiative conditions. This ignition propagates outwards in mass, followed by a single convective flash. The ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models have a more complex C-ignition that straddles the boundary between radiative and convective carbon burning (Timmes et al. 1996; Heger et al. 2000; Hirschi et al. 2004).

Figure 10 shows the differences in core C-ignition for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ model without mass-loss and with the mesa_160.net network. As the mass resolution increases, C-ignition transitions from radiative conditions at the coarsest resolution (${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$) to a mixture of convective and radiative conditions at intermediate resolutions (${\rm{\Delta }}{M}_{\max }\ =0.05$ ${M}_{\odot }$), and finally to purely convective conditions at higher resolutions (${\rm{\Delta }}{M}_{\max }=0.02$ ${M}_{\odot }$). Models with yet finer mass resolution are similar to the ${\rm{\Delta }}{M}_{\max }=0.02$ ${M}_{\odot }$ case.

Zoom In Zoom Out Reset image size

Carbon burning in all cases begins at the center. In the case of the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ model with the coarsest resolution shown in Figure 10, the central temperature is high enough (${T}_{9,{\rm{c}}}\approx 0.8$) to burn carbon even though the central regions are still dominated by thermal neutrino cooling (light purple). Only near $M\approx 0.3$ ${M}_{\odot }$ at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})$ ≈ 1.7 does the nuclear energy generation rate from burning overtake the thermal neutrino cooling rate, and thus the color in the Kippenhahn transitions to red.

As the C-ignition conditions transition from radiative to convective, the time spent in the core C-burning phase increases from ≈45 to ≈250 yr. After carbon burning ceases, the star has $\approx 1$ year until core collapse. As the mass resolution increases, the mass location of the innermost boundary of the ⁴He shell (see Figure 8) decreases, while the mass location of the outer boundary of the ¹²C shell increases by $\approx 0.03$ ${M}_{\odot }$. In addition, as the resolution increases, ${\mathrm{log}}_{10}({T}_{{\rm{c}}}/K)$ decreases from 8.90 to 8.88 and log₁₀ (${\rho }_{c}$/(g cm⁻³)) decreases from 5.36 to 5.20. The combustion at the top of the final carbon flash is predominantly neutron captures, with the neutrons provided by ²²Ne, onto the ashes of C-burning to form ²³Na,${}^{\mathrm{25,26}}$Mg and ²⁷Al.

Figure 11 shows the mass fraction profiles at core C-depletion for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ model with mass-loss, ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$, and the mesa_204.net reaction network. The composition of the core prior to Ne-ignition is the ashes of core C-burning, mainly ¹⁶O, ${}^{\mathrm{20,21}}$Ne, ²³Na, ${}^{\mathrm{24,25,26}}$Mg, ${}^{\mathrm{26,27}}$Al, to a smaller extent ${}^{\mathrm{29,30}}$Si and ³¹P, and other heavy elements present since the pre-MS from the initial composition and from any s-processing during helium burning (Arnett & Truran 1969; Arnett 1972a; Endal 1975; Lamb et al. 1976; Arnett & Thielemann 1985; Woosley & Weaver 1995; Woosley et al. 2002; Hirschi et al. 2004; Bennett et al. 2012).

Zoom In Zoom Out Reset image size

Carbon is not completely destroyed, as would be expected for complete combustion, because the final convective flash lasts only a short time. The initial flashes in Figure 10 last tens or hundreds of years depending on the resolution and number of isotopes in the reaction network. The final flash lasts only a few years, which is insufficient to burn all the carbon over the approximately few solar masses that the final convection zone covers. The flash does not survive longer because the convection is driven by the burning at the base of the convection zone. This burning front attempts to propagate toward the central regions, but encounters the ashes from the previous flash. This ash has insufficient ¹²C to sustain the proto-flame, and the convection dies. The situation has similarities to C-burning flames in super asymptotic giant branch models (Denissenkov et al. 2013; Chen et al. 2014; Farmer et al. 2015). The carbon left behind is thus predominantly concentrated at the top of the CO core, and may ignite again as a carbon shell (see Figure 18).

Neon is the next abundant nucleus to burn in balanced power via the ${}^{20}\mathrm{Ne}{(\gamma ,\alpha )}^{16}{\rm{O}}$ photodisintegration reaction (Arnett 1974) at a core temperature of T₉ ≈1.5 and core density of ${\mathrm{log}}_{10}\rho \approx 6.6$ g cm⁻³. The net result is that 2(²⁰Ne) $\to$ ¹⁶O + ²⁴Mg at a rate determined by how fast ²⁰Ne captures α-particles from the equilibrium set up between ¹⁶O and ²⁰Ne (e.g., Busso & Gallino 1985; Thielemann & Arnett 1985; Chieffi et al. 1998; Woosley et al. 2002).

For the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models, Ne-ignition occurs predominately at the center, which drives a $\approx 0.5$ ${M}_{\odot }$ convection zone for approximately one month. As ${M}_{\mathrm{ZAMS}}$ increases to 20 and 25 ${M}_{\odot }$, this initial Ne-flash can be followed by a subsequent Ne-flash that may or may not drive a convective region once O-burning has become vigorous. For the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models, the initial Ne-flash is followed by one or more additional Ne-flashes that propagate outwards in mass to $\approx 1$ ${M}_{\odot }$.

Figure 12 shows a variety of Ne-burning behaviors for the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models with mass-loss. The model using the approx22.net reaction network (upper left) shows a series of outward moving convective flashes starting at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})=-0.2$ and ending at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})=-0.8$. The model using the mesa_127.net reaction network (upper right) shows that neon ignites in a weak radiative flash that lasts $\approx 1$ month. For the mesa_160.net model (lower left), there is an extensive off-center radiative burning region that transitions into a convective flash. Finally, the mesa_204.net model contains a series of off-center flashes, where a pocket of convection persists from the first ignition of neon to the ignition of oxygen.

Zoom In Zoom Out Reset image size

Figure 13 shows the evolution of ${T}_{{\rm{c}}}$ and ${\rho }_{{\rm{c}}}$ for the models shown in Figure 12. The ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models start Ne-burning at ${\mathrm{log}}_{10}(T/K)\approx 9.17$ and log₁₀ (${\rho }_{c}$/(g cm⁻³)) ≈ 6.6. As Ne-burning progresses, the density and temperature increase until core O-burning begins with an accompanying creation of a central convection zone that lowers ${\rho }_{{\rm{c}}}$. The tracks are well converged. The maximum temperature difference is from ${T}_{9,{\rm{c}}}$ ≈ 1.65 to ≈1.50, which is an ≈8% difference, with the largest offset in the approx22.net. As the number of isotopes in the reaction network increases, the core is denser for a given temperature. The larger networks thus undergo increased neutrino cooling, which is density dependent, but not dependent on the isotopes in the network. This increased cooling rate prevents the neon from vigorously igniting at the center, and prevents the ignition from driving a central convection zone.

Zoom In Zoom Out Reset image size

These variations in burning structure that are due to changes in the number of isotopes in the nuclear reaction network lead to variations in the post Ne-burning abundance profiles. For the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models in Figure 12, the approx22.net model burns neon the longest amount of time and over the most mass, and thus shows the largest abundance changes. We find that X(²⁴Mg) is enhanced by a factor ≈6 compared to the pre-Ne-burning mass fraction. Models with the softwired mesa_127.net, mesa_160.net, and mesa_204.net networks show an X(²⁴Mg) enhancement factor that increases as the number of isotopes increases: from ≈1.5 for the mesa_127.net model to ≈2.0 for the mesa_204.net model, with a total of $\approx 0.05 \%$ X(²⁴Mg) left. The other main product of Ne-burning, ¹⁶O, increases from ≈0.7 to ≈0.8 in the inner 1.5 ${M}_{\odot }$, with the relative change increasing as the number of isotopes increases similar to the ²⁴Mg.

The large ¹⁶O mass fraction, coupled with ¹⁶O+¹⁶O being a true fusion reaction and not a photodisintegration-driven event like Ne-burning or Si-burning, ensures that O-burning is a key energetic and nucleosynthesis stage in the late phases of massive star evolution (Rakavy et al. 1967; Arnett 1972b; Woosley et al. 1972; Woosley & Weaver 1995; El Eid et al. 2004; Sukhbold et al. 2016).

In 1D stellar evolution instruments such as MESA , convective mixing and energy transport is modeled using MLT in a time-dependent manner (Vitense 1953; Böhm-Vitense 1958; Cox & Giuli 1968), which is usually tuned to reproduce solar properties (Asplund et al. 2009). However, thermal neutrino cooling speeds up O-burning to such an extent that the evolutionary timescales are close to the sound-crossing time, so that direct, multidimensional compressible numerical hydrodynamics must be applied for maximum fidelity to the underlying physics. Such studies show that nuclear burning tightly couples to turbulent convection so that fuel is consumed in chaotic episodes (Bazán & Arnett 1998; Asida & Arnett 2000; Meakin & Arnett 2007a). Core O-burning and shell O-burning are dominated by large-scale modes of fluid flow, which are of such low order that they do not cancel to a smooth spherical behavior (Couch et al. 2015; Chatzopoulos et al. 2016; Jones et al. 2016; Müller et al. 2016). Moreover, 3D simulations of O-burning suggest that MLT gives an incomplete representation of stellar convection (Arnett et al. 2015). Nevertheless, 3D simulations of core O-burning and beyond are resource intensive and to date have been run primarily to address hydrodynamic and transport aspects. Such 3D simulations have not yet been run to assess the detailed nucleosynthesis.

Figure 14 shows the O-core masses at the start of core oxygen convection. We define this as the location where $X{(}^{16}{\rm{O}})\gt 0.1$ and $X{(}^{12}{\rm{C}})\lt 0.01$ and measured when $X{(}^{16}{{\rm{O}}}_{{\rm{c}}})=0.7$. In general, as the mass resolution increases, the O-core mass decreases, with the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models having the smallest spread. This is because the 15 ${M}_{\odot }$ models have a more consistent behavior between the start of carbon and the end of neon burning. The O-core masses range from 1.35–1.4 ${M}_{\odot }$ for the 15 ${M}_{\odot }$ models, 1.5–2.5 ${M}_{\odot }$ for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models, 1.8–3.0 ${M}_{\odot }$ for the 25 ${M}_{\odot }$, and 2.25–3.15 ${M}_{\odot }$ for the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models. The spread in the O-core mass is smaller for the mass-losing models compared to models without mass-loss. Much of the fine structure seen in Figures 6 and 9 has been erased as a result of the variations in carbon and neon burning. Panels (e) and (f) still show a clear bifurcation, however, with the lowest resolution models having the highest O-core masses.

Zoom In Zoom Out Reset image size

Figure 15 shows the evolution during core O-burning for the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models without mass-loss, the mesa_204.net reaction network, and mass resolutions of ${\rm{\Delta }}{M}_{\max }\ =0.1$ ${M}_{\odot }$, 0.02 ${M}_{\odot }$, and 0.005 ${M}_{\odot }$. In each model O-ignition occurs at the center at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})\approx 0.4$ and O-burning continues until ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})\approx -0.5$. O-burning is concentrated over $\approx 0.1$ ${M}_{\odot }$ of the core, while a convective region is generated out to $\approx 1.0$ ${M}_{\odot }$. For ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$, the convective region undergoes a series of contractions and expansions (in mass coordinates), analogous to the CBP seen during core He-burning. During a core oxygen breathing pulse (OBP), when the convective region reaches its maximum extent (in mass coordinate) and ingests ¹⁶O, there is a brief increase in the extent of the central burning region and a corresponding increase in the energy released in the core. These OBPs have a smaller effect on the ¹⁶O abundance than CBPs because of the weaker density dependency of O-burning compared to He-burning. These OBPs are present for all four masses considered in this work, and unlike the helium counterparts, these OBPs are independent of the nuclear network used. As the mass resolution increases, the OBP behavior changes, the number of OBPs decreases, the maximum mass that the OBPs extend out to decreases, and they have a shorter lifetime. This suggests that careful treatment of the mass resolution is needed during this stage in the evolution of a star.

Zoom In Zoom Out Reset image size

Once core O-burning ceases, the core contracts until Si-ignition. During this contraction phase, the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models undergo a series of off-center flashes. These mimic the ignition at the center by having a series of small neon flashes that lead to a larger oxygen flash. The residual carbon left behind from carbon burning (Figure 11) can now also ignite as an additional flash near the edge of the CO core. For the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models in Figure 15, as the numerical resolution increases, the number of these flashes decreases, the ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ has three oxygen flashes and a number of smaller neon flashes, while the ${\rm{\Delta }}{M}_{\max }=0.02$ ${M}_{\odot }$ and ${\rm{\Delta }}{M}_{\max }=0.005$ ${M}_{\odot }$ has only one oxygen flash and a variable number of neon flashes. As the initial mass of the star increases, the number of neon/oxygen flashes decreases and the carbon flash ignites later, once silicon burning has commenced.

The ashes of O-burning are dominated by Si, S, Ar, and Ca in roughly solar proportions. Figure 16 shows that the isotopes ²⁸Si, ${}^{\mathrm{32,33,34}}$S, ${}^{\mathrm{35,37}}$Cl, ${}^{\mathrm{36,38}}$Ar, ${}^{\mathrm{39,41}}$ K, and ${}^{\mathrm{40,42}}$Ca are present in significant quantities (e.g., Woosley & Weaver 1995; Limongi et al. 2000; Limongi & Chieffi 2003; Sukhbold et al. 2016) for the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ model with mass-loss and ${\rm{\Delta }}{M}_{\max }\ =0.01$ ${M}_{\odot }$. During O-burning, several isotopes begin to be made as radioactive progenitors of stable isotopes, for example, ³⁷Cl as ³⁶Ar and ⁴¹K as ⁴¹Ca. The heaviest isotopes ($A\gg 50$) that were present in the initial composition at birth, along with those heavy isotopes made from the s-process during He-burning and C-burning (which we do not follow), begin to be destroyed by photodisintegration reactions; essentially melting them into the Fe-group. O-burning also features the first appearance of quasi-equilibrium clusters (Bodansky et al. 1968; Woosley et al. 1973; Woosley & Hoffman 1992; Hix & Thielemann 1996; Meyer et al. 1998). These clusters grow to encompass more isotopes as core O-burning proceeds. Weak interactions such as ³³S(e⁻, ${\nu }_{e}$)³³P, ³⁵35(e⁻, ${\nu }_{e}$)³⁵S, and ³⁷Ar(e⁻, ${\nu }_{e}$)³⁷Cl decrease the ${Y}_{{\rm{e}}}$ significantly during O-burning, especially at O-depletion (e.g., Heger et al. 2001). For example, we find ${Y}_{{\rm{e}},{\rm{c}}}$ ≈ 0.4778 at core O-depletion for the model shown in Figure 16, which has decreased from ${Y}_{{\rm{e}},{\rm{c}}}=0.499$ before O-burning.

Zoom In Zoom Out Reset image size

Silicon burning is the last exothermic burning stage and produces the Fe-peak nuclei. Owing to Coulomb repulsion, it is improbable that two ²⁸Si nuclei will fuse to ⁵⁶Ni. Instead, a photodisintegration-driven rearrangement of the abundances takes place, originating from equilibria established among individual reactions with their reverse reactions (Bodansky et al. 1968). When such equilibria occur among many reactions, the material reaches an equilibrium state where nuclei merge into clusters. Units of interaction are no longer nuclei, but the clusters themselves, which adapt their properties according to the local thermodynamic conditions (Woosley et al. 1973; Woosley & Hoffman 1992; Hix & Thielemann 1996; Meyer et al. 1998; The et al. 1998; Hix & Thielemann 1999; Magkotsios et al. 2010).

In general, not all reactions are in equilibrium. Consequently, this state is named QSE. The special case where all strong and electromagnetic reactions are balanced by their reverse reactions is called NSE, because all mass fractions may be described in terms of statistical properties of excited nuclear states (e.g., temperature-dependent partition functions) and nuclear structure variables (masses and Q values; Clifford & Tayler 1965; Hartmann et al. 1985; Jordan & Meyer 2004; Nadyozhin & Yudin 2004; Seitenzahl et al. 2008).

Weak interactions are excluded from these definitions since for conditions relevant to hadronic physics, they never attain equilibrium (e.g., Heger et al. 2001; Arcones et al. 2010). Hence, equilibrium notions are only related with strong and electromagnetic interactions. In practice, there is either one cluster in NSE or QSE, or two QSE clusters, one for the Si-group and one for the Fe-group nuclei. Even the main products of Si-burning depend quite sensitively on small changes in the electron fraction, temperature, and density. Although the physics of QSE/NSE is relevant for our models, we reiterate that no QSE or NSE approximations are made in our models.

During Si-burning, as for O-burning, the energetics of nuclear burning tightly couples to turbulent convection, and must be modeled with 3D simulations to assess the fidelity of the approximations made by 1D stellar evolution instruments (e.g., Arnett & Meakin 2011; Couch et al. 2015; Jones et al. 2016; Müller et al. 2016). For example, Couch et al. (2015) found that during the final minutes of Si-burning in a massive star that fluctuating, peak-convective speeds of ≈200 km s⁻¹ were common and that the speed of the convection increased to ≈500 km s⁻¹ as collapse approached and the core contracted. These speeds are not negligible relative to nominal infall speeds of 1000 km s⁻¹ for our core-collapse initial models. It remains a challgenge for future investigations to distill the essential features of 3D simulations into models suitable for 1D stellar evolution instruments.

Figure 17 shows the Si-core mass at Si-ignition. This is defined as the location where $X{(}^{28}\mathrm{Si})\gt 0.1$ and $X{(}^{16}{\rm{O}})\lt 0.01$, measured when ${\mathrm{log}}_{10}({T}_{{\rm{c}}}/K)=9.5$. As the mass resolution increases, there is a general trend, with substantial scatter, for the Si-core mass to decrease, with the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ and 20 ${M}_{\odot }$ models having the smallest spread. The Si-core masses range from 1.05 to 1.35 ${M}_{\odot }$ for the 15 ${M}_{\odot }$ models, 1.3–1.6 ${M}_{\odot }$ for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models, 1.0–1.7 ${M}_{\odot }$ for the 25 ${M}_{\odot }$ models, and 1.1–1.6 ${M}_{\odot }$ for the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models. The spread in the Si-core mass is about the same for the mass-losing models as for models without mass-loss. Nearly all of the fine structure seen in Figures 6 and 9, and even the coarser structure in Figure 14, has been largely erased before the onset of Si-burning. The bulk of the silicon core is built during the core oxygen-burning phase, but a number of short flashes, which can be seen in Figure 15 at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})\approx -1.5$, can occur before core Si-burning commences. These flashes introduce an additional fine structure into the core composition and hence the location of the silicon core boundary.

Zoom In Zoom Out Reset image size

Figure 18 shows the evolution during core Si-burning for ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ models with mass-loss and ${\rm{\Delta }}{M}_{\max }\ =0.005$ ${M}_{\odot }$ for the approx22.net, mesa_160.net and mesa_204.net reaction networks. Silicon ignites at the center within one day of core collapse. The initial phase of transforming ²⁸Si is sensitive to the thermodynamic conditions, and electron captures change the ${Y}_{{\rm{e}}}$ continuously during the buildup of the QSE clusters (e.g., Thielemann & Arnett 1985).

Zoom In Zoom Out Reset image size

For the approx22.net model in Figure 18, Si-burning propagates outwards to $\approx 1$ ${M}_{\odot }$, driving a convection region over this mass range. The Si-burning region then recedes to $\approx 0.5$ ${M}_{\odot }$ and leaves the convection zone behind. For the models using the mesa_160.net and mesa_204.net reaction networks, there is a central burning region occupying $\approx 0.2$ ${M}_{\odot }$. Outside this core region is a mixed region comprised of many small pockets of radiative/convective burning. The larger and more simply connected convective region for the approx22.net allows a larger fraction of the ²⁸Si to be converted into iron-group elements and to homogenize the spatial distribution of these iron-group isotopes. The final ${Y}_{{\rm{e}}}$ in the mesa_160.net and mesa_204.net models is determined to within ≈2% before the iron core begins to collapse, making shell Si-burning a key evolutionary state for determining the ${Y}_{{\rm{e}},{\rm{c}}}$ and the iron-core structure (e.g., Heger et al. 2001).

Figure 18 also shows a change in the qualitative behavior of the carbon shell. Carbon burns at the base of a convective region in the ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ models using the approx22.net reaction network, at $\approx 3.0$ ${M}_{\odot }$ and ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})\geqslant -2.25$, and has been since core O-burning. In contrast, the models using the mesa_160.net reaction network show C-burning in two distinct layers, one layer at $\approx 2.0$ ${M}_{\odot }$ and a second layer at $\approx 2.4$ ${M}_{\odot }$, once Si-burning becomes vigorous. The ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ models using the mesa_204.net reaction network do not ignite a carbon shell, which is chiefly due to an oxygen flash at $\approx 1.5$ ${M}_{\odot }$ and ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})=-2.25$, where a previous convection region extended out to $\approx 5$ ${M}_{\odot }$. This region partially burnt C, Ne, and O, leaving behind the mixed convection/semiconvection that is due to small-scale composition gradients, as seen in Figure 18 for mesa_204.net.

Figure 19 shows the composition of the cores at ${\mathrm{log}}_{10}(({\tau }_{\mathrm{cc}}-\tau )/\mathrm{yr})=-3.0$ for the ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$ models in Figure 18. The left panel corresponds to the stellar model using the approx22.net reaction network, and the right panel corresponds to the models using the mesa_204.net reaction network. First, we note that the approx22.net model has significantly smoother composition profiles in the inner $\approx 1.5$ ${M}_{\odot }$ because of the more simply connected convective region it exhibits (see Figure 18) relative to the multiple-connected convective regions exhibited by the mesa_204.net model. Second, the core in the approx22.net model is predominately ⁵⁴Fe, while at the same time, the mesa_204.net model is a mix of ${}^{\mathrm{54,56}}$Fe and ⁵²Cr until core collapse. At the center of the mesa_204.net model, the composition is ≈68% ⁵⁶Fe, ≈18% ⁵²Cr, and ≈2% ⁵⁴Fe, with the remainder in various iron-group elements. Third, there is a significant impact of the different shell O-burning behavior in the 2–5 ${M}_{\odot }$ range before silicon burning commences (see Figure 18). The model using the mesa_204.net reaction network has a smaller ¹⁶O mass fraction that extends deeper into the core than the model using the approx22.net reaction network. The O shell has also been polluted by the products of oxygen-burning and is also depleted in ¹²C.

Zoom In Zoom Out Reset image size

When the Fe core reaches its finite-temperature Chandrasekhar mass (Baron & Cooperstein 1990),

$$\begin{eqnarray}{M}_{\mathrm{Ch},\mathrm{eff}} & \simeq & 5.76\ {Y}_{{\rm{e}}}^{2}\\ & & \times \left[1-0.057+{\left(\displaystyle \frac{{s}_{e}}{\pi {Y}_{{\rm{e}}}}\right)}^{2}+1.21\left(\displaystyle \frac{{s}_{e}}{A}\right)\right]\,{M}_{\odot },\end{eqnarray} \tag{ 4 }$$

electron capture and photodisintegration drive the collapse of the Fe core with the largest infall speeds being reached near the outer edge of the Fe core. Here, ${s}_{e}={S}_{e}/({N}_{A}k)$ is a dimensionless average electron entropy and A is an average atomic weight. We terminate our MESA models when any mass coordinate within the Fe core exceeds an inward velocity of 1000 km s⁻¹. The structural and nucleosynthesis properties at this key evolutionary point (e.g., Woosley & Weaver 1995; Hirschi et al. 2004; Dessart et al. 2010; Chieffi & Limongi 2013) can have significant effects on the subsequent explosion (if achieved) and nucleosynthesis (e.g., Janka 2012; Dolence et al. 2013; Ott et al. 2013; Couch & O'Connor 2014; Couch et al. 2015; Pejcha & Thompson 2015; Bruenn et al. 2016; Jones et al. 2016; Müller et al. 2016).

Figure 20 shows the final Fe-core mass (which we define below) at core collapse, taken when the radial infall speed reaches $v\gt 1000$ $\mathrm{km}\,{{\rm{s}}}^{-1}$, as a function of the mass resolution ${\rm{\Delta }}{M}_{\max }$ and the number of isotopes in the reaction network. For a fixed ${\rm{\Delta }}{M}_{\max }$, in most cases, the final Fe-core mass for models using the mesa_127.net, mesa_160.net, and mesa_204.net reaction networks agree to within $\approx 0.05$ ${M}_{\odot }$. In some cases, the approx22.net and mesa_79.net reaction networks models produce a more massive Fe core (≈0.07–0.13 ${M}_{\odot }$). This is most noticeable for the 15 ${M}_{\odot }$ models without mass-loss, where mesa_160.net and mesa_204.net models yield an Fe core of ≈1.33 ${M}_{\odot }$ while the approx22.net model gives ${M}_{\mathrm{Fe}}$ ≈ 1.51 ${M}_{\odot }$, which is a difference of ≈14%.

Zoom In Zoom Out Reset image size

The definition of the final Fe-core mass is worth mentioning. We use the MESA definition, which is the maximum mass location where $X({}^{56}\mathrm{Fe})\gt 0.1$ and $X({}^{28}\mathrm{Si})\lt 0.01$. Other options for defining the Fe-core mass include the location of the maximum infall velocity, a fixed ${Y}_{{\rm{e}}}$ value, or the ${Y}_{{\rm{e}}}$ jump (Woosley & Weaver 1995; Heger et al. 2001). Comparing the different measures, we find that the median Fe-core mass using the MESA abundance definition is $\approx 0.01$ ${M}_{\odot }$ lower than that when using the peak velocity location. When we use the ${Y}_{{\rm{e}}}$ jump definition, the core mass is $\approx 0.03$ ${M}_{\odot }$ greater than the peak velocity location. A fixed ${Y}_{{\rm{e}}}$ definition was found to be unsuitable because of variations between the models in the ${Y}_{{\rm{e}}}$ value outside the Fe core, which could lead to changes as large as $\approx 0.5$ ${M}_{\odot }$ in the core mass location. Consideration needs to be given to finding a consistent and reliable definition for the final core mass to enable comparisons between different models. We do not claim that the MESA definition is superior.

For a fixed reaction network, say mesa_160.net, Figure 20 shows that mass resolution can have a significant impact on the final Fe-core mass. The largest spread occurs for the ${M}_{\mathrm{ZAMS}}=25$ ${M}_{\odot }$, where the Fe-core mass ranges from ≈1.38–1.8 ${M}_{\odot }$. For the 25 ${M}_{\odot }$ model the core masses split into two groups based on the resolution, a similar trend to that seen in its He-core mass (Figure 6), those models with lower He-core masses have grown smaller iron-core masses. However, for mass resolutions of ${\rm{\Delta }}{M}_{\max }\leqslant 0.01$, the final Fe-core mass is monotonic as the ZAMS mass increases, and agrees to within ±0.1 ${M}_{\odot }$ as the resolution increases. Whether the monotonic trend at the highest mass resolutions remains monotonic with a significantly finer grid of ZAMS masses is under consideration (I. Petermann et al. 2016, in preparation).

In comparing the mass-loss models against those without mass-loss, we find that the iron-core masses are slightly lower on average for models with mass-loss. However, the decrease in the mass of the cores is much smaller than the decrease in the final mass. For instance, the 30 ${M}_{\odot }$ models may lose $\approx 10$ ${M}_{\odot }$ of material over their lifetime (Table 1), but the iron cores will only be $\approx 0.1\mbox{--}0.2$ ${M}_{\odot }$ smaller. Mass loss rates are uncertain for the mass range considered here, but we have shown that the sizes of the iron cores are only mildly sensitive to the total amount of mass lost.

Table 1.  Median Pre-SN Model Properties, with Upper and Lower Bounds

Property	15 ${M}_{\odot }$		20 ${M}_{\odot }$		25 ${M}_{\odot }$		30 ${M}_{\odot }$
	$\dot{M}=0$	$\dot{M}\ne 0$	$\dot{M}=0$	$\dot{M}\ne 0$	$\dot{M}=0$	$\dot{M}\ne 0$	$\dot{M}=0$	$\dot{M}\ne 0$
${\mathrm{He}}_{\mathrm{core}}\,({M}_{\odot })$ a,b	${2.82}_{2.79}^{2.82}$	${2.77}_{2.72}^{2.78}$	${4.67}_{4.59}^{4.70}$	${4.57}_{4.52}^{4.59}$	${6.88}_{6.80}^{7.28}$	${6.56}_{6.52}^{6.67}$	${9.44}_{9.15}^{9.89}$	${8.67}_{8.62}^{8.78}$
${{\rm{C}}}_{\mathrm{core}}\,({M}_{\odot })$	${2.51}_{2.49}^{2.58}$	${2.44}_{2.43}^{2.53}$	${4.19}_{4.04}^{4.75}$	${4.07}_{3.69}^{4.08}$	${6.02}_{4.34}^{6.43}$	${5.75}_{5.53}^{5.92}$	${8.28}_{7.13}^{8.79}$	${7.62}_{7.22}^{7.82}$
${{\rm{O}}}_{\mathrm{core}}\,({M}_{\odot })$	${1.41}_{1.35}^{1.43}$	${1.40}_{1.32}^{1.42}$	${1.54}_{1.43}^{2.47}$	${1.57}_{1.41}^{2.05}$	${2.34}_{1.74}^{3.04}$	${1.81}_{1.76}^{2.47}$	${2.38}_{2.14}^{3.18}$	${2.39}_{2.26}^{3.06}$
${\mathrm{Si}}_{\mathrm{core}}\,({M}_{\odot })$	${1.15}_{1.02}^{1.38}$	${1.15}_{1.08}^{1.39}$	${1.38}_{1.30}^{1.65}$	${1.40}_{1.24}^{1.48}$	${1.19}_{0.91}^{1.61}$	${1.40}_{1.07}^{1.67}$	${1.16}_{1.08}^{1.64}$	${1.15}_{1.12}^{1.66}$
${\mathrm{Ye}}_{{\rm{c}},\mathrm{He}}$ b	${0.505}_{0.505}^{0.505}$	${0.505}_{0.505}^{0.505}$	${0.505}_{0.505}^{0.505}$	${0.505}_{0.505}^{0.505}$	${0.505}_{0.505}^{0.505}$	${0.505}_{0.505}^{0.505}$	${0.505}_{0.505}^{0.505}$	${0.505}_{0.505}^{0.505}$
${\mathrm{Ye}}_{{\rm{c}},{\rm{C}}}$	${0.499}_{0.499}^{0.500}$	${0.499}_{0.499}^{0.500}$	${0.499}_{0.499}^{0.500}$	${0.499}_{0.499}^{0.500}$	${0.499}_{0.499}^{0.500}$	${0.499}_{0.499}^{0.500}$	${0.499}_{0.499}^{0.500}$	${0.499}_{0.499}^{0.500}$
${\mathrm{Ye}}_{{\rm{c}},{\rm{O}}}$	${0.499}_{0.498}^{0.500}$	${0.499}_{0.498}^{0.500}$	${0.499}_{0.498}^{0.500}$	${0.499}_{0.498}^{0.500}$	${0.498}_{0.498}^{0.500}$	${0.499}_{0.484}^{0.500}$	${0.498}_{0.498}^{0.500}$	${0.498}_{0.490}^{0.500}$
${\mathrm{Ye}}_{{\rm{c}},\mathrm{Si}}$	${0.486}_{0.475}^{0.498}$	${0.486}_{0.475}^{0.498}$	${0.488}_{0.483}^{0.498}$	${0.489}_{0.483}^{0.498}$	${0.491}_{0.483}^{0.497}$	${0.490}_{0.482}^{0.499}$	${0.491}_{0.489}^{0.498}$	${0.490}_{0.488}^{0.497}$
${\mathrm{Ye}}_{{\rm{c}},\mathrm{Fe}}$	${0.431}_{0.419}^{0.461}$	${0.432}_{0.414}^{0.461}$	${0.438}_{0.425}^{0.462}$	${0.438}_{0.416}^{0.462}$	${0.438}_{0.414}^{0.462}$	${0.438}_{0.423}^{0.462}$	${0.444}_{0.437}^{0.462}$	${0.442}_{0.428}^{0.462}$
${\tau }_{{\rm{H}}}\,(\mathrm{Myr})$ c	${10.95}_{10.93}^{10.96}$	${10.99}_{10.94}^{11.00}$	${7.73}_{7.72}^{7.74}$	${7.78}_{7.76}^{7.79}$	${6.18}_{6.16}^{6.51}$	${6.24}_{6.22}^{6.38}$	${5.53}_{5.26}^{5.62}$	${5.34}_{5.33}^{5.43}$
${\tau }_{\mathrm{He}}\,(\mathrm{Myr})$	${1.69}_{1.48}^{1.71}$	${1.74}_{1.51}^{1.75}$	${1.11}_{1.05}^{1.72}$	${1.10}_{1.08}^{1.29}$	${0.81}_{0.77}^{1.19}$	${0.82}_{0.81}^{0.89}$	${0.65}_{0.63}^{0.73}$	${0.69}_{0.68}^{0.74}$
${\tau }_{{\rm{C}}}\,(\mathrm{yr})$	${78.37}_{75.65}^{82.62}$	${81.80}_{76.15}^{87.04}$	${111.36}_{25.90}^{115.61}$	${125.42}_{28.52}^{131.19}$	${27.06}_{14.81}^{28.42}$	${23.35}_{22.53}^{26.98}$	${16.56}_{7.97}^{22.45}$	${19.82}_{15.11}^{22.62}$
${\tau }_{{\rm{O}}}\,(\mathrm{yr})$	${3.83}_{3.28}^{5.05}$	${4.08}_{3.51}^{5.36}$	${1.44}_{0.16}^{2.76}$	${1.28}_{0.26}^{3.12}$	${0.14}_{0.08}^{0.29}$	${0.31}_{0.01}^{0.93}$	${0.13}_{0.02}^{0.19}$	${0.14}_{0.10}^{0.16}$
${\tau }_{\mathrm{Si}}\,(\mathrm{days})$	${3.43}_{0.74}^{5.05}$	${3.74}_{0.75}^{5.60}$	${0.81}_{0.32}^{1.56}$	${0.91}_{0.27}^{2.04}$	${0.61}_{0.18}^{6.58}$	${0.66}_{0.16}^{1.59}$	${0.35}_{0.15}^{1.92}$	${0.41}_{0.21}^{1.97}$
${\tau }_{\mathrm{Fe}}\,(\mathrm{hr})$	${33.16}_{11.21}^{57.11}$	${36.38}_{11.87}^{74.01}$	${11.11}_{3.66}^{35.22}$	${11.63}_{6.59}^{17.95}$	${6.57}_{4.04}^{25.25}$	${9.55}_{4.26}^{19.72}$	${4.74}_{3.23}^{24.33}$	${4.77}_{3.37}^{25.36}$
${\mathrm{He}}_{\mathrm{shell}}\,({M}_{\odot })$ d,e	${4.17}_{4.16}^{4.22}$	${4.09}_{4.08}^{4.17}$	${6.20}_{6.09}^{6.88}$	${6.05}_{5.84}^{6.07}$	${8.26}_{8.12}^{8.71}$	${7.95}_{7.71}^{8.10}$	${10.87}_{9.66}^{11.30}$	${9.99}_{9.58}^{10.24}$
${{\rm{C}}}_{\mathrm{shell}}\,({M}_{\odot })$	${2.51}_{2.49}^{2.58}$	${2.44}_{2.43}^{2.54}$	${4.19}_{4.04}^{4.75}$	${4.07}_{3.61}^{4.09}$	${6.03}_{5.84}^{6.44}$	${5.76}_{5.54}^{5.88}$	${8.26}_{7.09}^{8.73}$	${7.57}_{7.16}^{7.82}$
${{\rm{O}}}_{\mathrm{shell}}\,({M}_{\odot })$	${2.42}_{2.18}^{2.53}$	${2.33}_{2.13}^{2.49}$	${3.95}_{3.27}^{4.14}$	${3.82}_{2.61}^{3.93}$	${5.38}_{4.08}^{5.69}$	${5.25}_{2.80}^{5.63}$	${6.87}_{5.64}^{8.20}$	${6.47}_{4.68}^{7.20}$
${\mathrm{Si}}_{\mathrm{shell}}\,({M}_{\odot })$	${1.59}_{0.00}^{1.70}$	${1.61}_{1.24}^{1.65}$	${1.80}_{1.56}^{2.97}$	${1.81}_{0.00}^{2.11}$	${1.90}_{1.53}^{2.67}$	${1.88}_{1.62}^{2.22}$	${2.27}_{1.73}^{2.58}$	${2.30}_{1.50}^{3.14}$
${\mathrm{Fe}}_{\mathrm{core}}\,({M}_{\odot })$	${1.41}_{0.75}^{1.55}$	${1.42}_{0.77}^{1.53}$	${1.55}_{1.35}^{1.86}$	${1.57}_{1.38}^{1.74}$	${1.66}_{1.35}^{1.88}$	${1.59}_{1.46}^{1.83}$	${1.80}_{1.58}^{1.90}$	${1.79}_{1.39}^{1.90}$
${H}_{\mathrm{env}}\,({M}_{\odot })$ d,f	${6.96}_{6.94}^{7.21}$	${5.65}_{5.59}^{6.14}$	${8.41}_{7.58}^{8.55}$	${6.74}_{6.50}^{7.18}$	${9.76}_{9.17}^{9.96}$	${6.09}_{4.69}^{7.22}$	${10.47}_{9.99}^{11.02}$	${4.60}_{0.00}^{6.10}$
${M}_{\mathrm{final}}\,({M}_{\odot })$ d,g	⋯	${13.02}_{12.94}^{13.37}$	⋯	${17.26}_{16.91}^{18.11}$	⋯	${18.95}_{17.01}^{21.17}$	⋯	${20.27}_{17.38}^{22.06}$
${\mathrm{log}}_{10}(R/{R}_{\odot })$ d	${2.97}_{2.96}^{2.99}$	${2.98}_{2.97}^{3.00}$	${3.08}_{3.06}^{3.08}$	${3.07}_{3.05}^{3.07}$	${3.08}_{3.05}^{3.09}$	${3.07}_{3.05}^{3.08}$	${2.97}_{2.92}^{3.02}$	${2.99}_{2.95}^{3.03}$
${\mathrm{log}}_{10}({T}_{\mathrm{eff}}/K)$ d	${3.50}_{3.50}^{3.51}$	${3.50}_{3.49}^{3.52}$	${3.52}_{3.51}^{3.53}$	${3.52}_{3.51}^{3.53}$	${3.56}_{3.55}^{3.58}$	${3.56}_{3.54}^{3.73}$	${3.65}_{3.62}^{3.74}$	${3.63}_{3.60}^{3.67}$
${\mathrm{log}}_{10}(L/{L}_{\odot })$ d	${4.91}_{4.90}^{4.93}$	${4.90}_{4.88}^{5.01}$	${5.18}_{5.13}^{5.23}$	${5.15}_{5.09}^{5.20}$	${5.36}_{5.30}^{5.42}$	${5.33}_{5.25}^{6.03}$	${5.47}_{5.38}^{5.84}$	${5.49}_{5.40}^{5.56}$
${\xi }_{M=2.5}$ d,h	${0.10}_{0.04}^{0.15}$	${0.08}_{0.04}^{0.13}$	${0.24}_{0.10}^{0.63}$	${0.24}_{0.13}^{0.43}$	${0.32}_{0.16}^{0.66}$	${0.27}_{0.19}^{0.60}$	${0.60}_{0.31}^{0.69}$	${0.58}_{0.19}^{0.69}$

Notes.

^aCore mass values, see Figures 6, 9, 14, 17, and 20 for definitions. ^bMeasured at the corresponding ignition of each fuel, except for Fe, which is measured at core collapse. ^cApproximate time to transition to the next major fuel source. ^dMeasured at core collapse. ^eOuter mass coordinate where the element is the most abundant. ^fMass of H-rich envelope. ^gTotal mass of star. ^hCompactness parameter, with $M=2.5$ ${M}_{\odot }$.

Download table as:  ASCIITypeset image

Figure 20 shows that the choice of reaction network can directly affect the final value of ${M}_{\mathrm{Fe}}$. In turn, this determines the location of the maximum infall velocity. This is most notable in Figure 21 for the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models, where the infall velocity occurs at ≈1.51 ${M}_{\odot }$ and 1.55 ${M}_{\odot }$ for the approx22.net and mesa_79.net, respectively. For the larger networks, the core infall velocity occurs at arithmetic mean mass location of ${1.4}_{-0.04}^{+0.02}$ ${M}_{\odot }$. The $20\mbox{--}30$ ${M}_{\odot }$ stars in Figure 21 follow similar trends, with approx22.net and mesa_79.net producing the largest deviations about the arithmetic mean Fe-core value of ≈0.07–0.13 ${M}_{\odot }$. The 15 ${M}_{\odot }$ models, with ${\rm{\Delta }}{M}_{\max }=0.005$, mass-loss, and either the mesa_127.net or mesa_160.net have anomalously low ${M}_{\mathrm{Fe}}\lt 1.0$ ${M}_{\odot }$, these models achieve our definition of core collapse infall velocity in a few spatial zones, but they are not undergoing a collapse, and thus the evolution terminates early before the core could grow to its full extent.

Zoom In Zoom Out Reset image size

Figure 21 shows the ${Y}_{{\rm{e}}}$ and radial velocity profiles at core collapse for the ${M}_{\mathrm{ZAMS}}$ = 15, 20, 25, and 30 ${M}_{\odot }$ models with mass-loss and ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$ for different reaction networks. In most cases, we find that the models with the approx22.net reaction network underestimate the ${Y}_{{\rm{e}},{\rm{c}}}$ when compared to the larger reaction networks. The converse is found for models using the mesa_79.net network, which tend to overestimate ${Y}_{{\rm{e}},{\rm{c}}}$. For example, in the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ case, the models using the mesa_127.net, mesa_160.net, and mesa_204.net reaction networks converge to ${Y}_{{\rm{e}},{\rm{c}}}$ ≈ 0.43. However, models using the smaller approx22.net and mesa_79.net reaction networks give values of ${Y}_{{\rm{e}}}$ ≈ 0.41 and ${Y}_{{\rm{e}}}$ ≈ 0.46, respectively. A similar trend is found for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models. In the case of the ${M}_{\mathrm{ZAMS}}=25$ and 30 ${M}_{\odot }$ models, the approx22.net and mesa_79.net reaction networks show less disagreement in ${Y}_{{\rm{e}},{\rm{c}}}$ with the values found for the larger networks. The final ${Y}_{{\rm{e}},{\rm{c}}}$ is set by the final composition of the core, in the approx22.net, this is a mixture of ⁵⁶Fe and ⁶⁰Cr, the mesa_79.net are $\approx 95 \%$ ⁵⁶Fe (independent of mass), while in the larger networks the core becomes a mixture of iron-group elements. This is because the mesa_79.net lacks ⁵⁵Fe and ${}^{48-51}{\rm{V}}$ and ${}^{51-56}$Cr isotopes, which are present in the larger nuclear networks. These additional isotopes provide alternate decay routes for the ⁵⁶Fe in the core.

There is also considerable variation in the shape of the ${Y}_{{\rm{e}}}$ step at the edge of the iron core. This is not a well-defined jump in the ${Y}_{{\rm{e}}}$ value, but can show a number of substeps. For instance, for the 15 ${M}_{\odot }$ models in Figure 21, the approx22.net has a single large jump at 1.62 ${M}_{\odot }$, while with the larger nets the jump is at 1.3 ${M}_{\odot }$, but shows substructure in the ${Y}_{{\rm{e}}}$ values.

Figure 22 shows the oxygen mass fraction profile at core collapse for a subset of the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models. For each profile, all isotopes of oxygen in a particular network are summed to give the total oxygen mass fraction as a function of mass coordinate (see Figure 1 for the different isotopes of oxygen included for each network). The left panel plots the X(O) profile as a function of the number of isotopes in the network for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models without mass-loss at a mass resolution of ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$. The right panel plots the X(O) profile as a function of mass resolution for the ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models without mass-loss that use the mesa_160.net reaction network (this network follows ${}^{14-18}$O).

Zoom In Zoom Out Reset image size

The upper mass boundary of the O shell is well constrained at 4.4 ${M}_{\odot }$ for the different networks at fixed ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$ models. While the lower mass boundary varies over 0.5 ${M}_{\odot }$, the total mass of oxygen only varies between 1.9 and 2.0 ${M}_{\odot }$. For the models with a fixed reaction network and variable ${\rm{\Delta }}{M}_{\max },$ the upper mass boundary varies over $\approx 0.2$ ${M}_{\odot }$, similar to the spread in the lower mass boundary. Most models show a top-hat shape that is due to strong convection occurring above the shell, while the ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$, mesa_160.net shows a different behavior because an oxygen-burning wave propagates outwards from the base of the O shell. The almost step-like feature seen in the lower boundaries is due to a convective Si-burning shell, abutting the base of the O-burning shell and the top of the Fe core. This Si-burning shell only partially burns the silicon and the oxygen before core collapse.

Similarly, in Figure 23 we show mass fraction profiles for X(Si) for ${M}_{\mathrm{ZAMS}}=20$ ${M}_{\odot }$ models without mass-loss. This region shows considerably more variation than the O shells in Figure 22. The lower mass boundary of the silicon distribution denotes the edge of the iron core, while the upper mass boundary denotes the edge between the Si shell and the base of the O shell.

Zoom In Zoom Out Reset image size

For a fixed mass resolution of ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$, there are three groupings in the left panel of Figure 23: the approx22.net; the mesa_79.net, mesa_127.net, and mesa_160.net; and the mesa_204.net. The approx22.net model has the largest iron core because it has the most energetic shell Si-burning. This pushes the silicon shell outwards and allows the shell to burn over a larger oxygen-rich region. The larger networks are generally cooler, so that the Si-burning shell does not extend as far outwards. The sudden drop in the silicon mass fraction at the upper boundary, seen in the mesa_79.net, mesa_127.net, and mesa_160.net models, is caused by the convective lower mass boundary of the O shell, which mixes silicon outwards. The mesa_204.net shows a case when the iron core abuts the O shell, leaving no distinct Si shell.

For a fixed mesa_160.net reaction network, Figure 23 shows a trend for the Si shell to move inwards and concentrate at lower mass coordinates as the resolution increases. This is because the temperature and density of the Si shell decrease as the resolution increases. This slows the conversion of silicon into iron-group elements. For models with ${\rm{\Delta }}{M}_{\max }$ ≤ 0.02 ${M}_{\odot }$, the lower mass boundary of the silicon is sharply defined at $\approx 1.6$ ${M}_{\odot }$, with the core being composed primarily of ⁵⁴Fe up to the edge of the iron core. For the ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ and ${\rm{\Delta }}{M}_{\max }=0.5$ ${M}_{\odot }$ models, the lower mass boundary is more inclined because of the increased presence of ⁵⁶Ni. In the ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ case, this becomes a distinct ⁵⁶Ni shell between the iron core, which is mostly ⁵⁴Fe, and the silicon shell. In the ${\rm{\Delta }}{M}_{\max }=0.05$ ${M}_{\odot }$ case ⁵⁴Fe and ⁵⁶Ni are approximately equal components at ≈25% at the interface of the ⁵⁴Fe core and the Si shell.

Figure 24 shows the total number of mass zones at core collapse. As ${\rm{\Delta }}{M}_{\max }$ increases from ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$ to 0.005 ${M}_{\odot }$, there is an increase of a factor 4 in the number of zones, which is much less than the increase of a factor 20 in ${\rm{\Delta }}{M}_{\max }.$ This is because ${\rm{\Delta }}{M}_{\max }$ only provides an upper limit on the mass of a cell. For instance, a ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ model with ${\rm{\Delta }}{M}_{\max }=0.1$ ${M}_{\odot }$, assuming all cells have the same ${\rm{\Delta }}{M}_{\max },$ needs only 150 spatial zones—much fewer than the 1200 zones present in our models. On average, we find that the models with mass-loss need slightly fewer zones at core collapse. This is because they have a star with a lower mass. However, the zones are unevenly distributed because they are concentrated near the core and at the location where large changes in stellar properties occur, and therefore the scaling is not linear with stellar mass. There are a number of other MESA controls that affect the spatial resolution, like ${\delta }_{\mathrm{mesh}}$, which contributes to the total number of zones and to their relative distribution. We recommend taking a critical view of the spatial distribution of zones in stellar models and considering the effect of increased resolution on their results.

Zoom In Zoom Out Reset image size

Table 1 summarizes the pre-SN model properties measured either during the evolution of the model or at core collapse. We show the median value over all models with upper and lower bounds. The ${\mathrm{He}}_{\mathrm{core}},{{\rm{C}}}_{\mathrm{core}},{{\rm{O}}}_{\mathrm{core}},{\mathrm{Si}}_{\mathrm{core}}$ summarize the core masses at the point of the fuel ignition (see Figures 6, 9, 14, and 17, for an illustration of their spread). Similarly, the ${Y}_{{\rm{e}},{\rm{c}}}$ shows the electron fraction at the center of the star measured at the same time as the core masses. The electron fraction does not drop much bellow 0.5 until after oxygen-burning. Once silicon burning commences, the ${Y}_{{\rm{e}},{\rm{c}}}$ drops to $0.48\mbox{--}0.49$, with stars with higher initial mass having larger ${Y}_{{\rm{e}},{\rm{c}}}$. At core collapse, ${Y}_{{\rm{e}},{\rm{c}}}$ can drop to ${Y}_{{\rm{e}},{\rm{c}}}=0.43\mbox{--}0.44$, with higher initial masses having higher median values, but the upper/lower bounds overlap between all the masses.

Table 1 lists the final masses (${M}_{\mathrm{final}}$) at core collapse, but dominated by the mass at the TAMS, for models with mass-loss. Similarly to final ages, most of the variation in the final masses occurs as mass resolution increases (see Section 3.2). For the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ models, the total variation in the final mass is ≃2%. As the ZAMS mass increases, this variation increases to ≃10% for the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ models. In general, as the resolution increases, the final mass decreases, consistent with the findings from Section 3.1, where the higher resolution models live longer because of the change in behavior of the ICZ and SCZ.

The τ values in Table 1 denote the time between the depletion of the previous fuel to the core ignition of the next. As we progress between each major fuel source, the lifetime of the star decreases rapidly, as expected. As the initial mass increases, the lifetime for each fuel to burn decreases, and the fraction of the star lifetime that is spent burning each fuel (other than hydrogen) decreases. For instance, the MS lifetime of a 15 ${M}_{\odot }$ star is twice that of a 30 ${M}_{\odot }$, but eight times that during the final silicon shell burning. In addition, those stars that lose mass live longer than their companions that do not lose mass because of the slower evolution experienced by a lower mass star. ${\mathrm{He}}_{\mathrm{shell}},{{\rm{C}}}_{\mathrm{shell}},{{\rm{O}}}_{\mathrm{shell}},{\mathrm{Si}}_{\mathrm{shell}}$ denote the locations of each fuel shell at core collapse, while the ${\mathrm{Fe}}_{\mathrm{core}}$ is the iron-core mass at core collapse (see Figure 20). The shell locations increase as the initial mass increases, with mass-losing stars having lower shell masses. The ${\mathrm{Si}}_{\mathrm{shell}}$ lower bound of 0.0 does not denote that there is a distribution to 0.0, but that a few models have no silicon shell (see Figure 23) while the other models sit near the median value.

Last, we provide a summary of the surface of the star at core collapse. ${H}_{\mathrm{env}}$ is the mass of H-rich envelope, ${M}_{\mathrm{final}}$ is the final mass, and $R/{R}_{\odot },{T}_{\mathrm{eff}},L/{L}_{\odot }$ provide the surface radius, temperature, and luminosity. As the initial mass increases, the temperature and luminosity of the models increase, but the surface radius is maximum for the 20–25 ${M}_{\odot }$ models before decreasing for the 30 ${M}_{\odot }$ models. Mass loss does not significantly effect the final radius, temperature, or luminosity either, even though they can be significantly lower and have much lower ${H}_{\mathrm{env}}$ values. The final compactness parameter is given by

$$\begin{eqnarray}&&\xi =\displaystyle \frac{(M/\,{M}_{\odot })}{{R}_{M=2.5}/1000\,\mathrm{km}},\end{eqnarray} \tag{ 5 }$$

where $M=2.5$ ${M}_{\odot }$ and the radius is measured (in km) at the radius where $M=2.5$ ${M}_{\odot }$ (O'Connor & Ott 2011). In principle, it should be measured at the bounce, but Sukhbold & Woosley (2014) showed that no substantial accuracy is lost if it is measured in the pre-SN model. We see a steady increase in the final value as the initial mass increases, with mass-losing stars having slightly lower values. All the parameters in Table 1 show considerable spread in their values, which is due to changes in the spatial resolution and to the choice of nuclear network.

Figure 25 shows the stable isotopes from hydrogen to zinc for the pre-SN ${M}_{\mathrm{ZAMS}}=15$, 20, 25, and 30 ${M}_{\odot }$ models with mass-loss, ${\rm{\Delta }}{M}_{\max }=0.01$ ${M}_{\odot }$, and the mesa_204.net reaction network. These pre-SN yields are for the stable isotopes outside the Fe core. The x-axis is the atomic mass number. The y-axis is the logarithmic ratio of the model mass fraction to the Grevesse & Sauval (1998) solar mass fraction. The most abundant isotope of a given element is marked by an asterisk, and isotopes of the same element are connected by solid lines.

Zoom In Zoom Out Reset image size

Isotopes in Figure 25 heavier that about calcium will be strongly impacted by the explosion. Up to calcium, however, it generally makes little difference whether the pre-SN or the exploded yields are used (e.g., Timmes et al. 1995). In a forthcoming effort we anticipate exploding the pre-SN models examined in this paper. With this limitation of interpreting pre-SN yields in mind, the overall average production factor of $\approx 20$ for the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ model and rising to $\approx 100$ for the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ model are commensurate with existing pre-SN yields (e.g., Woosley & Weaver 1995; Limongi & Chieffi 2003; Chieffi & Limongi 2013) and the production factors needed by galactic chemical evolution models (e.g., Gibson et al. 2003) to reproduce the solar composition. Part of the spread below calcium is a consequence of uncertain physics in the MESA models: the treatment of convective, semiconvective, and overshoot mixing, the parameterization of mass-loss, nuclear reaction rates, and residual uncertainty in the measured solar abundances. Post-SN, the yields for isotopes with $A\gt 40$ will also sensitively depend on the choice of explosion dynamics (Paxton et al. 2015).

The light isotopes ${}^{\mathrm{6,7}}$Li, ⁹Be, and ${}^{\mathrm{10,11}}{\rm{B}}$ are not plotted in Figure 25 as the primary source of these isotopes is commonly taken to be Big-Bang nucleosynthesis (e.g., Coc et al. 2012; Cyburt et al. 2016), spallation of CNO nuclei in the interstellar medium by cosmic rays (e.g., Reeves et al. 1970; Olive & Schramm 1992; Fields et al. 2000; Ramaty et al. 2000; Prantzos 2012), or the ν-process during the explosion (e.g., Woosley et al. 1990; Balasi et al. 2015).

The isotopes ${}^{\mathrm{12,13}}{\rm{C}}$ and ¹⁴N are underproduced in Figure 25 by $\approx 3$ relative to the average overproduction factor for the ${M}_{\mathrm{ZAMS}}=15$ ${M}_{\odot }$ model and by $\approx 20$ relative to the average overproduction factor for the ${M}_{\mathrm{ZAMS}}=30$ ${M}_{\odot }$ model. This is consistent with the bulk of these isotopes being produced during CNO processing and dredged-up material from helium-shell flashes in intermediate- and low-mass stars (Iben & Truran 1978; Renzini & Voli 1981; Doherty et al. 2014; Karakas & Lugaro 2016).

The solar ²⁹Si abundance is $\approx 1.4$ times larger than the solar ³⁰Si abundance, in concordance with observations of some individual stars (Peng et al. 2013) and pre-solar silicon-carbide grains (e.g., Hoppe et al. 2010), but it disagrees with common models of 1D core-collapse supernova, SNe Ia, and AGB stars (e.g., Timmes & Clayton 1996; Lugaro et al. 1999; Lewis et al. 2013; Wasserburg et al. 2015). It is curious that the ³⁰Si abundance, normalized to solar in Figure 25, is larger than the ²⁹Si abundance in the ${M}_{\mathrm{ZAMS}}=20$ and 30 ${M}_{\odot }$ models, but traditionally smaller in the ${M}_{\mathrm{ZAMS}}=15$ and 25 ${M}_{\odot }$ models.

The isotope ⁴⁰K has the largest production factor in the ${M}_{\mathrm{ZAMS}}=25$ and 30 ${M}_{\odot }$ models of Figure 25. This isotope is also a special case as it is radioactive, but the half-life is long enough (1.248 × 10⁹ yr, Malonda & Carles 2002; Kossert & Günther 2004) that it is included in compilations of the solar composition. The overproduction factor may reflect an uncertain solar abundance.