Stellar Collapse Diversity and the Diffuse Supernova Neutrino Background

In this work, we use a combined set of 200 solar-metallicity progenitor models from Woosley & Heger (2007, 2015, "WH07" and "WH15") and Sukhbold & Woosley (2014, "SW14"), which was already applied in Sukhbold et al. (2016) and can be downloaded from the Garching Core-collapse Supernova Archive. ³ All models are nonrotating single stars, evolved with the Kepler code (Weaver et al. 1978) up to the onset of iron-core collapse. The resulting grid of progenitors, unevenly distributed over the zero-age main-sequence (ZAMS) mass interval of 9–120 M_⊙, spans the commonly assumed range of "conventional" iron-core collapse SNe (or BH-forming, failed SNe).

Below that, in the narrow band between 8.7 M_⊙ and 9 M_⊙, we additionally consider ECSNe of degenerate ONeMg cores as another channel for NS formation (Miyaji et al. 1980; Nomoto 1984, 1987); yet it should be stressed that the exact mass range of ECSNe in the local universe is not definitive according to current knowledge (see, e.g., Poelarends et al. 2008; Jones et al. 2013, 2016; Doherty et al. 2015; Kirsebom et al. 2019; Zha et al. 2019; Leung et al. 2020). We employ a simulation by Hüdepohl et al. (2010, "model Sf") for the neutrino signal of such core-collapse events. The upper-mass end of the ZAMS mass grid is similarly uncertain and depends strongly on the physics of mass loss. However, as will be detailed in Section 3.2, high-mass contributions are suppressed by the steeply declining initial mass function (IMF) and are therefore of subordinate importance for the DSNB. In Sections 5.2 and 5.3, we will further consider progenitors from binary systems. The helium-star models used in this context were published by Woosley (2019), and their explosions were investigated by Ertl et al. (2020).

Our stellar collapse and explosion simulations are performed with the Prometheus-HotB code (Janka & Müller 1996; Kifonidis et al. 2003; Scheck et al. 2006; Ertl et al. 2016, 2020). The innermost 1.1 M_⊙ of the nascent proto-NS (PNS) is excised and replaced by a contracting inner grid boundary and an analytic one-zone core-cooling model with tunable parameters (for details, see Ugliano et al. 2012). This "central neutrino engine" is calibrated to yield explosions in agreement with the well-studied cases of SN 1987A and the Crab SN (SN 1054). More specifically, for pre-SN stars with ZAMS masses above 12 M_⊙, which Sukhbold et al. (2016) termed "87A-like," a PNS core model is applied and adjusted such that a given progenitor in the range of 15–20 M_⊙, namely S19.8, N20, W18, W15, or W20 (as described in Sukhbold et al. 2016), reproduces the observed explosion energy ((1.2–1.5) × 10⁵¹ erg; Arnett et al. 1989; Utrobin et al. 2015), ⁵⁶Ni yield (∼0.07 M_⊙; Bouchet et al. 1991; Suntzeff et al. 1992), and basic neutrino emission features (Bionta et al. 1987; Hirata et al. 1987) of SN 1987A. The low-mass end (9–12 M_⊙) is connected to the 87A-like cases by an interpolation of the core model parameters. As a second anchor point, we use the progenitor z9.6 by A. Heger (2012, private communication), which explodes with low energy (∼10⁵⁰ erg; Janka et al. 2012; Melson et al. 2015) and a small ⁵⁶Ni yield (∼0.0025 M_⊙; Wanajo et al. 2018) in self-consistent simulations, in good agreement with the observational constraints for the Crab SN (Smith 2013; Tominaga et al. 2013; Yang & Chevalier 2015). For more details on our calibration procedure, the reader is referred to Sukhbold et al. (2016).

Depending on the engine model, we obtain more or less energetic or failed explosions over the range of considered pre-SN stars, as can be seen in Figure 1. While the S19.8 and N20 calibrations lead to the largest fraction of successful SNe (red), W20 is a rather weak engine, resulting in the largest fraction of BH-forming cases (black). W18 and W15 reside between these two extremes, as can also be seen in Table 1, which shows the IMF-weighted fractions of successful and failed explosions for the different neutrino engines. The outcome in the low-mass range (9–12 M_⊙) is the same for all five cases, since our interpolation toward z9.6 is independent of the high-mass calibration. The nonmonotonic pattern of successful SNe and BH-forming collapses in Figure 1 has been described in previous works (Ugliano et al. 2012; Pejcha & Thompson 2015; Ertl et al. 2016; Müller et al. 2016; Sukhbold et al. 2016; Ebinger et al. 2019). It is based on the progenitor structure, which strongly varies with the ZAMS mass (O'Connor & Ott 2011; Horiuchi et al. 2014; Nakamura et al. 2015; Sukhbold et al. 2018).

Zoom In Zoom Out Reset image size

Compared to the simulations of Ertl et al. (2016) and Sukhbold et al. (2016), the neutrino transport outside of the PNS core, which is treated by a gray approximation (Scheck et al. 2006; Arcones et al. 2007), is slightly improved such that we are able to follow cases of long-lasting mass accretion until late collapse to a BH. For numerical reasons, the neutrino–nucleon scattering rate (Equation (D.68) of Scheck et al. 2006) is now split into two separate source terms, one for absorption ($\propto \langle {\epsilon }_{\nu }^{4}\rangle$, with _ν denoting the neutrino energy) and one for emission ($\propto T\langle {\epsilon }_{\nu }^{3}\rangle$), to avoid sign fluctuations for high temperatures T (for details, see Appendix B of Stockinger et al. 2020). Furthermore, an adaptive grid is implemented to better resolve the steep density gradient at the PNS surface. Our new code is applied without recalibrating the core models, which leads to slightly increased explosion energies because of decreased neutrino luminosities as compared to the models reported by Sukhbold et al. (2016). Accordingly, a few scattered progenitors that fail to explode with the old code (see Sukhbold et al. 2016, Figure 13) yield successful SNe with our new treatment. (A detailed report of the code changes and the consequences for the model results is provided in the appendices of Ertl et al. 2020.) In the work at hand, we moreover neglect the neutrino emission from the late-time fallback in so-called fallback SNe, in which the fallback matter pushes the NS beyond the BH limit after a successful explosion is initiated. This is justified because such cases turn out to be rare in the considered set of solar-metallicity progenitors (Ertl et al. 2016; Sukhbold et al. 2016) and additionally reside in the IMF-suppressed high-mass regime. In the context of our paper we therefore consider fallback SNe as successful SN events with the corresponding neutrino emission from NS formation. BH-forming events are only those cases where the BH does not form by fallback but by continuous accretion, and we use the terms "BH formation" and "failed SN" equivalently.

For each progenitor, we obtain the total energy release in neutrinos through time integration of the time-dependent neutrino luminosities, ${L}_{{\nu }_{i}}(t)$, and the mean values of the energies of the radiated neutrinos by computing the (luminosity-weighted) time average of the time-dependent mean neutrino energies, $\langle {E}_{{\nu }_{i}}(t)\rangle$, for all three considered neutrino species ${\nu }_{i}={\nu }_{{\rm{e}}},{\bar{\nu }}_{{\rm{e}}},{\nu }_{x}$, where ν_x denotes a representative heavy-lepton neutrino (${\nu }_{\mu },{\bar{\nu }}_{\mu },{\nu }_{\tau },{\bar{\nu }}_{\tau }$). Successful SNe are simulated up to a post-bounce time of t = 15 s, when the neutrino luminosities from PNS cooling have already declined to an insignificant level (see Appendix A). In the cases of failed explosions, however, the continued accretion of infalling mass onto the PNS releases gravitational BE, leading to an ongoing accretion component of the neutrino luminosities. The signals of such cases are truncated only when the PNS is pushed beyond the (still unknown) limit of BH formation, for which we consider four different values of the baryonic mass, ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}$, namely 2.3, 2.7, 3.1, and 3.5 M_⊙, which are motivated as follows.

Assuming an NS radius of (11 ± 1) km for maximum-mass NSs ⁴ and utilizing Equation (36) of Lattimer & Prakash (2001), which we provide as Equation (B1) in Appendix B, a baryonic NS mass of 2.3 M_⊙ converts to a gravitating mass of ${1.95}_{-0.03}^{+0.02}\,{M}_{\odot }$. This is marginally below the largest currently measured pulsar masses of ∼2 M_⊙ (Demorest et al. 2010; Antoniadis et al. 2013; Özel & Freire 2016; Cromartie et al. 2020), setting a lower limit for the maximum NS mass.

From the first gravitational-wave observation of a binary NS merger (GW170817; Abbott et al. 2017a) and its electromagnetic counterparts (Abbott et al. 2017b), Margalit & Metzger (2017) placed a tentative upper bound on the maximum gravitational NS mass of 2.17 M_⊙ (at a 90% confidence level), which follows from their reasoning that the merger remnant was a relatively short-lived, differentially rotating hypermassive NS, disfavoring both prompt collapse to a BH and the formation of a long-lived, supermassive NS. Their mass limit is compatible with other recent publications (e.g., Shibata et al. 2017; Alsing et al. 2018; Rezzolla et al. 2018; Ruiz et al. 2018; Lim & Holt 2019; Essick et al. 2020). Consistently, we take our case of a baryonic mass of 2.7 M_⊙ (corresponding to a ${2.23}_{-0.04}^{+0.03}\,{M}_{\odot }$ gravitational mass), which is close to this bound, as our reference threshold for BH formation.

