Impact of Dark Photon Emission on Massive Star Evolution and Pre-supernova Neutrino Signal

Dark matter, which represents more than 80% (Tanabashi et al. 2018) of the matter density of the universe and whose nature remains one of the biggest mysteries in physics, could be part of a dark sector that weakly interacts with Standard Model (SM) particles. Such scenarios of dark sectors naturally appear in many extensions of the SM where dark matter particles only interact with the SM via a mediator. There is a rich experimental program searching for signatures of such mediators (Essig et al. 2013; Alexander et al. 2016; Aaboud et al. 2019). However, if the coupling to SM particles is too weak, these particles could evade the searches and remain hidden. Astrophysical probes can greatly extend the reach of the search for dark matter candidates, trading the precision associated with the controlled environment of a laboratory for the vast range of densities and temperatures of stars (Raffelt 1996). In one of the simplest extensions of the SM, the dark sector interacts with ordinary matter through the exchange of light vector bosons that couple to SM conserved currents (Holdom 1986; Rajpoot 1989; Nelson & Tetradis 1989; Batell et al. 2014a). Dark matter is charged under a local U(1) symmetry in which the mediator couples to the SM electric charge Q, and is described by the spin-one field ${A}_{\mu }^{D}$, called the dark photon, which mixes kinetically with the standard photon A_μ . Dark photons are characterized by two independent parameters, their mass m_A and the reduced coupling strength with normal matter ε. Several general dark matter searches have established experimental bounds on this parameter space (Essig et al. 2013; Batell et al. 2014b).

Complementary constraints are provided by core-collapse supernovae (SNe). For instance, parameters that would result in a noticeable reduction of the observed duration of the neutrino burst from SN 1987A can be excluded (Raffelt 1990; Rrapaj & Reddy 2016; Choplin et al. 2017; Chang et al. 2017; Sung et al. 2019). Additional constraints have been derived from observational signatures of our Sun and other low-mass stars (An et al. 2013). If particles beyond the SM can be produced at relatively low temperatures, they would affect stellar evolution during He burning and onward. This phase lasts millions of years and is sampled broadly by observations that can be used to derive stringent constraints (Raffelt & Dearborn 1987). Dark matter particles that are more massive, however, can only be produced at much higher temperatures and therefore only affect the short, advanced stages of the evolution of more massive stars.

This work is a continuation of Rrapaj et al. (2019) to study the potential impact of the dark photon on the final stages of a 15 M_⊙ star. During these stages temperatures in the star become high enough for dark photons to be emitted from electron–positron pair annihilation while the densities are still far below the regime where other processes, e.g., bremsstrahlung, become important. We therefore only include dark photons from pair annihilation in this work and consider values of m_A ranging from 2m_e to 10 MeV, where m_e is the electron mass. Due to energy–momentum conservation, dark photons with masses below 2m_e cannot be produced by electron–positron pair annihilation, and the temperatures needed for the emission of particles with masses more than 10 MeV are not reached before core collapse. We explore a wide range of coupling strengths ε, between 10⁻¹³ and 10⁻⁶, thereby probing a part of the parameter space unconstrained by neutrino observations of SN 1987A (Chang et al. 2017; Hardy & Lasenby 2017).

If the dark photon is lighter than all other particles in the dark sector it decays back into electron–positron pairs, either still inside the star or later in the interstellar medium. In the latter case, the γ-ray signature of the subsequent annihilation of positrons can be used to derive constraints (DeRocco et al. 2019). If the dark photons decay inside the star, this may act as an additional heating mechanism leading to further constraints from the observed explosion energies of SNe (Sung et al. 2019). In contrast to these scenarios, we assume that the dark sector contains lighter particles that the dark photon ultimately decays into, as also discussed in Rrapaj et al. (2019). In this case, any energy carried away by the dark photon is leaked into the dark sector, and the emission of the dark matter particles always acts as a cooling mechanism without additional signatures. Under this assumption the constraints discussed by DeRocco et al. (2019) and Sung et al. (2019) do not apply.

As a messenger of the stellar interior, we look at the emission of pre-SN neutrinos and find that, in a large part of our selected dark photon parameter space, the extra cooling leads to a speed-up of the final burning stages and systematically reduces the number of neutrinos emitted during the last day before core collapse. We also show that this reduction can be estimated with good accuracy just based on the baseline stellar model by adjusting the time-integration of the neutrino luminosity. Our results suggest that this effect may be used to constrain the dark photon parameter space if pre-SN neutrinos are detected in the future. In a very small region of the parameter space we also find that the effect of the extra cooling results in a slight increase in the neutrino emission. In addition, unstable and explosive O burning occurs for strong couplings. While the latter case may potentially also provide constraints on the dark photon parameters, improved simulations are required to better determine the final outcome and observational signatures.

Our paper is organized as follows. In Section 2 we describe the setup for our calculations. In Section 3 we discuss the details of our fiducial model of the evolution of a 15 M_⊙ star calculated without extra cooling. In Section 4 we provide an overview of the relevant parameter space and discuss details of the three types of effects that we find. In Section 5 we briefly outline the possibility of deriving constraints from future observations of pre-SN neutrinos and summarize our results.

