r-process Nucleosynthesis from Three-dimensional Magnetorotational Core-collapse Supernovae

We employ ideal GRMHD with adaptive mesh refinement (AMR) and spacetime evolution provided by the open-source Einstein Toolkit  (Löffler et al. 2012; Mösta et al. 2014a). GRMHD is implemented in a finite-volume fashion with WENO5 reconstruction (Tchekhovskoy et al. 2007; Reisswig et al. 2013) and the HLLE Riemann solver (Einfeldt 1988) and constrained transport (Tóth 2000) for maintaining $\mathrm{div}B=0$. We employ the K₀ = 220 MeV variant of the equation of state of Lattimer & Swesty (1991) and the neutrino leakage/heating approximations described in O'Connor & Ott (2010) and Ott et al. (2012). At the precollapse stage, we cover the inner ∼5700 km of the star with four AMR levels in a Cartesian grid and add five more during collapse. After bounce, the PNS is covered with a resolution of ∼370 m, and the AMR grid structures consists of boxes with extents 5674.0 km, 3026.1 km, 2435.1 km, 1560.3 km, 283.7 km, 212.8 km, 144.8 km, 59.1 km, and 17.7 km. The coarsest resolution is h = 94.6 km, and refined meshes differ in resolution by factors of two. We use adaptive shock tracking to ensure that the shocked region is always contained on the mesh refinement box with a resolution of h = 1.48 km.

We draw the 25 M_⊙ (at zero-age-main-sequence) presupernova model E25 from Heger et al. (2000). While this model includes rotation, we parameterize the initial rotation law to match the simulations in Mösta et al. (2014b). The rotation law is cylindrical and axisymmetric following Takiwaki & Kotake (2011, their Equation (1)) and as in Mösta et al. (2014b):

$$\begin{eqnarray}&&{\rm{\Omega }}(x,z)={{\rm{\Omega }}}_{0}\displaystyle \frac{{x}_{0}^{2}}{{x}^{2}+{x}_{0}^{2}}\,\displaystyle \frac{{z}_{0}^{4}}{{z}^{4}+{z}_{0}^{4}},\end{eqnarray} \tag{ 1 }$$

with an initial central angular velocity of Ω₀ = 2.8 rad s⁻¹. The fall-off of the angular velocity profile with the cylindrical radius and vertical position is controlled by parameters x₀ = 500 km and z₀ = 2000 km, respectively. The rotation parameter of our setup is β_rot = 0.1% where ${\beta }_{0}\equiv T/| W|$ is the ratio of kinetic and gravitational potential energy. This is consistent with typical GRB-oriented progenitors models (i.e., E25 in Heger et al. 2000 has β_rot ∼ 0.15%). We set up the initial magnetic field using a vector potential of the form

$$\begin{eqnarray*}&&{A}_{r}={A}_{\theta }=0;{A}_{\phi }={B}_{0}{({r}_{0}^{3})({r}^{3}+{r}_{0}^{3})}^{-1}\,r\sin \theta ,\end{eqnarray*}$$

where B₀ controls the strength of the field. In this way, we obtain a modified dipolar field structure that stays nearly uniform in strength within radius r₀ and falls off like a dipole at larger radii. We choose r₀ = 1000 km to match the initial conditions of model B12X5β0.1 of the 2D study of Takiwaki & Kotake (2011) and the 3D study of Mösta et al. (2014b). We perform simulations for two different initial magnetic field strengths, B₀ = 10¹³ G (B13 from here on) and B₀ = 10¹² G (B12). While specific constraints on core magnetic fields of massive stars are not available from observations or multi-D stellar evolution calculations, we deem model B13 with its ultra-strong poloidal precollapse field to be likely unrealistic. B12 and B12-sym have nearly equally high precollapse fields, but Mösta et al. (2015) have shown that the resulting 10¹⁵ G field after core bounce can be delivered by the MRI and dynamo action. We, therefore, consider these models to be more realistic.