Nonetheless, Margalit & Metzger (2017) pointed out several uncertainties related to their analysis. For instance, they neglected the effects of thermal pressure support on the stability of the compact merger remnant, which may change their conclusions. ⁵ Thermal effects are also important for the stability of hot PNSs on their way toward BH formation in the cases of failed SNe, possibly increasing the limiting mass as compared to the value for cold NSs (O'Connor & Ott 2011; Steiner et al. 2013; da Silva Schneider et al. 2020). For these reasons, we additionally explore two more extreme cases for the baryonic (gravitational) mass limit, namely 3.1 M_⊙ (${2.50}_{-0.05}^{+0.04}\,{M}_{\odot }$) and 3.5 M_⊙ (${2.75}_{-0.05}^{+0.05}\,{M}_{\odot }$). Eventually, further pulsar timing measurements (see Demorest et al. 2010; Antoniadis et al. 2013; Özel & Freire 2016; Cromartie et al. 2020) as well as an increased number of observed binary NS mergers (see, e.g., Abadie et al. 2010) should be able to shed more light on the maximum mass of NSs.

The most important results of the SN and BH formation simulations to be used in our DSNB calculations are the values of the time-integrated total energy release in neutrinos of all species and the time-averaged mean energies of the emitted electron antineutrinos. Figure 2 provides an overview of the corresponding values over the entire range of iron-core progenitors as a function of ZAMS mass for the exemplary case of the Z9.6 and W18 engine model, whose results will serve as a reference point in our later discussion (see Section 4). The three sets of pre-SN stars, WH15, SW14, and WH07, are separated by black vertical lines. Red bars indicate successful explosions and fallback SNe, whereas the outcomes of failed SNe are marked gray, dark blue, light blue, or cyan, depending on the different choices of critical baryonic mass for BH formation.

Zoom In Zoom Out Reset image size

The upper panel shows the explosion time, t_exp, for successful SNe, defined as the time when the shock passes 500 km (and not to be confused with the termination of our successful SN simulations and neutrino-signal calculations at 15 s, which is mentioned above). In cases of failed explosions, the time of BH formation, t_BH, is shown, which coincides with the sudden termination of the neutrino signal. Depending on the assumed NS mass limit and the progenitor-dependent mass-accretion rate, these times range from below 1 s up to 100 s in the most extreme cases (note the logarithmic scale). ⁶ This illustrates the need for a large set of long-time simulations to properly sample the neutrino contribution from the BH formation events.

The middle panel of Figure 2 displays the total radiated neutrino energies, ${E}_{\nu }^{\mathrm{tot}}$, computed as the time integrals of the summed-up neutrino luminosities of all species, ${L}_{\mathrm{tot}}(t)\,={L}_{{\nu }_{{\rm{e}}}}(t)+\,{L}_{{\bar{\nu }}_{{\rm{e}}}}(t)+\,4{L}_{{\nu }_{x}}(t)$, from core bounce (t = 0) to the end of the simulations (t = 15 s) or the termination of the signals (t = t_BH) for the SN or BH formation cases, respectively. Due to the aforementioned numerical improvements in the neutrino transport, these energies are slightly lower than those in Ertl et al. (2016) and Sukhbold et al. (2016). In Appendix B, we cross-check the values of ${E}_{\nu }^{\mathrm{tot}}$ by comparing them to the available budget of gravitational BE released during the cooling of the PNS, as estimated by using an analytic, radius-dependent approximate fit formula from Lattimer & Prakash (2001). We find good overall agreement, although our values might overestimate the neutrino energy loss by up to about 10%–20%, depending on the NS radius. In Section 7, we will discuss this and other uncertainties related to our DSNB predictions in more detail. In our work, we neglect contributions to the neutrino loss from fallback of matter after the successful launch of an explosion, since the amount of fallback is shown to be small (typically below 10⁻² M_⊙) for most progenitors (Ertl et al. 2016; Sukhbold et al. 2016) and since our values for the release of NS BE in neutrinos are on the high side anyway. In addition, fallback SNe with substantial late-time fallback (possibly turning NSs to BHs) are rare, as noted above.

The mean neutrino energies of the time-integrated energy emission are displayed in the bottom panel of Figure 2 for electron antineutrinos, which are most relevant for our study. Values around 15 MeV are the rather uniform outcome of successful SNe, in agreement with other publications (e.g., Mirizzi et al. 2016; Horiuchi et al. 2018). The mean energies from failed explosions, on the other hand, vary considerably between the progenitors and depend strongly on the NS mass limit. On the way to BH formation, the temperatures in the PNS's accretion mantle rise gradually, yielding increasingly harder neutrino spectra (see, e.g., Sumiyoshi et al. 2006, 2007, 2008; Fischer et al. 2009; Hüdepohl 2014; Mirizzi et al. 2016).

In a few cases, we need to extrapolate our neutrino signals, either because the simulations cannot be carried out to sufficiently late times or due to numerical problems, although such problems occur only after ∼10 s (see Appendix A and Figure A1 there). Furthermore, we should point out that the luminosities of heavy-lepton neutrinos are underestimated as compared to ν_e and ${\bar{\nu }}_{{\rm{e}}}$ in our simulations. This is a consequence of our approximate treatment of the microphysics (e.g., nucleon–nucleon bremsstrahlung is not included) and of the relatively modest contraction of the inner grid boundary and thus of underestimated temperatures in the accretion layer. To remedy this shortcoming as compared to more sophisticated transport calculations, we use information from such calculations with the Prometheus-Vertex code (Rampp & Janka 2002; Buras et al. 2006), which employs a state-of-the-art treatment of neutrino transport based on a Boltzmann moment closure scheme and a mixing-length treatment of PNS convection, to rescale the integrated energy loss in the different neutrino flavors, as detailed in Appendix C.

The neutrino signal of ECSNe from Hüdepohl et al. (2010, "model Sf") is followed for 8.9 s after core bounce and yields a total radiated neutrino energy of 1.63 × 10⁵³ erg, with a time-integrated ${\bar{\nu }}_{{\rm{e}}}$ mean energy of 11.6 MeV. These results are not shown in Figure 2, yet we use them for our DSNB estimates, which we describe in the next section.

For each progenitor, we compute the differential number spectrum $d{ \mathcal N }(t)/{dE}$ (in units of MeV⁻¹ s⁻¹) as a function of time t after core bounce from the time-dependent luminosity, $L(t)={L}_{{\nu }_{i}}(t)$, and mean energy, $\langle E(t)\rangle =\langle {E}_{{\nu }_{i}}(t)\rangle$, for all neutrino species (${\nu }_{i}={\nu }_{{\rm{e}}},{\bar{\nu }}_{{\rm{e}}},{\nu }_{x}$):

$$\begin{eqnarray}&&\displaystyle \frac{d{ \mathcal N }(t)}{{dE}}=\displaystyle \frac{L(t)}{ \langle E(t) \rangle }\displaystyle \frac{{f}_{\alpha }(E)}{{\int }_{0}^{\infty }{dE}\,{f}_{\alpha }(E)},\end{eqnarray} \tag{ 2 }$$

where we assume a spectral shape f_α (E) according to Keil et al. (2003),

$$\begin{eqnarray}&&{f}_{\alpha }(E)={\left(\displaystyle \frac{E}{\langle E(t)\rangle }\right)}^{\alpha }{e}^{-(\alpha +1)E/\langle E(t)\rangle }\,.\end{eqnarray} \tag{ 3 }$$

In our models, the shape parameter α of the spectrum is assumed to be constant over time. ⁸ Although this is a simplification, sophisticated simulations show that α does not change dramatically with time (e.g., Tamborra et al. 2012; Mirizzi et al. 2016), justifying this approximation. Instead, we vary α as a free parameter over a range of values (1 ≤ α ≤ 4), which we motivate in Appendix D.

For each progenitor and neutrino species, we then perform a time integration over the period of emission, from core bounce at t = 0 to a final time of t = t_f (with t_f = 15 s for successful explosions and t_f = t_BH for failed SNe):

$$\begin{eqnarray}&&\displaystyle \frac{{dN}}{{dE}}=\displaystyle \frac{\tilde{\xi }}{\xi }{\int }_{0}^{{t}_{{\rm{f}}}}{dt}\ \displaystyle \frac{d{ \mathcal N }(t)}{{dE}}\,.\end{eqnarray} \tag{ 4 }$$

Because the luminosities of heavy-lepton neutrinos ν_x are very approximate in our sets of simulations due to the incomplete microphysics and the relatively moderate core contraction mentioned in Section 2, we rescale each time-integrated spectrum with a factor $\tilde{\xi }/\xi$ (see Appendix C for details). By this procedure we adopt the total radiated neutrino energy (${E}_{\nu }^{\mathrm{tot}}$) from the simulated core-collapse models but redistribute them between the different neutrino species with weight factors obtained from SN and BH formation models with sophisticated neutrino treatment (see Table C1). $\tilde{\xi }={\tilde{\xi }}_{{\nu }_{i}}={\left({E}_{{\nu }_{i}}^{\mathrm{tot}}/{E}_{\nu }^{\mathrm{tot}}\right)}^{\mathrm{new}}$ thus constitutes the fraction of the total energy emitted in neutrino species ν_i . Correspondingly, $\xi ={\xi }_{{\nu }_{i}}={\left({E}_{{\nu }_{i}}^{\mathrm{tot}}/{E}_{\nu }^{\mathrm{tot}}\right)}^{\mathrm{old}}$ stands for the relative energy as originally computed in the core-collapse models considered in our study.

In Appendix D, we compare the shapes of our time-integrated spectra with results from sophisticated simulations (with detailed microphysics) by a few exemplary cases to examine the viability of our approximate treatment. We find good agreement with these simulations for values of the instantaneous shape parameter α of ∼3–3.5 for successful explosions and of ∼2 for failed SNe. In Appendix E, we provide, for a set of representative successful and failed SN models, the total radiated neutrino energies, the mean neutrino energies, and the spectral-shape parameters of the time-integrated neutrino (${\bar{\nu }}_{{\rm{e}}}$, ν_e, and ν_x ) spectra.

As mentioned in Section 2, our DSNB flux calculations also include the neutrino signal of the 8.8 M_⊙ ECSN simulated by Hüdepohl et al. (2010). The corresponding time-integrated spectra are computed according to Equations (2)–(4), but with the time-dependent shape parameters α = α(t) given by the simulation. We use the neutrino data of "model Sf," which takes into account the full set of neutrino interactions listed in Appendix A of Buras et al. (2006), including nucleon–nucleon bremsstrahlung, inelastic neutrino–nucleon scattering, and neutrino-pair conversions between different flavors, making rescaling of the spectra unnecessary, i.e., $\tilde{\xi }/\xi =1$ for all flavors.

The relative number of pre-SN stars depends on their birth masses. For our DSNB flux predictions, the time-integrated neutrino spectra ${dN}/{dE}$ for each core-collapse case therefore need to be weighted by an IMF (providing the number of stars formed per unit of mass as a function of the stellar ZAMS mass M). As in Hopkins & Beacom (2006), Horiuchi et al. (2011), and Mathews et al. (2014), we apply the modified Salpeter-A IMF of Baldry & Glazebrook (2003),

$$\begin{eqnarray}&&\phi (M)\propto {M}^{-\zeta }\,,\end{eqnarray} \tag{ 5 }$$

with ζ = 2.35 for birth masses M ≥ 0.5 M_⊙ and ζ = 1.5 for 0.1 M_⊙ ≤ M < 0.5 M_⊙. In our study, we consider masses up to 125 M_⊙. However, due to the steep decline of Equation (5), the high-mass end is suppressed and is thus of minor relevance for the DSNB.

The IMF-weighted neutrino spectrum ${{dN}}_{\mathrm{CC}}/{dE}$ of all core-collapse events can then be calculated as a sum over mass intervals ΔM_i associated with our discrete set of progenitor stars according to

$$\begin{eqnarray}&&\displaystyle \frac{{{dN}}_{{\rm{CC}}}}{{dE}}=\displaystyle \sum _{i}\displaystyle \frac{{\int }_{{\rm{\Delta }}{M}_{i}}{dM}\,\phi (M)}{{\int }_{8.7\,{M}_{\odot }}^{125\,{M}_{\odot }}{dM}\,\phi (M)}\displaystyle \frac{{{dN}}_{i}}{{dE}}\,,\end{eqnarray} \tag{ 6 }$$

where ΔM_i denotes the mass interval around the ZAMS mass M_i with the time-integrated spectrum ${{dN}}_{i}/{dE}$ of the corresponding SN, failed SN, or ECSN simulation. ⁹ Equation (6) is applied separately to the different neutrino species. As in Section 3.1, the indices ν_e, ${\bar{\nu }}_{{\rm{e}}}$, and ν_x are omitted here for the sake of clarity. In the following, we primarily focus on ${\bar{\nu }}_{{\rm{e}}}$, since the prospects for a first detection of the DSNB in upcoming detectors are the best for this species (see, e.g., Beacom & Vagins 2004; Yüksel et al. 2006; Horiuchi et al. 2009; An et al. 2016).

Nuclear burning proceeds fast in massive stars. As a consequence, the progenitors of core-collapse SNe (and failed SNe) have relatively "short" (<10⁸ yr) lives compared to cosmic timescales (see Kennicutt 1998). Therefore, the assumption is well justified that the cosmic core-collapse rate density R_CC(z) as a function of redshift equals the birth rate density of stars in the relevant ZAMS mass range (8.7 M_⊙ ≤ M ≤ 125 M_⊙), i.e.,

$$\begin{eqnarray}&&{R}_{{\rm{CC}}}(z)={\psi }_{\ast }(z)\displaystyle \frac{{\int }_{8.7\,{M}_{\odot }}^{125\,{M}_{\odot }}{dM}\,\phi (M)}{{\int }_{0.1\,{M}_{\odot }}^{125\,{M}_{\odot }}{dM}\,M\phi (M)}\simeq \displaystyle \frac{{\psi }_{\ast }(z)}{116\,{M}_{\odot }}.\end{eqnarray} \tag{ 7 }$$

Here, ψ_*(z) describes the cosmic star formation history (SFH) in terms of the star formation rate in units of M_⊙ Mpc⁻³ yr⁻¹, which can be deduced from observations (e.g., Hopkins & Beacom 2006; Reddy et al. 2008; Rujopakarn et al. 2010) and is thus independent of cosmological assumptions. In our study, we adopt the parameterized description by Yüksel et al. (2008),

$$\begin{eqnarray}&&{\psi }_{\ast }(z)={\dot{\rho }}_{0}\,{\left[{\left(1+z\right)}^{\alpha \eta }+{\left(\displaystyle \frac{1+z}{B}\right)}^{\beta \eta }+{\left(\displaystyle \frac{1+z}{C}\right)}^{\gamma \eta }\right]}^{\tfrac{1}{\eta }},\end{eqnarray} \tag{ 8 }$$

with the best-fit parameters from Mathews et al. (2014; see Table 1 therein). Note that the derivation of an SFH ψ_*(z) from observational data requires the use of an IMF, which should be consistent with the one employed in Equation (7). For this reason we use the Salpeter-A IMF (Baldry & Glazebrook 2003) to be consistent with the SFH data sample compiled by Mathews et al. (2014), which is based on the data sets by Hopkins & Beacom (2006) and Horiuchi et al. (2011). ¹⁰

Even though the cosmic core-collapse rate is not yet known to good accuracy (its impact on the DSNB flux is discussed, e.g., by Lien et al. 2010), our work is focused on variations of the neutrino source properties. To still account for the large uncertainty of R_CC, we additionally employ the ±1σ upper and lower limits to the SFH of Mathews et al. (2014), such that we obtain R_CC(0) = ${8.93}_{-3.01}^{+8.24}\times {10}^{-5}\,{\mathrm{Mpc}}^{-3}\,{\mathrm{yr}}^{-1}$ for the local universe. In Section 4, we further test parameterizations of the SFH by Madau & Dickinson (2014) and Fermi-LAT Collaboration et al. (2018). The cosmic metallicity evolution and its impact on the DSNB will be discussed briefly in Section 7.2.

For our DSNB calculations, we consider contributions up to a maximum redshift of ${z}_{\max }=5$. This limit is justified because, as pointed out in numerous previous works (Ando 2004; Keehn & Lunardini 2012; Mathews et al. 2014; Nakazato et al. 2015; Lunardini 2016), only sources at lower redshifts (z ≲ 1–2) noticeably add to the high-energy part of the DSNB, which is most relevant for detection (see Figure 3). Neutrinos from higher z are almost entirely shifted to energies below 10 MeV, where background sources dominate the flux and thus prevent a clear identification of the DSNB signal (see, e.g., Lunardini 2016).

Zoom In Zoom Out Reset image size

Throughout this work we assume standard ΛCDM cosmology with the present-day mass-energy density parameters Ω_m = 0.3 and Ω_Λ = 0.7 of matter and a cosmological constant, respectively, and the Hubble constant H₀ = 70 km s⁻¹ Mpc⁻¹. The expansion history of the universe is then given by dz/dt_c = $-{H}_{0}(1\,+\,z)\sqrt{{{\rm{\Omega }}}_{{\rm{m}}}{\left(1+z\right)}^{3}\,+\,{{\rm{\Omega }}}_{{\rm{\Lambda }}}}$. Using this together with Equation (1), we can write the DSNB flux spectrum (in units of MeV⁻¹ cm⁻² s⁻¹) as

$$\begin{eqnarray}&&\displaystyle \frac{d{\rm{\Phi }}}{{dE}}=\displaystyle \frac{c}{{H}_{0}}{\int }_{0}^{{z}_{\max }}\displaystyle \frac{{{dN}}_{\mathrm{CC}}}{{dE}^{\prime} }\displaystyle \frac{{R}_{\mathrm{CC}}(z)\ {dz}}{\sqrt{{{\rm{\Omega }}}_{{\rm{m}}}{\left(1+z\right)}^{3}+{{\rm{\Omega }}}_{{\rm{\Lambda }}}}}\,.\end{eqnarray} \tag{ 9 }$$

We do not vary the cosmological assumptions within our work, following most publications on the DSNB topic. For recent studies of the impact of different cosmological models on the DSNB, the reader is referred to Barranco et al. (2018) or Yang et al. (2019). Having described our computational model with all of its required inputs, we now proceed to a discussion of our results.

First, we study the impact of our engine model (as described in Section 2.2) on the DSNB flux spectrum. In the upper left panel of Figure 5, we show dΦ/dE for the various choices of central neutrino engines for our simulations. Sets with a higher percentage of failed explosions (see Figure 1 and Table 1) yield an enhanced DSNB flux, especially in the high-energy regime. This overall picture is in line with studies by Lunardini (2009), Lien et al. (2010), and Keehn & Lunardini (2012), who varied the fraction of BH-forming collapses while applying generic neutrino spectra and thus neglecting progenitor dependences. More recently, Priya & Lunardini (2017) and Møller et al. (2018) examined the fraction of failed SNe by assuming different ZAMS mass distributions, while Horiuchi et al. (2018), for the first time, employed a larger sample of simulations including seven BH formation cases, thus taking into account progenitor-dependent variations in the neutrino emission from failed explosions (by linearly interpolating the total energetics, mean energy, and shape parameter of their time-integrated neutrino spectra as a function of the compactness parameter of O'Connor & Ott 2011; see footnote 6). They explored relative fractions of BH formation cases between 0% and 45% by taking different threshold values for the compactness above which they assumed their progenitors formed BHs.

Zoom In Zoom Out Reset image size

Using our large sets of long-time simulations without predefined outcomes (also resulting in BH formation of less compact progenitors with low mass-accretion rates), we can confirm the common result of the previous studies: the larger the fraction of failed explosions, the stronger the enhancement of the DSNB at high energies. To better quantify this behavior, we follow Lunardini (2007) and fit the high-energy tail (20 MeV ≤ E ≤ 30 MeV) of our DSNB flux spectra with an exponential function:

$$\begin{eqnarray}&&\displaystyle \frac{d{\rm{\Phi }}}{{dE}}\simeq {\phi }_{0}\ {e}^{-E/{E}_{0}}.\end{eqnarray} \tag{ 10 }$$

Our model set Z9.6 and S19.8 with the lowest fraction of failed explosions (17.8%) features the steepest decline (i.e., E₀ = 4.5 MeV), while Z9.6 and W20 with 41.7% BH formation cases yields a flatter spectrum with E₀ = 5.1 MeV. The "normalization" ϕ₀, on the other hand, is hardly affected by our choice of engine model. Instead, it is determined by the uncertainty arising from the cosmic core-collapse rate, which shifts the entire flux spectrum vertically without changing the slope by more than ∼1%. ¹² The gray-shaded bands in Figure 5 indicate this severe normalization uncertainty (the +1σ upper limit to the SFH of Mathews et al. (2014) is taken for our highest-flux model, and the −1σ lower limit for our lowest-flux model). The aspect that the failed-SN fraction is likely to exhibit a dependence on metallicity (and thus on redshift) was pointed out by Nakazato et al. (2015) and Yüksel & Kistler (2015). We will come back to this point in Section 7.2.

The impact of the NS mass limit on the DSNB has been discussed in the literature to some extent (Lunardini 2009; Keehn & Lunardini 2012; Nakazato et al. 2015; Hidaka et al. 2016, 2018). Commonly, the spectra from exemplary simulations of BH formation with two different EoSs have been compared: the stiff Shen EoS (Shen et al. 1998, with incompressibility K = 281 MeV) and a softer EoS by Lattimer & Swesty (1991, "LS180" or "LS220," with K = 180 MeV or K = 220 MeV). Generally, a stiff EoS supports the transiently existing PNS of a failed SN against gravity up to a higher limiting mass than a soft EoS does. The final collapse to a BH therefore sets in after a longer period of mass accretion and neutrino emission with the consequence of higher spectral temperatures and an enhanced contribution to the DSNB flux.

Having a large compilation of long-time simulations at hand, we take a different (more rigorous) approach in our work: As described in Section 2, we directly vary the maximum baryonic NS mass, ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}$, without applying a certain EoS. Our neutrino signals from failed explosions are then truncated when the mass accretion from the collapsing progenitor star pushes the PNS mass beyond this critical threshold for BH formation. In the upper right panel of Figure 5, we show the DSNB flux spectra for our different choices of ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}$. Raising the NS mass limit from 2.3 M_⊙ to 3.5 M_⊙ drastically enhances the flux at higher energies, thus lifting the value of the slope parameter, E₀/MeV (see Equation (10)), from 4.4 to 5.6. This strong effect becomes immediately clear from Figure 2: a higher NS mass limit leads to enhanced time-integrated neutrino luminosities and generally hotter spectra, in line with studies by Lunardini (2009), Keehn & Lunardini (2012), Nakazato et al. (2015), and Hidaka et al. (2016, 2018).

