Nuclear Reaction Sensitivity in Magnetohydrodynamically Driven Supernovae

The MHD-SNe are an energetic class of SNe, also called hypernovae. Their explosion energy can reach about 10 B (1 B = 1 Bethe, 10⁵² erg), while being 10 times larger than the canonical CCSN. The MHD-SN is driven by its strong magnetic fields and fast rotation of the PNS; it occurs only when the energy of the magnetic field is roughly in equipartition with the kinetic energy. The MHD-SN is an alternative robust astronomical site for the r-process; the latest report on the element abundance pattern of the extremely metal-poor star SMSS J200322.54–114203.3 (Yong et al. 2021) was explained successfully by the MHD model. The explosion mechanism of MHD was proposed originally in the 1970s (Bisnovatyi-Kogan 1970; LeBlanc & Wilson 1970; Meier et al. 1976) and was considered as a source of r-process elements in the mid-1980s (Symbalisty 1984; Symbalisty et al. 1985).

A number of three-dimensional (3D) simulation efforts have been made to understand the mechanism of the MHD-SNe; the initial 3D model for a 15 M_⊙ progenitor, including sophisticated radiative transfer of neutrinos, is presented in Winteler et al. (2012). While the magnetic field amplifies from 5 × 10¹² G to 5 × 10¹⁵ G, corresponding to the ratio of magnetic to gravitational binding energy increasing 10,000 times, a successful full r-process up to the third peak A ∼ 195 is synthesized. A subsequential 3D model (Mösta et al. 2014) considers the microphysical finite-temperature equation of state (EOS) and a leakage scheme and lepton number exchange. It initially leads to strong bipolar jet formation, but afterward it is disrupted by a kink instability of the ultrastrong toroidal field near the rotation axis and then loses its jet feature. Further study (Mösta et al. 2015) has shown that the magnetorotational instability is a weak-field instability, while the scale toroidal field of the magnetar is small, and the MHD-driven explosions can be formed with the shape of a double lobe. But the r-process is not robust under this scenario (Mösta et al. 2018), and the elements beyond the second peak are reduced by a factor of 100 when the magnetorotational jets produced by the rapidly rotating core undergo a kink instability. The 3D simulations (Kuroda et al. 2020) of full general relativity with spectral neutrino transport show the shock wave is not revived after the bounce; expansion of the MHD-driven flows toward the polar directions is predominantly driven by the magnetic pressure. The most recent simulation (Obergaulinger & Aloy 2021) propagates the magnetically driven polar outflows with stronger magnetization, and the explosion reaches moderate energies of about 1 s. Generally speaking, the 3D simulations so far are still controversial, and the extremely high explosion energy of ∼10 B is hard to obtain by these models so far (Kuroda et al. 2020; Matsumoto et al. 2020); therefore, more studies are expected.

On the other hand, the 2D model continuously supports the MHD-SN as a conceivable site of r-process. Various numerical simulations are performed by investigating the features, including magnetic fields, rotation rates, and neutrino luminosity. The first 2D adiabatic model showed successful results in 2006 (Nishimura et al. 2006), although the effects of neutrino cooling and heating are ignored. Two types of the characteristic explosion are shown in Nishimura et al. (2015): prompt magnetic jets, with strong initial magnetic fields above 10¹⁵ G around the surface of the PNS and rapid rotation, promptly explode just after the bounce and launch a strong shock wave along with the rotational axis. The heavy nuclei including actinides can be produced here. The delayed magnetic jets, with weak initial magnetic fields and rotation, take about 10–100 ms after the bounce, while the magnetic fields increase owing to the wrapping. When the magnetic fields grow up to 10¹⁵ G, the stalled shock wave revives and ejects. The nuclei synthesized in this scenario cannot exceed the second peak A ∼ 130, where identified as the site of weak r-process. By considering the fission fragment distributions from the recycling of material (Shibagaki et al. 2016), the r-process abundance peaks are still well produced in the MHD-SNe. While introducing the magnetorotational instability into the MHD-SNe (Nishimura et al. 2017), not only the weak r-process elements but also the intermediate-mass nuclei are formed.