We implement a tabulation of the dark photon emission rates from Rrapaj et al. (2019) in the stellar evolution and hydrodynamics code KEPLER (Weaver et al. 1978; Woosley & Heger 2007) and calculate the evolution of a 15 M_⊙ star with an initial composition of solar metallicity (Lodders & Palme 2009). The mass-loss prescription is based on Nieuwenhuijzen & de Jager (1990) and the mixing length for convection is equal to the pressure scale height. Semi-convection is treated as described in Woosley & Weaver (1988), limiting the convective diffusion coefficient to 10% of the thermal value (see also Woosley et al. 2002). Overshoot and thermohaline mixing (Kippenhahn et al. 1980) are also included. Neutrino energy loss and the resulting total neutrino luminosity are based on Itoh et al. (1996).

We evolve the models until the onset of core collapse, which we define as the point when the infall velocity exceeds 5000 km s^–1. This limiting value is higher than the value of 1000 km s^–1 used in previous studies (Woosley et al. 2002; Woosley & Heger 2007) because the additional dark matter cooling tends to accelerate contraction.

In order to explore the dark photon parameter space we have calculated more than 470 stellar models, and to limit the computational cost we only use an approximate 19-isotope nuclear reaction network to calculate the energy generation rates. At a temperature of 3.5 GK the code switches to solving for the composition in quasi-statistical equilibrium (QSE) (Bodansky et al. 1968; Woosley et al. 1973). For even higher temperatures and when O is depleted full nuclear statistical equilibrium (NSE) is assumed. Furthermore, we limit the number of zones in our models to 2000.

Figure 1 shows the track of the fiducial stellar model without extra cooling in terms of central temperature T_C and density ρ_C. This track is in good agreement with the results from Woosley & Heger (2007) and also with those from different stellar evolution codes, such as MESA (Paxton et al. 2015). For the parameters that we are studying, dark photon production is not relevant for temperatures below ∼1 GK. The dotted contours in Figure 1 indicate where energy loss due to dark photon emission equals the neutrino loss, marking the boundary of the temperature and density regime where dark photon loss would dominate. This boundary depends on the values of m_A and ε. For instance, lighter dark photons with higher coupling strengths start to dominate at lower temperatures. The details of these contours were discussed by Rrapaj et al. (2019). As illustrated for m_A = 2 MeV and ε = 10⁻¹¹ in Figure 1, dark photons may already become important right after central C depletion, which occurs around 64 years before core collapse. Therefore, for context, we provide a short description of the evolution after C depletion.

Zoom In Zoom Out Reset image size

The convective burning phases over the last 60 yr are indicated in Figure 2, where hatched areas indicate convection and the color code shows the net energy generation or loss. After C core and shell burning, Ne burning ignites at the center at T_C = 1.35 GK, around 1500 days before collapse under slightly degenerate conditions with a value of the degeneracy parameter η = μ_e/k_B T ≈ 4, where μ_e is the chemical potential of electrons, T the temperature, and k_B the Boltzmann constant. The partial electron degeneracy leads to an initially rapid rise of the temperature, followed by core expansion at almost constant temperature. Once Ne is depleted, neutrino cooling leads to a decrease in temperature before the core continues to contract and eventually heats up again. This leads to the loop in the T_C–ρ_C diagram shown in Figure 1. Ne burning is visible as the small hatched peak in Figure 2, which also shows that this episode is very short.

Zoom In Zoom Out Reset image size

About 650 days before collapse, O burning ignites centrally at T_C = 1.7 GK, which leads to a rapid rise of the central temperature and subsequent expansion. In contrast to the Ne-burning loop, at the end of O burning in the core the temperature does not decrease and the core is stabilized by shell burning. Once the shell burning ceases, neutrino loss reduces T_C while the core contracts. At higher densities, neutrino loss is suppressed and the track resumes its upward climb, while a second O-burning shell is ignited. Between core O depletion and Si ignition, He and C shell burning are still active.

Convective Si burning ignites at the center around 3.5 days before collapse, after the second O-burning shell is extinguished. Around 5.6 hours before collapse, central Si burning finishes and leaves a hot 1.1 M_⊙ Fe core behind. This initial core grows further by shell burning that ignites about one hour before collapse and continues until the core mass exceeds its effective Chandrasekhar limit and collapses. Figure 3 shows the profiles of temperature, density, and mass fraction of the most important isotopes of our fiducial model at core collapse. The He core encompasses 4.26 M_⊙, the C/O core 3.02 M_⊙, the O/Ne core up to 2.70 M_⊙, and the Fe core 1.60 M_⊙ at a radius of 1336 km.

Zoom In Zoom Out Reset image size

Multidimensional hydrodynamics simulations of the last minutes before core collapse (Müller et al. 2016; Yoshida et al. 2019; Yadav et al. 2020) have shown deviations from spherically symmetric models using mixing length theory, but we do not expect our results to be qualitatively affected by these differences.

Because the evolution of the stellar core is mostly decoupled from the surface during the final stages, neutrinos are unique messengers that may provide detailed information about the processes and conditions in the core shortly before collapse (e.g., Guo & Qian 2016; Kato et al. 2020). Current and near-term neutrino detectors are expected to be able to detect the neutrinos from a nearby massive star only within the last day before collapse (e.g., Guo et al. 2019; Kato et al. 2020). The top panel of Figure 4 shows the neutrino luminosity during the last 10 days before core collapse for our fiducial model. In general, the luminosity increases as the star contracts and heats up, following the track in Figure 1. The two peaks visible in Figure 4 are caused by the ignition of nuclear burning that leads to expansion and cooling, temporarily delaying collapse. The peak at around 3.5 days before collapse corresponds to Si ignition at the center, and the peak at about one hour before collapse corresponds to the ignition of Si shell burning.

