Self-consistent microscopic calculations for electron captures on nuclei in core-collapse supernovae

Electron capture (EC) reactions on nuclei play important roles in the evolution of a variety of astrophysical phenomena [1], including intermediate-mass stars [2, 3], core-collapse supernovae (CCSNe) [4, 5, 6, 7], thermal processes in neutron-star crust [8, 9], and nucleosynthesis in thermonuclear supernovae [10, 11]. Although the results from this Letter are beneficial for a variety of astrophysical simulations, it is primarily focused on the collapse phase of CCSNe, during which ECs are responsible for regulating the electron fraction ($Y_{e}$) and the dynamics of the collapse. At the onset of this phase, the stellar temperature $T$, density $\rho$, and entropy $s$ are about 10 GK, 10¹⁰ g/cm³, and 1 $k_{B}$, respectively, and a nuclear statistical equilibrium (NSE) exists in the core [12]. ECs on nuclei reduce $Y_{e}$, and emitted electron neutrinos carry away energy and entropy, accelerating the collapse of the core [13]. It has previously been shown  [6, 14, 7] that the nuclei that contribute most strongly to the change in $Y_{e}$ are situated along the $N=50$ and $N=82$ shell closures near ⁷⁸Ni and ¹²⁸Pd, respectively. At $\rho\gtrsim$10¹² g/cm³, the dynamical timescale of the collapse becomes shorter than the electron-neutrino diffusion timescale and the electron neutrinos become trapped [15, 12]. The collapse proceeds up to $\rho\gtrsim$ n_sat $\approx$ $2.81\cdot 10^{14}$ g/cm³. An outward propagating shock wave is emitted at a radius where the velocity of the in-falling material equals the speed of sound. The mass of the remaining inner core is proportional to final ${Y_{e}}^{2}$ [15, 12].

Since the densities, and thus Fermi energies in the collapsing star, are high, ECs can occur to highly excited states in the daughter nuclei. As the temperatures are also high, excited states in the parent are thermally populated, and ECs can occur from these excited states [16]. The EC is induced by allowed (Gamow-Teller) and forbidden transitions  [17, 18, 19]. The relevant conditions cannot be reproduced in terrestrial laboratories, and theoretical modeling is required. The theoretical models must be developed based on and tested against experimental data at zero temperature. In the past decades, experimental data from charge-exchange reactions at intermediate energies, complemented with limited information from EC/$\beta^{+}$-decay data, have been used for this purpose [1].

The combined theoretical and experimental efforts led to the development of EC rate libraries, which employ a temperature-density grid first used by Fuller, Fowler, and Newman [21]. In Ref. [22], the EC rates have been calculated for a pool of 250 nuclei using a hybrid model of shell-model Monte Carlo (SMMC) and random phase approximation (RPA) with the occupation numbers determined from the SMMC [23] or from the Fermi-Dirac parametrization depending on the computational restrictions of the SMMC and sensitivity of the EC rates on detailed structure of the Gamow-Teller transitions. Further pool of EC rates for more than 2200 (medium-) heavy nuclei relevant for the collapse phase was calculated using the RPA with a Fermi-Dirac parameterization for fractional occupation numbers  [22]. More recently, an EC rate library [6, 24, 25] has been developed that uses various sets of EC rates based on microscopic nuclear structure calculations that were benchmarked by experimental data complemented with a single-state approximation for nuclei [23, 20] for which microscopic calculations were not available. Benefiting from the increased computational power now available, in this work, for the first time, EC rate calculations including allowed and forbidden transitions and the effects of finite temperature are calculated within a self-consistent finite-temperature covariant energy density functional (EDF) theory and the quasi-particle random-phase approximation (QRPA) for stable and unstable nuclei with atomic numbers $Z=20$ to $Z=52$. This makes it possible to evaluate the role of ECs in the collapse phase of CCSNe within a consistent and complete nuclear physics framework. By adding the new rates to the existing library [26], the impacts of using these new EC rates are evaluated in 1D spherically symmetric simulations of the collapse phase. Multi-dimensional simulations that do not assume spherical symmetry and that are computationally demanding are necessary to fully understand CCSNe [27, 28, 29, 30]. The detailed treatment of EC rates for all nuclear species in these models is not yet possible, and approximations are necessary. The library of Ref. [22] was previously used to develop realistic inputs. The results from the present work enable the development of inputs for these multi-dimensional simulations based on modern microscopic theory.