In the present work, the 2D model from Reichert et al. (2021a), developed by a series of studies (Obergaulinger & Aloy 2017, 2020; Aloy & Obergaulinger 2021; Obergaulinger & Aloy 2021), will be employed to run the reaction sensitivity study. This model includes parameterized rotation, magnetic field, and a self-consistent neutrino treatment. The simulation is based on the model 35OC (Woosley & Heger 2007), which has 35 M_⊙ for the initial mass of the star. While the collapse occurs, the mass is 28.1 M_⊙ with an iron core of mass 2.02 M_⊙. The angular velocity is 1.98 Hz at the center and 0.1 Hz at the surface of the iron core. The poloidal and toroidal components of the magnetic field at the core are 1.7 × 10¹⁰ G and 1.7 × 10¹¹ G, respectively.

Four models are calculated by varying the magnetic field, rotational velocity, timescale, explosion energy, mass ejection, and ratio of rotational/magnetic to gravitational energy of the PNS. A huge number of tracers of these models are generated to record the density, temperature, radius, electron fraction (Ye), neutrino luminosity, and energies during the evolution. The weak r-process, where the elements reach the second r-process peak around A ∼ 130, is successfully synthesized mainly by the models 35OC-Rs and 35OC-Rw. Ten tracers of the weak r-process group from model 35OC-Rs (Reichert et al. 2021b) will be implemented into the network calculation; the entropy versus Ye and maximum temperature versus density of these tracers are shown in Figure 1. By comparing these two plots with Figures 3 and 4 in Reichert et al. (2021b), the 10 tracers basically favor the low entropy of the whole sample, while the other three parameters are well representative. The mean value of the 10 tracers is a few kB/n less than the peak value of 35OC-Rs, which indicates that the present work generally reflects the character of MHD-SNe with relatively lower entropy conditions.

Zoom In Zoom Out Reset image size

Their abundances are shown in Figure 2, with the discrepancies mainly on the second peak A ∼ 130. Some tracers have a pretty strong second peak, similar to a solar abundance like tracer 2, while the amplitude of the second peak of some other tracer is an order of six less than the first peak, like tracer 6. Tracer 4 shows a special pattern with a feeble iron peak and elements beyond the second peak, even reaching to A ∼ 190.

Zoom In Zoom Out Reset image size

The nuclear reaction network code SkyNet (Lippuner & Roberts 2017), consisting of 7841 nuclei and 95,455 reactions, is employed in this literature. It includes a nuclear statistical equilibrium (NSE) solver, Helmholtz EOS, and electron screening corrections. SkyNet reads the density, electron fraction, and temperature as a function of evolution time from the MHD tracers, and the entropy is calculated by the EOS.

SkyNet uses nuclear reaction rates taken from the JINA-REACLIBV2.0 database (Cyburt et al. 2010), and the baseline (or reference) reaction rates in this work are also from it. A total of 1683 (α,n) reactions from Mn (Z = 25) to Sn (Z = 50) are extracted from the JINA-REACLIBV2.0 database. For each individual reaction, three calculations are run with all 10 trajectories: the baseline, reaction rate increased 100 times, and reaction rate decreased 100 times. In total, 50,490 network evolutions of (α,n) reactions are performed. The sensitivity of these reactions to a specific element is determined by the difference of abundance between the baseline and the varied rates.

Besides the (α,n) reaction, a couple of other types of reactions are supposed to play roles after NSE breakout (Woosley & Hoffman 1992) in the weak r-process. Here, the (n,γ), (p,γ), (α,p), (n,p) reactions of isotopes from Mn to Sn are investigated by their entire impact on the final abundance.