Zoom In Zoom Out Reset image size

Here we focus on the ${\bar{\nu }}_{{\rm{e}}}$ for two reasons. First, future scintillation detectors are likely to observe the pre-SN ${\bar{\nu }}_{{\rm{e}}}$, mainly by inverse beta decay (IBD, ${\bar{\nu }}_{{\rm{e}}}+{\rm{p}}\to {{\rm{e}}}^{+}+{\rm{n}}$). Second, towards the end of the life of a massive star ${\bar{\nu }}_{{\rm{e}}}$ are mostly produced by electron–positron pair annihilation. In contrast to the emission of ν_e, this thermal process is relatively insensitive to the details of the stellar composition and does not depend on the uncertainties related to electron capture on nuclei. We calculate the spectral neutrino flux from pair annihilation as in Guo & Qian (2016), using profiles of temperature, density, and electron fraction as functions of time from our stellar models. The ${\bar{\nu }}_{{\rm{e}}}$ luminosity due to pair annihilation alone is shown in the top panel of Figure 4 in comparison to the total neutrino luminosity. At one day before collapse, ${\bar{\nu }}_{{\rm{e}}}$ constitute almost 30% of the total luminosity. This fraction decreases towards collapse, as ν_e from electron captures become increasingly important. At 10 s before collapse, however, ${\bar{\nu }}_{{\rm{e}}}$ from pair annihilation still account for 10% of the total neutrino luminosity.

Based on the ${\bar{\nu }}_{{\rm{e}}}$ fluxes from the stellar evolution models, we calculate the expected pre-SN neutrino signal at the Jiangmen Underground Neutrino Observatory (JUNO) following Guo et al. (2019). With the spectral number luminosity ${\phi }_{{\bar{\nu }}_{{\rm{e}}}}({E}_{\nu },t)$ of the star, the expected energy-differential event rate at a distance d is

$$\begin{eqnarray}&&\displaystyle \frac{{d}^{2}N}{{{dE}}_{\nu }{dt}}=\displaystyle \frac{1}{4\pi {d}^{2}}{\epsilon }_{\mathrm{eff}}{N}_{p}{\sigma }_{\mathrm{IBD}}({E}_{\nu }){\phi }_{{\bar{\nu }}_{{\rm{e}}}}({E}_{\nu },t),\end{eqnarray} \tag{ 1 }$$

where _eff = 0.73 is the detector efficiency and N_p = 1.45 × 10³³ is the number of protons based on 20 kt of active detector material with a 12% proton fraction (An et al. 2016). The IBD cross section σ_IBD(E_ν ) is calculated as in Guo & Qian (2016). For the event rate we integrate Equation (1) over E_ν = 1.8–4 MeV, which is the optimal energy window for detection. Due to flavor transformations the detection rate depends on the neutrino mass ordering, and the expected number of events for the normal hierarchy is about 3.4 times higher than for the inverted hierarchy (e.g., Guo et al. 2019). Figure 4 shows the expected event rates from our fiducial model assuming a distance of 500 pc. We also show the main sources of background, geo-neutrinos and ${\bar{\nu }}_{{\rm{e}}}$ from nearby reactors. For the inverted hierarchy about 15 events are expected during the last day, which is approximately the same as the background. For the normal hierarchy, however, about 50 events are expected during the last day, which is a factor of ∼3 more than the background. Therefore, a pre-SN neutrino signal may be detected above the background during the final day before core collapse, especially in the case of the normal mass hierarchy.

We have calculated stellar models for a grid of the dark photon parameters m_A and ε laid out in Figure 5. For 3 MeV < m_A < 10 MeV we have looked at increments of 0.5 MeV and for 2m_e < m_A ≤ 3 MeV we have taken smaller steps of 0.1 MeV. For ε we looked at logarithmic intervals of 0.5 dex spanning values from 10⁻¹³ to 10⁻⁶. For reference, the region of the parameter space that is already excluded by the detection of the neutrinos from SN 1987A (Chang et al. 2017) is indicated as the top gray shaded area in Figure 5. Note also that the region corresponding to m_A < 2m_e is not studied here. The different symbols in Figure 5 indicate the outcomes of the models that we have calculated.

Zoom In Zoom Out Reset image size

We find three qualitatively different ways in which the extra cooling changes the late phases of stellar evolution and the last-day pre-SN neutrino signal. The green squares in Figure 5 mark models that exhibit only negligible deviations from the fiducial model. In contrast, the blue circles mark cases in which the pre-SN neutrino emission is reduced by more than 30% due to accelerated Si burning. Part of the parameter space with the blue circles is currently unconstrained, suggesting that future observations of pre-SN neutrinos may impose new constraints on the properties of dark photons.

Only three models with small dark photon masses, indicated by the red crosses in Figure 5, exhibit a slight increase in the pre-SN neutrino emission. As we will show in Section 4.2, this result is caused by an additional shell-burning episode that delays core collapse. For sufficiently small dark photon masses and intermediate coupling strengths, the dark photon emission predominantly originates not from the stellar center, but from the region of the final the O/Ne shell.

For the parameters in the red shaded region in the upper left corner of Figure 5, the nuclear timescale exceeds the hydrodynamic timescale during O burning, and thermonuclear runaway is to be expected, as we will explain in Section 4.3. We did not perform full calculations to cover this whole region because of the complications associated with this explosive scenario. While the final fate of the corresponding models requires further investigation beyond the scope of this paper, we can speculate that all these cases may even lead to thermonuclear SNe and the disruption of the whole star.

