Dark matter relic abundance and light sterile neutrinos

In this paper, we calculate the relic abundance of the dark matter particles when they can annihilate into sterile neutrinos with the mass ≲ 100 GeV in a simple model. Unlike the usual standard calculations, the sterile neutrino may fall out of the thermal equilibrium with the thermal bath before the dark matter freezes out. In such a case, if the Yukawa coupling yN between the Higgs and the sterile neutrino is small, this process gives rise to a larger ΩDMh2 so we need a larger coupling between the dark matter and the sterile neutrino for a correct relic abundance.


Introduction
The weakly interacting massive particles (WIMPs) are considered as the candidates of the dark matter (For a review, see ref. [1]). In this scenario, the dark matter particles are produced in the thermal bath of the early universe, then freeze out from the plasma as the temperature drops. It is well-known that the observed dark matter's relic abundance requires its thermally averaged annihilation cross section σv = 2-3 × 10 −26 cm 3 /s at the freezing-out temperature T ∼ mχ 20 , which is roughly the typical cross section of the weak interaction. This coincidence is called the "WIMP miracle".
Calculations of the relic abundance of the dark matter involve the Boltzmann equation (For derivation, see ref. [2]). In the case of the WIMP dark matter, some hypotheses are adopted in order to simplify the equation. One important hypothesis is that the annihilation products of the dark matter fall in thermal equilibrium with the thermal bath rapidly. This is true when the dark matter mainly annihilates into the standard model (SM) particles. However, in many new physics models, the dark matter might mainly annihilate into other beyond-SM particles. In this case, whether this hypothesis is valid needs to be carefully examined.
In the Type I see-saw model [3][4][5][6][7], the right-handed neutrinos (N ) couple with the left-handed neutrinos l ±,0 through the Higgs fields H. After the Higgs field acquires a vacuum expectation value (VEV), the majorana mass terms of the left-handed neutrino arise through the Type I See-saw Mechanisms. In the early universe, if there is no extra sector, the main processes that can generate the right-handed neutrinos are the decay and the inverse decay of the right-handed neutrinos and the Higgs bosons (For an example of calculations, see [8]). In the simplest Type I See-saw Mechanisms, if the mass of the righthanded neutrino is approximately 100 GeV, the Yukawa coupling constants of the N -l ±,0 -H couplings y N should be smaller than ∼ 10 −6 for the correct left-handed neutrino masses. However, it is too small for the right-handed neutrinos to reach in thermal equilibrium with JHEP01(2017)025 the thermal bath. Although there are some models [9][10][11][12][13][14][15], e.g., the inverse see-saw model, or the linear see-saw model, that can result in a larger y N ∼ 0.01, as the temperature drops, the thermal-averaged production rates of the sterile neutrinos Γ P ∝ e − m N T drops rapidly and then the sterile neutrinos decay out of equilibrium.
Ref. [16][17][18] calculate a general case of the secluded dark matter, in which the annihilation products fall out of thermal equilibrium from the thermal bath. Specifically, in the literature, there are some models that the dark matter can annihilate into light sterile neutrinos [19][20][21][22][23][24][25][26][27], unlike the models that the sterile neutrino itself plays the role of the dark matter [28][29][30]. However, as we have mentioned before, the sterile neutrinos might not be in thermal equilibrium with the SM particles in the early universe. Thus, the traditional calculations of the relic abundance might be unreliable and the standard Boltzman equation(s) should be modified. In this paper, in order to calculate these non-thermal effects, we rely on a simple model based on ref. [24]. We focus on the case that the masses of the dark matter and the right-handed neutrinos are less than ∼ 100 GeV and the dark matter particles only annihilate to right-handed neutrinos. In this case, the sterile neutrino mainly decays through the three-body final state channels. As we have calculated in ref. [31], this scenario can perfectly explain the gamma-ray excess from near the galactic center. We will also show that the nonthermal effects of the right-handed neutrinos can significantly modify the relic abundance when the Yukawa coupling constant y N lH becomes quite small.