After the NSE phase breakdown, various types of reactions play their roles in the evolution to move matter toward high Z. It is critical to know which ones among this huge number of reactions dominate the nucleosynthesis, in order to lead scientist focus on them in both theory and experiments. It would be more efficient to investigate the impact of the types rather than the tremendous reactions individually. The types of (α,n), (n,γ), (p,γ), (n,p) and (α,p) of isotopes from Mn (Z = 25) to Sn (Z = 50) are selected from the rate database. Each type includes about 1600 reactions; the rates of them and their reverse reactions are removed from the REACLIB database and entirely increased 10 times or decreased 10 times. Then, they were implemented into the network and run with the 10 trajectories. The relative discrepancies in percentage between them and the baseline on four representative tracers are shown in Figure 3. Tracers 2, 3, and 6 represent the weak r-process abundance pattern with a strong, medium, and weak second peak, respectively. The pattern of tracer 4 is between the classic weak r-process and main r-process, since it has lower electronic fraction, and thus more fuel is available for the nucleosynthesis to reach mass number A ∼ 190.

Zoom In Zoom Out Reset image size

The most active type is the (α,n) reaction; it changes the final abundances not only on the shown tracers in the figure but also on the remaining ones. The elements in the iron peak of tracer 2 are severely influenced by the (α,n) reaction; the decreased rates could enrich their abundance more than a couple of times, which is out of the range in the figure. This illustrates that the reaction flow toward high Z mainly relies on the (α,n) reaction under this condition, while the slower rate prevents the seeds' nuclei from going further and therefore stack; and vice versa while the reverse reaction rate is enhanced. The impacts of (α,n) to tracers 3 and 6 are mild compared to those of tracer 2, but still obvious over the whole range of mass number. This result means that the dominated reaction in the weak r-process under the current MHD-SN model is (α,n), which is in good agreement with the study of CCSNe. The timescales of different types of reactions show that the (α,n) reaction is fastest during the sensitive temperature and determine the CCSN nucleosynthesis evolution (Bliss et al. 2017). The dominated role of the (α,n) reaction in the weak r-process is thus supported by both timescale and sensitivity study of nuclear reactions, and in the scenario of both CCSNe and MHD-SNe.

The neutron capture reaction (n,γ), which dominates in the main r-process, has limited influence here. It favors the tracers with a medium second peak of the abundance pattern in Figure 2, like tracer 3. The decreased (n,γ) significantly enriched the elements located between 60 and 80 and the second peak. But the matter between them is barely changed, which means that other reactions push the synthesis from iron peak to the high-Z zone, while the (n,γ) only affects the elements' abundances around the second peak. This means that the (n,γ) is not the driving force at the key stage of the weak r-process. On the other tracers, the (n,γ) seldom plays a key role in the evolution.

The (n,p) reactions impact the elements around iron peak in a few tracers; this is consistent with the study from Bliss et al. (2017). They figured out that when the temperature is above 4.2 GK, the (n,p) reactions are main charge-particle reactions. This means that this type of reaction principally plays its role at the early stage of the weak r-process after NSE, which is why its influences concentrate around the seeds' nuclei but barely on the heavy nuclei zone, as shown in Figure 3.

Comparing to the key of (α,n), the (α,p) reactions have minor impact on the final abundance. Because the weak r-process takes place at the neutron-rich side of the stability valley. Once an α-particle and a nucleus fuse to a compound nucleus, it is much more easy to evaporate a neutron rather than a proton. Thus, we can see only a few percent difference, shown in Figure 3. The (p,γ) reactions play a key role in proton-rich environments like the νp-process; however, in these current tracers with Ye lower than 0.5, it has much less influence than others.

The impact factor R is defined as

$$\begin{eqnarray}&&R(i)=\displaystyle \frac{| {Y}_{\mathrm{varied}}(i)-{Y}_{\mathrm{base}}(i)| }{{Y}_{\mathrm{base}}(i)},\end{eqnarray} \tag{ 1 }$$

