Constraining the SelfInteracting Neutrino Interpretation of the Hubble Tension
Abstract
Large, nonstandard neutrino selfinteractions have been shown to resolve the tension in Hubble constant measurements and a milder tension in the amplitude of matter fluctuations. We demonstrate that interactions of the necessary size imply the existence of a forcecarrier with a large neutrino coupling () and mass in the keV – 100 MeV range. This mediator is subject to stringent cosmological and laboratory bounds, and we find that nearly all realizations of such a particle are excluded by existing data unless it carries spin 0 and couples almost exclusively to flavored neutrinos. Furthermore, we find that the light neutrinos must be Majorana, and that a UVcomplete model requires a nonminimal mechanism to simultaneously generate neutrino masses and appreciable selfinteractions.
I Introduction
The discrepancy between lowredshift and Cosmic Microwave Background (CMB) determinations of the presentday Hubble parameter, , has grown in significance to over the last several years Riess et al. (2016); Shanks et al. (2019); Riess et al. (2018); Aghanim et al. (2018); Riess et al. (2019). The standard cosmological model, CDM, may need to be augmented if this “ tension” is not resolved by observational systematics. Intriguingly, this tension cannot be addressed by modifying CDM at low redshift Vonlanthen et al. (2010); Verde et al. (2017); Evslin et al. (2018); Aylor et al. (2019), but adding new physics before recombination seems more promising Lesgourgues et al. (2016); Di Valentino et al. (2018); Poulin et al. (2018a); D’Eramo et al. (2018); Poulin et al. (2018b); Pandey et al. (2019); Escudero et al. (2019). Furthermore, low redshift measurements of the matter density fluctuation amplitude on 8 Mpc scales, , also appear to be lower than predicted by CDM from the CMB. This milder tension is not ameliorated in the models of Di Valentino et al. (2018); Poulin et al. (2018a); D’Eramo et al. (2018); Poulin et al. (2018b); Pandey et al. (2019); Escudero et al. (2019).
A particularly interesting resolution to these issues is a nonstandard neutrino selfinteraction of the form CyrRacine and Sigurdson (2014); Lancaster et al. (2017); Oldengott et al. (2017); Kreisch et al. (2019)
(1) 
where is a dimensionful coupling and flavor indices have been suppressed. If is much larger than the Standard Model (SM) Fermi constant, , neutrinos remain tightly coupled to each other until relatively late times. This inhibits their freestreaming and results in enhanced power on small scales and a shift in the acoustic peaks of the CMB spectrum relative to CDM Bashinsky and Seljak (2004). The effect of neutrino selfinteractions is degenerate with other parameters in the CMB fit, such as the angular scale of the sound horizon, the spectral index and amplitude of primordial fluctuations, and extra freestreaming radiation (as parametrized by ).
These approximate degeneracies have been shown to prefer in cosmological data CyrRacine and Sigurdson (2014); Lancaster et al. (2017); Oldengott et al. (2017); Kreisch et al. (2019) and relax the tension Oldengott et al. (2017); Lancaster et al. (2017); Kreisch et al. (2019). The latest results of Ref. Kreisch et al. (2019) extended previous analyses by allowing for finite neutrino masses and extra freestreaming radiation at the time of the CMB. They found that in the “strongly interacting” (SI) or “moderately interacting” (MI) regimes
(2) 
could simultaneously reduce the and tensions. Interestingly, the SI cosmology prefers a value of compatible with local measurements at the 1 level, even before including local data in the fit.
The favored range of in Eq. (2) vastly exceeds the strength of weak interactions, whose corresponding coupling is . As we discuss below, such an interaction can only arise from the virtual exchange of a force carrier (“mediator”) with mass and appreciable couplings to neutrinos. This mass range indicates that this scenario faces a variety of stringent cosmological and laboratory constraints.
We find that if the interaction in Eq. (1) resolves the tension, then:

SI is excluded: The SIrange of Eq. (2) cannot be realized in any of the consistent, lowenergy, weakly coupled models that we consider.

Vector forces are excluded: Constraints from the epoch of Big Bang Nucleosynthesis (BBN) exclude all selfconsistent vector mediators.

Dirac neutrinos are excluded: If neutrinos are Dirac particles, mediatorneutrino interactions efficiently thermalize the righthanded component of neutrinos, thereby significantly increasing the number of neutrino species at BBN.

Minimal seesaw models are excluded: Achieving the necessary interaction strength in Eq. (2) from a gaugeinvariant, UVcomplete model is challenging with minimal field content in a TypeI or II seesaw model; a nonminimal mechanism is required to generate neutrino masses.