We should mention that our study does not consider the possibility of a progenitor-dependent threshold mass for BH formation. O'Connor & Ott (2011) pointed out that thermal pressure support may be stronger for stars with high core compactness, lifting the maximum PNS mass to somewhat larger values. This might slightly reduce differences in the neutrino emission between individual progenitors. The results of Figure 5 should, however, remain essentially unchanged, because thermal stabilization of the PNS should be most relevant when the mass-accretion rate is high and the PNS becomes very hot. In such cases, however, the critical limit for BH formation is also reached quickly and the neutrino emission is not significantly extended. The wide range of values for ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}$ considered in our study should include the true NS mass limit, which depends on the still incompletely known high-density EoS of NS matter. Once the latter is better constrained by astrophysical observations and nuclear experiments and theory and thus once the maximum mass of cold NSs is better constrained, the question of a progenitor-dependent thermal effect on transient PNS stabilization can be addressed more thoroughly.

In our study, the spectral shape of the time-dependent neutrino emission is assumed to obey Equation (3) with a constant shape parameter α. Following our detailed analysis of the spectral shapes in Appendix D, we show in the lower left panel of Figure 5 how the DSNB flux spectrum changes when different values (between 1.0 and 3.0) of this instantaneous spectral-shape parameter α = α_BH are taken for the emission from failed explosions. For successful SNe, α is not varied but kept constant at the best-fit values of 3.0 and 3.5 (for M_NS,b > 1.6 M_⊙ and M_NS,b ≤ 1.6 M_⊙, respectively). Similar to its influence on individual failed-SN source spectra, a small value of α_BH broadens the shape of the DSNB such that its high-energy tail gets lifted relative to the peak (see Keil et al. 2003; Lunardini 2007, 2016). For α_BH = 1.0 (i.e., anti-pinched failed-SN source spectra), the exponential fit of Equation (10) yields E₀ = 5.4 MeV and ϕ₀ = 6.3 MeV⁻¹ cm⁻² s⁻¹ in the range of neutrino energies 20 MeV ≤ E ≤ 30 MeV. Choosing α_BH = 3.0, on the other hand, results in a more prominent peak at the cost of a suppressed flux at high energies (E₀ = 4.5 MeV; ϕ₀ = 12.3 MeV⁻¹ cm⁻² s⁻¹). Because the instantaneous spectral-shape parameter α is only varied for failed SNe while the contribution from successful SNe is unchanged, a slight "kink" becomes visible in the overall DSNB flux spectrum for cases of small α_BH, unveiling its "two-component" nature. Notice the crossings of the different curves at ∼3 MeV and ∼15 MeV. Accordingly, we construct the shaded band for the uncertainty of R_CC such that the lowest-flux and highest-flux models are considered in each segment.

In Table 3, we provide an overview of the two fit parameters ϕ₀ and E₀ for all models discussed in this section. We use the following naming convention for our DSNB models: "W18-BH2.7-α2.0" corresponds to our fiducial model with the Z9.6 and W18 neutrino engine ("W18"), with a baryonic NS mass limit of ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}=2.7\,{M}_{\odot }$ ("BH2.7"), and with the best-fit choice for the instantaneous spectral-shape parameter ("α2.0"; i.e., α_BH = 2.0). The two models "W20-BH3.5-α1.0" and "S19.8-BH2.3-α3.0," which employ the most extreme parameter combinations, yield the largest and smallest slope parameters E₀ among all our models and thus the highest and lowest fluxes at high energies, respectively.

Table 3. Exponential-fit Parameters of Equation (10) for a Subset of Our DSNB Models

