CONVECTIVE PROPERTIES OF ROTATING TWO-DIMENSIONAL CORE-COLLAPSE SUPERNOVA PROGENITORS

The final years in the lives of massive stars are characterized by vigorous convective shell burning, hydrodynamic and convective instabilities and, in many cases, episodic mass-loss events (Meakin & Arnett 2007; Quataert & Shiode 2012; Arnett et al. 2014; Shiode & Quataert 2014; Smith & Arnett 2014) that change their three-dimensional structure and the initial conditions (ICs) for the core-collapse supernova (CCSN) explosion. In addition, rotation and magnetic fields may further complicate the core-collapse process in a non-linear fashion.

The advanced stages of burning in massive stars have long been studied with 1D approximations (Clayton 1984; Woosley et al. 2002). For convection, the mixing-length theory (MLT; Böhm-Vitense 1958) remains the technique used and implemented in stellar evolution codes with options to choose the associated efficiency parameter, α_MLT, and the condition determining where convection becomes active by using the Schwarzschild or the Ledoux criterion accounting for suppression due to compositional gradients. Parametrizations for 1D models of convective overshoot, semiconvection, and thermohaline mixing are frequently adopted, as is the 1D shellular approximation to treat rotation (Zahn 1992; Meynet & Maeder 1997).

Advances in numerical algorithms, hydrodynamic software instruments, and computing power have allowed multi-dimensional studies of stellar convection that can assess the fidelity of 1D treatments of mixing. The properties of convective oxygen shell-burning in the progenitor of SN 1987A prior to collapse were studied by Bazán & Arnett (1998) in 2D hydrodynamics simulations. One of the implications of this study was the potential for post-explosion mixing of radioactive ⁵⁶Ni throughout the SN envelope. The first three-dimensional (3D) simulations of massive star oxygen shell convection were presented by Meakin (2006) and Meakin & Arnett (2007).This work showed that the boundaries between non-convective and convective regions are not stationary, as standard MLT theory predicted, but dynamical and the source of gravity waves. Arnett & Meakin (2011) followed up with 2D simulations of simultaneously active C, Ne, O, and Si burning shells for a 23 M_⊙ progenitor 1 hr prior to core-collapse, finding significant departures from spherical symmetry and strong dynamical interactions between shells. The 3D hydrodynamics simulations together with developments in mathematical methods indicate that MLT needs significant revision especially for late stages of stellar evolution (Arnett et al. 2015; Gilkis & Soker 2015).

The apparent departures from spherical symmetry that arise in the velocity fields of convective shells prior to iron core-collapse have a qualitative and quantitative impact on its susceptibility to explosion (Couch & Ott 2013, 2015; Müller & Janka 2015). Couch & Ott (2013, 2015) have shown that successful explosions occur for models whose ICs include velocity field perturbations due to convective burning, in contrast to the same models without these asphericities. This is due to non-radial motions in the accretion flow exciting stronger post-shock turbulence that aids shock revival. Full 3D simulations of the final minutes of iron core growth and collapse also suggest that non-spherical progenitor structure should not be ignored (Couch et al. 2015). One might characterize convection seen in 2D or 3D hydrodynamic simulations and then map realizations of convective velocity fields as ICs for simulations of CCSNe. Efforts to do so include the spherical Fourier–Bessel decomposition (Fernández et al. 2014) and nonlocal and time-dependent averaging (Meakin & Arnett 2007; Arnett et al. 2009; Viallet et al. 2013). In our analysis we use vector spherical harmonic (VSH) decomposition of the momentum density field as presented by Chatzopoulos et al. (2014).

The stochastic nature of 1D convection algorithms has a profound impact on the outcome of massive stellar evolution itself. For stars of very similar mass and metallicity the end points can be dramatically different. These give remnant masses and explosion properties which depend strongly on the pre-SN stellar structure and exhibit large variability even in narrow intervals of Zero Age Main Sequence (ZAMS) mass (Ugliano et al. 2012; Sukhbold & Woosley 2014).

Additional uncertainty arises from the dynamical boundaries between convective zones (Meakin & Arnett 2007), which are a source of gravity waves and might lead to episodic mass-loss events shortly before the SN explosion (Quataert & Shiode 2012; Shiode & Quataert 2014; Smith & Arnett 2014). The circumstellar (CS) environment that is formed around pre-SN stars due to this mass-loss history can have a significant impact on the radiative properties of the resulting explosion.