Interactions with are favored: Couplings to , with in the range of Eq. (2) are largely excluded by laboratory searches for rare decays and for neutrinoless doublebeta decay, except for a small island for the coupling. This means that the masseigenstate neutrinos only interact via their components.
This work is organized as follows: Sec. II demonstrates that a light new particle is required to generate the interaction in Eq. (1) with appropriate coupling strength; Sec. III presents the cosmological bounds on this scenario; Sec. IV discusses the corresponding laboratory constraints; Sec. V shows how Eq. (1) can arise in UV complete models; finally, Sec. VI offers some concluding remarks.
Ii The Necessity of a Light Mediator
The Boltzmann equations used in Refs. Oldengott et al. (2017); Lancaster et al. (2017); Kreisch et al. (2019); Barenboim et al. (2019) assume that lefthanded (LH), masseigenstate neutrinos participate in elastic scattering processes. They also assume that the interactions in Eq. (1) involve constant and flavoruniversal values of during all epochs relevant for the CMB. The largest CMB multipoles observed by Planck correspond to modes that entered the horizon when the universe had a temperature of .This temperature sets the characteristic energy scale of scattering reactions during this epoch, and it is important that the form of the Lagrangian shown in Eq. (1) is valid in this regime. At higher energies, however, this description can break down. In this section, we therefore emphasize the need to introduce new particle content to study laboratory and early universe processes that occur at energies .
As noted in Refs. Oldengott et al. (2017); CyrRacine and Sigurdson (2014); Lancaster et al. (2017); Kreisch et al. (2019), the operator in Eq. (1) is nonrenormalizable, and thus is necessarily replaced by a different interaction with new degree(s) of freedom at some energy scale higher than the energies probed by the CMB (see Ref. Cohen (2019) for a review). Since in Eq. (1) is momentumindependent, we will assume this interaction arises from “integrating out” a particle with mass and a perturbative coupling to neutrinos :
(3) 
where are twocomponent lefthanded neutrinos, and are flavor indices and the subscript “phen” indicates that we will use this Lagrangian for our phenomenological analysis. In Eq. (3) we have assumed that is a real scalar without loss of essential generality. In particular, our conclusions remain unaffected if is a CPodd or complex scalar. We focus on a scalar mediator here for clarity; introducing a new vector force instead, for example, follows the same reasoning but comes with additional, stronger constraints discussed in more detail below. We also explicitly allow for generic couplings between different neutrino species, and discuss implications of different choices of below.
Using Eq. (3), we see that the scattering amplitude is always proportional to two powers of multiplying the propagator, as shown schematically in Fig. 1. If the momentum transfer satisfies , we have
(4) 
where we have suppressed flavor indices and defined
(5) 
In the opposite limit, , then , leading to a qualitatively different energy and temperature dependence of neutrino selfinteractions; this regime was investigated in Refs. Forastieri et al. (2015, 2019), which found no improvement in the tension. Thus, for the remainder of this work we focus on models in which at energy scales relevant to the CMB. The intermediate regime, where the mediator mass is negligible for some CMB wavenumbers and not for others, is beyond the scope of this work. We therefore require a new degree of freedom with .
Throughout this cosmological epoch the neutrinos are relativistic, so the typical momentum transfer is . Therefore, the expansion indicated by the arrow in Eq. (4) is valid (and we may neglect the momentum and, hence, temperaturedependence in ) only if . Comparing the values in Eq. (2) to the expression for in Eq. (5), we see that
(6) 
so a new subGeV state is generically required to realize the selfinteractingneutrino solution to the tension. Since at horizon entry of the highest momentum modes relevant for CMB anisotropies, the validity of the effective interaction in Eq. (1) in the analyses of Refs. CyrRacine and Sigurdson (2014); Lancaster et al. (2017); Oldengott et al. (2017); Kreisch et al. (2019); Barenboim et al. (2019) requires (as already pointed out in Ref. Kreisch et al. (2019)). From Eq. (6), this condition translates to
(7) 
Finally, we note that the interaction in Eq. (3) is not gaugeinvariant at energies above the scale of electroweak symmetry breaking (EWSB). This is because is a component of the doublet of the group, whose gauge transformations rotate into the corresponding charged lepton . Thus, some particle with nontrivial electroweak quantum numbers must replace the interaction of Eq. (3) at energies above the Higgs vacuum expectation value . We defer a complete discussion of this issue to Sec V; for now we assume that the interaction of Eq. (1) is valid well below the electroweak scale.
Iii Cosmological Bounds
The successful predictions of BBN in standard cosmology (with only SM particle content) provide a powerful probe of additional light species. New particles in thermal equilibrium with neutrinos increase the Hubble expansion rate during BBN as extra relativistic degrees of freedom or by heating neutrinos relative to photons. Away from mass thresholds, both effects are captured by a constant shift in , the effective number of neutrinos. We find that the observed light element abundances constrain () at 95% CL for the SI (MI) preferred values of the baryon density, as explained in more detail in App. A.
Note that we do not use the Planck limit on Aghanim et al. (2018), since that result assumes freestreaming neutrinos. In fact, a large at the time of the formation of the CMB is crucial for reconciling the MI and SI cosmologies with the observed CMB power spectrum Kreisch et al. (2019). The concordance of the CMB and BBN measurements of within the selfinteracting neutrino framework thus seems to require an injection of energy between nucleosynthesis and recombination; we remain agnostic on this point. Given this discussion, we conservatively only apply a constraint on from considerations of BBN physics alone.
iii.1 Mediators and
The interaction in Eq. (3) induces decays and inverse decays; this process can bring into into thermal equilibrium with the neutrino bath before BBN.^{1}^{1}1Annihilation and scattering processes (e.g. and ) also contribute to equilibration, but the corresponding rates are suppressed by additional powers of the small coupling , so their contribution is always subdominant. From Eqs. (6) and (7), the mediator that generates sufficiently large neutrino selfinteractions must have mass in the range. Consequently, this particle is at least semirelativistic when neutrinos decouple from photons at MeV. We can estimate the resulting contribution to the expansion rate as parametrized by assuming is relativistic:
(8) 
where is the energy density of photons and is their temperature. If becomes nonrelativistic before the end of BBN, is increased further as decays raise the neutrino temperature (see, e.g., Ref. Boehm et al. (2012)). Here, we show that a mediator whose interactions realize the necessary in Eq. (2) necessarily equilibrates with neutrinos before .
Vector Mediators: If Eq. (1) arises from integrating out a vector particle with mass , the Lagrangian at energies above is
(9) 
where is the dimensionless gauge coupling. At low energies, the interaction in Eq. (9) yields
(10) 
The vector equilibrates before via if the corresponding thermally averaged rate exceeds Hubble when :
(11) 
where GeV and we have used Eqs. (10) and (7). Thus, this reaction is in equilibrium for all values of couplings and masses of interest. As a result, has a thermal number density at in both the MI and SI scenarios. Using Eq. (8), we find assuming remains relativistic throughout BBN. If, instead, becomes nonrelativistic between and keV at the end of BBN, then , so must become nonrelativistic well before . In Ref. Escudero et al. (2019) it was found that the Boltzmann suppression for massive vector particles is effective for at 95% CL. Using Eq. (10), this requires , which is excluded in all theoreticallyconsistent (or anomalyfree) vector models with neutrino couplings Bauer et al. (2018). Anomalyfree vectors, such as those coupled to currents, where is a lepton family number, would additionally introduce large interactions which would likely spoil the CMB fit.
Finally, we note that Eq. (9) can arise if active neutrinos mix with heavy sterileneutrinos, which couple directly to the vector Schmaltz and Weiner (2019); Bertuzzo et al. (2018); Ballett et al. (2019). The mediators in these models are subject to the same BBN constraint above, but because their coupling to active neutrinos is derived from the activesterile mixing angle, the resulting faces additional bounds from laboratory searches de Gouvêa and Kobach (2016) that exclude the interaction strengths in Eq. (2) required to resolve the tension.
Scalar Mediators: By a similar argument, any scalar mediator that realizes from the interaction in Eq. (3) with as required by Eq. (7) also has a thermal abundance at . Relativistic scalars in equilibrium with neutrinos contribute for a real (complex) , which has one (two) degree(s) of freedom. As for the vector case, the density must become Boltzmannsuppressed before neutrinophoton decoupling, leading to a lower limit on . We use AlterBBN 2.1 Arbey (2012); Arbey et al. (2018) as described in App. A to obtain the CL lower bounds
(12) 
for the SI preferred values of the baryon density (the corresponding MI bounds are somewhat weaker – see App. A). These bounds are are shown as red vertical lines labeled BBN in Fig. 2.
iii.2 Excluding Dirac Neutrinos
If neutrinos are Dirac particles, then all neutrino masses arise from the Higgsneutrino Yukawa interaction
(13) 
where is the Higgs doublet, is a lepton doublet, is a righthanded neutrino (RHN), and flavor indices have been suppressed. The Weyl fermions and become Dirac partners after EWSB and acquire identical masses. In the absence of additional interactions, the Yukawa coupling in Eq. (13) is insufficient to thermalize right handed helicity states, so the relic neutrinos consist exclusively of lefthanded neutrinos and righthanded antineutrinos Dolgov et al. (1995).
The light mediator interactions in Eq. (2) are much stronger than the weak force at late times, so and can both thermalize. Approximating the RHN production rate as , for eV we have
(14) 
where is the temperature at which RHN production is maximized relative to Hubble expansion. For a more detailed discussion, see App. B.
Observations of neutrino oscillations require that at least two of the light neutrinos are massive, with one heavier than eV and one heavier than eV Esteban et al. (2019). Thus, for all values of we consider in the SI range, at least one RHN will thermalize before BBN. For this reason, we assume that neutrinos are Majorana particles for the remainder of this work.
iii.3 Supernova 1987A
A light new particle with appreciable couplings to neutrinos can alter the qualitative behavior of neutrino emission that was observed in the aftermath of Supernova 1987A. We find that for the parameter space shown in Fig. 2 the particle is trapped and does not drain energy from the supernova. See App. C for further discussion.
Iv Laboratory Bounds
Because terrestrial experiments routinely reach energies above the MeV scale, there are many possible probes of the models we consider. Here we specialize to CPeven scalar mediators without loss of essential generality; we comment on CPodd scalars in Sec. V. The main constraints are as follows:
Double Beta Decay: If in Eq. (3) and is lighter than the value of a doublebetadecaying nucleus, the process may occur, contributing to measured rates. Measurements of this process constrain if MeV Agostini et al. (2015); Blum et al. (2018); Brune and Ps (2018). This constraint is shown in the top row of Fig. 2.
Meson Decays: Nonzero can allow for the meson decays if Blum et al. (2014); Berryman et al. (2018); Kelly and Zhang (2019); Krnjaic et al. (2019). Note that, unlike decays, this process is not helicity suppressed. The strongest constraints come from kaon decay. The measured ratio of constrains as shown in the top row of Fig. 2 Lessa and Peres (2007); Fiorini (2007). The limit on Tanabashi et al. (2018) constrains as shown by the purpleshaded region in bottomleft panel of Fig. 2. Nonzero still contributes when via offshell Kelly and Zhang (2019). Note that, if decays to visible matter via or through some other interaction (e.g. higgs mixing), the branching ratio must be subdominant to decays to avoid tight coupling below MeV, which would spoil the solution to the tension. Thus, the laboratory bounds presented here are unaffected by other possible couplings (see Ng and Beacom (2014) for a review).
Decays: Similar to meson decays, that couples to with allows for decays . For , this constrains Lessa and Peres (2007). This is depicted as a purple band in the bottomright panel of Fig. 2.
In Fig. 2 we summarize our main findings: values of from Eq. (2) favored by the tension are excluded if couples universally to all flavor and masseigenstates (top left panel), which was explicitly considered in Refs. CyrRacine and Sigurdson (2014); Lancaster et al. (2017); Oldengott et al. (2017); Kreisch et al. (2019); Barenboim et al. (2019), or if couples predominantly to or (topright and bottomleft panels, respectively). Similarly, we can rule out the possibility that couples to any single masseigenstate neutrino, since the  and composition of each mass eigenstate is similar. Moreover, in this case, the collisional Boltzmann equations would be much more complicated to solve (because different eigenstates will start to freestream at different times), and the results of Refs. CyrRacine and Sigurdson (2014); Lancaster et al. (2017); Oldengott et al. (2017); Kreisch et al. (2019); Barenboim et al. (2019) may not apply.
However, a flavorrestricted coupling leads to approximately the same neutrino masseigenstate interactions as in Refs. CyrRacine and Sigurdson (2014); Oldengott et al. (2017); Lancaster et al. (2017); Kreisch et al. (2019); Barenboim et al. (2019), since the flavor eigenstates are wellmixed in the mass basis. In particular, a only coupling, in which the matrix is zero except for the entry, is potentially viable since decays are less constraining than meson decays. Thus, we are unable to fully exclude an interaction . In this case, for defined in Eq. (5) and is a constant that accounts for the reduced scattering probability of each mass eigenstate. However, because mixed masseigenstate vertices are possible in this scenario, there are additional diagrams compared to the massdiagonal case. For this reason, we caution that the effect on the CMB anisotropies of flavorspecific neutrino selfinteractions can be mildly different than that considered in Refs. CyrRacine and Sigurdson (2014); Lancaster et al. (2017); Oldengott et al. (2017); Kreisch et al. (2019); Barenboim et al. (2019). Nonetheless, we expect that the preferred coupling range should shift slightly up relative to the flavoruniversal case. For this reason, the SI range is ruled out regardless of choice of flavor structure in Eq. (3). The MI range is still allowed in a flavoronly scenario, though a dedicated study is necessary and well motivated.
V Possible Ultraviolet Completions
In this section we quantify the difficulty of realizing the operator in Eq. (1) from renormalizable, gauge invariant interactions that also accommodate neutrino masses and mixings. In light of the discussion in Sec III, we only consider models of Majorana neutrinos. We focus on the typeI and typeII seesaw mechanisms of Majorana neutrino mass generation with an additional scalar particle . In both cases, we find the resulting Yukawa coupling is suppressed by factors of the light neutrino mass. In these minimal models, it is therefore impossible to simultaneously generate neutrino masses and a large enough in Eq. (2) to address the tension.
v.1 TypeI Seesaw
We consider the Lagrangian
(15) 
where we have suppressed flavor indices. After EWSB, the first term in Eq. (15) and the term generate light neutrino masses after is integrated out:
(16) 
These interactions also generate the coupling of to LH neutrinos in Eq. 3 with
(17) 
Thus, realizing requires . If this is true, the RHNs thermalize before BBN and spoil as discussed in analogy to the Dirac case discussed in Sec. III.2. Therefore, a minimal scenario where the same TypeI Seesaw generates the neutrino mass and the operator in Eq. (1) with the magnitude in Eq. (2) is not possible.
We also consider the case in which, instead of violating lepton number explicitly, we replace the term in parentheses in Eq. (15) with a complex scalar . This acquires a vacuum expectation value (VEV) of , dynamically generating a mass for Chikashige et al. (1981, 1980); Schechter and Valle (1982). The light neutrinos couple to the pseudoGoldstone boson , the Majoron, of this symmetry breaking. The coupling is of the form , so the Majoron is a pseudoscalar mediator. We note that our results from the previous sections apply to the Majoron, since all processes we have considered involve relativistic neutrinos for which the distinction between scalar and pseudoscalar is irrelevant. As with the above case, we see that the neutrino selfinteraction strength will be suppressed by the light neutrino masses.
v.2 TypeII Seesaw
We augment the Standard Model with a complex triplet . The typeII Seesaw Lagrangian is
(18) 
where the dimensionful coupling is a soft breaking of Lepton number. After integrating out and EWSB, neutrinos acquire masses of order
(19) 
and a coupling to of the form of Eq. 3 with
(20) 
The scale is bounded by nonobservation of rare leptonnumberviolating processes such as or Akeroyd et al. (2009); Dinh et al. (2012); Dev et al. (2017). Thus, for keV, the suppression is too severe to permit values of in Eq. (2), as required to resolve the tension.
As with the TypeI Seesaw, we may consider the Majoron case where is replaced by a field which acquires a VEV. The neutrinos couple to the pseudoGoldstone boson of this lepton number symmetry breaking, with coupling of the form . Again, this coupling is highly suppressed.
We conclude that the coupling to active neutrinos cannot arise from interactions with the same field that generates neutrino masses (the RH neutrino in TypeI seesaw and the triplet in TypeII), without being suppressed by neutrino masses. A large enough can therefore be obtained using separate seesaw mechanisms to generate the neutrino masses and the interaction. For instance, we can use the TypeI Seesaw with for the light neutrino masses and add the TypeII Seesaw with to allow for large . In this scenario, the size of is independent of the neutrino masses.
Vi Concluding Remarks
In this work we have shown that the selfinteracting neutrino explanation of the tension requires the existence of a light MeVscale mediator, which is subject to stringent cosmological and laboratory bounds. Consequently, most realizations of this scenario are robustly excluded by the BBNonly bounds on , accelerator constraints for rare decays, and precise results from neutrinoless doublebeta decay searches. If selfinteractions resolve the tension, they cannot viably arise in models with vectormediated forces or with Dirac neutrinos, and they cannot reside in the strongly interacting parameter space in Eq. (2).
Our analysis identifies viable phenomenological models with a 10 MeV scalar mediator and large couplings to flavored neutrinos only. However, realizing such interactions in UVcomplete, gaugeinvariant models is challenging. We find that sufficiently strong interactions cannot arise in models that generate neutrino masses via a single TypeI or II seesaw mechanism: the resulting neutrinoscalar coupling is suppressed by factors of where is the appropriate seesaw scale.
Looking forward, a dedicated exploration of the only scattering scenario favored by our analysis will be necessary to see if the moderately interactingscenario is indeed able to resolve the tension while not running afoul of laboratory measurements. Our results also motivate exploration of the “intermediate” mediatormass regime, where neutrino scattering is relevant for a partial range of redshifts explored by the CMB.
Acknowledgements.
We thank André de Gouvêa, FrancisYan CyrRacine, Joshua Isaacson, Martina Gerbino, Stefan Höche, Massimiliano Lattanzi, Kohta Murase, Jessica Turner, and Yue Zhang for helpful conversations. This manuscript has been authored by Fermi Research Alliance, LLC under Contract No. DEAC0207CH11359 with the U.S. Department of Energy, Office of Science, Office of High Energy Physics. The United States Government retains and the publisher, by accepting the article for publication, acknowledges that the United States Government retains a nonexclusive, paidup, irrevocable, worldwide license to publish or reproduce the published form of this manuscript, or allow others to do so, for United States Government purposes.Appendix A Calculating
The Hubble expansion rate at the time of BBN is sensitive to the energy density in neutrinos and other relativistic species when photon temperatures are below an MeV Serpico and Raffelt (2004); Boehm et al. (2013); Nollett and Steigman (2015). New relativistic particles or an injection of energy into the Standard Model neutrino bath (via the decays or annihilations of a coupled species) can increase the Hubble expansion rate at this time. A larger expansion rate at the time of BBN modifies the neutrontoproton ratio and the freezeout of deuteriumburning reactions, leading to larger yields of Helium4 and Deuterium. The primordial abundances of these elements are measured in pristine gas clouds to be Aver et al. (2015) and D/H Cooke et al. (2018), respectively. These observations are in good agreement with predictions assuming only Standard Model particle content at the time of BBN Cyburt et al. (2016).
We use the measurements of and D/H to constrain modifications to the expansion rate using the BBN Boltzmann code AlterBBN 2.1 Arbey (2012); Arbey et al. (2018). We follow the Monte Carlo procedure outlined in Ref. Berlin et al. (2019) to estimate theoretical uncertainties from nuclear reaction rates. In deriving limits we marginalize over a Gaussian prior on the baryon density corresponding to the bestfit points of CDM Aghanim et al. (2018), SI and MI Kreisch et al. (2019). Using the results of Refs. Steigman (2006); Cyburt et al. (2016) to convert to , we find that , , and in CDM, SI and MI, respectively. Note that the best fit values of in these cosmologies are all compatible with the range (95 % CL) extracted only from BBN data Tanabashi et al. (2018).
If the new particles remain relativistic throughout nucleosynthesis, their modification of the Hubble rate is specified by a constant shift of the number of relativistic degrees of freedom, . We find that the observed values of and D/H favor the values of shown in Tab. 1 for the different cosmologies. Note that the CDM result is higher than reported in, e.g., Ref. Cyburt et al. (2016), and it has a smaller uncertainty. This is because of updated D burning rates and observed abundances used in our analysis; for further discussion, see Ref. Berlin et al. (2019). These improvements actually weaken the BBN constraint compared to Ref. Cyburt et al. (2016) because of the slight preference for . The upper limits in SI and MI are weaker still, due to their larger central values and uncertainties of the baryon density: a larger baryon density reduces the Deuterium yield, which can be compensated by increasing Cyburt et al. (2016).
The upper bounds in Tab. 1 apply to light particles that are fully relativistic at the time of nucleosynthesis. If these particles have masses at the MeV scale, then their decays or annihilations can heat the neutrinos relative to photons; if they are much heavier, they transfer their entropy to neutrinos while those are still in equilibrium with photons. The resulting change in neutrino temperature and the corresponding can be estimated assuming instantaneous neutrinophoton decoupling at and by using entropy conservation Boehm et al. (2012). This crude estimate along with Tab. 1 suggests that neutrinocoupled scalars in the SI cosmology with MeV should be incompatible with the observed light element abundances, where the range corresponds to varying and the number of scalar degrees of freedom from 1 (real scalar) to 2 (complex scalar). The results of a full calculation using AlterBBN (which does not make these approximations), shown in Tab. 2, are compatible with this estimate.
Model  95% upper limit  