where i designates the index of atomic number. Either the increased or deceased rate of a reaction resulting in a R(i) less than 0.9 or more than 1.1 will be identified as an influential reaction. These reactions are tabulated in Table 1; they are mainly from the region Z from 29 to 34 and A from 80 to 90. The influences of the top several reactions spread over a pretty wide range, because they proceed in the most sensitive temperature range for the (α,n) reactions. As Bliss et al. (2017) figured out, at the early phase of weak r-process (α,n) reactions do not play the leading role but are replaced by other reactions. For instance, for tracer 2, the Ye is about 0.24 after the NSE breakdown, which promotes the neutron capture reactions driving the seed nuclei away from the stability valley quite vertically, so the (α,n) of iron peak isotopes even not sensitive to themselves. With decreasing temperature, the (α,n) gradually manages the synthesis and exhibits its most influence while the reaction flow reaches the zone of Z = 29–34. The increased reaction rates of these isotopes rapidly consume the seed nuclei and thus weaken the abundances of the iron peak. On the other hand, their impacts on the production of heavier elements are even much larger than those of neighboring nuclei. For example, for the element Mo, the (α,n) reactions of ¹⁰³Mo, ⁸⁹Nb, and a couple of other neighboring nuclei also have an impact on its abundance, but the R(i) are one order less than the top reactions on Table 1. Because both the temperature and density decrease as the synthesis goes on, the (α,n) reactions of heavier isotopes are not able to have universal influence, as in the early stage, and thus have limited impact to the nuclei of its neighbor zone. The isotopes on Table 1 extremely favor the decay mode of β delayed neutron emission, for instance, the ^79–81Zn uniquely have this mode among the zinc isotopic chain. The β delayed neutron emission (Roberts et al. 1939) affects the even–odd fluctuation in r-process abundance (Suzuki et al. 2018), which usually requires a high Q value of the neutron-rich isotope. Here the correlation between it and the (α,n) reaction is found in the weak r-process under the MHD-SN scenario.

The weak r-process explains the abundance of elements 36 < Z < 47 of UMP stars. Table 2 tabulates the impactful reactions for given elements from krypton to silver, including the value of reaction varied, the R(i), and the corresponding tracers. These reactions are selected by two categories: group I contains the ones with the largest R(i), and group II lists the neighbor isotopes of the target element. Group I basically coincides with the reactions in Table 1, which has universal impact over the whole range of elements. The ⁸⁰Zn, ⁸¹Ga, and ⁸²Ge among these isotopes are closed shell nuclei with N = 50, in particular ⁸⁰Zn, which holds the top position on every subtable. These isotopes with shell closure are bottlenecked at the early phase of the nucleosynthesis since the equilibrium of (n,γ)–(γ,n). The reaction flow is therefore slowed down here, so the enhanced reaction rates of (α,n) provide an additional driving force to revive it. The (α,n) reactions move along with the direction of slope = 2(Z+2,N+1) in the nuclei chart, so the heavier closed shell nuclei with N = 50 like ⁸³As, ⁸⁴Se, ⁸⁵Br, and ⁸⁶Kr are detoured and replaced by ^87,88As, ^88–90Se, etc. The reactions of ^78,79Zn, ^79,70Ga, and ^80,81Ge also have apparent influence over a wide range, as they already reach the gateway of the N = 50 shell closure. Most reactions in Table 2 expose their sensitivity, while the reaction rate increased 100 times, except for the ⁸²Ge. Both the increase and decrease of its rate impact the abundances on most elements at similar level besides Kr and Y. Unlike the other N = 50 closed shell nuclei, ⁸⁰Zn and ⁸¹Ga mainly show their influence under rate-enhancing conditions, the decreased rate of ⁸²Ge presents comparable effect indicates that the (α,n) reaction is irreplaceable at this point, once it is reduced the path for moving matter toward the high-mass region will be suppressed.

The R(i) of group II is generally dozens of times less than that of group I, and only has a limited effect on their adjacent elements. The (α,n) reactions of some isotopes have mild influence on three or four elements, such as ^93,94Kr and ⁹¹Br to Rb, Sr and Y; ¹¹⁷Rh, ¹²⁹Ag and ¹³⁴In to Ru, Rh, Pd and Ag; ¹¹⁰Sn to Rh, Pd and Ag; ¹⁰³Mo to Mo, pd and Ag. Few isotopes also have a wide range impact as the top items in group I, but with much less R(i), like ^119,121Tc, influence the abundances from Zr to Ag, while ⁸⁹Nb almost has the same coverage besides Tc. ¹²²Ru is sensitive to five elements, Tc, Ru, Rh, Pd, and Ag. The remaining reactions of group II just have a limited impact on one or two elements nearby.