Rapid rotation has been observed for many massive stars (Hunter et al. 2008; Vink et al. 2010; Dufton et al. 2011; Almeida et al. 2015) and can affect their pre-SN internal structure and composition via instabilities that alter the efficiency of angular momentum transport and the rate of chemical mixing (Brott et al. 2011a, 2011b; Ekström et al. 2012; see Maeder & Meynet 2012 for a review). The effects of rotation on convection in the extended envelopes of red giant stars has been studied with 3D simulations by Brun & Palacios (2009), who found the properties of turbulent convection to be sensitive to the rotation rate. Currently, the effects of rotation on the convective properties of massive, pre-SN stars and their implications for CCSNe have not been thoroughly investigated. In this paper we present 2D simulations of rotating pre-SN stars during convective C- and O-shell burning and Si-shell burning, and use the VSH method to quantify the behavior of convective velocity flows as a function of internal rotation rate.

The paper is organized as follows. In Section 2 we present the stellar evolution and 2D hydrodynamics simulations of a 20 M_⊙ pre-SN star during convective shell burning of C, O, and Si. In Section 3 we apply the method of VSH to decompose the convective velocity fields and obtain power spectra in order to characterize the convective properties of models of different rotation rates. Finally, in Section 4 we discuss our conclusions and implications for CCSNe.

Our analysis of 2D convection in rotation pre-SN progenitors has three distinct steps. First, we evolve models of a star with ZAMS mass of 20 M_⊙, solar metallicity, but different rotation rates using version 7503 of the stellar evolution code Modules for Experiments in Stellar Astrophysics (MESA; Paxton et al. 2011, 2013, 2015). Subsequently, MESA 1D profiles are extracted days prior to core-collapse, during core Si- and shell C- and O-burning, and ∼1 hr prior to collapse during Si-shell burning. These profiles are then mapped into 2D using the adaptive mesh refinement (AMR) multi-physics hydrodynamics code FLASH version 4.3 (Fryxell et al. 2000; Dubey et al. 2012) including rotation perpendicular to the plane of the simulation ("2.5D" approach). The 2D simulations are run for >3 convective turnover timescales, long enough to diminish the effects of the initial dynamical transient resulting from mapping 1D "convective" profiles to a multi-dimensional hydrodynamic grid. Lastly, FLASH simulation output at three different times is extracted and post-processed using the VSH method (Chatzopoulos et al. 2014) to obtain the power spectra of the convective motions and the energy cascade.

2.1. MESA Pre-SN Evolution

Our MESA models use initial rigid rotation rates of zero ("no-rot") and 50% of the critical Keplerial value ("rot-ST") at the equator. The initial rotation profiles impose rigid-body rotation (constant angular velocity) when the model first lands on the ZAMS. The rotating model assumes the transport of angular momentum and chemical mixing via the Spruit–Taylor mechanism (ST; Spruit 1999, 2002). Models without the inclusion of the ST mechanism are expected to produce different internal rotation structures for CCSNe progenitors since the coupling between the core and the envelope will be much weaker (Chieffi & Limongi 2013).

For MLT convection in MESA we adopt the Schwarzschild criterion and α_MLT = 1.6. We use an automatically extending nuclear reaction network starting from a basic 8-isotope network and reaching a 21-isotope network by the end of the run. We note that ∼100–150 isotopes are required to accurately represent core neutronization and neutrino cooling (Arnett & Meakin 2011; Arnett et al. 2015). The "Helmholtz" equation of state (EOS; Timmes & Swesty 2000) is used. Standard mass-loss prescriptions appropriate for massive stars are adopted (Vink et al. 2001; Glebbeek et al. 2009). After performing a resolution study, we choose a spatial resolution parameter (mesh_delta_coeff in MESA terminology) equal to 0.5 and a temporal resolution factor (varcontrol_target) of 10⁻³ where good convergence (≃10⁻² level) in terms of final carbon-oxygen core mass and iron core mass is achieved. The chosen grid resolution resulted in final output models with 1200–1600 Lagrangian zones.

