Revisiting Supernova 1987A Constraints on Dark Photons

We revisit constraints on dark photons with masses below ~ 100 MeV from the observations of Supernova 1987A. If dark photons are produced in sufficient quantity, they reduce the amount of energy emitted in the form of neutrinos, in conflict with observations. For the first time, we include the effects of finite temperature and density on the kinetic-mixing parameter, epsilon, in this environment. This causes the constraints on epsilon to weaken with the dark-photon mass below ~ 15 MeV. For large-enough values of epsilon, it is well known that dark photons can be reabsorbed within the supernova. Since the rates of reabsorption processes decrease as the dark-photon energy increases, we point out that dark photons with energies above the Wien peak can escape without scattering, contributing more to energy loss than is possible assuming a blackbody spectrum. Furthermore, we estimate the systematic uncertainties on the cooling bounds by deriving constraints assuming one analytic and four different simulated temperature and density profiles of the proto-neutron star. Finally, we estimate also the systematic uncertainty on the bound by varying the distance across which dark photons must propagate from their point of production to be able to affect the star. This work clarifies the bounds from SN1987A on the dark-photon parameter space.


Introduction
Dark sectors -containing particles that are not charged directly under the Standard Model forces -could provide an explanation for various shortcomings of the Standard Model, including dark matter candidates, neutrino masses, and the baryon asymmetry. Three "portals" allow for renormalizable interactions between the Standard Model and a dark sector: the neutrino portal, the Higgs portal, and the vector portal. In this paper we focus on the vector portal, in which (at low energies) the Standard Model photon interacts with a new gauge boson (the "dark photon") through kinetic mixing. Dark photons can easily remain hidden from conventional searches despite having masses well below the Weak scale. Over the last few years, much progress has been made to map out the viable parameter space and identify promising new search directions [1][2][3][4].
In this paper, we will focus on constraints on dark photons from the observation of a corecollapse supernova in 1987 in the Large Magellanic Cloud, called SN1987A. In the standard picture, the vast majority of energy liberated from the collapsing star, about 99% of the difference in the gravitational binding energy of the progenitor and the remnant, leaves the supernova in the form of neutrinos [5]. However, if new particles are (a) produced in large numbers and (b) able to travel macroscopic distances without transferring energy back to the stellar material, they provide a cooling mechanism that competes with the neutrinos. Since the neutrino burst observed in association with the explosion agrees qualitatively with predictions based on Standard Model-only simulations [6,7], such an anomalous energy loss is severely constrained. This means that dark-sector particles with masses below the characteristic supernova temperature of ∼ O(10's of MeV) must be either (a) so weakly coupled they are produced in negligible quantities or (b) coupled strongly enough that their mass and energy gets reprocessed by Standard Model particles. The purpose of this paper is to examine this scenario in detail, focusing on a realistic treatment of finite density effects on the Standard Model photon.
We will be interested in a minimal dark sector whose only light particle is a U (1) gauge boson. This gauge boson can mix with the Standard-Model photon, with field strength F µν , in the usual way, where F µν = ∂ µ A ν − ∂ ν A µ is the field strength for the dark photon field A and is the kinetic mixing parameter that allows the dark photon to couple to ordinary matter. The dark photon has a mass m . (We will assign a prime to any quantity associated with the dark photon.) We assume that the dark photon is the only new particle that could affect the supernova evolution. We leave to future work the possibility that there are additional stable, low-mass dark-sector scalars or fermions with U (1) charge [8] (see also [9][10][11]). Such particles, which could constitute the dark matter, could in principle provide another decay channel for dark photons and alter the supernova's evolution. In addition, we also imagine that a dark Higgs boson, which could be responsible for spontaneously breaking the U (1) gauge symmetry, is heavy enough not to effect the supernova evolution. We will eventually set bounds in the − m parameter space. Early studies indicate ∼ O(10 −10 ) is the order of magnitude of mixing angle that affects the supernova evolution [12][13][14][15].
Due to the mixing term F µν F µν /2 in Eq. (1.1), dark photons can be produced in any process that results in Standard Model photons. The primary processes of interest are depicted as the left two diagrams in Fig. 1. However, the kinetic mixing angle that dictates the production rate of dark photons can be substantially enhanced compared to its value in vacuum due to finite temperature and density effects. The Lagrangian parameter and the effective mixing angle are inequivalent when the Standard Model photon receives a nonzero "plasma mass." This is encoded by the polarization tensor Re Π ∼ ω 2 p = 4πα EM n i /E F,i [16], where the sum goes over charge carriers, n i is their number density, and E 2 F,i = m 2 i +(3π 2 n i ) 2/3 is their Fermi energy. At high density n i m 3 i , as encountered in the proto-neutron star, contributions to the plasma mass are dominated by the lightest charged particle present (i.e. electrons). The polarization tensor enters the gauge boson mass matrix as an additional mass term, which (dynamically) changes the diagonalization of the kinetic mixing parameter . The end result  Figure 1. Scattering and decay processes involving A particles: (left) bremsstrahlung in neutron-proton (n − p) scattering, (center) "semi-Compton" scattering, and (right) decay. The bremsstrahlung diagram is accompanied by two more diagrams at the same order (the A attached to the incoming p and a charge-exchange interaction). The grey blob in the bremsstrahlung diagram represents non-perturbative n − p scattering. The A can be produced inside the star by bremsstrahlung and semi-Compton scattering, but inverse versions of bremsstrahlung and semi-Compton contribute to absorption of A particles on their way out of the star. Decay to e + e − , possible only for m > 2m e , where m e is the electron mass, is not in equilibrium due to the low density of positrons in the protoneutron star.
is that the mixing angle in medium, m , depends on m and Π. On the mass shell the effective mixing angle is [17] where we reserve the symbol to always refer to the mixing angle in vacuum, which is the parameter that enters the Lagrangian. The transverse and longitudinal components need to be treated separately: each has its own polarization tensor and hence its own dispersion relation. We discuss how to calculate Π L,T in Sec. 3. The in-medium mixing angle from Eq. (1.2) controls both the production and absorption rates of dark photons, so that the amount of energy that flows out of the proto-neutron star is (schematically) L ∼ ωΓe −τ for a dark photon energy ω, production rate Γ, and "optical" depth τ . As we show below, τ has a strong dependence on ω, so there are a number of interesting limits for the production rate of the dark photon. We anticipate the results below by noting that at low masses and small mixing the bounds scale such that the product m is constant, m / MeV 2 × 10 −9 , independent of production rate and in contrast to prior results which indicated a flat asymptote of O(10 −10 ) at low mass. For large-enough mixing angles the bounds disappear: with m < 2m e and m / MeV 5 × 10 −5 , the absorption of dark photons through the inverse semi-Compton process is large enough to prevent the anomalous cooling of the supernova, while for m > 2m e and m / MeV 1 × 10 −6 , the dark photons rapidly decay to electron-positron pairs and are unable to carry away enough energy to affect the supernova evolution. We note that the high-limit we obtain in the range m < 2m e is up to two orders of magnitude stronger than it would be if we had taken a blackbody spectrum for escaping particles. Instead, we find that high-energy particles with ω few × T are capable of escaping the supernova interior without scattering, leading to relatively higher luminosity compared to that expected from a blackbody spectrum.
To establish these results, we will need to revisit several different aspects of supernova physics and finite-temperature field theory. In Sec. 2, we will discuss what is known and what remains uncertain about the explosion of SN1987A. We will focus on the interpretation of the luminosity limits and describe the requirements for transporting energy away from the zone of neutrino production. In Sec. 3, we will discuss how the A is produced and propagates in the supernova environment. As is evident from Eq. (1.2), a resonance is possible in the mixing angle; we will discuss the circumstances under which this resonance is attained, and what happens when it is suppressed. In Sec. 4, we will apply these lessons to draw excluded regions in the m − parameter space, and then we conclude. Several appendices provide technical details of our calculation and a comparison with previous work.

Supernova Constraints on New Particles
Despite great challenges in modeling the intermediate stages of core-collapse supernovae, the late stages of their evolution are adequately understood at a qualitative level. In particular, observations of the neutrino fluence and energy spectra at Kamioka [18], the IMB experiment [19], and Baksan [20] are in broad agreement with predictions based on simulations of core-collapse supernovae conducted by Lattimer and Burrows [6,7]. Despite O(1) uncertainties on the relevant astrophysics underlying the explosion, discussed at more length below, the agreement between theory and observation illustrates that there is only limited room for the possibility of a competing energy sink during the time of the core collapse.
This agreement has allowed constraints on a wide variety of models, including axions [21,22], sterile neutrinos [23,24], and compact extra dimensions [25,26]. Preliminary work on dark photons first appeared in [12], followed by calculations of the luminosity in the Born approximation [13,14] and then including nonperturbative effects from nucleon-nucleon scattering [15]. In this work, we revisit these constraints, focusing on the effects of finite temperature and density on the limits at low masses. In this section, we focus on how to apply these constraints, and discuss subtleties and uncertainties that hamper a complete analytic understanding of the bounds.

Using Luminosity to Obtain Constraints
Dark photons may be emitted as Standard-Model particles collide in the proto-neutron star. If dark photons are readily produced and are able to free-stream out, the total luminosity in dark photons may be so large that it drains the energy that powers electroweak interactions. This could potentially alter the luminosity of neutrinos that provide the terrestrial signal of the supernova explosion. Detailed simulations, including the backreaction of the new heat sink on the Standard-Model bath, are necessary to reliably assess the viability of any given point in the model parameter space. However, a convenient, approximate criterion was established in [5] Luminosity as a function of mixing angle for a fixed m . Above pr (m ), enough dark photons are produced that their luminosity may exceed L ν . For > tr (m ), dark photons produced far in the interior of the proto-neutron star are trapped. This means that enough of their energy will be reprocessed into neutrino energy, preventing serious energy depletion. The range pr (m ) ≤ ≤ tr (m ) is ruled out. If L(m , ) exceeds L ν for no value of at a given mass m , there is no constraint. If we neglected the optical depth entirely, we would find that the luminosity scales like ∼ 2 , shown as a dashed line labelled "No Trapping". Similarly, if we took the standard assumption that the emission is a blackbody spectrum from a surface at which the average dark photon has unit optical depth, we would find that the luminosity scales like ∼ e − 2 . This underestimates the contribution of very highenergy dark photons and leads to the bound shown with the dotted line labelled "Thermal Emission". The thermal emission bound is weaker than the real bound.
when the core reaches peak density ρ c ∼ 3 × 10 14 g / cm 3 and temperature T c ∼ 30 MeV, the duration of the neutrino burst from SN1987A is reduced by half, and the energy spectrum is inconsistent with observations. We will adopt this "Raffelt criterion" to find constraints. Stronger limits are potentially obtained by bounding the efficiency of the energy transport inside the core or the energy released over the entire neutrino burst, but these would also be subject to much larger systematic uncertainties from not knowing precisely the nature of the progenitor star. In Fig. 2 we plot a schematic representation of L for a fixed value of m . We see that at any given mass m there can be both upper and lower bounds on the mixing angle ; thus, there are two critical values of the mixing angle, which we call pr and tr . Below pr , dark photons are so weakly coupled that they cannot be produced readily enough to affect the evolution of the proto-neutron star. However, not all mixing angles > pr are ruled out: at sufficiently large mixing angles, above tr , dark photons are trapped before they exit the star, allowing energy to be efficiently reprocessed in the form of thermal neutrinos. Thus, only in the range pr ≤ ≤ tr can the dark photons cause a violation of the bound in Eq. (2.1).
The dark-photon luminosity is given by a volume integral of the differential power, dP/dωdV . Before focusing on the dark photon, we describe the luminosity for an arbitrary differential power. Of relevance is the luminosity in new particles that cannot be reprocessed efficiently as neutrino energy. We are thus interested in the luminosity from the zone of neutrino diffusion to a zone in which neutrinos are not produced efficiently. In other words, we want to find the energy that escapes from behind the "neutrinosphere" (the isotherm where Standard-Model neutrinos approximately transition from diffusion to free streaming) defined by some radius R ν . In general, this may be written where τ is the optical depth of a dark photon produced at radius r and the energy integral starts at 1 m . For a dark photon traveling radially outward, we have for an absorptive width Γ abs ; see App. B for the technical details of this calculation. As discussed below, dP/dV dω is proportional to the absorptive width, Γ abs using the principle of detailed balance, so we can simply calculate Γ abs for any given interaction. Two radii appear in Eq. (2.2). The radius R ν ∼ O(40 km) of the neutrinosphere is the radius outside of which most neutrinos free stream until arriving at Earth. For definiteness, we define R ν to be the radius where the temperature of the star has fallen to 3 MeV. This is roughly consistent with the condition for neutrino free streaming. The radius R far ∼ O(100 km −1000 km) is some "far radius" outside of which neutrinos no longer are produced efficiently. The reason these radii enter our calculations in such a fundamental way is that they provide a measure of the amount of energy that is diverted away from neutrinos in a manner that can alter the observed neutrino signal. If we were interested in the ability of terrestrial experiments to detect the new physics particles, R far would be the distance to the Earth, which takes into account interactions with the progenitor star material outside of R ν , interactions in the circumstellar medium, and the contribution from the column density of the interstellar medium across the remaining distance. However, the dark photons need not travel very far to have an effect on the evolution of the neutrino flux from the supernova explosion. Instead, the "Raffelt criterion" in Eq. (2.1) says that if energy is taken away from the core and deposited in a different region of the star, the energy is effectively lost because it becomes unavailable to neutrinos. This diversion depletes the fuel of the nuclear "engine" that allows the cooling timescale to be t cool ∼ O(10 sec).
Possible values for R far distance physical justification R gain 100 km for r > R gain , ν capture exceeds ν production R shock 1000 km for r > R shock , material is not yet shock heated Table 1. We consider two choices for the distant radius R far beyond which A particles must transport energy to affect the neutrino cooling phase.
There are a number of reasonable choices for R far ; the only strict requirement is R far > R ν . If R far is too close to R ν , we would erroneously conclude that arbitrarily well-mixed dark photons would carry away too much energy, because for dark photons produced in the final shell the integral in Eq. (2.3) would go over a very small range and τ → 0. However, this is unphysical: if the energy from the dark photons can be reprocessed by Standard-Model particles into neutrino energy, the dark photons do not provide an important energy sink. For this reason, we suggest that the lower bound on R far is the neutrino gain radius R gain ∼ O(100 km), outside of which neutrino production has a lower rate than neutrino absorption [28,29]. A reasonable upper limit on R far is the shock radius R s ∼ O(1000 km), outside of which matter is as yet uncompressed [5]. Here, we will parameterize our uncertainty on R far by using the gain radius and the shock radius, i.e. R far = 100 km and 1000 km, as representative values. We list these in Tab. 1 for reference.
We point out here that our method of calculating L differs from prior work: we do not split the calculation into free-streaming and trapped regimes. Instead, we allow τ to "speak for itself" and suppress the integration in parts of parameter space where the optical depth is large and the particles are mostly trapped. As we discuss in greater depth in Sec. 3.2, this has very important consequences at large mixing angles. At large mixing angles, we find that the energy spectrum is not thermal. For a given value of , assuming a thermal spectrum underestimates the luminosity and thus leads to weaker limits.

Uncertainties Regarding the Explosion of SN1987A
SN1987A is a promising environment for examining new physics because of the combination of the unique physical conditions attained in the star and the proximity of the explosion. However, constraints on new physics from the observation of SN1987A are inherently limited by difficulties in understanding the detailed process of the supernova even in the minimal case with no new physics. Many aspects of SN1987A remain poorly understood, from the nature of the progenitor to the primary driver of the "shock revival" required to sustain the supernova explosion. The mass of the progenitor star is only bracketed within a factor of two, and consequently the temperature and density profiles have large, qualitative uncertainties.
Given the uncertainties in modeling SN1987A, it is sufficient to use the conservative limit on the luminosity in Eq. (2.1) to derive bounds. To aid our analytic understanding, we will refer throughout the text to a "fiducial" model, a simple analytic supernova profile as  Figure 3. Profiles of temperature and density for the different profiles described in Sec. 2.2: the "fiducial model" described in Eq. (2.4), two obtained from the reference runs of [31], and one profile from [32,33]. We cut off each plot at R ν , which we define as the radius at which the temperature is T = 3 MeV for the profiles of Eq. (2.4) and [31] and the radius which equals the thermal neutrino mean free path for the profile of [33]. Note that the value of R ν differs across the simulations, and in the figures each profile is scaled by its corresponding value for R ν . advocated for in [5] This model is described by a core radius, R c , outside of which the temperature and density fall like power laws as functions of r. As fiducial parameters we will take k ρ = 0.2, k T = −0.5, ν = 5, and R c = 10 km with core temperature and density of T c = 30 MeV and ρ c = 3 × 10 14 g / cm 3 and a uniform proton fraction of Y = 0.3. Given these quantities, we also consistently solve for the electron chemical potential, as outlined in [30], which will be important for scattering and decay processes involving electrons. Our fiducial model, while convenient for analytic purposes, will inevitably mischaracterize the core collapse to some degree. To parameterize the uncertainties arising from the exact nature of the supernova explosion, we will also show bounds using numerical profiles from a variety of different sources [31][32][33]. The profile from [32,33] is the 13 solar mass progenitor at t = 1 sec post bounce and with a 100 msec shock revival time inserted by hand. This work includes a solution to a neutrino radiative hydrodynamical code before shock revival and a solution to the flux-limited diffusion equation after cooling has commenced [33]. We refer to this profile as "Nakazato 13M " in what follows. The profiles presented in [31] are of 11.8 or 18 solar-mass progenitors, taken at t = 1 sec after collapse. These solve the Boltzmann Table 2. Radii of interest for different profiles. We define the core radius R c as the radius of peak temperature. The maximum constrained mass m max in our final results in Sec. 4 is essentially determined by the Boltzmann tail of the core temperature m max ∼ 3T c . equation for neutrino transport with the AGILE-BOLTZTRAN code and an equation of state based on known nuclear isotopes and relativistic mean field models [34][35][36][37][38]. We refer to these as "Fischer 11M " and "Fischer 18M ", respectively. The fiducial model matches these profiles quite well in the interior of the star. We will emphasize the fiducial model in what follows because of its analytic nature and for ease of comparison with past work. However, this emphasis should not be interpreted as an endorsement of the fiducial profile over the numerical profiles we consider. We plot the temperature and density as functions of radius in Fig. 3 for the fiducial model and for the four numerical models, and we list some of their benchmark data in Tab. 2. As is apparent, these differ by as much as an order of magnitude in certain parts of the proto-neutron star. Such uncertainties imply that any given profile does not capture all of the subtleties of the behavior of the nuclear matter in this violent environment. Indeed, constraints derived with any given profile are necessarily approximate 2 : given the limited understanding we currently have of the core collapse mechanism and the nature of SN1987A, more observations of core collapse supernovae are necessary to determine whether any given point, especially one near the boundary of our excluded region, is in fact incompatible with observations.
Despite these caveats, we will ultimately see that these "boundary uncertainties" affect a subdominant part of the parameter space. Moreover, the bounds derived from these various profiles provide an estimate of the systematic uncertainties. We will refer to the region excluded for all profiles as the "robustly excluded" region, which we believe takes into account the dominant systematic uncertainties of supernova astrophysics.

Dark Photons in SN1987A
In this section, we discuss how to derive constraints from SN1987A on dark photons in the m − plane. Two new elements in our calculation, compared to previous work, are that we include finite-temperature and finite-density effects [17] and we do not require that the emitted dark photons follow a Planck distribution, which is relevant for m < 2m e and at large . We also include semi-Compton scattering, which is relatively important compared to inverse bremsstrahlung at high dark photon energies and in regions of low nuclear density. These changes have important effects on the resulting constraints.
Limits that neglect finite-temperature and finite-density effects are accurate only for dark photon masses m ω p , where ω p is the plasma mass of the Standard Model photon. For lower m , the effects of finite temperature and density on the gauge-boson mixing angle cause the constraints to be systematically weakened compared to the results in e.g. [15], which assumed a constant (vacuum) value of the mixing angle at all values of m and all radii. Instead, the mixing angle in a medium, which we denote m , is given by Eq. (1.2) for each polarization separately. The fact that m differs for the longitudinal and transverse polarizations of the dark-photon field has important phenomenological consequences, as emphasized in [17].
The main contribution to the polarization tensor of the photon, which determines the photon's plasma mass, comes from electrons, which are relativistic and degenerate for all radii within R ν , i.e. m e T µ e , where µ e is the chemical potential. In this region we can use Eq. (77) of [30] with the conventions of [17] to get the Standard-Model photon polarization tensor. As mentioned in Sec. 1, the real part of this tensor is proportional to the plasma mass of the Standard Model photon, In Eq. (3.1) we define a photon velocity, v, via | k| 2 = ω 2 v 2 , which can be small in medium. The imaginary part of the Standard-Model photon polarization function is proportional to the total width, Im Π L,T = −ω Γ abs|L,T − Γ prod|L,T , where the absorptive width Γ abs is the total rate at which Standard-Model photons of energy ω are absorbed, and Γ prod is the rate at which they are produced [17,44]. In local thermal equilibrium, these rates are related by Γ prod = e −ω/T Γ abs . Since Standard Model photons are in thermal equilibrium, the imaginary part of the photon polarization tensor is Together, Eqs. (3.1) and (3.2) allow one to calculate the mixing angle under all conditions. The differential power in photons is related to the absorptive width by [17] dP (r) where S L (S T ) = 1(2) is the spin degeneracy of a given polarization state of the Standard-Model photon. We are, of course, interested in the power carried by dark photons, which are exclusively produced through kinetic mixing with Standard-Model photons. The power dP of interest is obtained by making the substitution Given the absorptive width for the Standard-Model photon, we may then use Eq. (3.3) to calculate the differential power for dark photons of any mass and mixing angle. Since the optical depth τ is the integral of the absorptive width over all radii larger than r, as given in Eq. (2.3), it follows that we need only to calculate Γ abs to be able to find the luminosity in dark photons. We consider A production and absorption through bremsstrahlung in neutron-proton collisions, Compton-like scattering involving an e − , and decay to e + e − pairs (see Fig. 1). The absorptive width will generically be different for the longitudinal and transverse polarizations because the different polarizations have different dispersion relations in medium. We describe the calculations for these widths in detail in App. A; here we summarize the relevant results.
The dominant process in the core of the proto-neutron star is dark-photon bremsstrahlung, shown in the left panel of Fig. 1. The width for absorption through inverse bremsstrahlung is is the appropriately averaged neutron-proton dipole scattering cross section [15]. The subscript on the final fraction of Eq. (3.5) indicates that it should be included for Γ br|L only. As described in App. A, to the order of the soft radiation approximation that we can calculate [15,45], the rate for bremsstrahlung and inverse bremsstrahlung are equivalent. The reason is that this approximation is an expansion in ω/T CM , where T CM is the center-of-mass kinetic energy, so the dark photon energy should obey ω T , which formally corresponds to e −ω/T 1. However, if we do not enforce ω T , we discover large contributions from regions where the incoming nucleon energies are actually strongly Boltzmann suppressed. Thus, to consistently account for larger values of ω, we include a physical cutoff ω ≤ T CM , which allows a consistent reintroduction of a Boltzmann-like factor. With this substitution, the rate for bremsstrahlung production is no longer the same as the rate of inverse bremsstrahlung: we substitute σ We refer to emission of a dark photon in photon-electron scattering as "semi-Compton" scattering; this is shown in the center panel of Fig. 1. Semi-Compton scattering is important for m < 2m e despite the fact that it is subdominant to bremsstrahlung in the core. It has typically been neglected in these models, which is a reasonable approximation assuming a thermal spectrum. However, it has a different scaling with energy and density than bremsstrahlung, and we find that semi-Compton scattering can dominate absorption far out in the star for m < 2m e . Exact formulae for the different polarizations are cumbersome to write down due to the many mass and energy scales in the problem and the fact that the separate matrix elements are not Lorentz invariant. Nevertheless, we find an excellent fit to a complete calculation in the region ω 200 MeV and r ≤ R ν using the simple formula (in fact, we find that this approximation holds beyond R ν , but we only need it up to R ν ). This may underestimate Pauli blocking for low darkphoton energies, but qualitatively reproduces the energy dependence from a more complete calculation.
A final contribution to the dark photon optical depth comes from the decay to e + e − , which is allowed for the dark photon but not the Standard Model photon [30], and dominates the optical depth for m > 2m e . This has width where the limits of integration are , and the dimensionless matrix element squared m 2 L,T is given in Eq. (A.21). We emphasize that bremsstrahlung and inverse bremsstrahlung of photons are in equilibrium, np ↔ npγ, so that this process contributes both to dark-photon production and absorption. In contrast, positrons are very rare, so electron-positron coalescence is negligible, e + e −← →A , and this process only contributes to the absorptive width.

Luminosity in the Low-Mixing Regime
In this subsection, we consider small values of , for which the optical depth is small. In this region of parameter space, there is a "resonance" in m . From Eq. (1.2), we see that a necessary condition for such a resonance is Im Π Re Π. Using Eqs. (3.1) through (3.7), we verify that this condition is satisfied throughout the proto-neutron star. In the next section, we consider larger values of where there is no resonance.
On resonance we have Re Π = m 2 , and the differential power is As expected, the details 3 of the production process have fallen out. This follows from the scaling ∆ω * ∼ Γ prod ∼ ImΠ, where ∆ω * is the width of the Breit-Wigner function. The 3 We note that when the charge carriers are relativistic, like the electrons in the proto-neutron star, the resonant power has the same scaling with m for the different polarizations, unlike in the nonrelativistic limit [17]. Nonetheless, just as in the nonrelativistic case, most of the power for dark photons of low mass comes from the longitudinal mode: L modes are resonantly produced throughout the star where m < ωp, while T modes are resonantly produced only if radii exist at which ωp < m < We show ω * as a function of ω p /m in Fig. 4: for the longitudinal (transverse) mode, resonance occurs for ω p > m ( 2/3 m < ω p < m ). In the limit ω p → m , the resonance energy for each polarization approaches the mass, ω * → m , and the velocity becomes small.
Since the density decreases monotonically, we find that dark photons can be efficiently absorbed near the site of their production. Ther integral for the optical depth is resolved by a delta function just as the ω integral is for the power, giving Since e + e − coalescence is not in equilibrium we have Im Π ∝ Γ abs , and a ratio of widths appears in the numerator of Eq. (3.9). If the A has additional absorptive or decay channels, their widths will also appear in the numerator of this expression.
There can only be a resonance for m ≤ 3 2 ω p (r = 0) ∼ 15 MeV (for larger m , one can still produce enough dark photons off-resonance to violate the Raffelt criterion, see Sec. 3.2). We find that for ∼ 3 × 10 −9 /(m / MeV), the Raffelt criterion is violated, corresponding to pr (m ) in Fig. 2. At larger mixing angles, the resonant production is suppressed by resonant absorption. Using only the delta-function approximation for τ and P , we find that the bounds are lifted for mixing angles in excess of ∼ 10 −7 /(m / MeV), corresponding to tr (m ) in Fig. 2. With this simple approximation, we are able to rule out the region δ-function only: 3 × 10 −9 m / MeV 10 −7 , m 15 MeV . (3.10) This partial result already differs from earlier work using related techniques [12][13][14][15]. The difference in the magnitude of the lower bound above a few MeV compared to [14,15] is due to O(10) differences in our calculation of the differential power. The more consequential difference in the shape of our bound, which has a linear dependence on mass, is a consequence of including finite-mass effects on the mixing angle. The delta-function approximation, which gives our result in Eq. (3.10), fails in a dramatic but predictable way for large , which we discuss next.

Luminosity in the High-Mixing Regime
The delta-function approximation is not always suitable for obtaining the luminosity emitted in dark photons. The approximation fails because two delta functions appear in the integrals of Eq. (2.2): one arises in the optical depth τ , and another in the power dP/dV dω. Critically, the exponential of −τ * means that the delta function in the opacity competes against the delta function contained in dP/dV dω. For τ * 1, the differential luminosity on resonance may be suppressed compared to the differential luminosity at higher energies, even if the standard condition for a resonance, Im Π Re Π, holds. In other words, exponentiating τ at large values of may "depress" the peak of the would-be delta function in dP/dV dω while leaving contributions at other energies untouched. The end result is that we can put bounds above tr (m ) from Eq. (3.10). With a proper treatment of the absorption mechanism, we can even place bounds stronger than those from assuming that particles in the high-mixing regime are thermalized at some radius. We now discuss these effects in more depth.
At sufficiently small mixing angles where τ is small, dL dP is a good approximation, and the area under the peak of the resonance well approximates the full integral: even though this peak is very narrow, the height of the peak at resonance exactly compensates for the small width, as discussed around Eq. (3.8). This fails when the optical depth at the delta function peak is non-negligible. When is large, the dark photons from the delta function peak are trapped due to their large mixing angle, and the differential power at ω * is no longer large enough to overcome the small width of the resonance.
We demonstrate this behavior in Fig. 5, showing the case of dark photons with mass m = 1 and 1.1 MeV. These masses are chosen to be below and above the electron-positron pair production threshold, where the absorptive width scales differently. In this figure, we clearly see the breakdown of the delta function approximation in the large-mixing limit. At large mixing angles, the peak in dP/dV dω is entirely depressed and contributes negligibly to the luminosity. The energy spectrum is instead peaked somewhere above ω * . This location and its height relative to the delta-function peak depends on the absorptive modes available to the dark photon. The energy that maximizes the differential luminosity depends on both m and . At large and m < 2m e , the energy spectrum of dark photons can have a peak of We also compare to a Planck spectrum in the dot-dashed red line, shown at the radius R s for which the thermal luminosity L th (m ) = π 3 g * (m )R 2 s T (R s ) 4 /30 = L ν , following [15]. At large mixing and low mass, the energy spectrum is neither resonant nor thermal: the peak of luminosity is above the Planck spectrum and above the resonance peak for large when e + e − is not on shell. This implies stronger limits than when simply assuming a thermal spectrumsee Sec. 3.2 for details. higher magnitude and larger central value than a Planck or delta-function spectrum. Thus, assuming a resonant or thermal energy spectrum at large mixing angle will underestimate the luminosity emitted in dark photons.
This conclusion regarding the energy spectrum of dark photon emission has important differences from standard methods in the "trapping limit" familiar from supernova bounds on sterile-neutrino production [23,24]. Because sterile neutrinos scatter through weak processes with a rate that scales like Γ ν ∼ nsG 2 F and a mean optical depth that scales like τ ν ∼ ω 2 . At sufficiently large mixing angles the high-energy sterile neutrinos scatter many times to lower energies, and the ones that escape have approximately thermalized with the ambient matter at the radius at which the density and temperature are low enough for them to free stream. For this reason, it is accurate to say that the sterile neutrinos are emitted like a blackbody from some surface at large radius. However, if the dark photon scatters primarily through inverse bremsstrahlung, the optical depth goes like τ ∼ ω −5 for small masses, and the optical depth decreases with ω. Thus, high-energy particles are able to escape without scattering, providing a conduit for energy loss that is decoupled from local temperatures. The main consequence of this scaling is that low mass dark photons are emitted from within a sizable part of the star's volume, unlike the "surface emission" of sterile neutrinos.
For this reason, when the delta function approximation is inaccurate we need to integrate over the full range of dark photon energies to assess the luminosity in dark photons. This can be computationally expensive, especially if the resonance is subdominant but not negligible. However, because the optical depth decreases with energy, we find that a suitable approximation is to simply drop the low-energy tail: where ∆ω * = Im Π m 2 is the width of the delta function peak, and the factor of 4 gives > 99% containment of a Breit-Wigner peak of these approximate dimensions. The first term in Eq. (3.11) is a rescaled version of the delta function peak discussed in Sec. 3.1: the optical depth suppresses the height of the delta function, but otherwise the resonance operates as before. When is large τ * becomes considerable, and the second term will dominates the luminosity.
The consequence of this nonthermal energy spectrum at high is that we may constrain fairly large values of . Assuming that highly mixed A s are emitted with a thermal energy spectrum is bound to underestimate the luminosity in dark photons at masses below the e + e − threshold. This is because for a fixed 2 a fixed number of dark photons are produced. Assuming thermalization, we have L ∼ 2 T 4 R 2 s for some "A -sphere" radius R s , while with our approach we find that L ∼ 2 ω 3 s m 2 R 3 c , where ω s T is the maximum of the differential luminosity. We find that this energy ω 4 s ∼ n n n p σ (2) np may be quite high because it is not determined solely by the temperature, but rather by the energy scales of the absorptive processes. Since the proto-neutron-star core is close to nuclear saturation (so the interparticle spacing is small and n n and n p are large), the emission of high-energy A can overcome the Boltzmann suppression e −ωs/T and provide enough energy loss to violate the Raffelt criterion. While these arguments work well for m < 2m e , the decay A → e + e − invalidates them for m > 2m e and keeps L close to thermal. This can be seen in Fig. 5.

Results
By integrating the exact expressions for dL/dV dω, we can determine whether a point in the − m parameter space leads to a dark-photon luminosity that exceeds the luminosity in neutrinos. When the optical depth on resonance is τ * 1, this integral is well approximated by the resonant power. When the optical depth at the resonance energy is τ * O(1), we must integrate Eq. (3.11) over 0 ≤ r ≤ R ν to find the luminosity that is drained from inside the neutrinosphere.
Despite our quantitative understanding of this procedure as described in detail above, several major, qualitative sources of uncertainty linger in the calculation. We are unable to evade order-of-magnitude uncertainties regarding: the exact nature of the progenitor star, including its mass and peak temperature; the correct location for R far , beyond which energy is no longer available to be reprocessed for the neutrino cooling phase; and the form of the highenergy cutoff for Γ br , which, because we work in the soft-radiation approximation, is not given simply by a Boltzmann factor. Subdominant sources of uncertainty from the soft radiation approximation, such as Pauli blocking of initial state nucleons, relativistic corrections to the Boltzmann distribution, and low-energy effects like suppression of bremsstrahlung from the Landau-Pomeranchuk-Migdal mechanism [46] should also contribute corrections of order tens of percent to our calculation. More aggressive bounds, obtained in principle by limiting the amount of energy transport at the core radius, are subject to larger uncertainties due to the difficulty of modeling supernova explosions, so we do not pursue such bounds here. Irrespective of the uncertainties, we advocate for including finite-temperature effects with new particles in future simulations of core-collapse supernovae. Moreover, the assumption that the emitted particles (at large mixing) must follow a thermal distribution is not necessarily correct, as demonstrated by the dark-photon model, where the trapping probability decreases for larger energies. This would also be interesting to include in future studies. In Fig. 6, we estimate the size of the systematic uncertainties by plotting the excluded regions for all profiles and values of R far that were discussed in Sec. 2.2. We emphasize that the entire lower boundary is fairly robust against the uncertainties we consider. Uncertainties regarding the temperature profile of the progenitor contribute to significant uncertainty on the location of the maximum constrained mass, however, and also on the upper boundary where m < 2m e . The uncertainty in R far is the greatest contribution to uncertainties for the upper boundary where m > 2m e . This is because that boundary is mostly set by the decay length to e + e − , and the absorptive width scales linearly with R far . In contrast, for m < 2m e , the dark photon is stable against decays, 4 and the uncertainty is driven by the details of the stellar core rather than by the physics near the neutrinosphere. The combination of the uncertainties in profile and R far gives roughly up to an order of magnitude uncertainty at all points on the boundary. We can consider these uncertainties to be systematic error bars about the fiducial result. In Fig. 7, we display the fiducial result as a red line, and the envelope of other bounds provide the extent of the green region. The blue region is constrained regardless of the perturbations we have described above, and we refer to it as "robustly excluded".
Compared to other results on the dark photon model [13][14][15], our constrained region is somewhat smaller (see App. C for a direct comparison). The biggest change in our methodology is the inclusion of finite-temperature effects on the photon self-energy, which changes the asymptotic scaling at low mass to be constant in the product 2 m 2 (rather than mindependent). This has been accounted for in limits on dark photons from other stars [17], but our results are the first to include plasma effects on the constraints from SN1987A. Another important change is that we discard the assumption that A particles thermalize at high mixing angle. We give a sense for the difference this makes by also showing the upper boundary of the constrained region using the (incorrect) assumption of thermalization. For masses asymptotically smaller than 2m e , the difference is of order a few. The discrepancy increases around a few hundred keV as the thermal production falls off while the resonant production remains essentially fixed. Above 2m e , where decays to electron-positron pairs dominate the absorptive width (and consequently cut off the high-energy tail of the differential luminosity), the assumption of thermalization leads to an O(50%) error in the bound. Finally, due to ∼ O(10) differences in our calculation of the bremsstrahlung width, our lower boundary moves up by a factor of ∼ O(3) compared to the results in [15]. This lowers the maximum constrained mass because the highest constrained mass is determined by where the upper and lower boundaries (which are relatively steep functions of m ) intersect.
In Fig. 7 we also compare to other bounds for m ∼ 1 keV −1 GeV. Dark photons below O(100 keV) can disrupt the energetics of the Sun and other main sequence stars [17] (labelled "stars"). Decays of A → γγγ, which are not observed from the time of big bang nucleosynthesis until the present day, provide the constraints in the region labelled "late decays": at the present cosmological epoch these decays would be observed as an isotropic diffuse X-ray background, and an injection of energy at the time of Big-Bang Nucleosynthesis is prohibited by the success of standard predictions. Using the width calculated in the Euler-Heisenberg effective theory [47] and the model-independent bounds of diffuse X-ray background in [49] for a decay to two photons, which we multiply by 3/2 to account for the increased multiplicity in this model and rescale by the abundance from resonant production in the early Universe [48], we obtain the lower bound of this region. Requiring τ ≤ 1 sec so that decays do not interrupt nucleosynthesis gives the upper line of this gray region [48]. We caution that threshold effects near m  . Systematic uncertainties (green region) encompassing our "robustly excluded zone" (blue). The true boundary likely lies somewhere in the green region; we show our fiducial profile as the dotted black line. The blue area is excluded regardless of the perturbations we make to the physical inputs. We compare to bounds from other stars [17], decays to three photons on cosmological timescales [48,49], and beam dumps, meson decays, and other terrestrial experiments [4]. The comparison of the electron anomalous magnetic moment in two different systems is shown in the hatched region, which has not previously been shown in this mass range.
to two photons, should in general increase this partial width beyond the rate calculated in the effective theory. This would open up the high-range of this region near m = 2m e and deserves a more thorough investigation.

Conclusions
We have examined the effects of dark photons on SN1987A, incorporating for the first time the effects of the finite temperature and density in this environment on the kinetic mixing angle. Resonant production is important at low mixing. These effects cause the constraints to lift at small m in such a way that the product m is constant. Our result for the lower edge of excluded region is that we roughly require 2 (m / MeV) 2 3 × 10 −9 for m 15 MeV, and 10 −9 for 15 MeV m 120 MeV, when the production begins to be Boltzmann suppressed.
We have also carefully considered the non-thermal part of the dark-photon energy spectrum, discovering that the limits at very large mixing angles are somewhat more constraining than previously understood. The reason this is possible is that the energy spectrum at large mixing angle is not close to thermal. In fact, below the e + e − threshold the peak of the differential luminosity is up to an order of magnitude higher than the temperature, and our bounds are a factor of a few stronger than they would be if we assumed that the spectrum were thermal. This may have important implications for constraints on other dark-sector particles. We leave an investigation of this and variations of the simple dark-photon model to an upcoming publication [8].
where K µ = (ω, k) is the photon four-momentum, Ω is the photon scattering angle, is the collision operator for the dark photon in the proto-neutron star, P µ i = (E i , p i ) are the nucleon four-momenta, f n,p are the nucleon distribution functions, discussed in detail below, and |M| 2 γnp is the squared matrix element summed, not averaged, over all spins. The absorptive width defined in Eq. (A.1) enters the polarization tensor and the renormalized mixing parameter as discussed in Sec. 3. We now calculate this width for the processes shown in Fig. 1.

A.1 (Inverse) Bremsstrahlung
In this subsection, we detail the calculation of the rates for bremsstrahlung production and absorption of the Standard Model photon in the environment of the proto-neutron star. The extension to dark photons is addressed at the end.
We will use the soft radiation approximation [15,45] to calculate the rate for inverse bremsstrahlung. The soft radiation approximation is strictly appropriate only when photon energies are much below the other energy scales in the problem, ω T, T CM . However, this approximation is observed to reproduce data for energies in mild violation of the bound, ω ∼ T × few [15], so that we can use it for all of the parameter space.
With this approximation, the same calculation is formally appropriate for the production and absorption rates of dark photons: although production and absorption are in equilibrium as discussed in Sec. 3, and thus are exactly related by the Boltzmann factor e −ω/T , in the softradiation approximation we keep only the leading-order contribution, and this factor becomes unity. In contrast, errors incurred in constructing the matrix element are ∼ O(| p|/m N ), which scale like | p| ∼ ω/m N ω/m N . Therefore, corrections to the matrix element itself are of higher order than the Boltzmann factor. Nevertheless, we will find that a cutoff that parametrically scales like the Boltzmann factor necessarily arises with the requirement that ω not exceed the energy of the collision: without specifying a shape or a new energy scale, ensuring that collisions not violate energy conservation reintroduces a high energy cutoff.
In the soft-radiation approximation, using Lorentz invariance, gauge invariance, and the fact that the photon couples to the dipole current in n − p scattering, we find that the matrix element of interest 5 is [45] where 4πα EM = e 2 is the U (1) EM coupling; ε µ is the photon polarization; the dipole current is L µ = P µ 1 /P 1 · K − P µ 3 /P 3 · K; S L,T is the spin degeneracy for the polarizations of the Standard Model photon, with S L = 1 and S T = 2; and |M| np is the spin-averaged matrix element for neutron-proton scattering, discussed in [15]. We point out that the longitudinal and transverse modes of the dark photon have different matrix elements, since they have different dispersion relations in medium. The mean free path for bremsstrahlung using Eq. (A.1) is now written Because we work in the soft-radiation approximation, we drop K µ in the delta function in Eq. (A.4), which formally eliminates the distinction between bremsstrahlung and inverse bremsstrahlung. We make some simplifications to Eq. (A.4). First, we will assume that the nucleons are neither degenerate nor relativistic. This means that Pauli blocking is negligible: With this assumption, the second line of Eq. (A.4) is Lorentz invariant, so we may evaluate the second line in any frame that we choose. We choose the center-of-mass frame, E 3 = E 4 = m N + T CM /2 = E CM /2 (where m N is the nucleon mass), and we take this opportunity to define the center-of-mass kinetic energy T CM = ( p 1 − p 2 ) 2 /4m N = E CM − 2m N for future convenience. We also recall that for elastic scattering where dσ np /dΩ CM is measured in the center-of-mass frame. (We point out here that we differ from [15] in a factor of v rel .) In this frame the second line of Eq. (A.4) is now where the energy delta function is δ( We now use explicit forms of the distribution functions. For the non-degenerate nucleons, we take a rescaled Maxwell-Boltzmann distribution function, f n,p (E i ) = ce −E i /T . Here, c is defined such that g n,p d 3 p i (2π) 3 f n,p (E i ) = n n,p , where n n,p is the observed (or simulated) number density of the nucleon. Because there are no anti-nucleons in the proto-neutron star, g n,p = 2 and explicitly f n,p (p i ) = n n,p 2 Plugging Eq. (A.8) into Eq. (A.7) and simplifying, we have where we have done the trivial integration over dφ CM . We continue by changing coordinates This is not trivial to evaluate in generality, because E 1 and E 2 are functions of all six components of the incoming particle momenta. However, in the nonrelativistic limit we may take E 1 E 2 m 2 N . With this approximation, the width is which is only a function of | p + | and | p − |. Consequently, (ε · L) 2 = 8π 3 at leading order in T CM [15]. However, at finite temperature the longitudinal and transverse polarizations of the Standard-Model photon are not the same, as they have different dispersion relations, and thus Re Π L = Re Π T . Rather than using the vacuum polarization sum, we evaluate the polarization sums in a relativistic medium. Here we have where the photon velocity v 2 need not be 1. To leading order, we find that the projected integrals of L 2 over the dark-photon scattering angle dΩ are (A.14) Now we do the integral over θ CM by defining σ (2) np ≡ d cos θ CM (1 − cos θ CM ) dσnp d cos θ CM , and the integral | p + | 2 e −| p + | 2 /4m N T d| p + | = 2 π(m N T ) 3 is straightforward. With S T = 2 we have Γ SRA,T = 4 √ πα EM n n n p 15) and Γ SRA,L = 1 − v 2 Γ SRA,T . Finally, we change to a dimensionless variable We define the appropriately energy-and angle-averaged np dipole scattering cross section to be σ np (x), following [15], to further streamline the notation. To make this useful for calculating the absorptive width for the dark photon we include the correct mixing angles and the dark photon on-shell condition v 2 = 1 − m 2 /ω 2 . Thus, our final result for the absorptive width of the dark photon via inverse bremsstrahlung is which again is the same for production or absorption (although when considering total power produced we multiply by S L,T to account for the spin degeneracy of the transverse mode).
Assuming that we can compare to [15] by Γ For the production rate, the exact same calculation goes through as in the absorption rate because of the soft-radiation approximation. However, we stipulate that the dark-photon energy must not exceed the energy of the collision. This gives Γ br = Γ SRA Θ(T CM − ω), where the step function must consistently be included before doing the integral over 6 T CM . This takes ∞ 0 → ∞ ω/T in the definition of σ (2) np , so we define σ np (xT ) and subsequently Γ br,T = 32 3π Empirically, we find that this correction leads to ∼ O(1) higher luminosities (and thus ∼ 50% stronger bounds) compared to taking a simple rescaling Γ br|L,T = e −ω/T Γ ibr|L,T , as would be expected in an exact calculation where we should recover detailed balance.

A.2 Decay
For decay to an electron positron pair, we start with the relation Γ e + e − ;L,T = 1 2ω Assuming that the positrons are not Pauli blocked and that the electron distribution function is Fermi-Dirac, this simplifies immediately to Γ e + e − ;L,T = 1 32π 2 ω (A.20) Using the projection operators from Eq. (A.13) and S L = 1, S T = 2, we have where E 2 = ω−E 1 is understood everywhere. We verify as expected that |M| 2 T +|M| 2 L ∝ m 2 + 2m 2 e (with no independent spin average for the different polarizations in this sum). For ease of notation, we also define m(E 1 ) 2 L,T ≡ |M| 2 L,T /16πα EM m 2 2 L,T . Defining the common factor z( which should be ∼ O(1) for most energies.
Respecting the cylindrical symmetry of the process, m 2 L,T have no dependence on the azimuthal angle of the outgoing particles. Doing the dφ integral, our expression becomes where the extra factor under the measure is from changing variables in the delta function, and p * is the momentum that satisfies the delta function. Relativistic kinematics gives wherev = 1 − 4m 2 e /m 2 and β = 1 − m 2 /ω 2 = | k|/ω. We also observe that where we have used p * = E 2 1 − m 2 e . We plug back into Eq. (A.22) and after some cancellations the result is where the limits of integration are Eq. (A. 25) is not singular at ω → m because in this limit the bounds of the integration approach one another to cut off the division by zero. We are not able to integrate Eq. (A. 25) in generality due to the complicated form of the matrix element and the Pauli factor.

B Details for Calculating the Optical Depth
To calculate the optical depth for an outgoing dark photon, we must integrate its absorptive width along its path out of the proto-neutron star all the way to R far . In this appendix we describe some subtleties in this calculation. In general, the integral for the optical depth requires knowledge of the angular distribution of stellar material. This integral will furthermore be sensitive not only to the radius of production but also to the angle each dark photon makes with respect to the radial direction as it leaves the star. For simplicity, we will assume spherical symmetry of stellar material, and we will correct for the distribution of possible path length as follows. If the path of a dark photon produced at radius r is directed radially outward, the optical depth is simply τ radial, out (r) = R far r Γ abs (r)dr. (B.1) If a dark photon produced at radius r is directed radially inward towards the core, it traverses the distance to the center of the star once going inward and once going out. In this case, the optical depth is τ radial, in (r) = τ radial out (r) + 2 Dark photons that point inwards but do not cross through the core will not experience as large of a suppression. Thus averaging over the angles with which a dark photon may be emitted with respect to the radial direction, we have heuristically τ (r) = τ radial out (r) + a(r) where the correction factor a(r) accounts for the densities and temperatures along the chord that is traversed and takes into account the angles which encounter the core. We expect that a depends on r because, e.g., the fraction of solid angles that cross the core region (where the density is relatively high) will vary depending on how close to the core the dark photon is produced: for r < R c all paths go through the core, near r R c about half of the paths return through the core, but well outside R c this fraction is ∼ R c / R 2 c + r 2 . Given these estimates, we expect a to rise monotonically from 0 at r = 0 to something close to R c /R ν at r = R ν . In our numerical work, we will use the approximation τ (r) = 1 − r(r − R c ) 2R 2 ν τ radial out (r).

(B.4)
A similar geometric factor multiplying τ radial out (0) is found to equal 3/π 2 for an object of uniform density in [71], and determined to be close to 0.3 in the numerical work of [72] in the leakage approximation. The prefactor in Eq. (B.4) does not deviate far from one in practice. We also comment on the upper bound of the integration for the optical depth. Because we are interested in the luminosity that reaches regions outside of the neutrinosphere (possibly as far as the shock radius), we must integrate the partial widths at radii beyond the profiles of the simulations to get the total optical depth. To account for this subtlety, we approximate the optical depth as τ (r) = R far r Γ tot (r)dr = (R far − R ν ) Γ e + e − (R ν ) + Γ sC (R ν ) + 1 10 Γ ibr (R ν ) + Rν r Γ tot (r). (B.5) The factor of 1/10 in front of Γ ibr is present to account for the sharp dependence on r carried by Γ ibr . If Γ ibr ∝ r −a for a 1, we would expect the contribution at R far to be negligible compared to the region near R ν . However, the rate is unlikely to fall in this simple way out to arbitrarily high radii, so we adopt this approximation (which has a soft but nonzero scaling with R far ) as a very rough estimate. Fortunately, Γ ibr does not dominate the optical depth anywhere near the upper bounds of the − m parameter space, so our results are insensitive to this prescription.

C Comparison with Prior Work
As mentioned above, our constrained region goes to somewhat smaller values of compared to other results on the dark-photon model [12][13][14][15]. This is due to a number of partially offsetting factors, which we discuss in detail here.
The biggest change in our methodology, which leads to a qualitative change in the asymptotic scaling at masses below ∼ 10 MeV, is the inclusion of finite-temperature effects on the Standard Model photon self-energy. Thus, our limits are constant in the product 2 m 2 rather than m -independent. This has been included in limits on dark photons from other stars [17], but our results are the first to include plasma effects on the constraints from SN1987A. Another important change is that we do not assume that A particles thermalize at high mixing angle. Finally, our maximum constrained mass is somewhat lower in the fiducial model (analogous to the conditions assumed by previous work) because this mass is determined by where the upper and lower boundaries (which are relatively steep functions of m ) intersect.  Figure 8. Comparison of the bounds derived in this work with those from [12][13][14][15]. Note that the blue dotted line from [14] is not based on anomalous cooling of the supernova (see text for details).
The first (preliminary) bounds on this model space were presented in [12], using a simple rescaling of calculations in the axion case. This implicitly relies on a picture where a single pion exchange is used to model the nucleon interactions. The bounds in [13,14] also assume single pion exchange. As emphasized in [15], this is problematic across the entire range of energies considered in the proto-neutron star, and substantially overestimates the scattering cross section at the energy transfers expected in this environment. This leads to bounds that are too strong by almost an order of magnitude.
We also show a bound from [14] for the range 10 −13 10 −10 , denoted by a (d) on our plot. This bound relies on scaling arguments for the dark-photon production and assumptions about the outlying stellar material beyond the neutrinosphere. This estimate neglects finite-temperature effects, which we have shown are important in the core of the supernova below 10 MeV. We have not checked whether Pauli-blocking factors are important to include. Moreover, systematic uncertainties are at least as large as for the cooling bound.