Model	ϕ0 (MeV−1 cm−2 s−1)	E0 (MeV)
W18-BH2.7-α2.0 (fiducial)	${9.4}_{-3.3}^{+7.5}$ (${7.3}_{-2.6}^{+5.9}$, ${6.7}_{-1.4}^{+1.7}$, 4.4)	4.82+0.04 (4.84+0.04, ${5.1}_{-0.3}^{+0.2}$, 5.2)
W20-BH3.5-α1.0 (max.)	${6.6}_{-2.4}^{+5.3}$ (${4.8}_{-1.7}^{+3.8}$, ${4.9}_{-1.0}^{+1.3}$, 3.3)	6.79+0.05 (6.46+0.05, ${7.1}_{-0.3}^{+0.2}$, 7.1)
S19.8-BH2.3-α3.0 (min.)	${12.6}_{-4.3}^{+10.3}$ (${9.6}_{-3.3}^{+7.7}$, ${8.7}_{-1.7}^{+2.1}$, 5.6)	4.09+0.03 (4.32+0.04, ${4.4}_{-0.3}^{+0.2}$, 4.4)
S19.8-BH2.7-α2.0	${10.6}_{-3.7}^{+8.6}$ (${8.5}_{-3.0}^{+6.9}$, ${7.5}_{-1.6}^{+1.9}$, 4.9)	4.51+0.04 (4.60+0.04, ${4.8}_{-0.3}^{+0.2}$, 4.9)
N20-BH2.7-α2.0	${9.3}_{-3.2}^{+7.5}$ (${7.4}_{-2.6}^{+6.0}$, ${6.6}_{-1.4}^{+1.7}$, 4.3)	4.71+0.04 (4.75+0.04, ${5.0}_{-0.3}^{+0.2}$, 5.1)
W15-BH2.7-α2.0	${9.1}_{-3.2}^{+7.3}$ (${7.0}_{-2.4}^{+5.6}$, ${6.5}_{-1.3}^{+1.6}$, 4.2)	4.90+0.04 (4.90+0.05, ${5.2}_{-0.3}^{+0.2}$, 5.3)
W20-BH2.7-α2.0	${9.6}_{-3.4}^{+7.6}$ (${6.7}_{-2.3}^{+5.3}$, ${6.8}_{-1.3}^{+1.6}$, 4.4)	5.13+0.05 (5.14+0.05, ${5.5}_{-0.3}^{+0.3}$, 5.5)
W18-BH2.3-α2.0	${9.9}_{-3.4}^{+8.0}$ (${7.7}_{-2.7}^{+6.2}$, ${7.0}_{-1.4}^{+1.7}$, 4.5)	4.43+0.04 (4.57+0.04, ${4.7}_{-0.3}^{+0.2}$, 4.8)
W18-BH2.7-α2.0	${9.4}_{-3.3}^{+7.5}$ (${7.3}_{-2.6}^{+5.9}$, ${6.7}_{-1.4}^{+1.7}$, 4.4)	4.82+0.04 (4.84+0.04, ${5.1}_{-0.3}^{+0.2}$, 5.2)
W18-BH3.1-α2.0	${9.0}_{-3.1}^{+7.1}$ (${6.9}_{-2.4}^{+5.5}$, ${6.4}_{-1.3}^{+1.6}$, 4.3)	5.22+0.05 (5.14+0.05, ${5.5}_{-0.3}^{+0.2}$, 5.6)
W18-BH3.5-α2.0	${8.7}_{-3.1}^{+6.9}$ (${6.6}_{-2.3}^{+5.3}$, ${6.3}_{-1.3}^{+1.6}$, 4.2)	5.59+0.05 (5.44+0.05, ${5.9}_{-0.3}^{+0.2}$, 6.0)
W18-BH2.7-α1.0	${6.3}_{-2.2}^{+5.1}$ (${5.5}_{-1.9}^{+4.5}$, ${4.7}_{-1.1}^{+1.4}$, 3.2)	5.44+0.04 (5.25+0.04, ${5.7}_{-0.3}^{+0.2}$, 5.7)
W18-BH2.7-α1.5	${7.8}_{-2.7}^{+6.3}$ (${6.5}_{-2.3}^{+5.2}$, ${5.7}_{-1.3}^{+1.5}$, 3.8)	5.09+0.04 (5.02+0.04, ${5.4}_{-0.3}^{+0.2}$, 5.4)
W18-BH2.7-α2.0	${9.4}_{-3.3}^{+7.5}$ (${7.3}_{-2.6}^{+5.9}$, ${6.7}_{-1.4}^{+1.7}$, 4.4)	4.82+0.04 (4.84+0.04, ${5.1}_{-0.3}^{+0.2}$, 5.2)
W18-BH2.7-α2.5	${10.9}_{-3.8}^{+8.7}$ (${8.1}_{-2.8}^{+6.5}$, ${7.6}_{-1.5}^{+1.8}$, 5.0)	4.62+0.04 (4.70+0.04, ${4.9}_{-0.3}^{+0.2}$, 5.0)
W18-BH2.7-α3.0	${12.3}_{-4.3}^{+9.9}$ (${8.9}_{-3.1}^{+7.1}$, ${8.5}_{-1.6}^{+2.0}$, 5.5)	4.46+0.04 (4.58+0.04, ${4.8}_{-0.3}^{+0.2}$, 4.9)
W18-BH2.7-α2.0-He33	${8.0}_{-2.8}^{+6.4}$ (${6.3}_{-2.2}^{+5.1}$, ${5.7}_{-1.2}^{+1.4}$, 3.7)	4.76+0.04 (4.79+0.04, ${5.1}_{-0.3}^{+0.2}$, 5.1)
W18-BH2.7-α2.0-He100	${5.5}_{-1.9}^{+4.4}$ (${4.4}_{-1.5}^{+3.5}$, ${4.0}_{-0.8}^{+1.0}$, 2.5)	4.45+0.04 (4.59+0.04, ${4.7}_{-0.3}^{+0.2}$, 4.8)
S19.8-BH2.3-α2.0	${11.1}_{-3.8}^{+9.1}$ (${8.8}_{-3.0}^{+7.1}$, ${7.8}_{-1.6}^{+1.9}$, 5.0)	4.23+0.03 (4.42+0.04, ${4.5}_{-0.3}^{+0.2}$, 4.6)
W18-BH3.5-α1.0	${5.8}_{-2.1}^{+4.7}$ (${4.9}_{-1.7}^{+4.0}$, ${4.4}_{-1.1}^{+1.3}$, 3.0)	6.38+0.04 (5.96+0.04, ${6.6}_{-0.3}^{+0.2}$, 6.7)
W15-BH3.5-α1.0	${5.7}_{-2.1}^{+4.6}$ (${4.7}_{-1.7}^{+3.8}$, ${4.3}_{-1.0}^{+1.3}$, 2.9)	6.49+0.04 (6.08+0.04, ${6.8}_{-0.3}^{+0.2}$, 6.8)
W20-BH2.7-α1.0	${6.2}_{-2.2}^{+4.9}$ (${4.8}_{-1.7}^{+3.8}$, ${4.5}_{-1.0}^{+1.2}$, 3.1)	5.94+0.05 (5.73+0.04, ${6.2}_{-0.3}^{+0.2}$, 6.3)
W20-BH3.1-α1.0	${6.3}_{-2.3}^{+5.1}$ (${4.7}_{-1.7}^{+3.8}$, ${4.7}_{-1.0}^{+1.3}$, 3.2)	6.39+0.05 (6.12+0.05, ${6.7}_{-0.3}^{+0.2}$, 6.7)
W20-BH3.5-α2.0	${10.2}_{-3.6}^{+7.9}$ (${6.7}_{-2.4}^{+5.3}$, ${7.2}_{-1.3}^{+1.7}$, 4.8)	5.86+0.06 (5.77+0.06, ${6.2}_{-0.4}^{+0.3}$, 6.3)

Note. The fits are applied in the energy region 20 MeV ≤ E ≤ 30 MeV. The listed values correspond to the unoscillated ${\bar{\nu }}_{{\rm{e}}}$ DSNB flux spectra using the SFH from Mathews et al. (2014) with its associated ±1σ uncertainty. In parentheses, the values for the case of a complete flavor swap (${\bar{\nu }}_{{\rm{e}}}\leftrightarrow {\nu }_{x}$) are provided as well as the results for the SFH according to the EBL reconstruction model by Fermi-LAT Collaboration et al. (2018) and for the SFH of Madau & Dickinson (2014). The one-sided error intervals of E₀ in the cases with the SFH from Mathews et al. (2014) are caused by the fact that the functional fits to the SFH scale slightly differently with redshift (see footnote 12), with the best-fit case by Mathews et al. (2014) yielding the largest relative contribution from high-redshift regions and thus the smallest value of E₀ compared to both the +1σ and the −1σ limits.

Download table as:  ASCIITypeset image

As we mentioned in Section 2.1, the low-mass range of core-collapse SN progenitors is rather uncertain. It is widely believed that in degenerate ONeMg cores electron-capture reactions on ²⁰Ne and ²⁴Mg can counter the effects of oxygen deflagration, initiating collapse to an NS rather than thermonuclear runaway (Miyaji et al. 1980; Nomoto 1984, 1987). However, the conditions for the occurrence of such an ECSN in nature are still being debated (see, e.g., Jones et al. 2016; Kirsebom et al. 2019; Zha et al. 2019; Leung et al. 2020). Moreover, observations suggest that most massive stars are in binary systems (see, e.g., Mason et al. 2009; Sana et al. 2012), and evolution in binaries might lead to a larger population of degenerate ONeMg cores that produce ESCNe (Podsiadlowski et al. 2004).

In addition to these uncertain SN progenitors, three channels are being discussed that may lead to preferentially rather LM NSs, whose formation might contribute to the DSNB: Electron capture–initiated collapse may also occur when an ONeMg WD is pushed beyond the Chandrasekhar mass limit due to Roche-lobe overflow from a companion. Such an NS-forming event is referred to as AIC (see, e.g., Bailyn & Grindlay 1990; Nomoto & Kondo 1991; Ivanova & Taam 2004; Hurley et al. 2010; Jones et al. 2016; Wu & Wang 2018; Ruiter et al. 2019). Similarly, Saio & Nomoto (1985) suggested the MIC of two WDs as another possible scenario for the formation of a single NS (also see Ivanova et al. 2008; Schwab et al. 2016; Ruiter et al. 2019). Moreover, close-binary interaction might in some cases lead to the stripping of a star's hydrogen envelope and (most of its) helium envelope onto a companion NS, leaving behind a bare carbon–oxygen core (Nomoto et al. 1994; Dewi et al. 2002) and causing subsequent iron-core collapse. The explosion of such ultrastripped SNe (Tauris et al. 2013, 2015; Suwa et al. 2015; Müller et al. 2018) is discussed as the most likely evolutionary pathway leading to the formation of double-NS systems (Tauris et al. 2017; Mandel et al. 2021).

Previous works (e.g., Mathews et al. 2014; Horiuchi et al. 2018) have considered the contribution from ECSNe to the DSNB flux, and in a footnote, Lien et al. (2010) have already mentioned that, to a minor degree, neutrinos from the AIC of WDs might also add to the DSNB.

In our study we explore the consequences of additional formation channels of (rather) LM NSs on our DSNB predictions in a quantitative and systematic way, subsuming the possible contributions from ultrastripped SNe, AIC events, and MIC events in addition to the contribution from ECSNe that is included in our standard models. To this end, we employ a generic neutrino spectrum (${{dN}}_{\mathrm{LM}}/{dE}^{\prime}$) adopted from the ECSN calculations of Hüdepohl et al. (2010, "model Sf") since neutrino signals from sophisticated long-time simulations of AIC, MIC, and ultrastripped SNe are still lacking. We expect the neutrino emission properties of all three additional formation channels of LM NSs to be fairly similar to the case of ECSNe. Our approach is therefore meant to serve as an order-of-magnitude estimate, but it cannot capture any details connected to differences in the individual event rates and in the neutrino signals of the three channels of ultrastripped SNe, AIC events, and MIC events, which we combine to a single, additional LM NS formation component.

It should be mentioned here that the cosmic rates of such events are highly uncertain, because a large parameter space in the treatment of binary interaction (especially common-envelope physics) makes precise predictions difficult. Using population synthesis methods, Zapartas et al. (2017) found that core-collapse events in binary systems are generally delayed as compared to those of single stars. More particularly, Ruiter et al. (2019) showed that AIC and MIC can proceed in various evolutionary pathways, featuring a variety of delay times (from below 10² Myr up to over 10 Gyr) between starburst and eventual stellar collapse. For simplicity, we thus explore on the one hand different values of the comoving rate density, R_LM(z) = R_LM, which does not change with cosmic time ("LM_const"). On the other hand, we examine how our DSNB results differ in the case of an evolving rate for additional LM NS formation events ("LM_evolv"). The DSNB flux spectrum (Equation (9)) can be rewritten in the generalized form