The theoretical framework employed for calculating EC rates in this work is based on the relativistic EDF theory with the momentum-dependent D3C^∗ interaction [25, 31, 32]. The initial nuclear basis is determined within the finite-temperature Hartree-Bardeen-Cooper-Schrieffer (FT-HBCS) theory assuming spherical symmetry [19, 33]. Besides keeping the computational time manageable, such an approximation is reasonable since beyond $T\sim 2$ MeV, as in the collapse phase, most nuclei become spherical [34]. Unlike our previous work in Refs. [35, 25] which contained a simple monopole interaction, in the particle-particle channel, the pairing part of the more sophisticated Gogny D1S interaction is employed [36] with overall interaction strength $V_{pair}$ determined to reproduce pairing gaps of semi-magic isotopic chains. The final states are determined within the finite-temperature relativistic QRPA in the charge-exchange channel (FT-pnRQRPA) [19, 37, 38]. In order to obtain the expressions for the EC rates, we employ the current-current form of the weak-interaction Hamiltonian as in Refs. [39, 25]. Apart from the allowed Gamow-Teller (GT) transitions, we also consider the first-forbidden (FF) $J^{\pi}=0^{-},1^{-}$ and $2^{-}$ transitions, which largely determine the EC rate for neutron-rich nuclei. We assume that nuclei are fully ionized with the Fermi-Dirac distribution of electrons in the plasma, meaning that together with the nuclear properties, the rates are a function of the temperature $T$ and a product of stellar density with the electron-to-baryon ratio $\rho Y_{e}$. The EC rates can be expressed through the electron neutrino spectral functions $n(E_{\nu})$ [40]

where $E_{\nu}$ is the electron neutrino energy (neglecting a small rest mass). The neutrino spectrum can be subsequently utilized to calculate the neutrino energy loss rate, as well as the average neutrino energy required for CCSNe simulations [6, 22]. Model calculations include a large set of EC rates calculated from the proton to the neutron drip-line, starting from calcium $(Z=20)$ and terminating at tellurium $(Z=52)$, a total of 1652 individual nuclei. The neutron(proton) drip lines are determined from the neutron(proton) chemical potential condition $\lambda_{n(p)}>0$, based on the FT-HBCS calculation with the D3C^∗ interaction. The odd-$A$ and odd-odd systems are constrained with the particle number, similar to Refs. [25, 32]. The considered temperatures span a range from 0.01 to 50 GK, with a respective range in densities from $\rho Y_{e}=10^{1}$ to $10^{14}$ g/cm³, enough to consider the CCSNe trajectories where electron captures are of relevance.

The nuclear composition in the presupernova matter is obtained by considering the NSE [41]. Therefore, it is most relevant to check the rates of nuclei that have non-negligible NSE abundance. To get the NSE distribution, we employ the framework of Refs. [42, 22], with the Timmes and Arnett equation of state (EOS) from Ref. [43]. The NSE is determined by three parameters: stellar density $\rho$, electron-to-baryon ratio $Y_{e}$, and temperature $T$, along a CCSNe trajectory. Having the NSE abundances $Y_{i}$ for each nuclear species (which are normalized as $\sum\limits_{i}Y_{i}=1$), we can define the NSE-averaged rates as