In general, energy loss determines the timescales of the advanced burning stages in the evolution of massive stars (Woosley & Heger 2007). The balance between energy generation due to nuclear burning and energy loss determines the temperature of hydrostatic burning and thus the rates at which nuclear reactions occur. Furthermore, energy loss determines how fast the stellar core can contract because part of the resulting gain in gravitational binding energy needs to be radiated away. Thus, emission of dark matter particles that increases the energy loss rate generally leads to a speed-up of stellar processes. The consequences of this acceleration, however, depend on where the speed-up occurs and on whether the star can adjust sufficiently fast. Below we give more detailed explanations for the various effects that we observe.

4.1. Reduced Neutrino Emission

As long as the stellar path to core collapse is guided by energy loss, an additional cooling agent, such as dark photon emission, leads to a speed-up of the evolution. Once the core becomes hot enough, the bulk of the energy loss is taken over by the dark matter particles instead of the neutrinos, and the accelerated evolution leaves less time for neutrino emission.

As an example, Figure 6 shows the convective burning phases for a calculation with the parameters m_A = 2 MeV and ε = 10⁻¹². Comparison with Figure 2 shows that the stellar structure remains largely the same as in the fiducial model, but the evolution timescales change greatly. The period between the onset of convective Si burning and collapse is only around 0.9 day in Figure 6, compared to 3.5 days in the fiducial model (see Figure 2). That between the end of convective Si burning at the center and collapse is less than 3 hr, compared to 8 hr in the fiducial model. The radial extent of the convective regions, on the other hand, remains almost identical. The slight reduction of the extent of the convective zone due to the more efficient energy loss, along with the shorter time for shell burning due to the faster evolution of the core, results in a slightly smaller final Fe core mass of 1.58 M_⊙, compared with 1.60 M_⊙ for the fiducial model. With the remnant mass estimated by the point where the entropy decreases below 4 k_B/baryon (Woosley et al. 2002), dark photon emission slightly reduces the potential neutron star mass from 1.78 M_⊙ to 1.76 M_⊙.

Zoom In Zoom Out Reset image size

For the above values of m_A and ε, the overall changes to the stellar structure relative to the fiducial model are minor. In contrast, the integrated number of neutrinos emitted during the last day is already decreased by a factor of two. This result shows that the main effect of the additional cooling is a speed-up of Si burning and core contraction that scales with the energy loss rate. This speed-up of the evolutionary clock reduces the integrated neutrino emission, with dark photons also contributing to the energy loss. The effects on the ${\bar{\nu }}_{{\rm{e}}}$ luminosity from pair annihilation and the associated signal are illustrated in Figure 7 for models with m_A = 2 MeV and a range of ε values. In addition to a suppression of the luminosity, some features of the signal are shifted closer to the time of collapse. This shift is clearly visible for the peak associated with Si shell burning, which appears at about one hour before collapse for the fiducial model (see Section 3), but at only 30 minutes and 12 minutes before collapse for ε = 10⁻¹² and 10^−11.5, respectively. Note that for ε = 10⁻¹¹, the evolution of the ${\bar{\nu }}_{{\rm{e}}}$ luminosity becomes less smooth due to complications discussed in Section 4.3.

Zoom In Zoom Out Reset image size

For all dark photon masses covered in our study, we find the same type of reduction of the neutrino emission, and the necessary coupling strength to achieve the same reduction increases with mass. In the following we show that this reduction can be estimated semi-analytically, just based on the fiducial stellar model. Without dark photon cooling, the timescales during the final stages of evolution are determined mostly by the neutrino loss. We expect dt ∝ L_ν (t)⁻¹, where L_ν (t) is the total neutrino luminosity from the whole star at time t. With the additional dark photon cooling, we expect

$$\begin{eqnarray}&&{dt}^{\prime} =\displaystyle \frac{{L}_{\nu }(t)}{{L}_{\nu }(t)+{L}_{A}(t)}{dt},\end{eqnarray} \tag{ 2 }$$

where L_A (t) is the dark photon luminosity from the whole star. This change affects the results of any time-integrated quantity in two ways. First, the time measure itself is changed. Second, the effective integration limits are changed. The time of core collapse, t_CC, marks the final time with a well-defined physical condition. In order to accumulate a fixed time interval Δt, e.g., one day, up to t_CC, the starting point needs to be modified such that

$$\begin{eqnarray}&&{\rm{\Delta }}t={\int }_{{t}_{0}}^{{t}_{\mathrm{CC}}}{dt}={\int }_{{t}_{0}^{\prime} }^{{t}_{\mathrm{CC}}}{dt}^{\prime} .\end{eqnarray} \tag{ 3 }$$

With ${dt}^{\prime} \lt {dt}$ due to dark photon cooling, we have ${t}_{0}^{\prime} \lt {t}_{0}$, so Δt samples earlier stages of the evolution.

The above prescription allows us to estimate the energy emitted in ${\bar{\nu }}_{{\rm{e}}}$ from pair annihilation during the last day before collapse as

$$\begin{eqnarray}&&{E}_{1\mathrm{day}}^{\prime} ={\int }_{{t}_{0}^{\prime} }^{{t}_{\mathrm{CC}}}{L}_{{\bar{\nu }}_{{\rm{e}}},\mathrm{pair}}(t)\left(\displaystyle \frac{{L}_{\nu }(t)}{{L}_{\nu }(t)+{L}_{A}(t)}\right){dt},\end{eqnarray} \tag{ 4 }$$