Figure 1 shows the evolution of iron core mass (M_Fe) and central density (ρ_c), temperature (T_c), and average charge per baryon (Y_e) for both the "no-rot" and "rot-ST" models for the last 50,000 s of evolution and the formation of the iron core. The vertical dashed lines show the stage where we extracted the MESA models corresponding to Si-shell burning, at ∼1.3 hr prior to core-collapse, when the peak of _nuc profile reached the maximum value ((3–4) × 10¹⁸ erg g⁻¹ s⁻¹). The build-up of the Fe-core up to the Chandrasekhar mass is smoother for the "rot-ST" model predominantly due to the effects of ST and rotational mixing.

Zoom In Zoom Out Reset image size

Figures 2 and 3 show the distributions of nuclear energy generation rate (_nuc), composition (X_i), and temperature focused in the convective, shell-burning regions. The MESA models we calculate do not exhibit a phase of strongly evident shell Ne-burning, in contrast to the more massive (23 M_⊙) model used in Arnett & Meakin (2011). The vertical dashed lines indicate the inner and outer radial boundaries (${R}_{{\rm{in,sh}}}$ and ${R}_{{\rm{out,sh}}}$) that were chosen for the convective shells to be analyzed with the VSH method (see Section 3 for details on how their values were determined). The effects of enhanced mixing in the "rot-ST" model are clearly illustrated in the composition panels.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To further isolate the effects of rotation, we impose an artifical rotational velocity profile on the "no-rot" models:

$$\begin{eqnarray}&&{v}_{{\rm{rot}}}(r)=\displaystyle \frac{r{\rm{\Omega }}}{1+{(r/A)}^{2}},\end{eqnarray} \tag{ 1 }$$

where r is the radial (spherical) coordinate, Ω is the angular velocity, and A is the characteristic radius where the rotational velocity peaks as in Chatzopoulos et al. (2013). We denote models that use this rotation law as "rot-2." We carefully selected the values of A and Ω in order to capture a few full rotations (4–8) for the convective shell material within the simulated timescales in FLASH and assess the dynamical effects of the centrifugal force. For C- & O-shell burning we picked A = 3.46 × 10⁸ cm and Ω = 0.1 s⁻¹ and for Si-shell burning A = 2.2 × 10⁸ cm and Ω = 0.1 s⁻¹ respectively. Figure 4 shows the rotational velocity profiles for all models mapped in FLASH.

Zoom In Zoom Out Reset image size

2.2. 2D FLASH Hydrodynamics

The six 1D MESA models imported into FLASH are listed in Table 1. The C- & O-shell burning FLASH models are run for 1000 s (simulated time), and the Si-shell burning models for 500 s. The dynamical effects of rotation and angular momentum conservation are handled by the unsplit piecewise parabolic method (PPM) hydrodynamic solver in FLASH (Lee et al. 2009). We note that the newest implementation of the unsplit solver in FLASH handles species advection in a way that is nearly identical to consistent multifliud advection scheme methods (Plewa & Müller 1999). The 1D rotational velocity profiles were mapped on the 2D grid assuming "shellular" rotation as vectors with direction perpendicular to the R–z plane of the simulation and the rotation axis coincident with the polar axis. In this approximation the angular velocity, Ω(r), is constant for a particular spherical shell or on equipotential surfaces (Zahn 1992; Meynet & Maeder 1997). More specifically, each grid cell was given a rotational velocity ${v}_{{\rm{rot}}}={\rm{\Omega }}(r)R$, where ${\rm{\Omega }}(r)={v}_{{\rm{rot}}}(r)/r$, with ${v}_{{\rm{rot}}}(r)$ being MESA 1D rotational velocity and r the spherical radial coordinate ($r={({x}^{2}+{y}^{2})}^{0.5}$). We did not map and follow the evolution of the magnetic field in the FLASH simulations. We estimate that the plasma β parameter (ratio of the thermal to the magnetic pressure) in the convective regions simulated is much larger than 10⁸, implying that the effects of the magnetic pressure to the structure and size of the convective shells are minimal, particularly for the short timescales simulated.

Table 1.  Properties of MESA Models Mapped in FLASH

Model	Rin,sh (108 cm)	${R}_{{\rm{out,sh}}}$ (108 cm)	${v}_{{\rm{rot,sh}}}$ (km s−1)	${v}_{{\rm{conv}}}$ (km s−1)	τrot (s)	t1 (s)	t2 (s)	t3 (s)
			C- & O-shell Burning