$$\begin{eqnarray}&&\displaystyle \frac{d{\rm{\Phi }}}{{dE}}=\displaystyle \frac{c}{{H}_{0}}{\int }_{0}^{5}{dz}\,\displaystyle \frac{{R}_{\mathrm{CC}}(z)\tfrac{{{dN}}_{\mathrm{CC}}}{{dE}^{\prime} }+{R}_{\mathrm{LM}}(z)\tfrac{{{dN}}_{\mathrm{LM}}}{{dE}^{\prime} }}{\sqrt{{{\rm{\Omega }}}_{{\rm{m}}}{\left(1+z\right)}^{3}+{{\rm{\Omega }}}_{{\rm{\Lambda }}}}}\,.\end{eqnarray} \tag{ 11 }$$

In the lower right panel of Figure 5, we separately plot our fiducial DSNB prediction (dashed line; see Section 4) and the additional contribution from LM events for four different constant rate densities R_LM (solid lines), which we take as multiples of the local stellar core-collapse rate, R_CC(0) = 8.93 × 10⁻⁵ Mpc⁻³ yr⁻¹. However, since R_CC(z) varies strongly with redshift (it increases by over an order of magnitude from z = 0 to z = 1), we also consider the ratio of the comoving rate densities of LM NS formation events to those of "conventional" core-collapse SNe, both integrated over the cosmic history:

$$\begin{eqnarray}&&\chi =\displaystyle \frac{{\int }_{0}^{5}{dz}\,{R}_{\mathrm{LM}}(z)| {{dt}}_{{\rm{c}}}/{dz}| }{{\int }_{0}^{5}{dz}\,{R}_{\mathrm{CC}}(z)| {{dt}}_{{\rm{c}}}/{dz}| }\,.\end{eqnarray} \tag{ 12 }$$

This serves as a measure of the relative importance of how both types of neutrino sources contribute to the DSNB from the time of the highest considered redshifts (${z}_{\max }=5$) to the present day. In Table 4, we show the ratios

$$\begin{eqnarray}&&x=\displaystyle \frac{{\int }_{{E}_{1}}^{{E}_{2}}{dE}{\int }_{0}^{5}{dz}\,{R}_{\mathrm{LM}}(z)\tfrac{{{dN}}_{\mathrm{LM}}}{{dE}^{\prime} }| {{dt}}_{{\rm{c}}}/{dz}| }{{\int }_{{E}_{1}}^{{E}_{2}}{dE}{\int }_{0}^{5}{dz}\,{R}_{\mathrm{CC}}(z)\tfrac{{{dN}}_{\mathrm{CC}}}{{dE}^{\prime} }| {{dt}}_{{\rm{c}}}/{dz}| }\,,\end{eqnarray} \tag{ 13 }$$

i.e., the DSNB flux contributions from LM events relative to our fiducial model with R_LM = 0, integrated over different energy intervals [E₁,E₂] for our different choices of R_LM. To see an effect of at least 10% within the detection window (10–30 MeV), an additional (constant) LM rate R_LM = 1.55 × 10⁻⁴ Mpc⁻³ yr⁻¹ is required, which is nearly twice the local stellar core-collapse rate, R_CC(0), and corresponds to χ = 0.20. Such a fraction is well above present estimates for both AIC/MIC events (Metzger et al. 2009; Ruiter et al. 2019) and ultrastripped SNe (Tauris et al. 2013) of at most a few percent of the "conventional" core-collapse SN population. However, due to large uncertainties in the physics of binary interaction, the possibility of such a large population of LM NS formation events may not be ruled out completely.

Table 4. DSNB Contribution from Additional LM NS Formation Events

	(0–10) MeV	(10–20) MeV	(20–30) MeV	(30–40) MeV	(0–40) MeV
Fiducial DSNB Flux (${\bar{\nu }}_{{\rm{e}}}$), RLM = 0	22.7 cm−2 s−1	5.4 cm−2 s−1	0.6 cm−2 s−1	0.1 cm−2 s−1	28.8 cm−2 s−1
RLM = 1.0 × RCC(0), χ = 0.11	5.7% (6.1%)	5.7% (4.0%)	4.7% (2.5%)	2.6% (1.2%)	5.6% (5.6%)
RLM = 2.0 × RCC(0), χ = 0.23	11.3% (12.2%)	11.3% (7.9%)	9.3% (5.0%)	5.1% (2.4%)	11.2% (11.2%)
RLM = 3.0 × RCC(0), χ = 0.34	17.0% (18.4%)	17.0% (11.9%)	14.0% (7.5%)	7.7% (3.7%)	16.9% (16.9%)
RLM = 8.9 × RCC(0), χ = 1.00	50.0% (54.2%)	50.1% (35.0%)	41.3% (22.1%)	22.8% (10.8%)	49.8% (49.8%)

Note. First row: DSNB ${\bar{\nu }}_{{\rm{e}}}$-flux for our fiducial model (with R_LM = 0; see Table 2), integrated over different energy intervals. Rows 2–5: Flux contributions x (Equation (13)) from LM NS formation events (AIC, MIC, ultrastripped SNe) relative to the fiducial model for four different choices of the constant (LM_const) rate density R_LM. In parentheses, the values of x for an evolving LM NS formation rate (LM_evolv) with the same value of χ (Equation (12)) are given (see main text for details).

Download table as:  ASCIITypeset image

As a sensitivity check, we additionally consider a comoving rate density, R_LM(z), that linearly increases by a factor of 4 between z = 0 and z = 1 and stays constant at even larger redshifts, roughly following the observationally inferred rate of type Ia SNe (e.g., Graur et al. 2011). In the lower right panel of Figure 5, the LM flux contribution resulting from such an evolving rate is indicated by the pale red band (defined by 0.11 ≤ χ ≤ 1.00, as for the four cases of constant rates R_LM). The spectra are shifted toward lower energies, as expected from the relatively increased contribution from events at high redshifts. This can also be seen in Table 4 (values in parentheses). An enhancement of the DSNB flux by 10% at energies above 10 MeV would even require χ = 0.26 for an evolving LM rate, which would mean that, for example, roughly more than half of WD mergers lead to NS formation instead of a type Ia SN, if merging WDs explain the majority of SN Ia events. Although this seems to be disfavored on grounds of current observations and population synthesis models (e.g., Metzger et al. 2009; Ruiter et al. 2019), it might not be entirely impossible. Nevertheless, within the relevant detection window, the contribution from AIC/MIC events and ultrastripped SNe to the DSNB is likely to be hidden by the current uncertainty of the cosmic core-collapse rate (gray-shaded band in Figure 5). Only when this dominant uncertainty is reduced significantly may there be a chance to uncover a contribution to the neutrino background from such LM NS–forming events.

A large fraction of massive stars are expected to undergo binary interaction with a companion, possibly shedding their hydrogen envelopes (e.g., via Roche-lobe overflow or common-envelope ejection) and leaving behind bare helium stars (Sana et al. 2012). Taking this as a motivation, we explore how the inclusion of binary models affects our DSNB predictions. To this end, we employ a set of 132 helium stars with initial masses in the range of 2.5–40 M_⊙, originating from hydrogen burning in nonrotating, solar-metallicity stars, as reported in Woosley (2019). According to Equations (4) and (5) therein, this range of initial helium-core masses converts to ZAMS masses of 13.5–91.7 M_⊙. Stars with masses lower than those are assumed to form WDs, and thus do not contribute to the DSNB. ¹³ For the details of the pre-SN evolution (which includes wind mass loss), the reader is referred to Woosley (2019).

We use these progenitor models and perform SN simulations with the Prometheus-HotB code as done for single-star progenitors (see Section 2.2). A detailed and dedicated analysis of the explosions of these helium stars can be found in a recent paper by Ertl et al. (2020). In Figure 6, we show, for engine model Z9.6 and W18, the landscape of NS and BH formation events with the basic properties of the neutrino emission of relevance for our DSNB calculations. Compared to Figure 2, the range of stars experiencing core collapse is shifted toward higher ZAMS masses, starting only at 13.5 M_⊙. Moreover, there are no cases of BH formation below a ZAMS mass of 33 M_⊙. This can be understood as a consequence of mass loss by stellar winds during the pre-SN evolution of the helium stars, yielding less compact cores compared to stars that still possess their hydrogen envelopes (see Figures 1 and 10 in Woosley 2019).

Zoom In Zoom Out Reset image size

We note in passing that the values for the neutrino energy loss used in our present study differ in details from the numbers shown in Figure 5 of Ertl et al. (2020). First, we do not consider the additional neutrino energy loss from fallback accretion and consistently treat fallback SNe as NS formation events, whereas Ertl et al. (2020) took fallback into account in their estimates of the compact remnant masses and the associated release of gravitational BE through neutrinos. ¹⁴ Second, in our present study we extrapolate the neutrino emission of non-exploding cases until the accreting PNS in our Prometheus-HotB runs reaches the assumed and parametrically varied baryonic mass limit of stable, cold NSs, ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}$, and therefore collapses to a BH. In contrast, Ertl et al. (2020) employed for BH cases (with t_BH > 10 s) the radius-dependent fit formula of Lattimer & Prakash (2001) for the gravitational BE of an NS with the maximum mass assumed in their work. The energy release (${E}_{\nu }^{\mathrm{tot}}$) estimated that way is somewhat larger than our accretion-determined estimates (see Figure B2 in Appendix B).

Figure 7 illustrates how the inclusion of binary models impacts our DSNB predictions. In the left panel, we separately show the contributions from successful and failed explosions to the DSNB flux spectrum of electron antineutrinos for our fiducial model parameters (Z9.6 and W18; ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}=2.7\,{M}_{\odot };$ best-fit α; SFH from Mathews et al. 2014), assuming that all (100% of) progenitors evolve as helium stars ("W18-BH2.7-α2.0-He100"). Compared to single stars (Figure 3 and black dashed line in the right panel of Figure 7), the overall DSNB flux is reduced by a factor of ∼2 owing to the smaller fraction of stars experiencing core collapse. At the same time, the less frequent failed explosions produce a lower high-energy tail of the spectrum compared to our fiducial DSNB spectrum based on single stars (see also Table 3). If we assume that only 33% of all massive stars strip their hydrogen envelopes (as suggested by Sana et al. 2012; "W18-BH2.7-α2.0-He33"), the effects of helium stars on the DSNB spectrum are less dramatic, and the shifted spectrum lies within the uncertainty band associated with the SFH (gray-shaded band; see Figure 4).

Zoom In Zoom Out Reset image size

Applying the other neutrino engines considered in our work to the helium-star models, we obtain relative changes of the DSNB spectra similar to those in the case of our fiducial engine model Z9.6 and W18. We should stress at this point that if progenitors do not lose their entire hydrogen envelopes, end stages of stellar evolution more similar to those of single-star evolution can be expected (Woosley 2019).

So far we have not taken neutrino flavor oscillations into account but have identified the emission of electron antineutrinos (or neutrinos) by the considered astrophysical sources with the measurable DSNB flux of ${\bar{\nu }}_{{\rm{e}}}$ (or ν_e). However, on their way out of a collapsing star, neutrinos (and antineutrinos) undergo collective and matter-induced (MSW) flavor conversions (Wolfenstein 1978; Mikheyev & Smirnov 1985; Duan et al. 2010; Mirizzi et al. 2016). Hereafter, we discuss how such oscillations can affect our DSNB flux predictions.

Following Chakraborty et al. (2011) and Lunardini & Tamborra (2012), we write the DSNB flux spectrum of electron antineutrinos after including the effect of flavor conversions as

$$\begin{eqnarray}&&\frac{d{{\rm{\Phi }}}_{{\bar{\nu }}_{{\rm{e}}}}}{{dE}}=\bar{p}\,\frac{d{{\rm{\Phi }}}_{{\bar{\nu }}_{{\rm{e}}}}^{0}}{{dE}}+(1-\bar{p})\,\frac{d{{\rm{\Phi }}}_{{\nu }_{x}}^{0}}{{dE}},\end{eqnarray} \tag{ 14 }$$