Due to its simplicity, many CCSNe simulations employ analytic parameterizations for EC rates, assuming that rates depend strongly on EC $Q$-values through a series of Fermi integrals and fit coefficients to shell-model calculations in the $pf$-region [23, 20]. Such an approach omits detailed nuclear structure effects and has limited accuracy for neutron-rich nuclei where shell-model calculations are unavailable and extrapolations are necessary. Here, we compare our rates, $\langle\lambda_{ec}^{i,new}\rangle$, with the approximation from Ref. [20], $\langle\lambda_{ec}^{i,approx}\rangle$, averaged over the NSE distribution. The ratio between our NSE averaged rates and the approximation is shown in Fig. 1(a)–(c) for three points along the CCSNe trajectory. Starting from the FT-pnRQRPA rate table [26], the rates are interpolated using monotonic splines as discussed in Ref. [6]. We observe that with increasing stellar density $\rho$, the distribution of nuclei shifts to higher mass. At $\rho=10^{10}$ g/cm³, the NSE composition is dominated by nuclei up to $Z\sim 40$, and our calculations indicate higher rates than the approximation for more neutron-rich nuclei. At $\rho=10^{11}$ g/cm³, the ratio between the two rates for most nuclei is within one order of magnitude, with larger differences observed for more neutron-rich nuclei. Finally, at $\rho=10^{12}$ g/cm³, the NSE composition reaches $N=82$ region, and considerable differences between the FT-pnRQRPA rates and the approximation are obtained; for most of nuclei the new rates are larger, even up to 4 orders of magnitude for neutron-rich nuclei. The large discrepancy between the FT-pnRQRPA calculations and the approximation showcases a necessity to consider the rates within a microscopic framework, incorporating crucial effects of temperature on nuclear excitations. Another possible reason for the discrepancy is the significant contribution of FF transitions. In Fig. 1(d)–(f), we display the percentage of FF transitions in the total rate. At $\rho=10^{10}$ g/cm³, the EC rates are mostly dominated by allowed transitions. However, as the density increases to $10^{12}$ g/cm³, close to the neutrino trapping, the FF transitions dominate the rates for the whole NSE composition.

The CCSNe simulations are performed using the GR1D code, as in Refs. [6, 25, 24]. The GR1D considers early and post-bounce stages in spherical symmetry with general relativistic hydrodynamics, while the neutrino transport is handled through the NuLib library [44]. All simulations consider a 15-solar-mass, solar-metalicity, star progenitor (s15WW95 [45]) with SHFo EOS to determine the NSE [46]. We performed CCSNe simulations with the new FT-pnRQRPA EC rates for 1652 nuclei based on relativistic EDF theory, either including only the GT contribution [GT(EDF)] or GT with FF transitions [GT+FF(EDF)]. As a benchmark, we performed an additional simulation, but employing shell-model rates [13, 47, 48, 49] where available ($sd$- and $pf$-shell nuclei), and otherwise using the rates calculated in this work [GT(SM+EDF)]. Where available, the shell-model calculations generally achieve the highest accuracy [50] for nuclei up to the $pf$ shell but are not available for heavier nuclei and away from the valley of stability, and do not include contributions from the forbidden transitions.

First, in Fig. 2(a), we show the density evolution of the electron-to-baryon ratio $Y_{e}$. Up to $\rho=10^{11}$ g/cm³, we observe almost no impact of FF rates on the $Y_{e}$, as the pool of relevant nuclei consists of light and medium-mass nuclei for which mostly GT transitions are relevant. However, as the density increases beyond $10^{11}$ g/cm³ there is a significant impact of FF transitions, driving additional deleptonization and decreasing $Y_{e}$. As can be seen in Fig. 1(e), at $\rho=10^{11}$ g/cm³, the NSE averaged pool of nuclei extends up to $Z=45$, with ECs on nuclei beyond $Z=30$ having a significant contribution from the FF transitions. Ultimately, including FF transitions results in a 3% reduction of $Y_{e}$ close to the saturation density. By including the shell-model rates where available [GT(SM+EDF)], small differences in $Y_{e}$ can be noticed up to $10^{11}$ g/cm³. However, for $\rho>10^{11}$ g/cm³, the trajectory is almost independent of whether we use the SM rates or EDF rates for $pf$-shell nuclei. The density evolution of $Y_{e}$ directly affects the time evolution of the neutrino luminosity $L_{\nu}$, shown in Fig. 2(b). Before the bounce at $t\sim t_{bounce}$, the EC is dominated by $pf$-shell nuclei, and we observe almost no influence of FF transitions. However, due to faster deleptonization with FF transitions included, the peak neutrino luminosity is delayed by around 0.7 ms and decreases in amplitude by almost 5%, resulting in an overall considerably different luminosity profile. Finally, in Fig. 2(c), we display the central velocity as a function of the enclosed mass, which determines the effective mass of the underlying proto-neutron star (PNS). Including FF transitions within the new EC dataset leads to around a 7% reduction in the mass of the PNS.