no-rot	3.46	36.94	0.	136	⋯	545	749	1000
rot-ST	3.46	36.94	5.5	368	3957	545	749	1000
rot-2	3.46	36.94	185.6	142	126	545	749	1000
			Si-shell Burning
no-rot	1.64	36.00	0.	445	⋯	300	400	500
rot-ST	1.64	36.00	7.1	362	1505	300	400	500
rot-2	1.64	36.00	110.0	449	138	300	400	500

Note. "no-rot": non-rotating model. "rot-ST": rotating model that includes the magnetic field effects of the Spruit-Taylor dynamo. "rot-2": rotating model produced by introducing a rotational velocity profile to "no-rot" upon mapping to FLASH. The rotational timescale, τ_rot, corresponds to the rotational period (one full revolution around the rotation axis) for the given rotational speeds in the center of the convective shell.

Download table as:  ASCIITypeset image

The transition from 1D MESA to 2D FLASH was smooth in terms of the important physics involved; a hardcoded 21-isotope' network identical to the one used in MESA was employed (Couch et al. 2015) and the "Helmholtz" EOS was used. The main inconsistency is the treatment of convection in transition from 1D to 2D: in 1D the convective energy transport and cascade is treated via approximate MLT prescriptions, while in 2D and 3D convective flow naturally develops in the fluid within unstable regions. The mapping from 1D to 2D triggers an initial dynamical transient that lasts for  ∼100–300 s before a quasi-steady state is re-established. For our VSH analysis presented in Section 3 we use the FLASH output ("snapshots") at three different times, well after (>300 s for C- & O-shell burning and >200 s for Si-shell burning) the initial transient has transversed the computational domain. Table 1 also details the properties of the convective shells and simulation output.

All 2D FLASH simulations were run on the Texas Advanced Computing Center Stampede supercomputer. The size of all simulation domains was chosen to be 10¹⁰ cm including both the core and the convective shells of all models. For all of our simulations, we chose a reflecting inner boundary condition and a "diode" outer boundary condition allowing the flow of mass outside of the simulation box but prohibiting mass from entering it. The "diode" option in FLASH enforces a zero-gradient condition with a switch that depends on the sign of the velocity in the last zone interior to the simulation domain. As such, it can change with time. A more detailed discussion of the available boundary condition options in FLASH is provided in the FLASH user's guide. The maximum resolution chosen was 9 km, corresponding to convergence in total energy and mass at a ∼10⁻⁷ level over the course of the simulation. At that resolution we are able to resolve Eddy sizes covering a considerable range of the turbulent energy cascade (from ∼30,000 km down to ∼10 km). We should caution, however, that the turbulent energy cascade in 2D is inverted (Porter & Woodward 1994; Takiwaki et al. 2014; Couch & Ott 2015) and full 3D treatment is required to accuretly reproduce it. The morphology of the flow changes significantly from 2D to 3D, and the velocity scale is moderately higher in 2D (Meakin & Arnett 2007). We consider our simulations as an initial, exploratory step. Each simulation was run on 128 cores and the wallclock time ranged from 14 to 22 hr for a total of >12,000 core-hours used. Figures 5 and 6 show the O¹⁶ or Si²⁸ mass fraction and the 2D (x and y) velocity magnitude at the end of each simulation. The prevalence of convection is apparent in all cases with characteristic velocities reaching ∼150 km s⁻¹ for C- & O-shell burning and ∼450 km s⁻¹ for Si-shell burning in the "no-rot" and "rot-2" cases. It is noteworthy that by the end of all simulations the velocities in the "rot-ST" cases were much higher for C- & O-shell burning (∼370 km s⁻¹) but slightly lower for Si-shell burning (∼360 km s⁻¹).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The convective elements seen in the velocity magnitude panels span a range of sizes, with dominant large-scale motions that cover more than half of the size of the convective shells, as well as smaller scale vortices which are just visible (for example in the "rot-ST" C & O-burning shell). Convection is established ∼200–300 s after the start of each simulation, after the initial transient exits the computational domain. The evolution past 200 s shows large scale convective currents that cover the entire shell, breaking down to smaller structures and high velocity vortices interacting with each other. Mild mixing is also seen, with instabilities developing at the interfaces of the convective shells. This mixing can be due to the process of turbulent entrainment also seen by Meakin & Arnett (2007). In all simulations, convective elements interact with the inner (core) boundary and break into smaller structures that then subsequently reunite while rising upwards. Artifical flows are seen near the axis of the simulation (∼6°), a common issue of 2D cylindrical treatment. The radial distribution of angular momentum for the rotating models does not change significantly by the end of the FLASH simulations. We attribute this to the relatively short simulation timescale as compared to the time required to signficiantly redistribute angular momentum in the stellar interior (∼1 Myr for ST, which is the most efficient angular momentum transport mechanism in action).