We perform simulations in full, unconstrained 3D. For model B12, we add random perturbations with a magnitude of 1% of the velocity at the start of the simulation. In model B12-sym, we do not add perturbations, and the simulation, therefore, evolves identical to an octant symmetry 3D (90° rotational symmetry in the x–y plane and reflection symmetry across the x − y plane) simulation. In this way, we reproduce the dynamics of an axisymmetric simulation while keeping the tracer set up and distribution identical between the simulations. In model B13, the 10 times stronger initial poloidal magnetic field prevents the disruption of the jet by a kink instability (Mösta et al. 2014b), as the key quantity for instability is the ratio of toroidal over poloidal field (see Equation (2) in Mösta et al. 2014b). As a result, simulations with and without perturbations for model B13 are nearly identical, and we only present results for a simulation without added perturbations. Model B13 is closest in dynamics to the model presented in Winteler et al. (2012), while model B12-sym mimics the dynamics of the prompt (an explosion within 50 ms after core bounce) 2D jet explosions in Takiwaki & Kotake (2011) and Nishimura et al. (2015). We summarize the initial magnetic field strengths and perturbation setups used in the simulations in Table 1.

Table 1.  Initial Magnetic Field Strength and Perturbation Setup (in Velocity) for the Three Simulations Considered here

Simulation	B13	B12-sym	B12
Bpol (G)	1013	1012	1012
Perturbations	None	None	$0.01\times | {\boldsymbol{v}}|$

Download table as:  ASCIITypeset image

We extract the thermodynamic conditions of ejected material using Lagrangian tracer particles. We place 10⁵ tracer particles on each of the simulations. We limit particles to 30 km ≤ r ≤ 1000 km for all simulations to ensure high enough resolution in the particle mass but also guarantee that the infall time for the outermost shell of particles is longer than the simulated time. The tracers are uniformly spaced, so that they represent regions of a constant volume. Each tracer particle gets assigned a mass, taking into account the density at its location and the volume the particle covers. Tracer particles are advected passively with the fluid flow, and data from the 3D simulation grid are interpolated to the tracer particle positions. In this way, we record the thermodynamic conditions and neutrino luminosities the particles encounter as a function of time. To determine the ejected mass in the explosion, we only take dynamically unbound particles into account. We determine whether a particle is unbound by determining the total specific energy is positive (and define the specific energy as the sum of internal, kinetic, and magnetic energy).

For models B13 and B12-sym, we map the tracer particle distribution onto the simulation shortly after core bounce, as both of these models explode within the first 40 ms of postbounce evolution. Model B12 takes considerably longer to explode, and we map the tracer distribution onto the simulation shortly before the transition to an explosion at time t−t_bounce = 80 ms. This allows us to ensure that we have a sufficient number of tracer particles in the outflows along the rotation axis of the core. We postprocess the particles with the open-source nuclear reaction network SkyNet of Lippuner & Roberts (2017). The network includes 7843 isotopes up to isotope ³³⁷Cn. Forward strong rates are taken from the JINA REACLIB database (Cyburt et al. 2010), and inverse rates are computed assuming detailed balance. Weak rates are taken from Fuller et al. (1982), Oda et al. (1994), Langanke & Martínez-Pinedo (2000), or otherwise from REACLIB. REACLIB also provides nuclear masses and partition functions. SkyNet evolves the temperature via the computation of source terms due to the individual nuclear reactions and neutrino interactions.

Computations for each particle start from nuclear statistical equilibrium (NSE). We start the network as soon as the temperature drops below T = 25 GK. The initial conditions for the network calculation are taken at this time. The neutrino luminosity data from the particle trajectories are noisy due to interpolation effects and the very high time resolution at which the tracer particles record the neutrino luminosities. We, therefore, compute a moving-window time average of the neutrino luminosities of the form ${\nu }_{\mathrm{av},i}=\alpha \cdot {\bar{\nu }}_{i}+(1.0-\alpha )\cdot {\nu }_{\mathrm{av},i-1}$, where i denotes the current timestep data and i − 1 is the previous one. We choose a weight function for each data set in the moving average as $\alpha \,=2\cdot {(n+1.0)}^{-1}$, with n = 40, and keep the neutrino luminosities constant after the end of the particle data. In cases in which the particle data in the simulation does not reach temperatures low enough for the network calculation to start, we extrapolate the particle data assuming homologous expansion. We carry out the network calculations to 10⁹ s, which is sufficient to generate stable abundance patterns as a function of mass number A.