where $d{{\rm{\Phi }}}_{{\bar{\nu }}_{{\rm{e}}}}^{0}/{dE}$ and $d{{\rm{\Phi }}}_{{\nu }_{x}}^{0}/{dE}$ are the unoscillated spectra for electron antineutrinos (${\bar{\nu }}_{{\rm{e}}}$) and a representative heavy-lepton neutrino (ν_x ). $\bar{p}\,\simeq \,0.7$ or $\bar{p}\,\simeq \,0$ denotes the survival probability of ${\bar{\nu }}_{{\rm{e}}}$ in cases of normal (NH) or inverted (IH) mass hierarchy, respectively. ¹⁵ Recently, Møller et al. (2018) confirmed by numerically solving the neutrino kinetic equations of motion that (matter-induced) neutrino flavor conversions can be well approximated by the simplified analytic expression of Equation (14) for the small set of Prometheus-Vertex simulations that they used in their study and that we also employ in our work as reference cases to calibrate some degrees of freedom in our modeling approach (see Table C1 in Appendix C). We already mentioned earlier that the large sets of core-collapse simulations underlying our DSNB calculations do not provide reliable information on the heavy-lepton neutrino source emission, which is why we use Prometheus-Vertex SN and BH formation models to rescale the neutrino energy release in the different neutrino species (see Section 3.1). For the same reason, we also adjust the spectral parameters, $\langle {E}_{{\nu }_{x}}\rangle$ and ${\overline{\alpha }}_{{\nu }_{x}}$, of the time-integrated ν_x emission (the bar in the symbol ${\overline{\alpha }}_{{\nu }_{x}}$ indicates that the shape parameter refers to the time-integrated spectrum rather than the instantaneous spectrum), guided by the sophisticated Prometheus-Vertex models listed in Table C1, to get a useful representation of the unoscillated DSNB spectrum of heavy-lepton neutrinos, $d{{\rm{\Phi }}}_{{\nu }_{x}}^{0}/{dE}$ (see Appendix C for details).

In the left panel of Figure 9, we show our unoscillated, fiducial DSNB spectrum for ${\bar{\nu }}_{{\rm{e}}}$, $d{{\rm{\Phi }}}_{{\bar{\nu }}_{{\rm{e}}}}^{0}/{dE}$ (black dashed line), and the corresponding unoscillated DSNB spectrum for ν_x , $d{{\rm{\Phi }}}_{{\nu }_{x}}^{0}/{dE}$ (red solid line), for our fiducial model parameters (see Section 4). According to Equation (14), the latter represents the case of IH, where a complete flavor swap (${\bar{\nu }}_{{\rm{e}}}\,\leftrightarrow \,{\nu }_{x}$) takes place. If, instead, the case of NH is realized in nature, an outcome between the two plotted extremes can be expected. The uncertainty arising from the cosmic core-collapse rate (corresponding to the ±1σ interval of the SFH from Mathews et al. 2014) is indicated by shaded bands. In Table 5, we additionally provide the integrated ${\bar{\nu }}_{{\rm{e}}}$-flux for different energy intervals and a complete flavor swap (${\bar{\nu }}_{{\rm{e}}}\,\leftrightarrow \,{\nu }_{x}$) in analogy to what is given in Table 2 for the case of no flavor oscillations (${\bar{\nu }}_{{\rm{e}}}$). The most important difference is a reduced contribution from failed SNe. This can be understood by the small relative fraction of heavy-lepton neutrino emission, ${\tilde{\xi }}_{{\nu }_{x}}$, in our two Prometheus-Vertex reference models for BH formation, which we employ for our rescaling (Appendix C). At the same time, the contribution from successful explosions (including ECSNe) is largely unchanged, which reflects the approximate flavor equipartition in their neutrino emission.

Zoom In Zoom Out Reset image size

Table 5. DSNB Flux Components for the Case of a Complete Flavor Swap (${\bar{\nu }}_{{\rm{e}}}\,\leftrightarrow \,{\nu }_{x}$)

	(0–10) MeV	(10–20) MeV	(20–30) MeV	(30–40) MeV	(0–40) MeV
Total DSNB Flux (${\bar{\nu }}_{{\rm{e}}}$)	19.3 cm−2 s−1	4.3 cm−2 s−1	0.5 cm−2 s−1	0.1 cm−2 s−1	24.2 cm−2 s−1
ECSNe	3.0%	1.5%	0.8%	0.4%	2.7%
Iron-core SNe	68.4%	65.4%	52.3%	37.9%	67.5%
Failed SNe	28.6%	33.2%	47.0%	61.8%	29.9%

Note. First row: DSNB ${\bar{\nu }}_{{\rm{e}}}$-flux for the case of a complete flavor swap (${\bar{\nu }}_{{\rm{e}}}\,\leftrightarrow \,{\nu }_{x}$), integrated over different energy intervals. Rows 2–4: Relative contributions from the various source types (ECSNe/iron-core SNe/failed SNe). Our fiducial model parameters (Z9.6 and W18; ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}=2.7\,{M}_{\odot };$ best-fit α) are used. Compare with Table 2, where values for the unoscillated ${\bar{\nu }}_{{\rm{e}}}$-flux are provided.

Download table as:  ASCIITypeset image

Despite the less relevant contribution from failed SNe, the slope parameter E₀/MeV of the exponential fit of Equation (10) is increased marginally from 4.82 to 4.84 in the case of a complete flavor swap (see Table 3) because smaller values of the spectral-shape parameter ${\overline{\alpha }}_{{\nu }_{x}}$ for heavy-lepton neutrinos (see Table C1; ${\lambda }_{\overline{\alpha }}^{\mathrm{PV}}\,\lt \,1$) partly compensate for the reduced flux of ν_x in the high-energy region associated with BH cases. The mean energies of the time-integrated neutrino signals are fairly similar for ν_x and ${\bar{\nu }}_{{\rm{e}}}$ (see Table C1; ${\lambda }_{E}^{\mathrm{PV}}\,\sim \,1$), as suggested by state-of-the-art simulations (e.g., Marek et al. 2009; Müller & Janka 2014) and as a consequence of the inclusion of energy transfers (non-isoenergetic effects) in the neutrino–nucleon scattering reactions (see Keil et al. 2003; Hüdepohl 2014). In conflict with this result of modern SN models with state-of-the-art treatment of the neutrino transport, several previous DSNB studies have employed spectra with $\langle {E}_{{\nu }_{x}}\rangle$ being considerably higher than $\langle {E}_{{\bar{\nu }}_{{\rm{e}}}} \rangle$ (particularly for the emission from failed explosions).

In line with recent studies by Priya & Lunardini (2017) and Møller et al. (2018), we find that neutrino flavor conversions have a fairly moderate influence on the DSNB (for ${\bar{\nu }}_{{\rm{e}}}$), which is well dominated by other uncertainties. Nonetheless, for our highest-flux models (with a weak central engine and a high maximum NS mass), which possess a large DSNB contribution from BH-forming events, the oscillation effects become more pronounced. We will further comment on this in Section 8.1.

For the DSNB ν_e flux the effects of neutrino flavor oscillations can be described in an analog manner (see, e.g., Chakraborty et al. 2011; Lunardini & Tamborra 2012). In the most extreme case of NH (and purely MSW-induced flavor conversions), a complete flavor swap (ν_e ↔ ν_x ) can take place, whereas for IH a measurable DSNB ν_e-flux spectrum between the unoscillated spectra of ν_x and ν_e can be expected.

As we point out in Appendix B, the total radiated neutrino energies (${E}_{\nu }^{\mathrm{tot}}$) of our successful SNe might, on average, be overestimated by a few percent, whereas the neutrino emission from failed explosions could be slightly underestimated in our modeling approach for the neutrino signals. In the right panel of Figure 9, we therefore compare our fiducial DSNB prediction (black dashed line) with a spectrum where the ${E}_{\nu }^{\mathrm{tot}}$ of all exploding progenitors is reduced by 15% (lower edge of the red band). This choice of reduction is guided by a comparison of ${E}_{\nu }^{\mathrm{tot}}$ with the gravitational BEs (BE₁₂) of the corresponding NS remnants (Equation (B1) with R_NS = 12 km; see Figure B1 and Table B1), consistent with the cold-NS radius suggested by recent astrophysical observations and constraints from nuclear theory and experiments (see footnote 4). Analogously, the upper edge of the red band in Figure 9 indicates a model where the ${E}_{\nu }^{\mathrm{tot}}$ of all failed explosions is increased by 15%. This case is motivated by the circumstance wherein the maximum neutrino emission in our failed-SN models with late BH formation lies ∼10%–20% below the maximally available gravitational BE according to Equation (B1) of an NS at its mass limit (see Figure B2 and Table B2). Any mix of changes of the NS and BH energy release will lead to intermediate results. Note that the corresponding red uncertainty band is hardly visible on the logarithmic scale.

A somewhat stronger effect can be seen when we vary the mean energies, 〈E〉, of the time-integrated spectra by −10% for successful SNe or by +10% for failed explosions (lower or upper edge, respectively, of the blue shaded band). Particularly at high energies, the spectra fan out noticeably. Such an uncertainty range cannot be ruled out according to present knowledge. Again, changing 〈E〉 for both successful and failed SNe at the same time yields a result in between the given limits. In Appendix D, we show that the outcome of our simplified approach is in reasonable overall agreement with results from the sophisticated Prometheus-Vertex simulations; however, the mean energies of the time-integrated spectra do not match perfectly (they lie ∼1 MeV higher/lower than in the Vertex models for successful/failed SNe; see Figures D1 and D2). Besides this fact, we should emphasize that the neutrino emission characteristics depend considerably on the still incompletely known high-density EoS (e.g., Steiner et al. 2013; Schneider et al. 2019) and also on the effects of muons, which have been neglected in most previous stellar core-collapse models but can raise the mean energies of the radiated neutrinos (Bollig et al. 2017).

Despite these uncertainties associated with the neutrino source, the cosmic core-collapse rate R_CC still constitutes the largest uncertainty affecting the DSNB, especially at lower energies (see Figure 4). Accordingly, the gray-shaded band in the right panel of Figure 9 indicates the ±1σ variation of R_CC for the SFH from Mathews et al. (2014). Upcoming wide-field surveys such as LSST (Tyson 2002) should be able to pin down the visible SN rate (below redshifts of z ∼ 1) to good accuracy, opening up a chance for DSNB measurements to specifically probe the contribution from faint and failed explosions (Lien et al. 2010).

Finally, one should keep in mind that we only employ solar-metallicity progenitor models in our simulations. Obviously, this is a simplification, because the distribution of metals in the universe is spatially nonuniform (see, e.g., the low metallicities in the Magellanic Clouds) and evolves with cosmic time. Since the fraction of failed explosions depends on metallicity (e.g., Woosley et al. 2002; Heger et al. 2003; Langer 2012), Nakazato et al. (2015) and Yüksel & Kistler (2015) considered a failed-SN fraction that increases with redshift. On the other hand, Panter et al. (2008) suggested that the average metallicity does not decline dramatically up to z ∼ 2. Assuming solar metallicity should therefore be a sufficiently good approximation, in view of the fact that the DSNB flux in the energy window favorable to DSNB detection is produced almost entirely by sources at moderate redshifts (see Figure 3).

At this point we should also remind the reader that a core-collapse SN is an inherently multidimensional phenomenon (see, e.g., Müller 2016). While our simplified 1D approach should be able to capture the overall picture of the progenitor-dependent neutrino emission, the increasing number of fully self-consistent 3D simulations will have to validate our results eventually.

After discussing the dependence of the predicted DSNB spectrum on different inputs in Sections 5 and 7, we compare our results now with the most stringent ${\bar{\nu }}_{{\rm{e}}}$-flux limit set by the SK experiment (Bays et al. 2012): Φ_>17.3 ≡Φ(E >17.3 MeV) ≲ (2.8–3.1) cm⁻² s⁻¹.