Output from the FLASH simulation is taken at three instances, t₁, t₂, and t₃ (the end of the simulation); all are given in Table 1. A total of 18 snapshots for all 6 cases are post-processed using the VSH analysis implemented in the code. These timescales are chosen to represent different evolution phases during shell convection, well after the initial dynamic transient and several convective turnover timescales after that (≃10–100 for C- & O-shell burning and ≃10–50 for Si-shell burning).

The goal of the VSH analysis is to decompose the momentum density field within the selected shell into radial (A_nlm) and solenoidal (B_nlm and C_nlm) modes and then calculate the power spectrum for each of the modes to determine the global properties of 2D or 3D fluid motion. In our case we analyze a 2D momentum density field; therefore the solenoidal B_nlm modes are irrelevant because they cancel out and will not be discussed further. Also, because our simulations were conducted in a full domain possessing reflection symmetry, we expect the presence of an odd–even effect in the power spectra, a feature inherent to the 2D treatment. The first step is to determine the radial limits of the convective regions where VSH is applied. For this, volume-weighted radial momentum density profiles are calculated (see Equation (22) of Chatzopoulos et al. 2014). Locations where the radial component reaches a minimum are used as our final choice for R_in,sh and R_out,sh; these are also presented in Table 1.

The next step is to declare the maximum radial and angular resolution for the VSH components, both of which can be expressed as a length scale, λ_r. The chosen resolution scale then determines the number of radial (n) and angular (l and m) modes required for the expansion in the VSH components, which are given by the following formulae:

$$\begin{eqnarray}&&{n}_{\mathrm{max}}=\displaystyle \frac{2({R}_{2}-{R}_{1})}{{\lambda }_{{\rm{r}}}},\end{eqnarray} \tag{ 2 }$$

and

$$\begin{eqnarray}&&{l}_{\mathrm{max}}=\displaystyle \frac{\pi ({R}_{1}+{R}_{2})}{2{\lambda }_{{\rm{r}}}},\end{eqnarray} \tag{ 3 }$$

where the total number of modes in 2D momentum density field decomposition is

$$\begin{eqnarray}&&{N}_{\mathrm{total}}=({n}_{\mathrm{max}}+1)({l}_{\mathrm{max}}+1).\end{eqnarray} \tag{ 4 }$$

For the C- & O-shell burning models we choose ${\lambda }_{{\rm{r}}}=1.585\times {10}^{8}$ cm, corresponding to l_max = 40 and n_max = 42, while for the Si-shell burning models we choose λ_r = 1.495 × 10⁸ cm, l_max = 40, and n_max = 46. We choose to truncate radial n > 20 modes because we find that their contribution to the total power is minimal (10⁻⁵–10⁻⁴ level) and we can thus reduce computation time. With these choices, a total of 861 modes were calculated for each of the 18 FLASH snapshots. For the purposes of our study, we calculate reduced VSH spectra by first summing over all the "phase" (m) components. Then, we calculate the reduced angular (l) and radial n power spectra by summing over either n or l, respectively, for instance ${\alpha }_{l}^{\prime }\equiv {\sum }_{n}{\alpha }_{{nl}}$, ${\alpha }_{n}^{\prime }\equiv {\sum }_{l}{\alpha }_{{nl}}$, and similarly for the solenoidal modes. In our presentation of VSH power spectra later, we will simply refer to the radial and solenoidal modes as A and C, respectively (for the exact definitions of A and C consult Equations (6)–(8) of Chatzopoulos et al. 2014).

3.1. Shell O- and C-burning

