title: Einstein-Cartan Portal to Dark Matter
authors: Mikhail Shaposhnikov Andrey Shkerin Inar Timiryasov Sebastian Zell
text: 
  
  
  
  
  
  
  Institute of Physics, Laboratory for Particle Physics and Cosmology,
  École Polytechnique Fédérale de Lausanne, CH-1015 Lausanne, 
  Switzerland
  
  
  
  It is well-known since the works of Utiyama and Kibble that the gravitational force can be obtained by gauging the Poincaré group, which puts gravity on the same footing as the Standard Model fields. The resulting theory – Einstein-Cartan gravity – inevitably contains four-fermion and scalar-fermion interactions that originate from torsion associated with spin degrees of freedom. We show that these interactions lead to a novel  mechanism for producing singlet fermions in the Early Universe. These fermions can play the role of dark matter particles. The mechanism is operative in a large range of dark matter particle masses: from a few keV up to ∼ 10^8GeV. We discuss potential observational consequences of keV-scale dark matter produced this way, in particular for right-handed neutrinos. We conclude that a determination of the primordial dark matter momentum distribution might be able to shed light on the gravity-induced fermionic interactions.
  
  
  Einstein-Cartan Portal to Dark Matter
      Sebastian Zell
      December 30, 2023
  =====================================
  
  
  
  
  
  Introduction.—Gravity is a universal force, and it inevitably reveals some properties of yet-to-be-discovered constituents of Nature – dark matter (DM) and dark energy (see, e.g., <cit.> for reviews). Indeed, all evidence for DM in galaxies, galaxy clusters and at large scales is due to the gravitational interaction. Thus, we know well about how DM gravitates, but little about what it is. The list of DM particle candidates includes, e.g., weakly interacting massive particles (see <cit.> for important early contributions), sterile neutrinos <cit.>, axions <cit.> and axion-like particles <cit.>. Regardless the nature of DM, any model of it should explain how it was produced in the Early Universe and how its abundance is maintained.
  
  In this paper, we argue that gravity might be able to tell us not only about the distribution of DM in the Universe, but also about the mechanism of its production. The two key ingredients of the framework allowing for such a statement
  are the Einstein-Cartan (EC) formulation of gravity <cit.> and the assumption about the fermion nature of a DM candidate.
   Regarding the second one, we take a DM particle to be a Dirac or Majorana fermion which is a singlet under the gauge group of the SM. 
  
  
  We begin with a brief discussion of EC gravity.[See, e.g., <cit.> for reviews.] There is no doubt that below the Planck scale, General Relativity (GR) provides an elegant and accurate description of gravity. Nevertheless, this leaves still unanswered the question about which formulation of GR one should employ. An important alternative to the most widely used metric gravity is the EC formulation. In this theory, the role of fundamental fields is played by the tetrad and the spin connection, in terms of which the metric and the Christoffel symbols are introduced. The latter are, in general, not symmetric in the lower indices, hence EC gravity contains torsion. Still, the number of propagating degrees of freedom – two of the graviton – is the same as in the metric formulation. A conceptual advantage of EC gravity is that it can be viewed as a gauge theory of the Poincaré group, thus allowing for a similar treatment of all fundamental forces.[Moreover, we note that in EC gravity there is no need for the Gibbons-Hawking-York boundary term <cit.> for the variational problem to be well-posed.]
  
  In the absence of matter, the EC and metric formulations of gravity are equivalent. This changes once scalar fields coupled non-minimally to the Ricci scalar are introduced. The resulting theory is then equivalent to the Palatini formulation of gravity <cit.> (see also <cit.>). In this version of gravity, the metric and the Christoffel symbols are viewed as independent fundamental variables, and the latter are assumed to be symmetric in the lower indices. Another way to break the equivalence of the EC and metric formulations is to take into account fermions since the latter source torsion. The theory admits an exact solution for the torsionful part of the connection. Plugging this solution back in the action, one arrives at an equivalent theory in the metric formulation. The difference of the two theories of gravity is then manifested in the appearance of dimension-six terms representing an interaction of fermionic axial currents <cit.>.[Of course, one could have started from the beginning in a metric theory with additional higher-dimensional operators. In this case, however, a consistent effective field theory approach would dictate that all possible higher-dimensional interactions (consistent with relevant symmetries) are taken into account.] Their strength is fixed and suppressed by 1/M_P^2.
  
  In the statements thus far, we have considered in EC gravity the same action as in the metric theory. Due to the presence of torsion, however, one can form more terms of mass dimension not bigger than four than in the metric case. Specifically, the fermion kinetic terms can be generalized by introducing non-minimal fermion couplings <cit.>. In an equivalent metric theory, they lead to vector-vector, axial-vector and axial-axial fermion current interactions. Furthermore, the gravitational action of the EC theory can be extended by adding the Holst term <cit.>, which modifies the four-fermion interaction. Note that the additional terms come with a priori unknown coupling constants and the strength of current-current interactions depends on these couplings.
  
  The present paper uses the results of <cit.>, where a theory of scalar and fermion fields coupled to gravity was studied. There, we included additional terms of mass dimension not bigger than four that are specific to the EC formulation and derived the equivalent metric theory. When applied to cosmology and experiment, the scalar field can be associated with the SM Higgs field, and the fermions with the SM quarks and leptons as well as, possibly, additional species such as right-handed neutrinos. The phenomenology of the scalar-gravity part of the EC theory has already been investigated in <cit.>. There, we considered the Higgs field as an inflaton. Our study was motivated by the well-known fact that the Higgs field can be responsible for inflation provided that it couples non-minimally to the Ricci scalar <cit.>. The models of Higgs inflation in the metric <cit.> and Palatini <cit.> formulations of gravity find their natural generalization within the EC theory <cit.> (see also <cit.>).
  
  The goal of this paper is to study phenomenology of the fermionic sector of the EC theory. Namely, we show that the four-fermion interaction originating from EC gravity can be responsible for DM production. We compute the abundance and the spectrum of produced (Dirac or Majorana) particles and show that the right amount of DM can be generated for a wide range of fermion masses. We also discuss an interesting case of warm DM where the primordial  momentum distribution characteristic for EC gravity can potentially be observable. 
  
  
  Einstein-Cartan gravity and fermions.—In this work, we focus on the fermion-gravity part of the general theory studied in <cit.>. To simplify the presentation, we only keep the Einstein-Hilbert term and the non-minimal couplings of fermions in the action.  We comment on the inclusion of other terms later; see also appendix A.
  Then, for each (Dirac or Majorana) fermion species Ψ the relevant part of the action reads as follows <cit.>:[We work in natural units ħ=c=1 and use the metric signature (-1,+1,+1,+1). The matrix γ^5 is defined as γ^5=-iγ^0γ^1γ^2γ^3.]
  
      Ł    =1/2M_P^2R + i/2Ψ̅(1-iα -iβγ^5)γ^μ D_μΨ
          -i/2D_μΨ(1+iα+iβγ^5)γ^μΨ ,
  
  where M_P = 2.435× 10^18 GeV is the Planck mass and D_μ is the covariant derivative containing the connection field. The real couplings α, β are chosen to be the same for all generations of fermions which implies the universality of gravity in the fermionic sector. Allowing for the couplings to depend on a generation index yields qualitatively the same results. In metric gravity, the non-minimal terms sum up to a total derivative, but in the torsionful case they contribute to the dynamics of the theory. 
  
  The theory (<ref>) can be resolved for torsion explicitly <cit.>. Upon substituting the solution for torsion back to the action, one obtains an equivalent metric theory with extra higher-dimensional fermion interaction terms. They read:
  
      0.92[1]Ł_4f= 3α^216  M_P^2 V^μ V_μ + 3  α β8 M_P^2 V^μ A_μ
      -3-3β^216 M_P^2 A^μ A_μ ,
  
  where V^μ =N̅γ^μ N+∑_XX̅γ^μ X is the vector fermion current and A^μ is the analogous axial current (with γ^μ replaced by γ^5 γ^μ). The sum is performed over all SM fermion species  X. For convenience, we wrote separately the terms containing N which plays the role of DM  and can be Dirac or right-handed Majorana fermion. The interaction (<ref>) vanishes only if α=0, β=± 1, and in what follows we do not consider this particular choice of the couplings. 
  
  Adding the Holst term to Eq. (<ref>) modifies the couplings in Eq. (<ref>). Furthermore, including a scalar field coupled non-minimally to gravity results in additional scalar-fermion interactions. Below we focus on the four-fermion interaction (<ref>) whose contribution to the DM production dominates for small masses of N. We discuss the general case in appendix A, which includes Refs. <cit.>.
  
  
  
  Thermal production of singlet fermions.—The four-fermion interaction (<ref>) opens up the production channel of N-particles through the annihilation of the SM fermions X, via the reaction X + X̅→ N + N̅.[The production of singlet fermions due to some higher-dimensional operators was considered in <cit.>. However, the four-fermion interaction which appears in EC gravity was not accounted for.] The kinetic equation corresponding to this reaction takes the form
  
      (∂/∂ t-H q_i∂/∂ q_i) f_N(t, q⃗)=R(q⃗, T)  ,
  
  where f_N is the phase-space density of N, H is the Hubble rate and R is the collision integral, also referred to as a production rate. In an isotropic background, both f_N and R depend only on the absolute value of the spatial momentum |q|≡ |q⃗|.
  
  In what follows, we assume that all SM particles, including fermions, are in thermal equilibrium at the moment of DM production. To check the validity of this assumption, one would need a careful examination of dynamics of bosonic and fermionic SM species at and after preheating, which goes beyond the scope of the present paper. We expect, however, that deviations from thermality do not change qualitatively our results.
  
  As long as the concentration of N remains small and we can neglect the inverse processes, the collision integral in <ref> reads 
  
      R=1/2 |q|∑_X ∫Π  |ℳ_X|^2  f_X( p_1) f_X̅(p_2)  ,
  
  where the sum runs over the SM species (24 left- and 21 right-handed fermions), ℳ_X is the amplitude of the process summed over all spinor indices, and f_X is the distribution of X which we assume to be the Fermi-Dirac one. The differential Π is the phase space volume element accounting for energy conservation; see appendix B for its definition, where also details on the derivations of subsequent res_rategeneral_result are given. 
  The typical values of the momenta in <ref> are large compared to the mass of N, so we neglect all masses when computing ℳ_X.
  Introducing the dimensionless variable y = E/T, where T is the temperature of the cosmic plasma   and E = |q|,
  we arrive at 
  
      R(E,T) = T^5  C_f(α, β)/M_P^4 r(y)  .
   
  Here C_f(α, β) is a combination of the non-minimal fermion couplings whose precise form depends on whether N is Majorana (f=M) or Dirac (f=D): 
  
      C_M = 94(24a^4 + 21 b_+^4 )  ,   C_D = 94(45a^4+ 21 b_+^4 + 24 b_-^4)  ,
  
  where a^2=1+α^2-β^2, b_±^2=1-(α±β)^2. Next, r(y) is a function that we compute numerically (following <cit.>). It is accurately approximated by[This expression is exact if instead of the Fermi-Dirac distribution, one uses the Boltzmann distribution for f_X.]
  
      r(y) ≃1/24π^3 y f_X  .
  
  <ref> can now be easily integrated, leading to
   
      f_N(y) = C_f   T^3_prod M_0(T_prod)/3   M_P^4 r(y)  ,
  
  where T_prod is the temperature at which the DM production begins,
  M_0(T)=M_P √(90/π^2 g_eff(T)), and g_eff(T_ prod)=106.75 is the number of effectively massless degrees of freedom at high temperature. 
  Plugging in the numbers, we obtain for the abundance of N-particles:[In deriving this result, we assume that all other possible interactions of N-particles with particles of the SM are not essential for the DM production.]
  
      Ω_N/Ω_DM≃ 3.6 · 10^-2  C_f ( M_N/10 keV)  ( T_prod/M_P)^3  ,
  
  where Ω_DM is the observed DM abundance and the coefficient C_f is defined in <ref>. <ref> shows that, depending on the value of C_f, the right amount of DM can be generated in a broad range of fermion masses M_N. 
  
  To proceed further, we need an estimate for the production temperature. We obtain it within the framework of Higgs inflation in the Palatini formulation of gravity <cit.>. In this model, preheating is almost instantaneous <cit.>, and one can take T_ prod∼ T_ reh where
  
      T_ reh≃(15λ/2π^2 g_ eff)^1/4M_P/√(ξ)
  
  is the preheating temperature, λ is the Higgs field self-coupling and ξ is the non-minimal coupling of the Higgs field to the Ricci scalar. Both λ and ξ are taken at a high energy scale. Using ξ=10^7 and λ=10^-3 <cit.>, we obtain T_ prod≃ 4× 10^13 GeV.
  
  Now we can investigate two particularly interesting cases. The first one is the limit of vanishing non-minimal couplings, α = β = 0. Then, from <ref> we obtain that Ω_N≃Ω_DM if M_N≃ 6× 10^8 GeV for the Majorana fermion and M_N≃ 3× 10^8 GeV for the Dirac fermion. We conclude that heavy fermion DM can be produced in EC gravity even if the action of the EC theory is identical to that of the metric theory.[Such a heavy DM can also be produced due to graviton exchange between fermionic currents.] Interestingly, the given bounds on M_N are close to the bound M_N≲ 10^9 GeV above which N-particles are overproduced due to the varying geometry at the radiation-dominated stage of the Universe <cit.>.
  
  The second case corresponds to setting α∼β∼√(ξ). With this choice, the scale of suppression of the interaction (<ref>) coincides with the inflationary cutoff scale which in Palatini Higgs inflation is of the order of M_P/√(ξ) <cit.>. 
  For both the Majorana and Dirac cases, <ref> becomes
  
      0.92[1]Ω_NΩ_DM≃ 1.4   √(ξ)λ^3/4g_ eff^3/4 (α + β)^4ξ^2( M_N10 keV)  ( T_prodT_ reh)^3 .
  
  Thus, the right amount of DM is generated for a keV-scale M_N.
  
  
  Einstein-Cartan portal to warm dark matter.—Let us discuss in more detail the second choice of the non-minimal couplings. Then, N is an example of warm DM, and its free-streaming length affects structure formation. Since DM is produced at very high temperatures, its spectrum is redshifted. Consequently, the average momentum is only ≃ 0.61 of the equilibrium momentum at T=1MeV. Depending on the history of reionization, such colder DM candidate can provide a better fit to the Lyman-α data than pure cold DM <cit.>.
  
  
  In <ref> we show the spectrum f_N of DM produced in EC gravity. For comparison, we also show the spectrum of resonantly produced sterile neutrino taken from <cit.>[This result depends on the assumptions about the lepton asymmetry and the values of Yukawa coupling. The curve in <ref> corresponds to the “case a”, sin ^2(2 θ)=20 × 10^-11 in terms of <cit.>.] and the thermal spectrum at 1MeV. We see that the DM momentum distribution (<ref>) has a unique form and  differs from both non-resonantly <cit.> and resonantly produced sterile neutrinos <cit.>. This potentially allows one to distinguish DM produced via the EC portal from other DM candidates.
  We conclude that, at least for certain values of M_N, the distribution of DM in the Universe can bear the information about the gravity-induced fermionic interactions.
  
  
  νMSM.—The lack of a DM candidate is a famous shortcoming of the SM, but it is not the only one. It also cannot explain neutrino oscillations and the baryon asymmetry of the Universe. A possible way to address these three issues at once is the Neutrino Minimal Standard Model (νMSM), which extends the particle content of the SM by three right-handed neutrinos N_1,2,3 <cit.>. One of them, N_1, can be the DM candidate. A lower bound on its mass, M_N_1≳ 1keV, comes from small-scale structures in the matter power spectrum, as inferred, e.g., from Lyman-α measurements. The mixing of N_1 with active neutrinos is bounded from above by X-ray constraints. The other two right-handed neutrinos N_2,3 have nearly-degenerate GeV-scale masses and are responsible for generating the baryon asymmetry. Moreover, the parameters of the model can be chosen such that the observed pattern of neutrino oscillations is obtained
  (for reviews see <cit.>).
  
  So far, the production of DM sterile neutrino was associated with the mixing of N_1 with ordinary neutrinos and typical temperatures 100-300 MeV <cit.>. The so-called non-resonant production of N_1 <cit.> has already been excluded by X-ray searches of radiatively decaying DM <cit.>. The resonant production of DM sterile neutrino <cit.> requires large lepton asymmetries which can be produced in interactions of N_2,3 <cit.>. It can be successful albeit a fine-tuning of parameters is required <cit.>. For resonantly produced DM, the non-observation of X-ray decays of DM in galactic halos leads to an upper bound on the mass, M_N_1≲ 50 keV <cit.>, thereby constraining M_N_1 in a quite narrow range.
  
  The present work supplies the νMSM with a different mechanism to produce the N_1-particles. Sterile neutrino DM produced this way may be absolutely stable and, therefore, is not subject to any X-ray constraints. This opens up a new interesting mass window up to M_N_1∼ 10^8 GeV. The prediction of the νMSM that the lightest active neutrino is effectively massless remains in force also for large values of M_N_1; see appendix C.[Note that this prediction can potentially be tested by the Euclid space mission <cit.>.]
  
  It is very intriguing that if the perturbative cutoff is universal for both the fermion and scalar-gravity sectors of the EC theory (∼ M_P/√(ξ)), then the mass of the DM sterile neutrino is required to be in the keV range. This is the domain where the warm nature of the DM particle is most visible and in which
  the most intensive searches of the radiatively decaying DM are being carried out.
  
  We note that Majorana fermions N_2,3  will be also produced by the same mechanism. Their number density is the same as that of N_1, n_2,3≃  n_1, and the latter is fixed by the DM abundance n_1≃  10^-2n_eq (10  keV/M_N_1). However, this density is too small to affect the analysis of the νMSM leptogenesis carried out previously (e.g. <cit.>), see appendix D.
   
  
  Discussion and outlook.—Once gravity is coupled to matter, such as fermions or a non-minimally coupled scalar field, its different formulations are no longer equivalent. As long as they are consistent, only observations can help us to distinguish between them. In the present work, we have shown that properties of fermionic dark matter may be able to discern the Einstein-Cartan (EC) theory of gravity from the most commonly used metric formulation. In particular, the universal dimension-six interactions of fermionic currents in EC theory can cause the production of the observed amount of dark matter for a wide range of fermion masses. Moreover, they lead to a characteristic momentum distribution of dark matter, which can serve to confirm or exclude our proposed production mechanism. 
  
  On the one hand, these findings are relevant for any extension of the Standard Model (SM) by sufficiently heavy fermions. On the other hand, an exciting unified picture of gravity and the SM emerges. It relies on the EC formulation of gravity and the extension of SM by three right-handed neutrinos, i.e. the Neutrino Minimal Standard Model (νMSM) <cit.>. When non-minimally coupled to the Ricci scalar, the Higgs field can assume the role of the inflaton <cit.>. As discussed in <cit.>, the resulting inflationary scenario generalizes the model of Palatini Higgs inflation <cit.> and is fully compatible with observations. On its own, the νMSM is able to explain neutrino oscillations, baryogenesis and provides a dark matter candidate in the form of a right-handed neutrino. The result of the present work is that the EC formulation of gravity can lead to the production of this sterile neutrino in an amount that matches the observed abundance of dark matter. Finally, we remark that the Palatini formulation of gravity and, by generalization, EC gravity is convenient for addressing the question of the big difference between the Electroweak and the Planck scales <cit.> (see also <cit.>).
  
  
  
  Acknowledgments.—We thank Alexey Boyarsky and Oleg Ruchayskiy for useful discussions. The work was supported by ERC-AdG-2015 grant 694896 and by the Swiss National Science
  Foundation Excellence grant 200020B 182864.
  
  
  
  
  
  
  
  
  
  Appendix A.—<ref> is not the most general dimension-four Lagrangian allowed in EC theory. Due to the presence of torsion, its gravitational part can be extended by adding the Holst term <cit.>:
  
      Δℒ_ grav.=M_P^2/4γ̅ϵ^μνρσR_μνρσ .
  
  The coupling γ̅ is called the Barbero-Immirzi parameter <cit.>. The presence of the Holst term modifies the coefficients in the four-fermion interaction terms: <ref> becomes <cit.>
  
      ℒ_4f =    3/16M_P^2α^2γ̅^2/γ̅^2 +1V^μ V_μ
         - 3/8M_P^2αγ̅^2/γ̅^2 +1(1/γ̅-β) V^μ A_μ
         - 3/16M_P^2γ̅^2/γ̅^2 +1(1+2 β/γ̅-β^2)A^μ A_μ .
  
  In the limit γ̅→∞, <ref> is restored. In the limit γ̅→ 0, all coefficients vanish. This agrees with the fact that the latter limit corresponds to the case of vanishing torsion (see, e.g., <cit.>). 
  
  The action of EC theory can also be extended by including non-minimal coupling of a scalar field h (such as the Higgs field) to gravity. This coupling does not contribute to the four-fermion interaction in the metric description. Instead, it gives rise to scalar-fermion interaction terms <cit.>.
  The scalar-fermion part of the action reads
  
      S_s f   = ∫^4x√(-g)3α/4( ∂_μΩ^2/Ω^2+γ/γ^2+1( ∂_μη̅/Ω^2+∂_μγ) ) V^μ
          + ∫^4x√(-g)3/4( β∂_μΩ^2/Ω^2+1+γβ/γ^2+1( ∂_μη̅/Ω^2+∂_μγ) ) A^μ
  
  with Ω^2=1+ξ h^2/M_P^2 and
  
      γ(h) =1/γ̅Ω^2(1+ξ_γ h^2/M_P^2)  ,    η(h) =ξ_η h^2/M_P^2 .
  
  In the expressions above h is the Higgs field in the unitary gauge, whereas ξ_γ and ξ_η are the non-minimal couplings of h to the Holst term and Nieh-Yan invariant respectively, see <cit.> for details.
  We study this action in the limit h≪ M_P/√(ξ) and h≪ M_P/√(ξ_γ).
  Upon integrating by parts, we obtain the following interaction Lagrangian:
  
      ℒ_s f    = -3α/4(ξ+γ̅ξ_η+ξ_γ-ξ/γ̅^2+1)(h/M_P)^2∂_μ V^μ
          - 3/4(βξ+(γ̅+β)(γ̅ξ_η+ξ_γ-ξ)/γ̅^2+1)(h/M_P)^2∂_μ A^μ .
  
  The vector current is conserved, ∂^μ V_μ = 0, whereas the divergence of the axial current is proportional to the fermion mass M_N. As a result, in the presence of the additional terms non-minimally coupled to the Higgs field, a new channel h + h → N + N̅ contributes to the production of singlet fermions. The amount of DM produced via this channel can be calculated analogously to the four-fermion case. The abundance of N generated via the Higgs annihilation is given by
  
      Ω_NΩ_DM≃    3.2  (βξ + (γ̅+β)(γ̅ξ_η+ξ_γ-ξ)/γ̅^2+1)^2
         × ( T_ prodM_P)^3 ( M_NT_ prod)^2 ( M_N10  keV)  ,
  
  where we specialized to a Dirac fermion N; the Majorana case would lead to an additional factor of 1/2. 
  
  We see that, as compared to the four-fermion channel <ref>, the DM production due to the scalar-fermion interaction is suppressed by a factor of (M_N/T_ prod)^2. This means that for small masses of the DM particle, the four-fermion channel provides the dominant contribution to the DM abundance. Moreover, it is possible to choose the non-minimal couplings such that the scalar-fermion channel vanishes regardless the mass of N. For example, this is the case for
  
      1/γ̅=ξ_η=β=0  ,
  
  or when the coupling β is given by the following combination of other parameters:
  
      β=-ξ_ηγ̅^2+ξγ̅-ξ_γγ̅/ξγ̅^2+ξ_ηγ̅+ξ_γ .
  
  In a more general situation, the observed amount of DM can be generated mainly in the scalar-fermion channel in certain regions of the parameter space. As an example, let us take γ̅ =1, ξ≃ 10^7, ξ_γ=ξ_η=α=0 corresponding to a particular realization of Holst inflation <cit.>. If we then take β = √(ξ), the contribution of the process h + h → N + N̅ is suppressed compared to the four-fermion one by several orders of magnitudes. If, on the other hand, we take β = 0, the scalar-fermion channel dominates over the four-fermion one, and the observed amount of DM is generated for m_N=3.6 × 10^7 GeV.
  
  We note that the value of Ω_N depends on T_ prod, which in turn is sensitive to the details of preheating. These may generically depend on the non-minimal couplings. Nevertheless, the main conclusions about the possibility of gravitational production of DM via the four-fermion or scalar-fermion interactions and about the relative strength of different production channels are independent of the details of preheating.
  
  
      
  Appendix B.—Here we provide some details about the calculation of the DM abundance.
  We start from computing the amplitude entering <ref>.
  It is convenient to split the production channels as
  
      X_L(p_1) + X̅_L(p_2)    → N(q) + N̅(p_3) ,
      
      X_R(p_1) + X̅_R(p_2)    → N(q) + N̅(p_3) .
  
  As mentioned in the main text, we consider both Majorana and Dirac masses for N. In the Majorana case, we assume that N is a right-handed fermion, and then the amplitudes corresponding to <ref> and <ref> read
  
      |ℳ_X_L|^2    = 64/M_P^4[ C_A A-C_V V]^2( p_1 · p_3 )  ( p_2 · q ) ,
      
      |ℳ_X_R|^2    = 64/M_P^4[ C_A A+C_V A + C_V V]^2( p_1 · q )  ( p_2 · p_3 ) ,
  
  where by C_VV, C_VA and C_AA we denote the coefficients of the vector-vector, vector-axial and axial-axial currents, respectively. In the Dirac case, the amplitudes are
  
      |ℳ_X_L|^2    = 64/M_P^4{[ C_A A-C_V V]^2( p_1 · p_3 )  ( p_2 · q ) . 
         +. [ C_A A-C_V A + C_V V]^2( p_1 · q )  ( p_2 · p_3 ) } , 
      
      
      |ℳ_X_R|^2    = 64/M_P^4{[ C_A A-C_V V]^2( p_1 · p_3 )  ( p_2 · q ) .
         +. [ C_A A+C_V A + C_V V]^2( p_1 · q )  ( p_2 · p_3 ) } .
  
  If there is no asymmetry between particles (X) and anti-particles (X̅), the integral in <ref> is symmetric under the replacement p_1 ↔ p_2. As a consequence, both structures yield the same function r(y) (cf. <ref>).
  The dimensionless function is defined simply as r(y) = r(y T, T) / T^5, where y = E/T and 
  
      r(E,T) =   1/2 E∑_X ∫3 p⃗_1/(2π)^3   2 E_13 p⃗_2/(2π)^3   2 E_23 p⃗_3/(2π)^3   2 E_3
         × (2π)^4δ^(4)( p_1+p_2-q-p_3 )  ( p_1 · q )  ( p_2 · p_3 ) 
         × f_X( p_1) f_X̅(p_2)  .
  
  The integration in <ref> can be performed numerically following the recipe in <cit.>.
  We assume that all SM fermions thermalize rapidly after inflation and, hence, we take the Fermi-Dirac distribution f_X = f_X̅ = (exp(y)+1)^-1.
  Note that the computation of the exact fermion spectrum right after inflation is a highly non-trivial task, and we do not attempt to address it here.
  The result of the numerical integration is shown in <ref>. In the same figure we also show the analytic approximation given in <ref>.
  
  
  We can rewrite <ref> in the form
  
      T  H  ∂ f_N/∂ T=-C_f/M_P l^4T^5  r(y) ,
  
  where the coefficient C_f (f = M,  D) is given by <ref>. In the more general case of non-vanishing 1/γ̅ the coefficient can be read from <ref> and M_MLM_DR. Now <ref> can be easily integrated leading to <ref>.
  
  In order to compute the DM abundance, it is convenient to introduce the variable Y_N=n_N/s which remains constant in the expanding Universe. Here s corresponds to the entropy density and n_N=2· 1/(2π)^3 ∫^3 q⃗ f_N is the number density. In the latter, the factor 2 accounts for N̅. Then we get
  
      Y_N = 1/2π^2/45 g_eff T^3  2· T^3 ∫_0^∞ y  y^2/2 π^2 f_N(y) .
  
  The abundance of N is given by
  
      Ω_N≡M_N n_N/ρ_cr = M_N Y_N/ρ_cr/s ,
  
  where ρ_cr is the critical energy density of the Universe. Plugging in the numbers, we arrive at <ref>.
  
  
  
  Appendix C.—As was mentioned in the main text, the lightest active neutrino is practically massless in the νMSM regardless the mass of the DM sterile neutrino N_1.[The statement about keV-scale DM has been made in <cit.>.] Indeed, in the seesaw mechanism the contribution of N_1 to the mass matrix m_ν of the active neutrinos can be estimated as
  
      δ m_ν∼ M_N_1θ^2,   θ^2≡∑_α |θ_α|^2,
  
  where θ_α is the mixing between N_1 and the neutrino of flavor α = e, μ, τ.
  The mixing θ^2 implies that N_1 is not absolutely stable and processes like N_1 →ννν, N_1 →νγ are allowed. The stability of DM and the non-observation of a monochromatic signal from the decay N_1→νγ constrain θ^2 and, by virtue of <ref>, the contribution of N_1 to the active neutrino mass matrix. 
  The other two right-handed neutrinos N_2, 3 of the νMSM can provide two contributions to m_ν matching the two observed mass differences (the solar and the atmospheric ones). Therefore, the contribution of N_1 can be identified with the mass of the lightest active neutrino.
  
  In the mass range 1-50 keV where the resonant production is effective, the combination θ^2 has been bounded from above by the X-ray observations <cit.> (for a recent review, see <cit.>). One can expect that for larger masses of N_1, the X-ray (and gamma) constraints also lead to the strongest bounds. However, for our purposes the requirement of the DM stability is already enough. Indeed, the width of the three-body decay of N_1 reads <cit.>
  
      Γ_N_1 → 3 ν=G_F^2 M_N_1^5/96 π^3θ^2,
  
  where G_F is the Fermi constant. Requiring that 1/Γ_N_1 → 3 ν exceeds the age of the Universe, we arrive at an upper bound on θ^2. Using <ref> we get 
  
      δ m_ν≲ 3.3  eV ×( 10 keV/M_N_1)^4.
  
  We see that already for M_N_1=100 keV, the requirement of the stability of N_1 yields δ m_ν being ∼ 25 times less than the solar mass difference √(Δ m_sol)≃ 8.6 × 10^-3 eV. It is also clear that this bound rapidly gets stronger when M_N_1 increases.  
  
  
  
  Appendix D.—The possibility of successful leptogenesis with two Majorana fermions N_2,3 has first been shown in Ref. <cit.>. Inspecting the perturbative solution of Ref. <cit.> one can see that the produced Baryon Asymmetry of the Universe (BAU)
  is proportional to ( n_2,3-n_eq)/s, where n_2,3 is the number density of  Majorana fermions N_2,3 and n_eq is the equilibrium number density of relativistic fermions, and s is the entropy density. The fermions N_2,3 are produced by the same processes as N_1. Since the masses can be safely neglected, the yields Y_X=n_X/s will be the same for all right-handed neutrinos. Using <ref> and assuming that N_1 forms all DM, we can find that
  
      Y_N_2,3 =  Y_N_1≃ 4.4× 10^-510 keV/M_N_1 .
  
  This value should be compared with the equilibrium yield Y_eq=n_eq/s≃ 1.9× 10^-3. As one can see, for a 10 keV N_1, the baryogenesis results are affected at a percent level. 
  We have also verified numerically (using the kinetic equations of ref. <cit.>), that the same conclusion holds in the regions of the parameter space where the perturbative solutions is not applicable. 
  
  
  
  
  
  