where ${L}_{{\bar{\nu }}_{{\rm{e}}},\mathrm{pair}}(t)$ is the ${\bar{\nu }}_{{\rm{e}}}$ luminosity due to pair annihilation from the whole star. Note that L_ν (t) and L_A (t) are based on the temperature and density profiles of the fiducial model without dark photon emission. This approximation is valid as long as the effects of extra cooling on stellar structure are negligible. In this case, the change in the ${\bar{\nu }}_{{\rm{e}}}$ emission can be estimated semi-analytically without computing new stellar models with the additional energy loss.

Figure 8 displays the energy emitted in ${\bar{\nu }}_{{\rm{e}}}$ from pair annihilation over the last day before core collapse, normalized to the fiducial model, as a function of the dark photon coupling strength for a range of m_A . The symbols show the results from full stellar evolution calculations and demonstrate that increasing dark photon emission leads to a reduction in the neutrino emission. The solid curves show the estimates based on Equation (4), which successfully describe the onset of the effects of dark photons as well as the overall trend for many of the calculated models, especially when the reduction of the neutrino emission is ≲50%. Deviations from the actual calculations are expected when the dark photon cooling starts to change the temperature and density profiles for sufficiently large values of ε. For a large part of the dark photon parameter space, however, the main effect is the relatively straightforward reduction of the neutrino emission without significant changes to the stellar structure as discussed above. Those cases are indicated by the blue circles in Figure 5 and cover an unconstrained region of the parameter space.

Zoom In Zoom Out Reset image size

Figure 8 also shows large deviations of the estimates from the actual calculations for m_A = 1.2 MeV. As we will explain in Section 4.2, these deviations are due to changes in the shell-burning evolution.

4.2. Increased Neutrino Emission

As discussed in Section 4.1, the reduction in the pre-SN neutrino emission results almost entirely from the modification of the timescales of core burning and contraction. Below we show that the deviations from this explanation for low m_A (see Figure 8) are due to the appearance of shell-burning phases with a stabilizing effect, which leads to the cases of the increased neutrino emission marked by the red crosses in Figure 5.

For the dark photon parameters under our consideration, only three models with m_A = 2m_e, 1.1 MeV, and 1.2 MeV, all with ε = 10⁻¹², have an increase of ≲50% in the pre-SN neutrino emission relative to the fiducial model. The comparison of the neutrino emission with the fiducial model is shown for the case of m_A = 1.2 MeV in Figure 8. The Kippenhahn diagram for this case is shown in Figure 9, which reveals an additional episode of convective O shell burning between an enclosed mass of 1.6 M_⊙ and 2 M_⊙. This episode starts 0.25 days before collapse, following the end of core Si burning, but before the onset of Si shell burning. Such an episode occurs neither in the fiducial model nor in the model with the reduced neutrino emission shown in Figure 6. Heating by this additional O-shell-burning phase relieves the pressure on the Si core and delays the ignition of Si shell burning and hence core collapse. This delay gives the star more time, thereby increasing the total number of pre-SN neutrinos emitted.

Zoom In Zoom Out Reset image size

For m_A = 1.2 MeV, the additional O-shell-burning phase also occurs for larger values of the dark photon coupling, but it cannot compensate for the effects of the accelerated burning and core contraction. Consequently, the integrated neutrino emission for such cases is still reduced (see Figure 8), although not by as much as estimated from Equation (4), which does not take into account the additional O shell burning. We can only find a net increase in the pre-SN neutrino emission for a very narrow range of ε, which is strong enough to cause an additional burning episode but weak enough to limit the acceleration of the evolution.

In order to understand the apparent association of the additional burning episodes with the lightest dark photons in our study, we look at the time-integrated energy loss as a function of the enclosed mass M_r inside radius r,

$$\begin{eqnarray}&&{E}_{\mathrm{loss},x}({M}_{r})={\int }_{{t}_{0}}^{{t}_{{CC}}}{\dot{q}}_{x}({M}_{r},t){dt},\end{eqnarray} \tag{ 5 }$$

where x is either ν for the neutrino with the energy loss rate ${\dot{q}}_{\nu }({M}_{r},t)$ or A for the dark photon with the energy loss rate ${\dot{q}}_{A}({M}_{r},t)$, and t₀ is the time of C ignition when the energy loss due to neutrinos and dark photons becomes relevant. Figure 10 shows the ratio R_A/ν (M_r ) = E_loss,A (M_r )/E_loss,ν (M_r ) relative to the value for the central zone R_A/ν (M_r,0) based on the fiducial model. Due to the scaling of R_A/ν (M_r ) with ε, the quantity R_A/ν (M_r )/R_A/ν (M_r,0) is independent of ε. It is, however, quite sensitive to m_A .

Zoom In Zoom Out Reset image size

For masses above 2 MeV the emission of dark photons dominates that of neutrinos mostly in the final Fe core because the temperature is generally higher at smaller radii. For masses below 2 MeV, however, the dark photon emission dominates further outside the core. This different behavior for the lower values of m_A can be understood as follows. Because the number of e^± pairs is reduced in degenerate conditions, dark photon production by pair annihilation is suppressed at high densities in the core (see Rrapaj et al. 2019 for details). This suppression does not occur in the region at lower densities outside the core, but the temperature of this region can only facilitate significant emission of the lower-mass dark photons. Therefore, shell-burning phases are accelerated and additional burning episodes with stabilizing effects as described above tend to occur for m_A < 2 MeV. For heavier dark photons, their production requires energetic e^± pairs that can only be provided by the final core-burning stages. Therefore, the assumption that the effect of dark photon cooling is limited to the core-burning regions is justified for m_A > 2 MeV, but the evolution of the outer regions is noticeably affected for smaller dark photon masses.