Lastly, we aim to quantify the impact of individual nuclei on the CCSNe dynamics. The study in Ref. [6] determined two regions of the nuclear chart near $N=50$ and $N=82$, which have the highest contribution to the deleptonization process $\Delta Y_{e}^{i}$. Unfortunately, due to a lack of theoretical calculations in these regions, the study utilized analytic approximations for the rates. Therefore, we implement the new set of FT-pnRQRPA EC rates in the analysis of the impact of individual nuclei on the deleptonization $\Delta Y_{e}^{i}$, calculated by summing $|\dot{Y}_{e}|$ for each individual nuclear species and integrating it from the onset of collapse phase up to the neutrino trapping $t=t_{trapping}$. Results are shown in Fig. 3(a) for the EC rate set calculated with the FT-pnRQRPA, including GT and FF transitions. Simulations employing the new rate set confirm the results obtained in Ref. [6], that regions around $N=50$ and $N=82$ have the most significant impact driving the deleptonization. Including FF transitions plays almost no role for nuclei around $N=50$ shell closure. However, for the $N\sim 82$ nuclei, we notice a shift towards lower-$Z$ nuclei as indicated in Fig. 3(c) which contains FF contribution, compared to Fig. 3(b) with only GT. In particular, the deleptonization in ¹³²Sn changes by an order of magnitude once the FF transitions are included.

In conclusion, we have developed a set of stellar density- and temperature-dependent EC rates for nuclei ranging from $Z=20$–52 within the self-consistent finite-temperature covariant EDF theory based on the FT-pnRQRPA, which includes contributions from allowed and first-forbidden transitions from ground and thermally-populated excited states. Presented calculations provide a fundamental improvement over previous EC rate estimates by not using approximations and extrapolations that do not consider detailed nuclear structure effects. The implementation of the new rates in astrophysical simulations of CCSNe, compared to the previous version of weak-rate library [26], indicates that the EC rate is enhanced, resulting in a reduction of the electron fraction at saturation, the peak neutrino luminosity, and the enclosed mass at the bounce. By employing a consistent microscopic description of the EC rates, we have demonstrated that nuclei most critical to the deleptonization process near $N=82$ shell closure become less proton rich, which could provide guidance for future experimental efforts. The EC rates established in this work provide an essential foundation for advancing the accuracy of the future multi-dimensional supernova simulations, and better understanding the evolution of stellar explosions.

Acknowledgement. We thank Gabriel Martínez-Pinedo for insightful discussions and providing the NSE code. We also thank Toshio Suzuki for his insights and support. This work was supported by the U.S. Department of Energy under Award No. DOE-DE-NA0004074 (NNSA, the Stewardship Science Academic Alliances program), the U.S. National Science Foundation under award number NSF-2209429, Windows on the Universe: Nuclear Astrophysics at FRIB, and by the Croatian Science Foundation under the project number IP-2022-10-7773. This work was supported in part through computational resources and services provided by the Institute for Cyber-Enabled Research at Michigan State University. We also acknowledge support by the U.S. National Science Foundation under Grant No. PHY-1927130 (AccelNet-WOU: International Research Network for Nuclear Astrophysics [IReNA]).