Collapse and early postbounce evolution proceed identically in B12 and B12-sym. Core bounce occurs ∼350 ms after the onset of collapse for model B12 and ∼450 ms for model B13. The delay in core bounce for model B13 is due to the additional support by the 100 times higher magnetic pressure. Shortly after core bounce, the poloidal and toroidal B-field components reach B_pol, B_tor ∼ 10¹⁵ G for model B12 and B_pol, B_tor ∼ 10¹⁶ G for B13.

In model B13, the hydrodynamic shock launched at bounce, while still approximately spherical, never stalls, but continues to propagate into a jet explosion along the rotation axis (see left panels of Figures 1 and 2). The jet is powered by the extra pressure and stress from the strong magnetic field. The shock propagates at mildly relativistic speeds (v_jet ∼ 0.1–0.2 c) and reaches 1000 km at around 35 ms after core bounce. The jet is stabilized against the MHD kink instability by its large poloidal magnetic field (see stability condition Equation (2) in Mösta et al. 2014b). A mild m = 0 deformation is visible in the outflow in Figure 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In models B12 and B12-sym, the bounce shock stalls after ∼10 ms at a radius of ∼110 km. At this time, there is strong differential rotation in the region between the PNS core and the shock. This differential rotation powers rotational winding of the magnetic field and amplifies its toroidal component to 10¹⁶ G near the rotation axis within 20 ms of bounce. The strong polar magnetic pressure gradient, in combination with hoop stresses exerted by the toroidal field, then launches a bipolar outflow. As in  Mösta et al. (2014b), B12-sym now continues into a jet explosion and reaches ∼900 km after ∼100 ms. The expansion speed at that point is mildly relativistic (v_jet ≃ 0.1–0.15 c).

The 3D simulation with perturbations B12 starts to diverge from its symmetric counterpart B12-sym around ∼15 ms after bounce due to the non-axisymmetric spiral MHD kink instability (Mösta et al. 2014b). The subsequent 3D evolution is fundamentally different from both the effectively axisymmetric B12-sym model and model B13. The jet is strongly disrupted by the kink instability, which causes the outflow to cover a larger solid angle. The shock also propagates at a lower velocity than in simulations B13 and B12-sym (v_jet ≃ 0.03–0.05 c).

The developed jet structures of models B13, B12-sym, and B12 are depicted in Figure 1 (meridional slices) and Figure 2 (volume renderings). In both B13 (left) and B12-sym (center), a clean jet emerges and propagates at mildly relativistic speeds into the outer layers of the core and star. For simulation B12, the explosion propagates in a dual-lobe fashion, as in Mösta et al. (2014b), at nonrelativistic speeds. Material in the outflows of all three simulations is highly magnetized (β = P_gas/P_mag ≪ 1), ranges between $10\,{k}_{{\rm{B}}}\,{\mathrm{baryon}}^{-1}\leqslant s\leqslant 20\,{k}_{{\rm{B}}}\,{\mathrm{baryon}}^{-1}$ in specific entropy, and is neutron-rich (0.1 ≤ Y_e ≤ 0.4).