Model descriptions
The model discussed in this paper contains a majorana fermion χ and a real-scalar boson φ. Both these two fields are SM-singlets and are odd under a dark Z 2 discrete symmetry. The sterile neutrino together with the SM-fields are all even under this Z 2 symmetry. In this paper, we discuss two cases. In one case there is one majorana right-handed neutrino N , and in the other case there are a pair of pseudo dirac sterile neutrino Weyl-fields N 1,2 .
In the majorana right-handed neutrino case, the general Lagrangian is given by where σ 2 is the second Pauli-matrix. In eq. (2.1), all the y χ,χ5,N i , λ φ,φH and m χ,φ,N are real numbers. In fact, y χ and y χ5 are respectively the real part and the imaginary part of a single coupling constant JHEP01(2017)025 (y χ + iy χ5 )χ w · N w φ + h.c., where χ w and N w are the Weyl components of the χ and N fields. We should note that all the complex phases in the m χ,N and y N i can be rotated away by redefining the fields and the y χ + iy χ5 . However, for simplicity, in the numerical calculations we just omit the y χ5 and set it to be zero.
In the pseudo-Dirac sterile neutrino case, the general Lagrangian is given by is a Dirac four-spinor composed of two Majorana fields N 1 and N 2 , m D is the Dirac mass term of the sterile neutrinos, and µ 1,2 are the Majorana mass terms, y χD , y χD5 , y N i , y N Ci are the Yukawa coupling constants which is different from the eq. (2.1). This time, we can redefine all the fields in order for the m χ , m N D and y N i to be real numbers, while the y χD , y χD5 , y N Ci , µ 1 , µ 2 might be complex numbers and in general their phases cannot be rotated away. In this paper, for simplicity, we only discuss the case that all these numbers are real and y χD5 = 0.
In general, µ 1,2 and y N Ci terms violate the lepton number and cause the mass splitting of the two components of the N D . If y N Ci = 0 and µ 1,2 m N D , this is called the "inverse see-saw", and if µ 1,2 = 0 and y N Ci = 0, this is called the "linear see-saw". Thus, the masses of the left-handed neutrinos are mainly decided by the strength of the lepton number violating terms and y N i can be much larger than the simplest see-saw models. Usually the smallness of the left-handed neutrino masses requires the smallness of the lepton number violating terms µ 1,2 and y N Ci . Therefore, in the early universe, the effects of the µ 1,2 and y N Ci terms are negligible. Therefore, these parameters are set to zero in our study and the N D is regarded as a Dirac fermion during the calculation processes.

Calculations of the dark matter's relic abundance
The calculations of the relic abundance of the dark matter are based on the Boltzmann equations (We derive the equations in this paper according to ref. [8,32]). In principle, a full solution to the Boltzmann equations should involve the evolutions to the distribution functions of the particles. However, we assume that the elastic scatterings are fast enough as the particles in the thermal bath can maintain kenitic equilibrium as usual. For simplicity, we only consider the case that m φ > m χ + m N (D) . In the Majorana right-handed neutrino case, the Boltzmann equations are given by where the Y A = n A s is the actual number of the constituent A per-comoving-volume, and the Y Aeq = n Aeq s is the equilibrium number of the constituent A per-comoving-volume, n A(eq) is the (equilibrium) number density of the constituent A, s is the entropy density, z = mχ T , and T is the temperature, H is the Hubble constant. σv AB→CD is the thermally averaged cross section times velocity 2) where δ CD = 1(0) if C and D are identical(different) particles, g A and g B are the degrees of freedoms of particle A and B, K 1 is a Bessel function, σ AB→CD (s ) is the cross section of the process AB → CD with the total energy in the center of mass frame is √ s . The definition of theΓ φ→χN is given bȳ where m t φ is the thermal mass of φ which will be defined later. TheΓ N is a little bit complicated. We need to consider the decay/inverse-decay processes N → h ±0 l or h ±0 → lN . However, as the temperature drops below the electroweak symmetry breaking (EWSB) critical temperature T c , we should considers the processes N ↔ W ± /Z/h. In this paper, we adopt the approximation method described in ref. [33] to calculate the N → h ±0 l or h ±0 → lN with all the four states of the Higgs doublets having the Higgs boson mass m h (T ) below T c . If m N > m h (T ), at the zero temperature, we set in order to let the right-handed neutrino decay. We calculate the thermal masses in the following procedures. The thermal effects on the fermions are neglected. As for the scalar bosons, the effective potential in a finite temperature T is [34][35][36][37] where λ is the self-interacting coupling constant of the Higgs boson, µ 2 < 0 is the mass term at zero temperature of the Higgs potential. The definition of the c and c φ is given by The critical temperature T c of the EWSB is Then the temperature dependent masses of the Higgs boson and the scalar φ are given by where As for the pseudo-Dirac sterile neutrino case, all the N 's in the above formulas should be replaced with N D andN D . One needs to note that the condition of YN D = Y N D always holds and a summation over the particle and anti-particle should be considered. That is to say, in (3.1), σv AA→N N should be replaced with σv AA→N D N D + σv AA→N D N D + σv AA→N DND ,Γ φ→χN should be replaced withΓ φ→χN D +Γ φ→χN D , σv N N →AA should be replaced with σv N D N D →AA + σv N DND →AA . In (3.8), the 2y 2 χ + 2y 2 χ5 should also be replaced with 4y 2 χD + 4y