Figure 7 shows the reduced VSH power spectra for all models undergoing C- & O-shell burning. The evolution of the total power, summed over all components, is shown in the top row of Figure 9. For all cases the bulk of convective power is concentrated in large scales (l < 10, n < 5), with the peak values implying angular scales of ∼1.6–3.2 × 10⁹ cm and radial scales of approximately the same range, revealing a nearly circular characteristic shape for the convective eddies. These eddies can be comparable in size to the shell itself. For some spectra, a secondary peak of power is observed at smaller values of l (∼15–20), indicating that a small fraction of the total power is possessed by smaller scales (∼4.3 × 10⁸ cm). The overall slope of the spectra remains nearly fixed over time. As expected for a convective velocity field confined in a shell, the power in the solenoidal modes is clearly dominant over that in the radial modes.

Zoom In Zoom Out Reset image size

Regardless of the degree of rotation, the total power in dominant solenoidal modes declines over the course of the simulations as the convective energy cascade settles. A comparison between the "no-rot" and "rot-2" case shows that the addition of rotation in an otherwise equal stellar structure has little effect on the global properties of convection. A small reduction of total power is apparent, possibly related to the effects of centrifugal forces and the existence of an extra degree of freedom (movement perpendicular to the simulation domain). On the other hand, the algorithm for calculating the effects of rotation during stellar evolution leads to a pre-SN star with clearly stronger convection during the C- & O-shell burning phase ("rot-ST" model, upper middle panel in Figure 9). Indeed, the final peak 2D-velocity magnitude in the "rot-ST" model is more than double those of the "no-rot" and "rot-2" models. We discuss this in more detail in Section 4.

3.2. Shell Si-burning

The VSH power spectra for the cases of shell Si-burning are shown in Figure 8. The evolution of the total power is shown in the bottom row of Figure 9. The reduced spectra reveal that large scales also dominate over smaller scales during shell Si-burning, with l peaking in the range 4–6 throughout the evolution. The corresponding angular length scales range from 1.5 × 10⁹ cm to $2.0\times {10}^{9}$ cm. The radial scales are also in the same range implying a nearly circular shape for the convective elements, which are at about half the size of the convective shell. Secondary peaks occur at smaller scales throughout the VSH spectral evolution (∼10⁹ cm). As in the case of C- & O-shell burning, the spectral slopes remain consistent over time and solenoidal models dominate radial models by a factor of ∼10,000 for l < 5.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The evolution of the total power in the solenoidal (C) modes shows that convection during Si-shell burning is about 10 times stronger than convection during C- & O-shell burning, as expected from the higher rates of local energy generation. The total power does not seem to vary significantly over time for all rotation rates. A small reduction of total power in the C-modes is seen by the end of the simulations for the "rot-ST" and "rot-2" cases. For models initially identical, modulo the inclusion of rotation ("no-rot" versus "rot-2"), this effect may be due to the dynamical impact of centrifugal forces and the extra degree of freedom, as argued for the case of C- & O-shell convection. We do, however observe a qualitative difference between the two stages of shell burning: during Si-shell burning, the "rot-ST" model exhibits nearly identical and, at late times, somewhat lower convective power as compared to the "no-rot" model, while the opposite behavior was observed for C- & O-shell burning. We return to this issue in the following section.

In this paper we explored the effects of rotation on the convective properties of a 20 M_⊙ pre-SN progenitor star, more specifically during C- & O-shell burning and Si-shell burning in the hours to months prior to core collapse. We studied rotation both by including it self-consistently during the evolution of the star but also by imposing it on an otherwise identical, non-rotating model in order to better isolate its effects. We run 2D FLASH simulations with the rotational velocity field mapped and pointing inwards perpendicular to the plane of the simulation domain. We employed VSH decomposition of the momentum density field to characterize convection at different times during the hydrodynamic simulations. Our initial exploration of the influence of rotation on the properties of convection deep in the core of a massive star can (and should) be improved by using 3D simulations, where the turbulent cascade goes from large scales to small scales.

We find that, regardless of the degree of rotation, the characteristic size of the convective elements is somewhat larger during the C- & O-shell burning phase compared to the Si-shell burning phase, while the characteristic convective velocities are more than twice as large during Si-shell burning as compared to C- & O-shell burning. The dominant scales implied by the resulting VSH power spectra span more than 50% of the size of the convective shells, with secondary typical scales as small as (0.5–1) × 10⁹ cm. The slope of the VSH power spectra stays consistent over time regardless of the presence or absence of rotation. In all cases we find nearly circular shapes for the characteristic convective elements. In terms of the 2D convective flow properties, our results are in agreement with those of Meakin & Arnett (2007; see also Meakin et al. 2015; Cristini et al. 2016). In particular, we find a highly intermittent flow and mixing due to turbulent entrainment.