The stabilizing effect of additional shell-burning phases due to the accelerated evolution outside the core cannot be easily captured by a semi-analytical prescription like the one presented in Section 4.1. Detailed stellar evolution calculations are essential to the discovery and understanding of this effect.

4.3. Runaway O Burning

The stellar models with those dark photon parameters indicated by the red stars in Figure 5 exhibit a thermonuclear runaway during O burning. This result may be understood as follows. With the dark photon cooling, the accelerated burning timescales get shorter than the convective and eventually also the hydrodynamic timescales. In addition, the increased electron degeneracy further prevents the star from adjusting its structure in response to the nuclear energy release. Consequently, the central temperature rises rapidly to ∼5 GK and O is exhausted in a convective region covering less than the innermost 0.1 M_⊙ of material. Afterwards, O burning propagates outward in a thin shell as a deflagration front. Figure 11 shows the evolutionary tracks of T_C and ρ_C for m_A = 2 MeV and ε = 10⁻¹², 10⁻¹¹, 10^−10.5, and 10⁻¹⁰, in comparison with the fiducial model. The degeneracy parameter η is also shown in the bottom panels of this figure. The onset of the runaway can be seen clearly for the largest two values of ε: T_C increased almost vertically to ∼5 GK when we stopped the calculations due to the required very small time steps.

Zoom In Zoom Out Reset image size

Compared to the fiducial model, the tracks for ε = 10⁻¹² in Figure 11 do not show much deviation, but the density tends to be higher for the same temperature, thereby increasing the degeneracy. While the track of T_C versus ρ_C for the fiducial model already exhibits a loop, which arises from burning under degenerate conditions (see also Figure 1), the track for ε = 10⁻¹¹ shows two such loops and several backward kinks resulting from shell flashes, signifying ignition of shell burning under degenerate conditions. For ε = 10^−10.5, O burning becomes unstable and the central temperature explosively increases up to ∼5 GK. It can be seen from Figure 11 that in this case the core expands somewhat initially, but the decrease in density fails to stop the runaway heating. As T_C increases rapidly, convection also quickly becomes inefficient at distributing the energy released from the nuclear reactions. The detailed evolution in this case is discussed in the Appendix. We find the same type of behavior as in this case for even larger values of ε. For example, for ε = 10⁻¹⁰, the density response to the rise in temperature due to O burning is barely visible in Figure 11, indicating a further acceleration of the runaway process.

While the runaway burning described above results from the interplay of several factors, including the competition of burning timescales with convective and hydrodynamic timescales, as well as partial electron degeneracy, below we present a simple argument for the association of a thermonuclear runaway with large values of ε and provide a conservative estimate for the dark photon luminosity (per unit mass) at which hydrostatic burning is no longer possible. In configurations of stable nuclear burning, the energy release from nuclear reactions is compensated by energy losses, leading to a stable consumption of the nuclear fuel. For advanced burning phases without dark photons, the temperature for hydrostatic burning can be estimated by equating the neutrino energy loss and nuclear energy release (Woosley et al. 2002). Assuming

$$\begin{eqnarray}&&\rho ={10}^{6}{\left(T/{\rm{G}}{\rm{K}}\right)}^{3}\,{{\rm{g}}{\rm{c}}{\rm{m}}}^{-3}\end{eqnarray} \tag{ 6 }$$

for O burning, the intersection of the O-burning energy generation and the neutrino loss rate gives the typical O-burning temperature of ≈1.8 GK (see blue curves in Figure 12).

Zoom In Zoom Out Reset image size

The additional energy loss due to dark photons changes the temperature for hydrostatic O burning. Due to the very steep temperature dependence of the nuclear reaction rate, a small change in temperature is sufficient to compensate for a large increase in the energy loss. As shown in Figure 12, a factor of ≈100 increase in the energy loss for m_A = 2 MeV and ε = 10⁻¹¹ only slightly increases the nominal O-burning temperature to ≈2 GK. Nevertheless, this temperature increase is responsible for the speed-up of the evolution discussed in Section 4.1. Of course, the burning temperature must increase substantially for strong dark photon couplings. For example, it reaches ≈3 GK for ε = 10⁻⁹ (see Figure 12).

Complications arise if the timescale τ_n of nuclear burning becomes shorter than the hydrodynamic timescale on which stellar structure can adjust to the temperature change. We take the latter timescale to be the free-fall timescale

$$\begin{eqnarray}&&{\tau }_{\mathrm{ff}}=\displaystyle \frac{1}{\sqrt{G\rho }},\end{eqnarray} \tag{ 7 }$$

where the density ρ can be estimated from the temperature assuming Equation (6). The timescale of nuclear burning can be estimated as

$$\begin{eqnarray}&&{\tau }_{{\rm{n}}}={XQ}/{\dot{q}}_{\mathrm{nuc}},\end{eqnarray} \tag{ 8 }$$

where Q = 5 × 10¹⁷ erg g^–1 (Woosley et al. 2002) is the effective energy release, X = 0.7 is the mass fraction of the ¹⁶O fuel, and ${\dot{q}}_{\mathrm{nuc}}$ is the specific energy generation rate per unit mass. For O burning, τ_n = τ_ff occurs at a critical specific energy generation rate of