CDM  
SI  
MI 
Model  95% CL lower bound on (MeV)  

real  complex  
CDM  
SI  
MI 
We point out that the contribution of a scalar to at the 95% CL limit from our BBN analysis is also compatible with the extra radiation density at the best fit point in Ref. Kreisch et al. (2019) at only the level. In order to be in better agreement with for the SI mode, new semirelativistic degrees of freedom may need to come into equilibrium with the neutrino bath after BBN is complete Berlin et al. (2019).
Appendix B Dirac Neutrino Thermalization
If neutrinos are Dirac fermions, their righthanded (RH) components must not come into equilibrium with the Standard Model plasma before , since they would give during nucleosynthesis, which is clearly incompatible with the results of App. A. In this section we derive the conditions for RH neutrino thermalization.
In the Standard Model, RH neutrinos (neutrinos with the “wrong” helicity) can be created in any interaction that produces LH neutrinos, i.e. weak interactions. The characteristic production rate of RH neutrinos in scattering reactions is Dolgov et al. (1995). If weak interactions are the only interactions producing neutrinos, then RHNs do not thermalize as long as , which is comfortably satisfied in the Standard Model. However, the interaction introduced in Eqs. (1) and (2) is many orders of magnitude stronger than its weak counterpart, and the resulting RH neutrino production rate is much larger.
Thermalization is described by a Boltzmann equation of the form
(21) 
where is the energy density in “wrong helicity” neutrinos and is the collision term encoding the processes that produce these neutrinos. Thermalization/equilibration occurs when . The dominant process responsible for RH neutrino production is the decay , where () is a helicity label corresponding to RH (LH) neutrinos; scattering reactions such as are suppressed by an additional factor of . We therefore evaluate the collision term for :
(22) 
where is the Lorentzinvariant phasespace (including the momentum conservation delta function) and , are the phasespace distributions of and LH neutrinos. We have neglected the inverse decay contribution and the RH neutrino Pauliblocking factors. These approximations are adequate for estimating the onset of equilibrium, assuming the initial abundance of RH neutrinos is negligible, i.e. . We have also multiplied the righthand side by two to account for decays which can also produce wrong helicity neutrinos.
We use the Lagrangian in Eq. (3) (with a complex for a Dirac neutrino) to evaluate the matrix element squared in the plasma frame, without summing over one of the helicities Dreiner et al. (2010)
(23) 
where is the angle between the and direction of motion, and we have kept only the leading terms in . Note that when , the “wrong” helicity amplitude vanishes at this order as a result of angular momentum conservation (the next term is ). Using this result in Eq. (22), we find that
(24) 
where the factor in parentheses is the number density of relativistic and the function is for and becomes Boltzmannsuppressed for .
We can apply the nucleosynthesis bound as computed in Appendix A if thermalization occurs before . By using the approximate thermalization criterion below Eq. (21), we find that RH neutrinos thermalize before BBN if
(25) 
This bound (evaluated using the full numerical ) is shown in Fig. 2 as a black dashed line.
Appendix C Supernova 1987A
A new weakly coupled particle can change the behavior of the neutrino emission that was observed from the explosion of Supernova 1987A. The protoneutron star cooling phase that was observed in large water Cherenkov detectors was qualitatively similar to the Standard Modelonly expectation Burrows and Lattimer (1987). If a new particle species carried away too much energy during the protoneutron star cooling phase, the time over which neutrinos arrived would have been unacceptably reduced Burrows et al. (1989, 1990). A semianalytic criterion that the luminosity of this particle should obey is at times of order 1 second after the core bounce Raffelt (1996). Following the procedure described in more detail in Chang et al. (2017), we have
(26) 
where: the has fourmomentum ; is the absorptive width; is the production rate, which is related to the absorptive width in equilibrium by ; is the radius of the neutrinosphere, outside of which neutrinos freestream; and is the radius inside of which neutrinos gain energy on average in elastic scattering events. We calculate including decay and annihilation to neutrino pairs, both of which are important for the masses of interest. We have not included contributions from the neutrino effective potentials, which may be significant at small Brune and Ps (2018). We have also neglected neutrino Pauli blocking in , which is important near the protoneutron star core, since this will become unimportant between and .
We find that given by Eq. (26) exceeds if is roughly in the range
(27) 
The sharp change in the shape of the bound at is due to the fact that rate of decay and inverse decay, , becomes subdominant to the annihilation rate, , for masses , where is the core temperature. These bounds are approximately compatible with those shown in Brune and Ps (2018) at masses above 10 keV. From Eq. (27), we see that bounds arising from the luminosity of particles from Supernova 1987A are generally below the coupling range of interest in this work.
It is also interesting to understand the constraints on arising from deleptonization of the core, which have been obtained in the mass range in Brune and Ps (2018) and which approximately overlap the range in Eq. (27). At lower masses, this likely has an effect on the early phases of collapse and the collapse progenitor. Such a possibility was suggested in Kolb et al. (1982) and was studied in the aftermath of Supernova 1987A by Fuller et al. (1988). This latter study found a constraint , neglecting any dependence. Because these bounds are determined by physics at the beginning of the core collapse, when temperatures are , the bound likely cuts off at , similar to the bounds cited above. A detailed study is of interest, but beyond the scope of this work.
References
 Riess et al. (2016) A. G. Riess et al., Astrophys. J. 826, 56 (2016), arXiv:1604.01424 [astroph.CO] .
 Shanks et al. (2019) T. Shanks, L. Hogarth, and N. Metcalfe, Mon. Not. Roy. Astron. Soc. 484, L64 (2019), arXiv:1810.02595 [astroph.CO] .
 Riess et al. (2018) A. G. Riess, S. Casertano, D. Kenworthy, D. Scolnic, and L. Macri, (2018), arXiv:1810.03526 [astroph.CO] .
 Aghanim et al. (2018) N. Aghanim et al. (Planck), (2018), arXiv:1807.06209 [astroph.CO] .
 Riess et al. (2019) A. G. Riess, S. Casertano, W. Yuan, L. M. Macri, and D. Scolnic, (2019), arXiv:1903.07603 [astroph.CO] .
 Vonlanthen et al. (2010) M. Vonlanthen, S. Rsnen, and R. Durrer, JCAP 1008, 023 (2010), arXiv:1003.0810 [astroph.CO] .
 Verde et al. (2017) L. Verde, E. Bellini, C. Pigozzo, A. F. Heavens, and R. Jimenez, JCAP 1704, 023 (2017), arXiv:1611.00376 [astroph.CO] .
 Evslin et al. (2018) J. Evslin, A. A. Sen, and Ruchika, Phys. Rev. D97, 103511 (2018), arXiv:1711.01051 [astroph.CO] .
 Aylor et al. (2019) K. Aylor, M. Joy, L. Knox, M. Millea, S. Raghunathan, and W. L. K. Wu, Astrophys. J. 874, 4 (2019), arXiv:1811.00537 [astroph.CO] .
 Lesgourgues et al. (2016) J. Lesgourgues, G. MarquesTavares, and M. Schmaltz, JCAP 1602, 037 (2016), arXiv:1507.04351 [astroph.CO] .
 Di Valentino et al. (2018) E. Di Valentino, C. Behm, E. Hivon, and F. R. Bouchet, Phys. Rev. D97, 043513 (2018), arXiv:1710.02559 [astroph.CO] .
 Poulin et al. (2018a) V. Poulin, T. L. Smith, D. Grin, T. Karwal, and M. Kamionkowski, Phys. Rev. D98, 083525 (2018a), arXiv:1806.10608 [astroph.CO] .
 D’Eramo et al. (2018) F. D’Eramo, R. Z. Ferreira, A. Notari, and J. L. Bernal, JCAP 1811, 014 (2018), arXiv:1808.07430 [hepph] .
 Poulin et al. (2018b) V. Poulin, T. L. Smith, T. Karwal, and M. Kamionkowski, (2018b), arXiv:1811.04083 [astroph.CO] .
 Pandey et al. (2019) K. L. Pandey, T. Karwal, and S. Das, (2019), arXiv:1902.10636 [astroph.CO] .
 Escudero et al. (2019) M. Escudero, D. Hooper, G. Krnjaic, and M. Pierre, JHEP 03, 071 (2019), arXiv:1901.02010 [hepph] .
 CyrRacine and Sigurdson (2014) F.Y. CyrRacine and K. Sigurdson, Phys. Rev. D90, 123533 (2014), arXiv:1306.1536 [astroph.CO] .
 Lancaster et al. (2017) L. Lancaster, F.Y. CyrRacine, L. Knox, and Z. Pan, JCAP 1707, 033 (2017), arXiv:1704.06657 [astroph.CO] .
 Oldengott et al. (2017) I. M. Oldengott, T. Tram, C. Rampf, and Y. Y. Y. Wong, JCAP 1711, 027 (2017), arXiv:1706.02123 [astroph.CO] .
 Kreisch et al. (2019) C. D. Kreisch, F.Y. CyrRacine, and O. Dor, (2019), arXiv:1902.00534 [astroph.CO] .
 Bashinsky and Seljak (2004) S. Bashinsky and U. Seljak, Phys. Rev. D69, 083002 (2004), arXiv:astroph/0310198 [astroph] .
 Barenboim et al. (2019) G. A. Barenboim, P. B. Denton, and I. M. Oldengott, Phys. Rev. D99, 083515 (2019), arXiv:1903.02036 [astroph.CO] .
 Cohen (2019) T. Cohen, (2019), arXiv:1903.03622 [hepph] .
 Forastieri et al. (2015) F. Forastieri, M. Lattanzi, and P. Natoli, JCAP 1507, 014 (2015), arXiv:1504.04999 [astroph.CO] .
 Forastieri et al. (2019) F. Forastieri, M. Lattanzi, and P. Natoli, (2019), arXiv:1904.07810 [astroph.CO] .
 Blum et al. (2014) K. Blum, A. Hook, and K. Murase, (2014), arXiv:1408.3799 [hepph] .
 Berryman et al. (2018) J. M. Berryman, A. De Gouvêa, K. J. Kelly, and Y. Zhang, Phys. Rev. D97, 075030 (2018), arXiv:1802.00009 [hepph] .
 Kelly and Zhang (2019) K. J. Kelly and Y. Zhang, Phys. Rev. D99, 055034 (2019), arXiv:1901.01259 [hepph] .
 Krnjaic et al. (2019) G. Krnjaic, G. MarquesTavares, D. Redigolo, and K. Tobioka, (2019), arXiv:1902.07715 [hepph] .
 Agostini et al. (2015) M. Agostini et al., Eur. Phys. J. C75, 416 (2015), arXiv:1501.02345 [nuclex] .
 Blum et al. (2018) K. Blum, Y. Nir, and M. Shavit, Phys. Lett. B785, 354 (2018), arXiv:1802.08019 [hepph] .
 Brune and Ps (2018) T. Brune and H. Ps, (2018), arXiv:1808.08158 [hepph] .
 Boehm et al. (2012) C. Boehm, M. J. Dolan, and C. McCabe, JCAP 1212, 027 (2012), arXiv:1207.0497 [astroph.CO] .
 Bauer et al. (2018) M. Bauer, P. Foldenauer, and J. Jaeckel, JHEP 07, 094 (2018), arXiv:1803.05466 [hepph] .
 Schmaltz and Weiner (2019) M. Schmaltz and N. Weiner, JHEP 02, 105 (2019), arXiv:1709.09164 [hepph] .
 Bertuzzo et al. (2018) E. Bertuzzo, S. Jana, P. A. N. Machado, and R. Zukanovich Funchal, Phys. Rev. Lett. 121, 241801 (2018), arXiv:1807.09877 [hepph] .
 Ballett et al. (2019) P. Ballett, S. Pascoli, and M. RossLonergan, Phys. Rev. D99, 071701 (2019), arXiv:1808.02915 [hepph] .
 de Gouvêa and Kobach (2016) A. de Gouvêa and A. Kobach, Phys. Rev. D93, 033005 (2016), arXiv:1511.00683 [hepph] .
 Arbey (2012) A. Arbey, Comput. Phys. Commun. 183, 1822 (2012), arXiv:1106.1363 [astroph.CO] .
 Arbey et al. (2018) A. Arbey, J. Auffinger, K. P. Hickerson, and E. S. Jenssen, (2018), arXiv:1806.11095 [astroph.CO] .
 Dolgov et al. (1995) A. D. Dolgov, K. Kainulainen, and I. Z. Rothstein, Phys. Rev. D51, 4129 (1995), arXiv:hepph/9407395 [hepph] .
 Esteban et al. (2019) I. Esteban, M. C. GonzalezGarcia, A. HernandezCabezudo, M. Maltoni, and T. Schwetz, JHEP 01, 106 (2019), arXiv:1811.05487 [hepph] .
 Lessa and Peres (2007) A. P. Lessa and O. L. G. Peres, Phys. Rev. D75, 094001 (2007), arXiv:hepph/0701068 [hepph] .
 Fiorini (2007) L. Fiorini (NA48/2), Tau lepton physics. Proceedings, 9th International Workshop, TAU’06, Pisa, Italy, September 1922, 2006, Nucl. Phys. Proc. Suppl. 169, 205 (2007), [,205(2007)].
 Tanabashi et al. (2018) M. Tanabashi et al. (Particle Data Group), Phys. Rev. D98, 030001 (2018).
 Ng and Beacom (2014) K. C. Y. Ng and J. F. Beacom, Phys. Rev. D90, 065035 (2014), [Erratum: Phys. Rev.D90,no.8,089904(2014)], arXiv:1404.2288 [astroph.HE] .
 Chikashige et al. (1981) Y. Chikashige, R. N. Mohapatra, and R. D. Peccei, Phys. Lett. 98B, 265 (1981).
 Chikashige et al. (1980) Y. Chikashige, R. N. Mohapatra, and R. D. Peccei, Phys. Rev. Lett. 45, 1926 (1980), [,921(1980)].
 Schechter and Valle (1982) J. Schechter and J. W. F. Valle, Phys. Rev. D25, 774 (1982).
 Akeroyd et al. (2009) A. G. Akeroyd, M. Aoki, and H. Sugiyama, Phys. Rev. D79, 113010 (2009), arXiv:0904.3640 [hepph] .
 Dinh et al. (2012) D. N. Dinh, A. Ibarra, E. Molinaro, and S. T. Petcov, JHEP 08, 125 (2012), [Erratum: JHEP09,023(2013)], arXiv:1205.4671 [hepph] .
 Dev et al. (2017) P. S. B. Dev, C. M. Vila, and W. Rodejohann, Nucl. Phys. B921, 436 (2017), arXiv:1703.00828 [hepph] .
 Serpico and Raffelt (2004) P. D. Serpico and G. G. Raffelt, Phys. Rev. D70, 043526 (2004), arXiv:astroph/0403417 [astroph] .
 Boehm et al. (2013) C. Boehm, M. J. Dolan, and C. McCabe, JCAP 1308, 041 (2013), arXiv:1303.6270 [hepph] .
 Nollett and Steigman (2015) K. M. Nollett and G. Steigman, Phys. Rev. D91, 083505 (2015), arXiv:1411.6005 [astroph.CO] .
 Aver et al. (2015) E. Aver, K. A. Olive, and E. D. Skillman, JCAP 1507, 011 (2015), arXiv:1503.08146 [astroph.CO] .
 Cooke et al. (2018) R. J. Cooke, M. Pettini, and C. C. Steidel, Astrophys. J. 855, 102 (2018), arXiv:1710.11129 [astroph.CO] .
 Cyburt et al. (2016) R. H. Cyburt, B. D. Fields, K. A. Olive, and T.H. Yeh, Rev. Mod. Phys. 88, 015004 (2016), arXiv:1505.01076 [astroph.CO] .
 Berlin et al. (2019) A. Berlin, N. Blinov, and S. W. Li, (2019), arXiv:1904.04256 [hepph] .
 Steigman (2006) G. Steigman, JCAP 0610, 016 (2006), arXiv:astroph/0606206 [astroph] .
 Dreiner et al. (2010) H. K. Dreiner, H. E. Haber, and S. P. Martin, Phys. Rept. 494, 1 (2010), arXiv:0812.1594 [hepph] .
 Burrows and Lattimer (1987) A. Burrows and J. M. Lattimer, Astrophys. J. 318, L63 (1987).
 Burrows et al. (1989) A. Burrows, M. S. Turner, and R. P. Brinkmann, Phys. Rev. D39, 1020 (1989).
 Burrows et al. (1990) A. Burrows, M. T. Ressell, and M. S. Turner, Phys. Rev. D42, 3297 (1990).
 Raffelt (1996) G. G. Raffelt, Stars as laboratories for fundamental physics (1996).
 Chang et al. (2017) J. H. Chang, R. Essig, and S. D. McDermott, JHEP 01, 107 (2017), arXiv:1611.03864 [hepph] .
 Kolb et al. (1982) E. W. Kolb, D. L. Tubbs, and D. A. Dicus, The Astrophysical Journal 255, L57 (1982).
 Fuller et al. (1988) G. M. Fuller, R. Mayle, and J. R. Wilson, The Astrophysical Journal 332, 826 (1988).