The various parameter combinations considered in our study lead to a wide spread between the DSNB flux spectra, as can be seen in the left panel of Figure 10. At high energies, the spectral tails of our different models fan out over more than an order of magnitude, with our most extreme cases yielding an integrated flux Φ_>17.3 that clearly exceeds the SK limit. To guide the eye, we roughly mark the region of such disfavored models (with Φ_>17.3 ≳ 3.1 cm⁻² s⁻¹) by a red shaded band, while flux spectra with Φ_>17.3 ≲ 3.1 cm⁻² s⁻¹, including our fiducial prediction (dashed line; see Section 4), lie in the gray band. We take the specific model "W20-BH3.5-α2.0" (i.e., Z9.6 and W20 neutrino engine, ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}=3.5\,{M}_{\odot }$, α_BH = 2.0) with the best-fit parameters taken for the SFH from Mathews et al. (2014) as a bounding case; it yields an integrated flux Φ_>17.3 = 3.09 cm⁻² s⁻¹, just within the uncertainty range of the SK limit (2.8–3.1 cm⁻² s⁻¹). We should emphasize, however, that this does not define a rigorous borderline, since spectra with quite different values of the slope parameter E₀ (Equation (10)) can yield similar integrated fluxes in the energy range above 17.3 MeV.

Zoom In Zoom Out Reset image size

In the right panel of Figure 10, we therefore plot Φ_>17.3 as a function of the fit parameter E₀ for a selection of models reaching close to (or beyond) the SK bound, which is marked by the red shaded region (with its uncertainty indicated by two slightly shifted lines). This plot is also intended for facilitating comparison with other works (see, e.g., Table 1 in Lunardini & Peres 2008 and Figure 19 in Bays et al. 2012). The tendency of greater integrated fluxes Φ_>17.3 for higher values of E₀ is obvious, yet there is significant scatter. In particular, the large uncertainty connected to the cosmic core-collapse rate (±1σ interval from Mathews et al. 2014, indicated by error bars) impedes definite conclusions. Nonetheless, models of our study with the most extreme combinations of parameters, such as the different cases of W20-BH3.5, which possess a strong contribution from failed SNe and thus large values of E₀ (see Section 5.1 and Table 3), are already disfavored, because their fluxes Φ_>17.3 reach beyond the SK limit (unless a minimal R_CC is taken). Also a less extreme value of the NS mass limit or a neutrino engine with a lower fraction of BH formation events can lead to an integrated flux close to the SK bound: models W20-BH3.1-α1.0 and W15-BH3.5-α1.0 (not shown in Figure 10) yield ${{\rm{\Phi }}}_{\gt 17.3}={2.7}_{-0.9}^{+2.3}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$ and ${{\rm{\Phi }}}_{\gt 17.3}\,={2.6}_{-0.9}^{+2.2}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$, respectively, with a dominant fraction (85% and 81%, respectively) of ${\bar{\nu }}_{{\rm{e}}}$ above 17.3 MeV originating from BH formation events. In Table 6, we provide the total integrated fluxes (Φ_tot), the fluxes within the observational window of 10–30 MeV (Φ₁₀₋₃₀), and the flux integrals above 17.3 MeV (Φ_>17.3) for a subset of our DSNB models.

Table 6. Total Integrated DSNB Flux (Φ_tot), Flux within the Observational Window of 10–30 MeV (Φ₁₀₋₃₀), and Flux above 17.3 MeV (Φ_>17.3) for the Same Subset of Our DSNB Models as Listed in Table 3

Model	Φtot (cm−2 s−1)	Φ10−30 (cm−2 s−1)	Φ>17.3 (cm−2 s−1)
W18-BH2.7-α2.0 (fiducial)	${28.8}_{-10.9}^{+24.6}$ (${24.2}_{-9.2}^{+20.7}$, ${20.8}_{-5.3}^{+6.6}$, 15.2)	${6.0}_{-2.1}^{+5.1}$ (${4.8}_{-1.7}^{+4.0}$, ${5.0}_{-1.6}^{+2.0}$, 3.4)	${1.3}_{-0.4}^{+1.1}$ (${1.0}_{-0.3}^{+0.9}$, ${1.2}_{-0.5}^{+0.6}$, 0.8)
W20-BH3.5-α1.0 (max.)	${41.7}_{-15.8}^{+35.7}$ (${30.4}_{-11.5}^{+26.0}$, ${30.1}_{-7.7}^{+9.5}$, 22.0)	${10.8}_{-3.8}^{+8.9}$ (${7.0}_{-2.5}^{+5.8}$, ${8.6}_{-2.5}^{+3.1}$, 6.0)	${3.5}_{-1.2}^{+2.9}$ (${2.1}_{-0.7}^{+1.8}$, ${3.0}_{-1.0}^{+1.2}$, 2.0)
S19.8-BH2.3-α3.0 (min.)	${24.4}_{-9.2}^{+20.9}$ (${22.8}_{-8.6}^{+19.5}$, ${17.6}_{-4.5}^{+5.6}$, 12.9)	${4.5}_{-1.6}^{+3.8}$ (${4.2}_{-1.5}^{+3.5}$, ${3.8}_{-1.3}^{+1.6}$, 2.6)	${0.7}_{-0.3}^{+0.7}$ (${0.8}_{-0.3}^{+0.7}$, ${0.7}_{-0.3}^{+0.4}$, 0.5)
S19.8-BH2.7-α2.0	${27.7}_{-10.5}^{+23.7}$ (${24.7}_{-9.4}^{+21.2}$, ${20.0}_{-5.1}^{+6.3}$, 14.6)	${5.5}_{-1.9}^{+4.6}$ (${4.7}_{-1.6}^{+3.9}$, ${4.6}_{-1.5}^{+1.8}$, 3.1)	${1.0}_{-0.4}^{+0.9}$ (${0.9}_{-0.3}^{+0.8}$, ${1.0}_{-0.4}^{+0.5}$, 0.7)
N20-BH2.7-α2.0	${27.4}_{-10.4}^{+23.4}$ (${23.6}_{-8.9}^{+20.2}$, ${19.8}_{-5.1}^{+6.2}$, 14.4)	${5.6}_{-1.9}^{+4.7}$ (${4.5}_{-1.6}^{+3.8}$, ${4.6}_{-1.5}^{+1.8}$, 3.2)	${1.1}_{-0.4}^{+1.0}$ (${0.9}_{-0.3}^{+0.8}$, ${1.0}_{-0.4}^{+0.5}$, 0.7)
W15-BH2.7-α2.0	${28.7}_{-10.9}^{+24.6}$ (${23.7}_{-9.0}^{+20.3}$, ${20.7}_{-5.3}^{+6.5}$, 15.2)	${6.1}_{-2.1}^{+5.1}$ (${4.7}_{-1.6}^{+3.9}$, ${5.1}_{-1.6}^{+2.0}$, 3.5)	${1.3}_{-0.4}^{+1.1}$ (${1.0}_{-0.3}^{+0.9}$, ${1.2}_{-0.5}^{+0.6}$, 0.8)
W20-BH2.7-α2.0	${32.6}_{-12.3}^{+27.9}$ (${24.9}_{-9.4}^{+21.3}$, ${23.5}_{-6.0}^{+7.4}$, 17.2)	${7.4}_{-2.6}^{+6.1}$ (${5.2}_{-1.8}^{+4.3}$, ${6.1}_{-1.9}^{+2.3}$, 4.2)	${1.7}_{-0.6}^{+1.5}$ (${1.2}_{-0.4}^{+1.0}$, ${1.5}_{-0.6}^{+0.7}$, 1.1)
W18-BH2.3-α2.0	${24.8}_{-9.4}^{+21.2}$ (${21.7}_{-8.2}^{+18.6}$, ${17.9}_{-4.6}^{+5.7}$, 13.1)	${4.8}_{-1.7}^{+4.0}$ (${4.1}_{-1.4}^{+3.4}$, ${4.0}_{-1.3}^{+1.6}$, 2.8)	${0.9}_{-0.3}^{+0.8}$ (${0.8}_{-0.3}^{+0.7}$, ${0.8}_{-0.3}^{+0.4}$, 0.6)
W18-BH2.7-α2.0	${28.8}_{-10.9}^{+24.6}$ (${24.2}_{-9.2}^{+20.7}$, ${20.8}_{-5.3}^{+6.6}$, 15.2)	${6.0}_{-2.1}^{+5.1}$ (${4.8}_{-1.7}^{+4.0}$, ${5.0}_{-1.6}^{+2.0}$, 3.4)	${1.3}_{-0.4}^{+1.1}$ (${1.0}_{-0.3}^{+0.9}$, ${1.2}_{-0.5}^{+0.6}$, 0.8)
W18-BH3.1-α2.0	${32.3}_{-12.2}^{+27.6}$ (${26.2}_{-9.9}^{+22.4}$, ${23.3}_{-6.0}^{+7.3}$, 17.0)	${7.3}_{-2.6}^{+6.1}$ (${5.4}_{-1.9}^{+4.5}$, ${6.0}_{-1.9}^{+2.3}$, 4.1)	${1.7}_{-0.6}^{+1.5}$ (${1.2}_{-0.4}^{+1.1}$, ${1.5}_{-0.6}^{+0.7}$, 1.1)
W18-BH3.5-α2.0	${35.4}_{-13.4}^{+30.3}$ (${28.1}_{-10.7}^{+24.0}$, ${25.5}_{-6.5}^{+8.1}$, 18.7)	${8.6}_{-3.0}^{+7.2}$ (${6.1}_{-2.1}^{+5.1}$, ${7.1}_{-2.1}^{+2.6}$, 4.9)	${2.2}_{-0.8}^{+1.9}$ (${1.5}_{-0.5}^{+1.3}$, ${2.0}_{-0.7}^{+0.9}$, 1.4)
W18-BH2.7-α1.0	${28.8}_{-10.9}^{+24.6}$ (${24.1}_{-9.1}^{+20.6}$, ${20.8}_{-5.3}^{+6.6}$, 15.2)	${6.0}_{-2.1}^{+5.0}$ (${4.8}_{-1.7}^{+4.0}$, ${5.0}_{-1.5}^{+1.9}$, 3.4)	${1.4}_{-0.5}^{+1.2}$ (${1.1}_{-0.4}^{+0.9}$, ${1.3}_{-0.5}^{+0.6}$, 0.9)
W18-BH2.7-α1.5	${28.8}_{-10.9}^{+24.6}$ (${24.2}_{-9.2}^{+20.7}$, ${20.8}_{-5.3}^{+6.6}$, 15.2)	${6.1}_{-2.1}^{+5.1}$ (${4.8}_{-1.7}^{+4.0}$, ${5.0}_{-1.6}^{+1.9}$, 3.4)	${1.3}_{-0.5}^{+1.2}$ (${1.0}_{-0.4}^{+0.9}$, ${1.2}_{-0.5}^{+0.6}$, 0.9)
W18-BH2.7-α2.0	${28.8}_{-10.9}^{+24.6}$ (${24.2}_{-9.2}^{+20.7}$, ${20.8}_{-5.3}^{+6.6}$, 15.2)	${6.0}_{-2.1}^{+5.1}$ (${4.8}_{-1.7}^{+4.0}$, ${5.0}_{-1.6}^{+2.0}$, 3.4)	${1.3}_{-0.4}^{+1.1}$ (${1.0}_{-0.3}^{+0.9}$, ${1.2}_{-0.5}^{+0.6}$, 0.8)
W18-BH2.7-α2.5	${28.8}_{-10.9}^{+24.6}$ (${24.2}_{-9.2}^{+20.7}$, ${20.8}_{-5.3}^{+6.6}$, 15.2)	${6.0}_{-2.1}^{+5.1}$ (${4.7}_{-1.6}^{+4.0}$, ${5.0}_{-1.6}^{+2.0}$, 3.5)	${1.2}_{-0.4}^{+1.1}$ (${1.0}_{-0.3}^{+0.8}$, ${1.1}_{-0.4}^{+0.6}$, 0.8)
W18-BH2.7-α3.0	${28.8}_{-10.9}^{+24.6}$ (${24.2}_{-9.2}^{+20.7}$, ${20.8}_{-5.3}^{+6.6}$, 15.2)	${6.0}_{-2.1}^{+5.0}$ (${4.7}_{-1.6}^{+4.0}$, ${5.0}_{-1.6}^{+2.0}$, 3.4)	${1.1}_{-0.4}^{+1.0}$ (${0.9}_{-0.3}^{+0.8}$, ${1.1}_{-0.4}^{+0.6}$, 0.8)
W18-BH2.7-α2.0-He33	${23.7}_{-9.0}^{+20.3}$ (${20.2}_{-7.7}^{+17.3}$, ${17.2}_{-4.4}^{+5.4}$, 12.5)	${4.9}_{-1.7}^{+4.1}$ (${4.0}_{-1.4}^{+3.3}$, ${4.1}_{-1.3}^{+1.6}$, 2.8)	${1.0}_{-0.3}^{+0.9}$ (${0.8}_{-0.3}^{+0.7}$, ${0.9}_{-0.4}^{+0.5}$, 0.7)
W18-BH2.7-α2.0-He100	${13.6}_{-5.2}^{+11.6}$ (${12.4}_{-4.7}^{+10.6}$, ${10.0}_{-2.6}^{+3.2}$, 7.2)	${2.7}_{-0.9}^{+2.3}$ (${2.4}_{-0.8}^{+2.0}$, ${2.3}_{-0.7}^{+0.9}$, 1.6)	${0.5}_{-0.2}^{+0.4}$ (${0.5}_{-0.2}^{+0.4}$, ${0.5}_{-0.2}^{+0.3}$, 0.3)
S19.8-BH2.3-α2.0	${24.4}_{-9.2}^{+20.9}$ (${22.8}_{-8.6}^{+19.5}$, ${17.6}_{-4.5}^{+5.6}$, 12.9)	${4.6}_{-1.6}^{+3.8}$ (${4.2}_{-1.5}^{+3.5}$, ${3.9}_{-1.3}^{+1.6}$, 2.6)	${0.8}_{-0.3}^{+0.7}$ (${0.8}_{-0.3}^{+0.7}$, ${0.8}_{-0.3}^{+0.4}$, 0.5)
W18-BH3.5-α1.0	${35.3}_{-13.4}^{+30.2}$ (${28.1}_{-10.6}^{+24.0}$, ${25.5}_{-6.5}^{+8.0}$, 18.6)	${8.4}_{-3.0}^{+7.0}$ (${6.0}_{-2.1}^{+5.0}$, ${6.8}_{-2.0}^{+2.5}$, 4.7)	${2.5}_{-0.8}^{+2.1}$ (${1.6}_{-0.6}^{+1.4}$, ${2.1}_{-0.7}^{+0.9}$, 1.5)
W15-BH3.5-α1.0	${35.4}_{-13.4}^{+30.3}$ (${27.8}_{-10.5}^{+23.8}$, ${25.6}_{-6.5}^{+8.0}$, 18.7)	${8.6}_{-3.0}^{+7.2}$ (${6.1}_{-2.1}^{+5.1}$, ${7.0}_{-2.1}^{+2.5}$, 4.8)	${2.6}_{-0.9}^{+2.2}$ (${1.7}_{-0.6}^{+1.4}$, ${2.2}_{-0.8}^{+1.0}$, 1.5)
W20-BH2.7-α1.0	${32.5}_{-12.3}^{+27.8}$ (${24.8}_{-9.4}^{+21.2}$, ${23.5}_{-6.0}^{+7.4}$, 17.2)	${7.3}_{-2.6}^{+6.1}$ (${5.2}_{-1.8}^{+4.3}$, ${6.0}_{-1.8}^{+2.2}$, 4.1)	${2.0}_{-0.7}^{+1.7}$ (${1.3}_{-0.5}^{+1.1}$, ${1.7}_{-0.6}^{+0.8}$, 1.2)
W20-BH3.1-α1.0	${37.2}_{-14.1}^{+31.8}$ (${27.7}_{-10.5}^{+23.7}$, ${26.9}_{-6.9}^{+8.5}$, 19.6)	${9.0}_{-3.2}^{+7.5}$ (${6.1}_{-2.1}^{+5.0}$, ${7.3}_{-2.1}^{+2.6}$, 5.0)	${2.7}_{-0.9}^{+2.3}$ (${1.7}_{-0.6}^{+1.5}$, ${2.3}_{-0.8}^{+1.0}$, 1.6)
W20-BH3.5-α2.0	${41.8}_{-15.9}^{+35.8}$ (${30.4}_{-11.5}^{+26.0}$, ${30.2}_{-7.7}^{+9.5}$, 22.1)	${11.1}_{-3.9}^{+9.2}$ (${7.2}_{-2.5}^{+5.9}$, ${9.0}_{-2.7}^{+3.3}$, 6.2)	${3.1}_{-1.1}^{+2.6}$ (${1.9}_{-0.7}^{+1.7}$, ${2.7}_{-1.0}^{+1.2}$, 1.9)