Material that gets accreted across the shock into the postshock flow is pushed to higher densities and temperatures as it is advected toward the PNS. Eventually, some of this infalling material is entrained in the outflow and ejected. We show the evolution of the temperature and density for a typical tracer particle trajectory from simulation B13 in Figure 3. Additionally, we show neutrino luminosities for both electron and electron antineutrinos in Figure 4, as recorded by three representative particles from simulations B13, B12-sym, and B12. After the initial neutronization burst that is visible for simulations B13 and B12-sym in the first 20 ms of postbounce evolution, the neutrino luminosities for both electron neutrinos and antielectron neutrinos converge toward values of $\simeq 5\,\times {10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ and stay approximately constant for the duration of the simulations.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Material that is ejected often reaches conditions where weak reactions proceed rapidly enough for weak equilibrium (or β equilibrium) to nearly take hold. In weak or β equilibrium, the rate of neutron destruction balances the rate of proton destruction. For a closed, thermalized system, β equilibrium is characterized by the condition ${\mu }_{{\nu }_{e}}+{\mu }_{n}={\mu }_{e}+{\mu }_{p}$, where μ_i is the chemical potentials of electron neutrinos, neutrons, electrons, and protons, respectively. For a fixed lepton fraction ${Y}_{L}={Y}_{{\rm{e}}}+{Y}_{{\nu }_{{\rm{e}}}}$, the condition of β equilibrium determines Y_e and the net electron neutrino fraction, ${Y}_{{\nu }_{{\rm{e}}}}$. When neutrinos are not trapped, the material moves toward dynamic β equilibrium, which is not determined by chemical potential equality but rather by rate balance (e.g., Arcones et al. 2010). Considering only captures on neutrons and protons, dynamic β equilibrium is given by the condition

$$\begin{eqnarray}{\dot{Y}}_{{\rm{e}}} & = & [{\lambda }_{{e}^{+}}(\rho ,T,{Y}_{e})+{\lambda }_{{\bar{\nu }}_{e}}]{Y}_{n}\\ & & -\ [{\lambda }_{{e}^{-}}(\rho ,T,{Y}_{e})+{\lambda }_{{\nu }_{e}}]{Y}_{p}=0,\end{eqnarray} \tag{ 2 }$$

where the free proton and neutron fractions are set by NSE, and λ_e−, λ_e+, ${\lambda }_{{\nu }_{{\rm{e}}}}$, and ${\lambda }_{{\bar{\nu }}_{{\rm{e}}}}$ are the rates of electron, positron, electron neutrino, and electron antineutrino capture, respectively. Equation (2) is an implicit equation for Y_e in β equilibrium, ${Y}_{{\rm{e}},\beta }={Y}_{{\rm{e}},\beta }(\rho ,T)$. When material is out of β equilibrium, weak interactions will push its Y_e toward ${Y}_{{\rm{e}},\beta }$ on a timescale given by ${\tau }_{\mathrm{weak}}=({\lambda }_{{e}^{-}}(\rho ,T,{Y}_{e})$ + ${\lambda }_{{e}^{+}}{(\rho ,T,{Y}_{e})+{\lambda }_{{\nu }_{e}}+{\lambda }_{{\bar{\nu }}_{e}})}^{-1}$. If ${\tau }_{\mathrm{weak}}$ is shorter than the dynamical timescale, ${\tau }_{d}=\rho /\dot{\rho }$, then the material should relax to a composition determined by ${Y}_{e,\beta }$.

The value of ${Y}_{{\rm{e}},\beta }$ depends both on the imposed neutrino fluxes and the thermodynamic state of the material, which determines the lepton capture rates. Generally, the electron capture rate dominates at high densities where electrons are degenerate, which pushes ${Y}_{{\rm{e}},\beta }$ to values less than ∼0.3. For fixed entropy at densities between ${10}^{9}\mathrm{and}{10}^{12}\,{\rm{g}}\,{\mathrm{cm}}^{-3}$, the β-equilibrium electron fraction goes down with the increasing density. Therefore, rapidly ejected material that has reached a higher density will often have a lower Y_e at the time that r-process nucleosynthesis begins because it experienced more electron captures. On the other hand, when neutrino captures dominate (and when there are no nuclei present), one finds ${Y}_{{\rm{e}},\beta }\approx {\lambda }_{{\nu }_{{\rm{e}}}}{({\lambda }_{{\bar{\nu }}_{{\rm{e}}}}+{\lambda }_{{\nu }_{{\rm{e}}}})}^{-1}$ (Qian & Woosley 1996). Since the neutrino emission in CCSNe is fairly similar in all flavors, β equilibrium driven by neutrino captures generally predicts ${Y}_{{\rm{e}},\beta }\gt 0.4$. Material ejected in the jet goes through density regimes where the electron captures dominate and then later through regimes where the neutrino captures dominate. In many cases, the material is never able to fully attain β equilibrium.

In the jet-driven SNe considered here, weak reactions at high density play a dominant role in setting Y_e just before nucleosynthesis starts. This in turn greatly impacts nucleosynthesis in the ejecta, since the electron fraction is the determining factor in whether or not a robust r-process occurs. In particular, the dynamics of the MHD explosion influence the conditions under which the electron fraction of ejected material is set, since the dynamics determine when τ_weak becomes close to the dynamical timescale and when weak interactions freeze out. In Figure 5, we show the evolution of Y_e, ${Y}_{{\rm{e}},\beta }$, τ_weak, and τ_d for representative particles from the three simulations. The four colored lines indicate the Y_e evolution as obtained from the nuclear reaction network calculation using four constant neutrino luminosities. In these calculations, we assume ${L}_{{\nu }_{e}}={L}_{{\nu }_{\bar{e}}}$ and constant mean neutrino energies of $\langle {\epsilon }_{{\nu }_{e}}\rangle =10\,\mathrm{MeV}$ and $\langle {\epsilon }_{\bar{{\nu }_{e}}}\rangle =14\,\mathrm{MeV}$. The solid black lines show the Y_e evolution using the neutrino luminosities, as obtained from the tracer particles. There is an initial decrease in Y_e as material is advected inward to higher density, causing ${Y}_{e,{eq}}$ to go down and τ_weak to decrease. Except in a limited number of cases, ${\tau }_{d}\lt {\tau }_{\mathrm{weak}}$ and β equilibrium is never obtained, although Y_e is always moving toward ${Y}_{e,\beta }$. Then, as the particle moves outward with the jet and evolves toward lower densities, the electron capture rate is reduced, and neutrino captures begin to dominate the weak reaction rates, causing Y_e to increase (this occurs before 0.1 s in all of the plots). Finally, once the temperature reaches T ≈ 5 GK, r-process nucleosynthesis begins and β⁻-decays of heavy nuclei cause Y_e to increase. The initial value of Y_e that is relevant at the start of r-process nucleosynthesis is seen as the plateau near 1 s in Figure 5.

Zoom In Zoom Out Reset image size

In the most energetic model, B13, the jet is formed at high densities very soon after bounce. As the jet propagates out, it entrained collapsing material that has not fully deleptonized, trapped neutrinos, and reached weak equilibrium. Therefore, after the ejecta begins to move out to larger radii and smaller densities, its Y_e still goes down due to electron captures trying to move the material toward the lower Y_e predicted by neutrino-free β equilibrium. The particles in B13 do not reach their minimum Y_e at their maximum density, but rather continue to experience electron captures that drive Y_e down toward its neutrino-free β-equilibrium value. This generally pushes the ejected material in B13 to small Y_e values (Y_e ≃ 0.15).

In model B12-sym, the evolution differs since it takes $\simeq 20\mbox{--}30\,\mathrm{ms}$ before a jet explosion is launched, and the propagation of the jet is slower than in model B13. Tracer particles that accrete toward the PNS reach slightly lower maximum densities than in simulation B13, but these particles reach their lowest Y_e at their maximum density. As particles get ejected in the outflow, they evolve to higher Y_e due to neutrino interactions. This effect is more pronounced for higher neutrino luminosities and is enhanced compared to simulation B13 due to the longer dwell time of particles in the vicinity of the PNS.

In model B12, the explosion dynamics are drastically different than in models B13 and B12-sym. The shock starts to expand only after ≃80 ms. Additionally, particles do not accrete to minimum radii as low as in simulations B13 and B12-sym, hence they reach smaller maximum densities before being ejected. The eventually ejected material comes close to reaching β equilibrium when it reaches its maximum density. Since this maximum density is lower than the maximum densities encountered in B12-sym, the minimum Y_e reached by the ejected B12 material is systematically higher. The propagation speed of the explosion is slower than in models B13 and B12-sym. However, the dwell time of particles in the vicinity of the PNS before being ejected is similar to that of model B12-sym. This is due to a similar ejection speed in the initially forming outflow near the PNS. It is only the shock surface itself that propagates at slower expansion speeds as the outflow material spirals away from the rotation axis. As the ejected material interacts with neutrinos, it evolves to higher Y_e values. This rise in Y_e is similar to the evolution in simulation B12-sym.

We show selected particles from simulations B13, B12-sym, and B12 as a scatter plot in Y_e at T = 5 GK and specific entropy in Figure 6. This figure illustrates the behavior described above for the individual tracer particles mentioned above. The symbols for each particle are color coded with the maximum density reached. For simulation B13, particles reach the highest densities as they reach the smallest minimum radii. The Y_e values at the time when the particles last exceed a temperature of 5 GK (approximately the temperature threshold for r-process nucleosynthesis) are peaked at low Y_e ≃ 0.2. The entropy values for the particles are similar to those for simulation B12-sym but are lower than for simulation B12. The low Y_e values for simulation B13 are almost exclusively set by β equilibrium since neutrino irradiation has less of an effect on the Y_e distribution, even for high neutrino luminosities, since material gets ejected very rapidly and efficiently. In model B12-sym, particles are at similar entropy but at lower maximum densities and higher Y_e values at T = 5 GK compared to simulation B13. The Y_e values at T = 5 GK for the particles from simulation B12-sym are set by two effects. First, the β-equilibrium Y_e values for these particles at lower densities are higher than for the particles at higher densities in simulation B13. Second, the dwell time for particles in the vicinity of the PNS is an order of magnitude longer in simulation B12-sym than in B13. This causes the Y_e values of the particles to shift to higher values as they cool to T = 5 GK. In the full 3D simulation B12, particles are at the lowest maximum densities and highest entropies. Their Y_e values in β equilibrium are, therefore, higher than for both B13 and B12-sym. The shift in the Y_e distribution as the particles evolve toward T ≃ 5 GK is similar to the evolution in simulation B12-sym. This is caused by the similar dwell time of material at small radii before being ejected in the outflow and hence a similar amount of neutrino interactions. The distribution of Y_e values for particles from simulation B12 is considerably wider than for simulations B13 and B12-sym.

Zoom In Zoom Out Reset image size

In addition to the high-density lepton captures, neutrino captures at lower densities can also impact Y_e at the beginning of nucleosynthesis. We parameterize the neutrino luminosities for the network calculation to determine how much of an impact uncertainties in our neutrino transport approximation have on the nuclear network calculation. In all simulations, higher neutrino luminosities push the particle Y_e values more quickly toward the higher end. This is particularly pronounced for neutrino luminosities ${L}_{\nu }={10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ and ${L}_{\nu }={10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$. The neutrino luminosities recorded from the tracer particles peak at a few ${L}_{\nu }={10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ (see Figure 4) and are bracketed by the ${L}_{\nu }={10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ and ${L}_{\nu }={10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ constant luminosity cases.

The ejecta properties vary significantly between the simulations. For the r-process nucleosynthetic signature of the explosion, the most important factor is how neutron-rich the ejected material is. In Figure 7, we show the distribution of the electron fraction Y_e for all particles in the ejected material when the temperature for the particles is last above 5 GK. This is representative of Y_e at the beginning of neutron-capture nucleosynthesis and leads to different ejecta properties between jet explosions (simulations B13 and B12-sym) and the 3D dual-lobe explosion (B12). We show results for both the leakage neutrino luminosities and our assumed constant neutrino luminosities. In the case of the leakage neutrino luminosities, the luminosities are also assumed to be constant in the network calculation after the end of the tracer particle data.

Zoom In Zoom Out Reset image size

For zero neutrino luminosities, the distributions for all simulations are peaked at Y_e ≲ 0.2. B12-sym is peaked at lower Y_e ≃ 0.15 than B12 at Y_e ≃ 0.21. The distribution for B12 is significantly broader than for B12-sym and B13. There are more particles at low Y_e values for simulation B12-sym than for B12. This is caused by particles reaching higher densities before they get turned around and swept up in the outflow (see Figure 6). In the ${10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ luminosity case, neutrino interactions shift the distributions to higher Y_e for all simulations. For model B13, where the dwell time of particles in the neutrino field is a factor of ≃10 shorter than in models B12-sym and B12, this shift is not large, but for simulations B12-sym and B12, the distributions are shifted by almost ≃0.1. As a result, there is effectively no material at Y_e ≲ 0.2 for simulation B12-sym and no material below Y_e ≃ 0.22 for simulation B12. The results obtained with the neutrino luminosities from the tracer particles show this effect even more clearly. Here, the effect of neutrino interactions is large enough that even the Y_e distribution for simulation B13 is shifted to values of Y_e ≳ 0.2, the distribution for B12-sym is now centered at Y_e ≃ 0.34, and the Y_e distribution for simulation B12 is shifted to Y_e ≃ 0.36.

The variations in the distribution of Y_e have consequences for the eventual nucleosynthesis, since one must have Y_e ≲ 0.25 to make the third r-process peak (Lippuner & Roberts 2015). Figure 8 shows abundance patterns for all three simulations, B13, B12-sym, and B12. We show the fractional abundance pattern averaged over all particles in the ejecta as a function of mass number A.

Zoom In Zoom Out Reset image size

If no neutrino luminosities are taken into account in the nucleosynthesis calculation, we find a robust r-process pattern in all three simulations. This is also true for a constant neutrino luminosity of ${L}_{{\nu }_{e}}={L}_{\bar{{\nu }_{e}}}={10}^{51}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$. For neutrino luminosity ${L}_{{\nu }_{e}}={L}_{\bar{{\nu }_{e}}}={10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$, all simulations still show a robust second r-process peak. B13 still has robust third-peak abundances, while B12-sym and B12 have reduced abundances in their third peaks (with B12 seeing the larger reduction). For a neutrino luminosity of ${L}_{{\nu }_{e}}={L}_{\bar{{\nu }_{e}}}={10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$, none of the simulations show significant amounts of material synthesized beyond A = 135. In all simulations, the reduction in the fractional abundance beyond A = 135 is accompanied by an overproduction of nuclei with A < 135 compared to the lower neutrino luminosity cases.

The abundance patterns calculated with the neutrino luminosities, as recorded from the tracer particles, fall in between the ${L}_{{\nu }_{e}}={L}_{\bar{{\nu }_{e}}}={10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ and ${L}_{{\nu }_{e}}={L}_{\bar{{\nu }_{e}}}={10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ constant luminosity cases. For simulation B13, material beyond A > 135 is reduced by a factor of 10 relative to the L_ν = 0 case, but for simulations B12-sym and B12, the results with luminosities as recorded from the tracer particles follow the ${L}_{{\nu }_{e}}={L}_{\bar{{\nu }_{e}}}\,={10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ case closely. There is no or very little second- or third-peak r-process material synthesized.

For a more direct comparison, we show nucleosynthesis calculations with constant neutrino luminosities of ${L}_{{\nu }_{e}}={L}_{\bar{{\nu }_{e}}}\,={10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$ for simulations B13, B12-sym, and B12 in Figure 9. While model B13 matches the solar abundance pattern well, model B12-sym falls short in the amount of material synthesized beyond A ≃ 170 by a factor of a few. For third-peak r-process material, the reduction in abundance between models B13 and B12-sym is slightly more than a factor of 10. For model B12, the reduction in material beyond the second peak is even more severe. Material beyond A = 135 is underproduced by two orders of magnitude with respect to simulation B13 and the solar abundance pattern. This underproduction is accompanied by an overproduction of nuclei with mass numbers of 50 ≤ A ≤ 80.

Zoom In Zoom Out Reset image size

We can also compare the elemental abundance patterns produced by our models to the elemental abundance patterns observed in low-metallicity halo stars. In Figure 10, we show the elemental abundances of B12, B12-sym, and B13 along with the observed abundances of the low-metallicity halo stars CS22892-052 (Sneden et al. 2000) and HD122563 (Honda et al. 2006), similar to Nishimura et al. (2017a). CS22892-052 has an abundance pattern that is consistent with the solar r-process abundances and, therefore, provides a good match to the abundances of B12-sym and B13. Although the nucleosynthetic pattern of B12 does not match the solar r-process abundance pattern or the pattern of CS22892-052, it is reasonably consistent with the incomplete r-process abundance pattern of HD122563.

Zoom In Zoom Out Reset image size

Table 2 summarizes the mass of the total and r-process ejecta material for models B13, B12-sym, and B12. The ejecta mass for simulation B13 is an order of magnitude larger than for simulations B12-sym and B12. This is due to the immediate jet launch after core bounce and the propagation speed of v ≃ 0.15c. All of the ejected mass measurements are only lower bounds on the total ejecta mass, since it is still increasing at the end of each of the simulations. Our lowest neutrino luminosity scenario that is still within the uncertainty of the Leakage luminosities from the tracer particles is the constant neutrino luminosity of ${10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$. This acknowledges a factor of up to a few uncertainty in the neutrino luminosities and average energies from the Leakage scheme (E. O'Connor 2018, private communication). For this luminosity, the r-process ejecta mass in model B13 is comparable to that found in Winteler et al. (2012) and Nishimura et al. (2015). For neutrino luminosities taken from the tracer particles and for constant neutrino luminosities of ${10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$, the r-process ejecta mass is reduced by an order of magnitude. In simulations B12-sym and B12, the r-process ejecta mass for our most optimistic scenario is already an order of magnitude smaller than for the same neutrino luminosity in simulation B13. For neutrino luminosities taken from the tracer particles and for constant neutrino luminosities of ${10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$, the r-process ejecta mass is effectively zero.

Table 2.  Total and r-process Ejecta Masses (Material with 120 ≤ A ≤ 249) for the Three Simulations, B13, B12-sym, and B12, for the Four Constant Neutrino Luminosities and the Neutrino Luminosities as Obtained from the Tracer Particles

Simulation	B13	B12-sym	B12
${M}_{\mathrm{ej},\mathrm{tot}}\,({M}_{\odot })$	0.0356	0.0043	0.0048
${M}_{\mathrm{ej},{\rm{r}}}$ ${L}_{\nu }=0\,\mathrm{erg}\,{{\rm{s}}}^{-1}\,({M}_{\odot })$	0.0337	0.0042	0.0038
${M}_{\mathrm{ej},{\rm{r}}}$ ${L}_{\nu }={10}^{51}\,\mathrm{erg}\,{{\rm{s}}}^{-1}\,({M}_{\odot })$	0.0336	0.0042	0.0037
${M}_{\mathrm{ej},{\rm{r}}}$ ${L}_{\nu }={10}^{52}\,\mathrm{erg}\,{{\rm{s}}}^{-1}\,({M}_{\odot })$	0.0320	0.0034	0.0018
${M}_{\mathrm{ej},{\rm{r}}}$ Lν from tracer $\,({M}_{\odot })$	0.0038	$5.4\times {10}^{-5}$	4.0 × 10−7
${M}_{\mathrm{ej},{\rm{r}}}$ ${L}_{\nu }={10}^{53}\,\mathrm{erg}\,{{\rm{s}}}^{-1}\,({M}_{\odot })$	0.0012	0.0	0.0

Note. Masses are in solar masses, M_⊙.

Download table as:  ASCIITypeset image