$$\begin{eqnarray}&&{\dot{q}}_{{\rm{n}}{\rm{u}}{\rm{c}},{\rm{c}}{\rm{r}}{\rm{i}}{\rm{t}}}\approx 9\times {10}^{16}{(T/{\rm{G}}{\rm{K}})}^{3/2}\,{\rm{e}}{\rm{r}}{\rm{g}}\,{{\rm{g}}}^{-1}\,{{\rm{s}}}^{-1}.\end{eqnarray} \tag{ 9 }$$

The contour for τ_n = τ_ff is shown in Figure 12. It intersects the curve for the energy release due to O burning at T ≈ 3 GK and ${\dot{q}}_{{\rm{nuc}},{\rm{crit}}}\approx 4.6\times {10}^{17}\,{\rm{erg}}\,{{\rm{g}}}^{-1}\,{{\rm{s}}}^{-1}$. At this temperature and the corresponding density from Equation (6), if the specific dark photon energy loss rate exceeds ${\dot{q}}_{\mathrm{nuc},\mathrm{crit}}$, i.e.,

$$\begin{eqnarray}&&{\dot{q}}_{A}(T=3\,{\rm{G}}{\rm{K}})\gt 4.6\times {10}^{17}\,{\rm{e}}{\rm{r}}{\rm{g}}\,{{\rm{g}}}^{-1}\,{{\rm{s}}}^{-1},\end{eqnarray} \tag{ 10 }$$

a thermonuclear runaway is expected. As shown in Figure 12, this result occurs for m_A = 2 MeV and ε = 10⁻⁹.

The red shaded region in Figure 5 indicates the general parameter space of m_A and ε for which the condition in Equation (10) is fulfilled. While the boundary of the region roughly reflects the appearance of the red stars indicating runaway O burning, the stellar models already exhibit such burning for much weaker couplings for m_A ≲ 2 MeV. This result may be explained by the complications due to electron degeneracy, which tends to prevent the star from reacting to nuclear energy release, and by the competition of nuclear burning with convection. Both these complications are ignored in the simple picture assumed for deriving Equation (10). At O ignition, our fiducial model already exhibits some degree of degeneracy with η ∼ 2 (see Figure 11). When dark photon cooling is included, the temperature of stable burning increases, which requires the stellar core to contract to a higher density. Because degeneracy rises more steeply with density than the temperature, values up to η ∼ 8 are reached in the core at O-burning temperatures when O burning becomes unstable. As the adjustment of the core is slowed down, the O-burning luminosity exceeds the maximum luminosity that can be transported by convection (see the Appendix), which also favors a local runaway.

The boundary of the red shaded region in Figure 5 reflects the dependence of the specific dark photon energy loss rate on m_A and ε. It also tends to predict too small threshold values of ε for runaway burning for higher values of m_A . Specifically, for (m_A /MeV, ε) = (4, 10⁻⁷), (4.5, 10^−6.5), (4.5, 10⁻⁶), and (5, 10⁻⁶), which are above the boundary, although O burning leads to a rapid rise in temperature, the dark photon energy loss increases fast enough to avoid a thermonuclear runaway.

We have studied the impact of extra cooling due to particles beyond the SM on the evolution of a 15 M_⊙ star. Specifically, we have implemented the dark photon emission from e^± pair annihilation in the stellar evolution code KEPLER, assuming that these particles decay into other dark sector components, thereby representing an additional mechanism of energy loss for the star. We have considered dark photon masses of m_A = 2m_e to 10 MeV and couplings of ε = 10⁻¹³ to 10⁻⁶, and found that the dark photon can affect the O- and Si-burning phases. There are three types of potentially observable effects, which are summarized in Figure 5.

For broad ranges of m_A and ε, the extra cooling allows the star to contract faster and increases the temperature at which nuclear burning proceeds in equilibrium with energy losses. These effects speed up the burning processes and reduce the number of neutrinos emitted during the last day before core collapse. We have developed a semi-analytical approach that describes the reduction of the pre-SN neutrino emission for the relevant ranges of m_A and ε (see Section 4.1).

We further find that, due to the density dependence of the emission process, dark photons with m_A < 2 MeV are also produced outside the core, which causes additional shell-burning episodes, thereby increasing the pre-SN neutrino emission for a very narrow range of parameters (see Section 4.2).

Because extra cooling increases the nominal nuclear burning temperatures, sufficiently strong dark photon couplings are expected to produce a thermonuclear runaway (see Section 4.3). We have found many cases of runaway O burning in our models and speculate that thermonuclear SNe may be the outcome in these cases (see the Appendix). Finding the actual outcome and possible observables in these cases, however, requires more refined simulations that can properly follow the propagation of narrow nuclear burning fronts. Were dark photon cooling to cause complete disruption of a star, the relevant parameter space may be constrained by the observed inventory of neutron stars and stellar black holes, which requires most massive stars to undergo core-collapse SNe.

The first type of effect, i.e., the suppression of the pre-SN neutrino emission by dark photon cooling, may be constrained by the detection of such neutrinos. In the most conservative approach, this suppression means that part of the parameter space is only consistent with a nondetection of pre-SN neutrinos. Therefore, any positive detection of a pre-SN neutrino signal from a well-understood progenitor excludes this part of the dark photon parameter space. For illustration, we take our 15 M_⊙ model to explode as a core-collapse SN at a distance of 500 pc. Assuming the normal neutrino mass hierarchy, we can find the values of (m_A , ε) that reduce the pre-SN neutrino emission sufficiently to give an expected number of events less than the background count of about 17 events, in contrast to about 50 events without dark photon cooling. The corresponding parameter space is bounded by the red curve in Figure 13. The detection of a number of events significantly above the background level would exclude this part of the parameter space, if the progenitor is similar to the model we have studied here. Note that a small fraction of this parameter space may lead to a thermonuclear runaway, and therefore could be constrained by the observation of a core-collapse rather than thermonuclear SN.