In terms of the effects of rotation on the strength and properties of convection alone we find our 2D simulations suggest minimal impact and generally lead to an overall small reduction of the total convective power stored in solenoidal motions ("no-rot" versus "rot-2" models) regardless of the nature of convective burning (C, O, or Si shell). We suggest this is due to the effects of the centrifugal forces mildly expanding the star, and thus changing the locations of the convective shells over time, subjecting the shell material to lower temperatures that, in turn, trigger lower nuclear burning rates and weaker convection. That alone may lead to initial velocity perturbations of smaller amplitude than in the case of zero rotation which may have a small impact on the susceptibility to a successful SN explosion following iron core collapse (Couch & Ott 2013).

We emphasize that, to answer the question of whether the inclusion of rotation significantly changes the ICs to the core-collapse SN mechanism, self-consistent evolution with the effects of rotation included is the proper approach to take. In our analysis we have followed the evolution of a 20 M_⊙ SN progenitor with both the effects of rotation and magnetic fields included in the transport of angular momentum and chemical mixing using the MESA code ("rot-ST" models). We find differences in the sign of the effect depending on the nature of convective shell-burning: during C- & O-shell burning there is more power stored in the solenoidal components than in the case of no rotation while during Si-shell burning the effects are very small, with hints of even reduction of the solenoidal mode power by the end of the simulation. We attribute this effect to differences in the initial MESA models for the two distinct shell-burning stages. A careful look in the upper panel of Figure 2 shows that the "rot-ST" _nuc profile during C- & O-shell burning has a secondary peak (due to O-burning) that is nearly an order of magnitude greater than the corresponding one for the "no-rot" model. On the contrary, in Figure 3 we see that during Si-shell burning the peaks in the _nuc profile for the "no-rot" model in the region (4–6) × 10⁸ cm are greater than the corresponding ones for "rot-ST." This is due to the enhanced chemical mixing by the ST mechanism during the C- & O-shell burning phase that effectively recycles fresh fuel from outer layers to deeper and hotter regions enabling faster specific nuclear energy generation rates. In contrast, during the later and more short-lived Si-shell burning phase the ST mechanism does not have the same radial extent and efficiency to instigate similar effects. This result illustrates that the presence of efficient mixing mechanisms need to be studied self-consistently and in more detail since they can alter the convective properties and structure of massive stars prior to CCSNe quantitatively.

If the Spruit-Taylor mechanism is even roughly correct, our simulations suggest that the cores of most massive stars do not rotate rapidly enough for rotation to be dynamically relevant to the CCSN mechanism (see also Heger et al. 2005; Gilkis et al. 2015). The inclusion of rotation and attendant angular momentum-transporting instabilities in the stellar evolution calculation, however, does significantly impact the nature of the convection surrounding the pre-collapse iron core. This could have important implications for the CCSN mechanism itself following core collapse. Recently, Mösta et al. (2015) suggested that above a certain rotational thershold, the magneto-rotational instability (MRI) can drive an inverse cascade of the magnetic energy, generating a large-scale magnetic field that can provide the conditions for gamma-ray burst jets and explain the origins of Type Ib/c SNe as well as some superluminous supernovae (SLSNe) powered by the spin-down of newly born magnetars. In the limit of the slow and typical rotation rates explored here we do not expect MRI-induced turbulence to have an important effect on the progenitor properties for single-star evolution. Binary evolution seems to offer an alternative channel, with the possibility of rapidly rotating core collapse in which the MRI may be effective as an explosive and jet-forming mechanism.

We thank J. Craig Wheeler for useful conversations and the anonymous referee for their comments. EC thanks the Enrico Fermi Institute for its support via the Enrico Fermi Fellowship. FXT thanks NASA support under the Theoretical and Computational Astrophysics Networks grant NNX14AB53G. The authors acknowledge the Texas Advanced Computing Center (TACC) at The University of Texas at Austin for providing HPC, visualization, and storage resources that have contributed to the research results reported within this paper. An award of computer time was provided by the Innovative and Novel Computational Impact on Theory and Experiment (INCITE) program. This research used resources of the Argonne Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract DE-AC02-06CH11357.