Numerical solutions and results
To solve the differential equations (3.1), we use the ready-made function [38,39] embeded in the micrOMEGAs [40] for computing the stiff equations eq. (3.1). We also use the CalcHEP [41] embeded in the micrOMEGAs to calculate the σv (s) and the widths of the particles. The model file is implemented and output by FeynRules [42]. We adopt the g * and g * S implemented in the micrOMEGAs in order to calculate the Hubble constant M pl , and s = 2π 2 45 g * S T 3 . Here M pl = 1.22 × 10 19 GeV is the Planck energy. If the λ φh is too small, and when T m φ , φ becomes in thermal equilibrium with the SM particles through the φφ ↔ HH interactions. The χ and N (D) then fall into thermal equilibrium through the φ-portal. However, once the temperature drops below the mass of the φ, the number density of φ rapidly becomes so small that the χ and N (D) decouple from the thermal bath altogether. Finally, χ, N (D) decouple from each other.
During the calculations, we simplify the eq. (3.1) by eliminating all the terms involving Y φ once Y φ < 0.01Y χ . In order to present our result, we fix m φ = 180 GeV, λ φ = 0.5 and λ φH = 0.45, y N 2 = y N 3 = 0. We plot our results on the m N D -y χ plane in the different combinations of the values of y N = y N 1 = 10 (−7),(−6),(−5),(−4),(−3),(−2) (most of the values are far beyond the current collider bounds, for the related discussions, see ref. [43][44][45][46][47][48][49][50][51][52][53]), and m N ( D) = 25, 52, 76 GeV. For each m N (D) , we find one y χ that results in 0.117 < Ω DM h 2 < 0.120 [54]. For comparison, we also present the result using the traditional Boltzmann equation The results are shown in figure 1, 2 and 3. We should note that in these figures, some lines are left over since they are nearly identical with the drawn lines. When y N 10 −3 , the deviation from the standard calculation by the old Boltzmann equation is quite small. Note that when y N   in this case. Therefore, some lines may disappear before the most right-handed side in figure 1-3.
To study the decoupling processes in details, we plot the z-evolution of In the region where z > 45, although Y N Y N eq still arises as the z accumulates, y N actually drops as y N eq decreases much faster. When χ and N (D) decouples with each other while y N is small, Y N (D) is usually larger than Y N (D) eq , which gives rise to the Y χ at the freeze-out point. That is the reason we need a larger y χ to suppress the Y χ∞ ∝ Ω DM h 2 in the case of small y N 10 −3 as shown in figure 1-3.

Discussions
In the case that the fermionic χ to be the dark matter candidate, there is no tree-level diagrams contributing to the direct detection processes. However, As has been mentioned in ref. [23,24], one-loop diagram will result in χ(I, or iγ 5 )χH † H operators, which give rise to not only the direct detection processes through exchanging a Higgs boson with the target nucleon, but also lead to the Higgs invisible decays. These effects are all suppressed by the loop factor and are proportional to λ 2 φH y 4 χ . The ref. [24] calculated these constraints in a similar model, and the result was shown in its figure 2. From the left panel we can see that we do not need to worry about these contraints once λ φH 1. Although in this paper, we need a larger y 4 χ than the usual standard calculations in order to get an appropritate dark matter relic abundance, we can at least escape the constraints by assigning a smaller λ φH accordingly. Such an assignment usually does not affect the relic abundance of the dark matter significantly, because the main annihilation channels do not involve λ φH .
Ref. [31] has calculated the galactic center gamma-ray excess in such kind of scenario. In ref. [31] we have pointed out that an approximately 10-60 GeV sterile neutrino together with a heavier dark matter particle can perfectly explain the the observed spectrum. The annihilation cross section σv is within the range 0.5-4 × 10 −26 cm 3 /s. Considering the uncertainties of the parameters of the dark matter profile, e.g., the local dark matter density ρ which varies from 0.2-0.6 GeV/cm 3 , σv real ∝ ρ −2 can differ by one order of magnitude. In this paper, we need a larger y 4 χ , which amplifies the σv real ∝ y 4 χ from the standard WIMP cross section ∼ 3 × 10 −26 cm 3 /s in order for a correct dark matter relic abundance. Such an amplification factor is typically ∼ 2-100 depending on the m χ,N when y χ 10 −4 , leaving us enough room to adjust the dark matter profile parameters to fit the gamma-ray data from the galactic center.

JHEP01(2017)025
Finally, we should note that even in the extreme case that the y N is as small as 10 −7 which lead to a long-life right-handed neutrino, the Big-Bang neucleosynthesis (BBN) is not affected (For a review of the BBN, see the section 24 in the ref. [55], and for references therein). For example, the width Γ N ∼ 10 −17 GeV when m N ∼ 50 GeV, but this is still much larger compared with the Hubble constant H ∼ 10 −22 GeV at the BBN temperature T ∼ 10 MeV, so nearly all of the out-of-equilibrium sterile neutrino decay before they may have an impact on the BBN.

Conclusions
We have calculated the relic abundance of the dark matter particles when they annihilate into sterile neutrinos with the mass m N < m χ 100 GeV. In the model we have relied on, the sterile neutrino will become in thermal equilibrium with the thermal bath when T m N (D) , however it will decouple from the thermal bath before the dark matter freezes out if y N is small. This gives rise to a larger Ω DM h 2 so we need a larger coupling between dark matter and the sterile neutrino for a correct relic abundance. In the future, we will continue to dedicate ourselves in some more detailed research in such kind of scenarios.