Zoom In Zoom Out Reset image size

Pre-SN neutrinos have not been detected yet. Our results, however, indicate that their observation in the future may offer a unique probe of not only the dark photon but also other dark matter particles that are efficiently produced in the temperature regime of ∼1–10 GK during stellar evolution.

This work was supported in part by the US Department of Energy [DE-FG02-87ER40328 (UM)]. Calculations were carried out at the Minnesota Supercomputing Institute. E.R. was supported by the NSF (PHY-1630782) and the Heising–Simons Foundation (2017-228). G.G. acknowledges support from the Academia Sinica by Grant No. AS-CDA-109-M11. We thank Alexander Heger for providing access to the KEPLER code.

In order to understand the phenomenon of runaway O burning, we investigate the model for m_A = 2 MeV and ε = 10^−10.5 in more detail. Figure 14 shows a sequence of snapshots from the evolution of this model. The top row shows the temperature and density profiles, the second row shows the nuclear energy release and combined energy loss due to neutrinos and dark photons, the third row shows the mass fractions of the major isotopes to illustrate the progress of nuclear burning, and the bottom row shows the velocity profile, where positive (outward) velocities are shown as solid curves and negative (infall) velocities are represented by dashed curves. The red shaded regions are convective zones. Time runs from panel (a) on the left to panel (c) on the right, spanning 357 s. Panel (a) shows the onset of O burning, which is confined to a convective region reaching up to 1000 km and encompassing about 0.1 M_⊙ of material. At T ≈ 2 GK and ρ ≈ 8 × 10⁷ g cm^–3, the energy release from nuclear reactions exceeds the combined neutrino and dark photon loss. In response to the excess heat, positive velocities develop, indicating the onset of core expansion. Panel (b), however, shows that O burning proceeds much faster than the adjustment of stellar structure. The fuel in the center is almost depleted and the temperature has risen to 4.7 GK while the density has only been reduced by about 40% compared to panel (a).

Zoom In Zoom Out Reset image size

At a temperature above 4 GK, the luminosity due to O burning exceeds the maximum value that can be transported by convection (Woosley & Heger 2015),

$$\begin{eqnarray}&&{L}_{\max }\approx 4\pi {r}^{2}\rho {v}_{\mathrm{conv}}{{fC}}_{{\rm{P}}}T,\end{eqnarray} \tag{ A1 }$$

where v_conv is the convective velocity, C_P T is the thermal energy content, and f ≪ 1 indicates the efficiency of convection to remove the internal energy. With f = 0.1 we find ${L}_{max}\lt {10}^{50}\,{\rm{erg}}\,{{\rm{s}}}^{-1}$, corresponding to an energy generation rate of ${\dot{q}}_{max}\approx 5\times {10}^{17}\,{\rm{erg}}\,{{\rm{g}}}^{-1}\,{{\rm{s}}}^{-1}$ for the conditions of O burning assuming an initial core mass of 0.1 M_⊙. The nuclear energy release in panel (b) exceeds this luminosity already by four orders of magnitude and convection cannot distribute the heat and the fuel throughout the convective region on the burning timescale. Therefore, O is only depleted in the innermost 300 km although the convective zone reaches out to 2000 km.

Without efficient convection, the very rapid burning gets confined into a narrow burning front at the bottom of the nominally convective layer, as shown in panel (c) of Figure 14. The change in composition and continuing energy generation also show that Si burning immediately follows because of the high temperature of the central region. In our model, the burning region consists only of a few zones and the fuel in one zone is consumed before the zone above it ignites. The propagation of the burning depends critically on the heat transport and would need to be described by a model for flame propagation. For a very narrow burning front, multidimensional effects are important for the propagation of such a deflagration (Fryxell & Woosley 1982) and the outcome is highly sensitive to turbulence and mixing at the boundary of the burning front (Jones et al. 2013, 2016). The bottom row of Figure 14 also shows that significant outward velocities develop but remain subsonic. The above situation cannot be adequately simulated by our model due to the lack of resolution and the employed small reaction network and mixing length treatment of convection. Therefore, our results on the final outcome are speculative.

We have followed the calculation of our model further and find that the shell burning steepens into a supersonic shock once it reaches the steep density gradient at the upper edge of the convective zone around the radius of 2000 km, where Ne burning also provides additional energy. The outward velocity reaches several thousand km s^–1, which takes the pressure off the core, eventually facilitating a rapid expansion and cooling. T_C drops below 0.1 GK and ρ_C below 10⁵ g cm^–3. When the shock reaches the H envelope we find that the kinetic energy of 3.6 × 10⁵⁰ erg exceeds the gravitational binding energy of 1.4 × 10⁵⁰ erg for the star. The shock may thus disrupt the whole star or at least unbind a significant fraction of the H envelope and lead in either case to an optical transient. Even if the core eventually collapses, the runaway O burning would delay the collapse by many hours up to days. Here we caution again that our model is probably inadequate to describe this situation accurately and a better treatment is needed to predict potentially observable signatures.