Note. The given values correspond to the unoscillated ${\bar{\nu }}_{{\rm{e}}}$ DSNB flux spectra using the SFH from Mathews et al. (2014) with its associated ±1σ uncertainty. In parentheses, the values for the case of a complete flavor swap (${\bar{\nu }}_{{\rm{e}}}\leftrightarrow {\nu }_{x}$) are provided as well as the results for an SFH according to the EBL reconstruction model by Fermi-LAT Collaboration et al. (2018) and for the SFH of Madau & Dickinson (2014).

Download table as:  ASCIITypeset image

Unlike the experimental DSNB flux limits of Malek et al. (2003), those provided by Bays et al. (2012) depend on the DSNB model employed. However, for an energy threshold close to ∼20 MeV, the flux limits are rather insensitive to the shape of the DSNB spectrum as pointed out by Lunardini & Peres (2008). In any case, the (Fermi–Dirac) spectral temperatures (3 MeV ≤ T_ν ≤ 8 MeV) that Bays et al. (2012) considered for their modeling of a "typical" SN source spectrum lead to DSNB spectra with slope parameters E₀ that cover the range of values obtained in our work. ¹⁶ Repeating their analysis of computing upper DSNB flux limits with our DSNB models should therefore lead to comparable bounds. Instead, we simply compare the experimental flux limit of (2.8–3.1) cm⁻² s⁻¹ to a subset of our model predictions in Figure 10. Naturally, this cannot replace a sophisticated statistical analysis, which is beyond the scope of this work.

Our fiducial model (W18-BH2.7-α2.0) yields an integrated flux of ${{\rm{\Phi }}}_{\gt 17.3}={1.3}_{-0.4}^{+1.1}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$, which is just below the SK bound, possibly not even by a factor of 2. Intriguingly, Bays et al. (2012) pointed out that there might already be a hint of a signal in the SK-II and SK-III data, giving hope that the first detection of the DSNB is within close reach now (see Beacom & Vagins 2004; Yüksel et al. 2006; Horiuchi et al. 2009; Keehn & Lunardini 2012; An et al. 2016; Priya & Lunardini 2017).

Since the SN neutrino emission is different for heavy-lepton neutrinos compared to electron antineutrinos (see Section 7.1), a complete (or partial) flavor swap (${\bar{\nu }}_{{\rm{e}}}\,\leftrightarrow \,{\nu }_{x}$) would affect our previous conclusions: In the case of IH (${\bar{\nu }}_{{\rm{e}}}$ survival probability $\bar{p}\,\simeq \,0$), the integrated flux above 17.3 MeV of our most extreme model (W20-BH3.5-α1.0) decreases by 39% from ${3.5}_{-1.2}^{+2.9}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$ to ${2.1}_{-0.7}^{+1.8}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$, which is just below the SK bound (however, note the large uncertainties due to the cosmic core-collapse rate). For the case of NH ($\bar{p}\,\simeq \,0.7$), we obtain ${{\rm{\Phi }}}_{\gt 17.3}={3.1}_{-1.1}^{+2.6}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$, which is still somewhat above the SK limit. Applying neutrino flavor conversions to our fiducial model, the effects are much reduced, as described in Section 7.1 (see left panel of Figure 9 and Table 5): Φ_>17.3 decreases by only 7% (21%) from ${1.3}_{-0.4}^{+1.1}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$ to ${1.2}_{-0.4}^{+1.0}({1.0}_{-0.3}^{+0.9})\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$ for NH (IH), still reaching close to the SK bound. In Table 6, the first values in parentheses correspond to the case of a complete flavor swap.

Taking alternative SFHs such as the ones from Madau & Dickinson (2014) or Fermi-LAT Collaboration et al. (2018), which we discussed in Section 4, leads to lower predictions of the DSNB flux compared to our fiducial model, which employs the SFH from Mathews et al. (2014). This implies weaker constraints by the experimental limit, as can be seen in Table 6, where the second and third values in parentheses show the results for an SFH according to the EBL reconstruction by Fermi-LAT Collaboration et al. (2018) and for the SFH from Madau & Dickinson (2014), respectively. Independently of the chosen model parameters, the integrated fluxes are reduced as compared to the cases with the SFH from Mathews et al. (2014). The flux values of Φ_>17.3 for the case of the SFH from Madau & Dickinson (2014) lie about one-third below the ones when taking the best-fit SFH from Mathews et al. (2014; roughly corresponding to their −1σ lower limit case), whereas there is still significant overlap between the flux values for the SFH of Fermi-LAT Collaboration et al. (2018) and our fiducial flux values (also note the large uncertainty ranges). At the same time, the values of the slope parameter E₀ are increased (i.e., the spectral tails are lifted) when taking the SFH of Fermi-LAT Collaboration et al. (2018) or the one from Madau & Dickinson (2014) (see Table 3 and Figure 4). Apparently, the large degeneracy between the parameters entering the flux calculations impedes both precise predictions and the exclusion of models.

Finally, we compare our DSNB flux predictions with the results of other recent works. For instance, Priya & Lunardini (2017) found a ${\bar{\nu }}_{{\rm{e}}}$-flux above 11 MeV in the range of (1.4–3.7) cm⁻² s⁻¹, with their highest-flux model being a factor of ∼3 below the SK limit of Bays et al. (2012). In contrast, our fiducial model yields flux values of ${4.6}_{-1.6}^{+3.9}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$ (${3.9}_{-1.3}^{+3.3}\,\,{\mathrm{cm}}^{-2}\,{{\rm{s}}}^{-1}$) above 11 MeV in the case where we consider neutrino oscillations for NH (IH) (to follow Priya & Lunardini 2017), reaching very close to the SK bound (see Section 8.1). Likewise, a recent study by Møller et al. (2018) suggested a clearly lower DSNB flux compared to our result (see their Figures 3 and 10). These differences between our DSNB estimates and previous results can be understood by the large variations of the neutrino outputs between the different core-collapse events in our sets of SN and BH formation models, as shown in Figure 2. While progenitors at the low end of the considered ZAMS mass range radiate ${E}_{\nu }^{\mathrm{tot}}\,\simeq \,2\,\times \,{10}^{53}\,\mathrm{erg}$, the emission increases to values of (3–4) × 10⁵³ erg for progenitors above about (11–12) M_⊙. On the other hand, Priya & Lunardini (2017) and Møller et al. (2018) applied the low-energy neutrino signals (${E}_{\nu }^{\mathrm{tot}}\,\simeq \,2\,\times \,{10}^{53}\,\mathrm{erg}$) of the s11.2c and z9.6co models they considered for the entire mass interval between ∼8 M_⊙ and ∼15 M_⊙, which received a high weight by the IMF in the integration over all core-collapse events. Moreover, both studies made use of failed-SN models that formed BHs relatively quickly (within ≲2 s after bounce) and therefore radiated less energy (≲3.7 × 10⁵³ erg) than most of our failed explosions. Each of these two aspects accounts for the reduction of the integral flux by several tens of percent as compared to our work.

Horiuchi et al. (2018) for the first time employed a larger set of neutrino signals in their DSNB study, including seven models of BH-forming, failed explosions. However, the total neutrino energies ${E}_{\nu }^{\mathrm{tot}}$ radiated from their failed SNe are in general below ∼3.5 × 10⁵³ erg (see their Figure 5). In contrast, we find total neutrino energies in the cases of failed explosions of up to 5.2 × 10⁵³ erg for an NS mass limit of ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}=2.3\,{M}_{\odot }$ and of up to 6.7 × 10⁵³ erg when using our fiducial value, ${M}_{\mathrm{NS},{\rm{b}}}^{\mathrm{lim}}=2.7\,{M}_{\odot }$ (see Figures 2 and B2), enhancing the integral flux by a few tens of percent as compared to Horiuchi et al. (2018). Accordingly, our study suggests that in particular the inclusion of slowly accreting progenitors that lead to late BH formation (not considered in previous works) is responsible for a significant contribution to the DSNB.