Radiative Production of Non-thermal Dark Matter

TL;DR

This paper investigates dark matter production mechanisms in the early universe, highlighting that Standard Model loop processes can generate the correct relic abundance even without direct inflaton coupling.

Contribution

It demonstrates that Standard Model loop processes can produce the correct dark matter relic abundance without direct inflaton coupling, expanding understanding of non-thermal dark matter origins.

Findings

Standard Model loops can produce sufficient dark matter relic abundance.

Non-thermal production via inflaton decay is comparable to thermal processes.

Loop processes may dominate dark matter production even without direct coupling.

Abstract

We compare dark matter production from the thermal bath in the early universe with its direct production through the decay of the inflaton. We show that even if dark matter does not possess a direct coupling with the inflaton, Standard Model loop processes may be sufficient to generate the correct relic abundance.

Figures

Equations

$$\frac{d{n}_{\chi }}{dt}+3H{n}_{\chi }=R(T)$$

$$R(T)=\int f_1{f}_2\frac{E_1{E}_{2}d{E}_{1}d{E}_{2}d\cos {\theta }_{12}}{1024{\pi }^6}\int |{\mathcal{M}}_{fi}{|}^{2}d{\Omega }_{13},$$

$$\frac{d{\rho }_{\varphi }}{dt}+3H{\rho }_{\varphi }=-{\Gamma }_{\varphi }{\rho }_{\varphi }$$

$$\frac{d{\rho }_R}{dt}+4H{\rho }_R={\Gamma }_{\varphi }{\rho }_{\varphi }.$$

$$H^2=\frac{{\rho }_{\varphi }}{3M_P^2}+\frac{{\rho }_R}{3M_P^2},$$

$$T\propto a^{-3/8},$$

$$Y=\frac{n}{T^8},$$

$$\frac{dY}{dT}=-\frac{8}{3}\frac{R(T)}{H{T}^9}.$$

$$H(T)\approx \frac{5}{6}\frac{\alpha }{{\Gamma }_{\varphi }{M}_P^2}{T}^4.$$

$$\frac{dY}{dT}=-\frac{24c}{5\sqrt{3\alpha }}\frac{R(T)M_P{T}_{RH}^2}{T^{13}}$$

$$\frac{3}{2}cH(T_{\mathrm{RH}})=\frac{3}{2}c\sqrt{\frac{\alpha }{3}}\frac{T_{\mathrm{RH}}^2}{M_P}={\Gamma }_{\varphi }$$

$$R(T)=R_{annihilation}(T)+R_{decay}(T).$$

$$R_{annihilation}(T)=\frac{T^{n+6}}{{\Lambda }^{n+2}},$$

$$\mathrm{n}<6:n(T_{RH})=\frac{24c}{5}\frac{T_{RH}^{n+4}{M}_P}{(6-n)\sqrt{3\alpha }{\Lambda }^{n+2}}$$

$$\mathrm{n}=6:n(T_{RH})=\frac{24c}{5\sqrt{3\alpha }}\frac{M_P{T}_{RH}^{10}}{{\Lambda }^8}\ln \left(\frac{T_{\max }}{T_{RH}}\right)$$

$$\mathrm{n}>6:n(T_{RH})=\frac{24c}{5}\frac{T_{\max }^{n-6}{M}_P{T}_{RH}^{10}}{(n-6)\sqrt{3\alpha }{\Lambda }^{n+2}}$$

$$\Omega =\frac{\rho }{{\rho }_c}\approx \frac{n(T_0)\times m_{\chi }}{10^{-5}{h}^2\mathrm{GeVc}{\mathrm{m}}^{-3}}$$

$$\Rightarrow \Omega h^2\approx 10^5\frac{n(T_{RH})}{\mathrm{c}{\mathrm{m}}^{-3}}\left(\frac{g(T_0)}{g(T_{RH})}\right){\left(\frac{T_0}{T_{RH}}\right)}^3\left(\frac{m_{\chi }}{1\mathrm{GeV}}\right)$$

$$\simeq 5.9\times 10^6\left[\frac{n(T_{RH})}{T_{RH}^3}\right]\left(\frac{m_{\chi }}{1\mathrm{GeV}}\right)$$

$$R_{decay}(T)={\Gamma }_{\varphi }\frac{{\rho }_{\varphi }}{M_{\varphi }}{B}_R=\frac{25}{12}\frac{{\alpha }^2{T}^8}{{\Gamma }_{\varphi }{M}_P^2{M}_{\varphi }}{B}_R.$$

$$\frac{dY}{dT}=-\frac{20}{3}\frac{\alpha }{M_{\varphi }{T}^5}{B}_R,$$

$$Y(T_{RH})\simeq \frac{5}{3}\frac{\alpha }{M_{\varphi }{T}_{RH}^4}{B}_R=\frac{g(T_{RH}){\pi }^2}{18M_{\varphi }{T}_{RH}^4}{B}_R$$

$$\Rightarrow n(T_{RH})=\frac{g(T_{RH}){\pi }^2}{18M_{\varphi }}{T}_{RH}^4{B}_R.$$

$$\Omega h^2{|}_{\mathrm{decay}}=\left(\frac{B_R}{2.9\times 10^{-9}}\right)\left(\frac{T_{RH}}{M_{\varphi }}\right)\left(\frac{m}{1\mathrm{GeV}}\right)$$

$$\simeq 0.1\left(\frac{B_R}{8.7\times 10^{-7}}\right)\left(\frac{3\times 10^{13}\mathrm{GeV}}{M_{\varphi }}\right)\left(\frac{T_{RH}}{10^{10}\mathrm{GeV}}\right)\left(\frac{m}{1\mathrm{GeV}}\right)$$

$$\Omega h^2{|}_{\mathrm{tot}}\simeq \Omega h^2{|}_{annihilation}+\Omega h^2{|}_{\mathrm{decay}}$$

$$\Omega h^2{|}_{n<6}\simeq 0.1\left[\frac{2.8\times 10^{8}c}{(6-n)}\frac{T_{RH}^{n+1}{M}_P}{\sqrt{3\alpha }{\Lambda }^{n+2}}+\frac{B_R}{2.9\times 10^{-10}}\times \frac{T_{RH}}{M_{\varphi }}\right]\left(\frac{m}{1\mathrm{GeV}}\right)$$

$$\Omega h^2{|}_{n=6}\simeq 0.1\left[2.8\times 10^{8}c\frac{T_{RH}^7{M}_P}{\sqrt{3\alpha }{\Lambda }^8}\ln \left(\frac{T_{max}}{T_{RH}}\right)+\frac{B_R}{2.9\times 10^{-10}}\times \frac{T_{RH}}{M_{\varphi }}\right]\left(\frac{m}{1\mathrm{GeV}}\right)$$

$$\Omega h^2{|}_{n>6}\simeq 0.1\left[\frac{2.8\times 10^{8}c}{(n-6)}\frac{T_{RH}^7{T}_{max}^{n-6}{M}_P}{\sqrt{3\alpha }{\Lambda }^{n+2}}+\frac{B_R}{2.9\times 10^{-10}}\times \frac{T_{RH}}{M_{\varphi }}\right]\left(\frac{m}{1\mathrm{GeV}}\right)$$

$$\Omega h^2{|}_{n=6}^{tot}\simeq 0.1\left[{\left(\frac{T_{\mathrm{RH}}}{10^{11}\mathrm{GeV}}\right)}^7{\left(\frac{7.3\times 10^{12}\mathrm{GeV}}{\Lambda }\right)}^8\ln \left(\frac{T_{\max }}{T_{\mathrm{RH}}}\right)+\left(\frac{B_R}{8.7\times 10^{-8}}\right)\left(\frac{3\times 10^{13}\mathrm{GeV}}{M_{\varphi }}\right)\left(\frac{T_{\mathrm{RH}}}{10^{11}\mathrm{GeV}}\right)\right]\left(\frac{m}{1\mathrm{GeV}}\right).$$

$$\mathcal{L}=y_f\varphi \bar{f}f+y_{\chi }\varphi \bar{\chi }\chi +\frac{1}{{\Lambda }^2}\bar{f}f\bar{\chi }\chi$$

$${\Gamma }_{\varphi }^{loop}\simeq N_f\frac{y_f{y}_{\chi }}{16{\pi }^3}\frac{M_{\varphi }^3}{{\Lambda }^2}+N_f\frac{y_f^2}{128{\pi }^5}\left(1+\frac{{\pi }^2}{4}\right)\frac{M_{\varphi }^5}{{\Lambda }^4}$$

$${\Gamma }_{\varphi }^f=\frac{N_f{y}_f^2}{8\pi }{M}_{\varphi };{\Gamma }_{\varphi }^{\chi }=\frac{y_{\chi }^2}{8\pi }{M}_{\varphi }$$

$$B_R$$

$$W=\sqrt{3}{m}_{3/2}(z+\nu )$$

Full text

††preprint: LPT–Orsay 19-01††preprint: UMN–TH–3812/19, FTPI–MINN–19/03

I I. Introduction

Despite indirect but clear evidence planck of the presence of large amount of dark matter in our Universe, its nature still remains elusive. The absence of any signal in direct detection experiments XENON XENON , LUX LUX and PandaXPANDAX question the weakly coupled dark matter paradigm. Simple extensions of the Standard Model such as the Higgs-portal hp ; Higgsportal , $Z$-portal Zportal , or even $Z^{\prime}$-portal Zpportal as well as more complex extensions such the minimal supersymmetric standard model mssm ; Go1983 ; ehnos have large part of their parameter space (if not all) excluded when combining direct, indirect and accelerator searches (for a review on WIMP searches and models, see Arcadi:2017kky ). As a consequence, it is useful to look for different scenarios, including those with ultra-weak couplings fimp (see Bernal:2017kxu for a review), or the possibility that dark matter production occurred at very early stages of reheating after inflation as in the case of gravitino production nos ; ehnos ; kl .

In this context, it has been confirmed that dark matter production is naturally feasible in a wide variety of models such as SO(10) grand unification Mambrini:2013iaa , anomaly free U(1)’ models Bhattacharyya:2018evo , spin-2 portal Bernal:2018qlk high scale supergravity Benakli:2017whb ; grav2 ; grav3 ; gravitino ; highsc or moduli portal Chowdhury:2018tzw . In all these models, it has been shown that effects of non-instantaneous reheating Giudice:2000ex ; Garcia:2017tuj and non-instantaneous thermalization Garcia:2018wtq are non-negligible.

As in the case of gravitino production during inflationary reheating, in many models, dark matter is produced from the annihilation of thermal Standard Model particles, and one neglects the direct production of dark matter $\chi$ from the decay of the inflaton, $\phi$. It was shown egnop ; grav2 ; grav3 that the branching ratio $BR(\phi\rightarrow\chi\chi)$, is constrained to be very small (below $10^{-8}$ for a 100 GeV dark matter and a reheating temperature $T_{RH}\gtrsim 10^{10}$ GeV). In this paper, we compare the dark matter production rate from both the thermal bath and direct decay. We show that if dark matter is produced from the thermal bath, independent of the tree level branching ratio of inflaton decay to dark matter, one cannot ignore the radiative decay of the inflaton into dark matter. We further show that the radiative decay may well dominate the production rate thus providing the main source for relic abundance in the Universe. As a particular example, we also consider the radiative contribution for gravitino production in models of high scale supersymmetry.

The paper is organized as follows. We first compute the dark matter abundance is section II, then apply it to a generic microscopic model in section III before looking to consequences in supersymmetric scenarios in section IV. We conclude in section V.

II II. Dark matter production

II.1 Generalities

We presume that the dark matter is not produced in thermal equilibrium during reheating and as a consequence, the dark matter abundance, $n_{\chi}$, in the early universe is much lower than the thermal density $n_{\gamma}$. If we further assume that the dark matter annihilation process, $n_{\chi}~{}n_{\chi}\rightarrow n_{\gamma}~{}n_{\gamma}$, is also out of equilibrium, since $n_{\chi}^{2}\ll n_{\gamma}^{2}$, we can write the Boltzmann equation as

[TABLE]

with

[TABLE]

for a process 1 2 $\rightarrow$ 3 4 with 1 and 2 representing standard model particles in the thermal bath and 3, 4 representing the dark matter candidate, with $f_{1}$ and $f_{2}$ being the thermal distribution functions of the incoming particles 1 and 2 and $d\Omega_{13}$ is the solid angle between the particle 1 and 3 in the center of mass frame.

To compute the relic abundance, the strategy is straightforward. To solve Eq.(1), one needs to know the relation between $T$ and the cosmological scale parameter, $a$. We can then use the Hubble parameter $H=\frac{\dot{a}}{a}$ to express $t$ as function of $T$ before integrating Eq. (1). The dependence of $T$ on $a$ is obvious in a pure radiation dominated or a pure matter dominated era ($T\propto a^{-1}$) due to entropy conservation. However, during reheating, the temperature grows from effectively 0 (when the Universe is still dominated by inflaton oscillations) to a maximum temperature, $T_{\rm max}$, Garcia:2017tuj ; Giudice:2000ex . Subsequently, the temperature decreases to $T_{\rm RH}$ when the inflaton is (nearly) fully decayed. The maximum temperature is determined by the inflaton decay width, $\Gamma_{\phi}$ and is approximated by $T_{\rm max}\simeq 0.5(M_{\phi}/\Gamma_{\phi})^{1/4}T_{\rm RH}$, where $M_{\phi}$ is the inflaton mass. During this period, the density of Universe is a mixture between the inflaton and the growing radiation density. It has been shown in Garcia:2017tuj ; Giudice:2000ex that solving the set of conservation of energy conditions for the inflaton density $\rho_{\phi}$, the radiation density $\rho_{R}$ combined with the Friedmann equation :

[TABLE]

where $\Gamma_{\phi}$ is the width of the inflaton111A more precise computation should distinguish $\Gamma_{\phi}$ into the part contributing to the thermal bath $\Gamma^{\gamma}_{\phi}$ and the part contributing to the dark sector $\Gamma^{\chi}_{\phi}$. However, we will always consider (and justify) $\Gamma_{\phi}^{\chi}\ll\Gamma^{\gamma}_{\phi}$ throughout our study., and $M_{P}$ the reduced Planck mass222$M_{P}=\frac{1}{\sqrt{8\pi}}M_{Pl}=2.4\times 10^{18}$ GeV. one obtains

[TABLE]

as $T$ decreases from $T_{\rm max}$ to $T_{\rm RH}$. Defining the comoving yield $Y$

[TABLE]

the Eq. (1) becomes

[TABLE]

If we further make the approximation that the inflaton dominates the total energy density between $T_{\rm max}$ and $T_{\rm RH}$, we can express the Hubble parameter in terms of $\rho_{\phi}$ which falls as $a^{-3}\sim T^{8}$ with a constant of proportionality depending on the inflaton decay rate,

[TABLE]

In Eq. (7), $\alpha=\frac{g(T)\pi^{2}}{30}$, where $g(T)$ counts the relativistic degrees of freedom of the thermal bath at the temperature $T$ ($g_{*}=106.75$ in the Standard Model) so that $\rho_{R}=\alpha T^{4}$. In this case, Eq. (6) becomes

[TABLE]

where we have taken

[TABLE]

to define the reheating temperature, where $c\approx 1.2$ is a constant obtained from a numerical integration ps2 ; egnop . The particular choice of a dark matter candidate and its interactions will determine $R(T)$ and allow for the integration of Eq. (8).

II.2 Dark matter production from the thermal bath

It has been shown in Mambrini:2013iaa that dark matter can be produced in the very early stages of the Universe, even if it is not directly coupled to the Standard Model, through the exchange of a massive mediator. Indeed, thermal gravitino production nos ; ehnos ; kl was an early example of this type of process. In Bhattacharyya:2018evo , Bernal:2018qlk gravitino ; grav2 ; grav3 ; highsc and Chowdhury:2018tzw this type of dark matter production has been extended to Chern-Simons type couplings, spin-2 mediators, supergravity, and moduli–portal scenarios respectively taking into account the effects of non instantaneous reheating we discussed above Garcia:2017tuj ; Garcia:2018wtq .

We can distinguish between the annihilation production processes and the decay rate by defining

[TABLE]

We parametrize the production rate

[TABLE]

where $\Lambda$ is some beyond the Standard Model mass scale333This parametrization corresponds to a cross section $\sigma\propto\frac{T^{n}}{\Lambda^{n+2}}$.. If one looks at specific models, $n=2$ could correspond to the exchange of a massive particle with mass $\Lambda>T_{\rm max}$ or two non-renormalizable mass suppressed couplings. $n=6$ appears in processes invoking two mass suppression couplings and the exchange of a massive particle, which is the case in high scale supergravity gravitino ; grav2 ; grav3 ; highsc or moduli–portal scenarios Chowdhury:2018tzw .

Inserting Eq. (11) in Eq. (8) we obtain, after integration, from $T_{\rm max}$ down to $T_{\rm RH}$, for the relic abundance at $T=T_{RH}$ 444These expressions agree (using $c=1$) with the more exact derivation in Garcia:2017tuj with the substitution of $\Lambda^{n+2}=(\pi/\lambda)(\pi^{2}/g_{\chi}\zeta(3))^{2}M^{n+2}$ where $g_{\chi}$ accounts for the internal number of degrees of freedom of the dark matter $\chi$. We note that the cross section was assumed to be $\sigma v=\lambda T^{n}/\pi M^{n+2}$ yielding a rate $R=n_{\chi}^{2}\sigma v$, with $n_{\chi}=g_{\chi}\zeta(3)T^{3}/\pi^{2}$. :

[TABLE]

from which we can deduce the present relic abundance at the temperature $T_{0}$

[TABLE]

with $g(T)$ accounting for the relativistic degrees of freedom, where we have considered only Standard Model degrees of freedom ($g(T_{RH})=106.75$ and $g(T_{0})=3.91$).

II.3 Dark matter production from inflaton decay

We can also compute the dark matter density produced directly from the decay of the inflaton. Indeed, in Eq. (12), we neglected any direct couplings of the dark matter to the inflaton. But we can easily compute the relic abundance obtained if we allow a branching fraction $B_{R}$ of the inflaton decay into dark matter. Here $B_{R}$ is defined as the number of dark matter quanta produced per inflaton decay.

We will assume non-instantaneous thermalization and solve the Boltzmann equation Eq. (8) defining the production rate (number of dark matter particles produced per unit of time and volume) as

[TABLE]

Plugging the rate Eq. (14) into Eq. (8), we obtain

[TABLE]

which gives after integration between $T_{\rm max}$ and $T_{\rm RH}$,

[TABLE]

Interestingly, this approximate result is remarkably close to the exact integration done in egnop where the coefficient of $\alpha T_{\rm RH}^{4}B_{R}/M_{\phi}$ in $n$ is $(8/3)(c/3)^{1/2}$ which is almost exactly 5/3 when $c=1.2$ 555Notice that if we considered instantaneous reheating, one would have $n(T_{RH})=B_{R}\times\frac{\rho_{\phi}}{M_{\phi}}=B_{R}\times\frac{g(T_{RH})\pi^{2}}{30M_{\phi}}T_{RH}^{4}$..

We can then combine Eq. (13) with Eq. (16) to obtain

[TABLE]

As we indicated earlier, we observe that one needs a very low branching ratio to avoid the overabundance of the dark matter. The total relic abundance is then obtained by adding the annihilation and decay contributions. Combining Eq. (12) and Eq. (17), we obtain

[TABLE]

or

[TABLE]

which gives, for $n=6$ for example

[TABLE]

We show in Fig. 1 the result on the scan of the parameter space (in $B_{R}$, $\Lambda$, and the dark matter mass, $m$) in the case $n=6$ using the relic density constraint $\Omega h^{2}\simeq 0.11$, for $T_{RH}=10^{11}$ GeV and $T_{max}=10^{13}$ GeV. For each value of $\Lambda,B_{R}$, the mass needed to obtain the correct relic density is color coded by the scale at the right of the figure. We notice that the range of dark matter masses allowed is very large (from the MeV to the EeV scale). This is a direct consequence of the strong power dependence of the relic abundance on the scale $\Lambda$. We also remark that a branching ratio of $\mathcal{O}(1)$ is possible for dark matter masses of order 100 keV, whereas EeV scale dark matter requires very tiny branching fraction of the order of $B_{R}\simeq 10^{-16}$ to avoid the over-closure of the Universe. Note also, that over most of the parameter space with large values of $\Lambda$, the thermal production from annihilations contributes negligibly to the relic density.

III III. Generic microscopic models

We next apply the generic analysis we did using an effective field theory formulation. We first consider the effective Lagrangian666Without loss of generalities, we will work with a fermionic dark matter candidate. The extension to a scalar candidate is straightforward. between the inflaton $\phi$, standard model fields, $f$, and the dark matter candidate, $\chi$:

[TABLE]

This potential can be viewed as an effective interaction between standard model fermions and the dark sector through the exchange of a massive field of mass $M$ $\simeq\Lambda\gg T_{max}$.

One of our key points is that there will be direct production of dark matter through inflaton decay, even if the dark matter does not couple to the inflaton at tree level. That is, even if $y_{\chi}=0$, dark matter will be produced radiative via the diagram shown in Fig. 2. Indeed, from Eq. (20) one can deduce

[TABLE]

via the loop of the $N_{f}$ families of the standard model fermions $f$, where the first term is the interference between the tree and loop diagrams, and we have taken the massless limit for $f$. The other contributions to the inflaton width $\Gamma_{\phi}$ are given by

[TABLE]

We can express $B_{R}$ as function of the microscopic parameters,

[TABLE]

where we assumed that $N_{f}y_{f}\gg y_{\chi}$ and $\Lambda\gg M_{\phi}$. These conditions are necessary to avoid the overabundance of dark matter. Also plotted in Fig. 1 is the corresponding branching ratio in the plane ($\Lambda$, $B_{R}$) for two different values of $y_{\chi}$ ($y_{\chi}=0$ and $y_{\chi}=10^{-4}y_{f}$). First of all, we notice that the lines lie in the region where the decay of the inflaton is responsible of the total amount of dark matter. The thermal production does not contribute in this area. Secondly, the plateau seen by the dashed curve corresponds to the points where the dark matter is completely produced by the tree level decay, uniquely determined by $y_{\chi}$ and thus independent of $\Lambda$. Moreover, even if $y_{\chi}=10^{-4}y_{f}$, for low values of $\Lambda\lesssim 5\times 10^{14}$ GeV, the loop contribution dominates over the tree-level production.

IV IV. Supergravity

We next consider a more concrete supergravity model where the inflaton decays to only a pair of Higgs bosons at tree level, while a pair of gravitinos are produced via the loops of Higgses and Higgsinos. As a specific application, we would like to reconsider the high scale supergravity model described in grav3 . The model is based on no-scale supergravity no-scale .The inflationary and supersymmetry breaking sector contains three chiral fields, $T,\phi$ driving inflation ENO6 , and a Polonyi-like field, $z$ pol whose superpotential we take

[TABLE]

with $\nu$ being a constant. The Polonyi field is assumed to be twisted and strongly stabilized dine ; ego ; dlmmo so that the Kähler potential takes the form

[TABLE]

where one of the $\phi_{i}$ is related to inflation, and the rest are matter fields. Choosing Cecotti ,

[TABLE]

leads to Starobinsky-like inflation Staro and is consistent with Planck observations planck when the inflaton mass, which we now designate as $M_{T}\simeq 3\times 10^{13}$ GeV.

If $\phi$ in (26) is fixed so that $\langle\phi\rangle=0$, the inflaton is associated with $T$ and its coupling to Standard Model fields leads to reheating, where the canonically normalized inflaton field $t$ is defined as $T\simeq(1/2)(1+\sqrt{2/3}t)$ about the minimum of $\langle T\rangle=1/2$. In a high scale supersymmetry model, all of the superpartners, except the gravitino, are assumed to be more massive than the inflaton. In particular, the MSSM $\mu$-parameter is also large, $\mu>M_{T}$. In this case, the dominant decay mode for the inflaton is $t\rightarrow H_{u,d}{H^{*}}^{u,d}$ which ultimately corresponds to a decay of $t\rightarrow hh$ where $h$ is the SM Higgs boson. The decay rate to two Higgs bosons is EGNO4

[TABLE]

with $y^{2}\equiv\mu^{4}/(6m_{T}^{2}M_{P}^{2})$. Then the reheating temperature can be expressed as $T_{\rm RH}\simeq 0.5(y/2\pi)\sqrt{M_{T}M_{P}}$ grav2 ; egnop ; Pradler:2006hh and $T_{\rm max}\simeq 0.5(8\pi/y^{2})^{1/4}T_{\rm RH}$.

Gravitinos are produced by the pair annihilation of SM particles, ${\rm SM}~{}{\rm SM}\rightarrow$ two gravitinos, where the resultant gravitino abundance depends on the reheating process as discussed above. The reaction rate of the gravitino production is given by Benakli:2017whb

[TABLE]

with $\Lambda^{8}\equiv 9m_{3/2}^{4}M_{\rm P}^{4}/21.65$. Then, from Eqs. (12) and (13), we obtain

[TABLE]

where the instantaneous thermalization of SM particles is implicitly assumed. On the other hand, when the non-instantaneous thermalization effect in the gravitino production is incorporated, we have Garcia:2018wtq

[TABLE]

where only gluon pair annihilation was assumed. We evaluate the SU(3)C gauge coupling $\alpha_{3}=g_{s}^{2}/4\pi$ at $T_{RH}$ by solving renormalization group equations at two-loop level.

In this model, gravitinos are also produced by the tree-level decay of inflaton. Although the tree level coupling between gravitino and inflaton vanishes when $\langle\phi\rangle=0$, supersymmetry breaking shifts $\langle\phi\rangle$, giving rise to the tree level decay given by777By comparing to the result of Ref. EGNO4 , the given expression is corrected by taking into account higher order terms of vacuum expectation values of $T,~{}\phi,$ and $z$, and the mixing between $t$ and $z$.

[TABLE]

The radiative decay of inflaton to a pair of gravitinos is induced by the interactions given by

[TABLE]

where $C_{tSS}\equiv\sqrt{2/3}\mu^{2}$ and $C_{tFF}\equiv\mu^{2}/2\sqrt{6}$ EGNO4 .888Note that those interaction terms are obtained by using equations of motion. A more rigorous calculation would not use equations of motion when the interaction terms are relevant for loop diagrams. However, a more rigorous treatment of this particular process is not our main focus, and is expected to have only a small quantitative effect on our result. The other relevant terms in the supergravity Lagrangian are given by

[TABLE]

where we denote a chiral multiplet $(\varphi,\chi_{L})$. The relevant diagrams for the radiative decay of inflaton $t$ into a pair of gravitinos are shown in Fig. 3, where the dominant contribution is coming from the upper two diagrams (A and B). While a detailed discussion of the decay width is given in Appendix VI, in the case of $m_{3/2}\ll M_{T}\ll\mu$ we obtain an approximate expression given by

[TABLE]

where $\mu_{\rm ren}$ is the renormalization scale, and we take $\mu_{\rm ren}=M_{T}$ in our analysis. Then, by using

[TABLE]

we can evaluate $\Omega h^{2}|_{\rm decay}$ given in Eq. (17).

Figure 4 shows the required relation between the $\mu$-parameter and the gravitino mass when total gravitino abundance given by $\Omega h^{2}|_{\rm tot}=\Omega h^{2}|_{\rm ann}^{\rm inst/non-inst}+\Omega h^{2}|_{\rm decay}=0.12$, where we take $M_{T}=3\times 10^{13}$ GeV. The top panel shows the case assuming instantaneous thermalization using Eq. (29), and the bottom panel takes the non-instantaneous thermalization effect into account using Eq. (30). In the both panels, the black dotted line shows the contribution from thermal production alone, using only $\Omega h^{2}|_{\rm ann}$. The short spaced dotted line shows the relation with only the loop contribution included. As one can see, the contribution from the 1-loop decay diagrams dominate over the thermal annihilation processes when $m_{3/2}\lesssim 10^{-1}M_{T}$. The solid blue contour corresponds the sum of the annihilation and loop decay contributions in $\Omega h^{2}|_{\rm tot}=0.1197$ planck evaluated without the approximation on the Passarino-Veltman functions Passarino:1978jh , while the dashed blue line uses the approximate formula given in Eq. (34). To evaluate the Passarino-Veltman functions, we have utilized LoopTools Hahn:1998yk . For both blue lines, we take the tree-level decay contribution to be negligibly small by assuming $\Lambda_{z}\ll M_{P}$. The green dot-dashed lines are the total abundance with the tree level decay. The region below the solid blue line is the allowed parameter space.

V V. Conclusion

The origin of the dark matter relic density is unknown. The thermal production mechanism hlw , works quite well for weakly interacting massive particles. However, this mechanism does not work when the dark matter candidates are super-weakly interacting and decouple very early (before inflation) in the history of the Universe. Gravitinos are a well known example of dark matter candidates for which their relic abundance is not determined by thermal annihilations. Instead, these particles are produced during the reheating process after inflation.

Some fraction of the dark matter may also be produced directly from inflaton decay if the inflaton couples to the dark matter candidate. However, even if the inflaton is not directly coupled to the dark matter, we have shown that direct production necessarily still occurs through radiative processes as in Fig. 2. We further showed in this work that direct production of dark matter through the radiative decay of the inflaton can in fact dominate the relic abundance in the Universe.

We highlighted this result with a specific example from no-scale supergravity. The model makes use of no-scale Starobinsky-like inflation, with supersymmetry breaking through a strongly stabilized Polonyi mechanism. We considered the case of high-scale supersymmetry in which all superpartners except the gravitino are more massive than the inflaton. Therefore, the thermal content of the universe is just that of the Standard Model. The dominant decay of the inflaton in this is model is to a pair of Standard Model Higgs bosons. The branching fraction of direct decays to gravitinos is suppressed (as given in Eq. (36)), though the degree of suppression depends on the strong stabilization parameter $\Lambda_{z}$. However the branching fraction to gravitinos through the loop process shown in Fig. 3 and given in Eq. (37) dominates over the thermal production when $m_{3/2}<10^{-1}M_{T}\sim 3\times 10^{12}$ GeV.

Acknowledgments: The authors want to thank especially E. Dudas and M. Garcia for very insightful discussions. This research has been supported by the (Indo-French) CEFIPRA/IFCPAR Project No. 5404-2. Support from CNRS LIA-THEP and the INFRE-HEPNET of CEFIPRA/IFCPAR is also acknowledged. This work was also supported by the France-US PICS MicroDark. Y. Mambrini acknowledges partial support from the European Union Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie: RISE InvisiblesPlus (grant agreement No 690575) and the ITN Elusives (grant agreement No 674896). The work of K. Kaneta and K. A. Olive was supported in part by the DOE grant DE-SC0011842 at the University of Minnesota.

Appendix

VI Calculation of the radiative gravitino decay

The relevant gravitino interactions in the supergravity Lagrangian999See, for instance, Ref. Lsugra that involves Higgs and Higgsino are given by

[TABLE]

where $G\equiv K+\ln|W|^{2}$, $G_{i}(G_{ij^{*}})$ represents derivative of $G$ with respect to $\varphi^{i}$ ($\varphi^{i}$ and $\varphi^{*j}$), and $\Phi_{\rho}\equiv G_{i}\partial_{\rho}\varphi^{i}-G_{i^{*}}\partial_{\rho}\varphi^{*i}$ which vanishes for $i=H_{u}$ and $H_{d}$. We also have $G_{H_{u}H_{u}^{*}}=G_{H_{d}H_{d}^{*}}\simeq 1$. The second term of Eq. (39) gives Higgs-Higgs-gravitino-gravitino coupling. In the case of $m_{3/2}\ll\mu$, the dominant terms are induced by

[TABLE]

Thus, we have the four diagrams, shown in Fig. 3, that induce the radiative inflaton decay into a pair of gravitinos.

In calculating the decay amplitudes, the gravitino equation of motion (EOM) may be used to simplify the expressions, which is given as

[TABLE]

Using the EOM, we also obtain the following relations:

[TABLE]

In the following discussion, we will use these equations to simplify the amplitudes whenever possible.

First, let us consider diagrams A and B. For a generic mass spectrum, the amplitudes are given by

[TABLE]

where we have omitted the factors arising from changing the basis of Higgs and Higgsinos into neutralinos and charginos, and the propagator in the denominators are defined as $D_{0}=q^{2}-m_{0}^{2}$ and $D_{i}=(q+p_{i})^{2}-m_{i}^{2}~{}(i=1,2)$. The gravitino wave functions are denoted as $u_{\mu}$ and $v_{\nu}$ which satisfy the Majorana condition $u_{\mu}=v_{\mu}^{c}$. The momenta in the propagators are taken to be $p_{1}=p$ and $p_{2}=-k$, while $m_{0}$ and $m_{1}(=m_{2})$ are the Higgsino and Higgs boson masses for diagram A, and are the Higgs and Higgsino masses for diagram B, respectively. Then, the resultant amplitudes are given by

[TABLE]

where $C_{ij}^{I}$ and $C_{ijk}^{I}$ are the Passarino-Veltman functions Passarino:1978jh for the diagram $I=A,B$, which are defined as

[TABLE]

These functions may in general be written as

[TABLE]

The terms denoted by I, II, III, … are the amplitudes simplified by the EOM and the relevant relations, and are given by

1. I.

$\overline{u}_{\mu}P_{L}\gamma^{\alpha}\gamma^{\mu}\gamma_{\alpha}\gamma^{\nu}\cancel{p}_{i}v_{\nu}=0$, 2. II.

$\overline{u}_{\mu}P_{L}\cancel{p}_{i}\gamma^{\mu}\gamma^{\beta}\gamma^{\nu}\gamma_{\beta}v_{\nu}=0$, 3. III.

$\overline{u}_{\mu}P_{L}\gamma^{\alpha}\gamma^{\mu}\cancel{p}_{1}\gamma^{\nu}\gamma_{\alpha}u_{\nu}=4m_{3/2}g^{\mu\nu}\overline{u}_{\mu}P_{R}v_{\nu}$, 4. IV.

$\overline{u}_{\mu}P_{L}\gamma^{\alpha}\gamma^{\mu}\cancel{p}_{2}\gamma^{\nu}\gamma_{\alpha}u_{\nu}=4m_{3/2}g^{\mu\nu}\overline{u}_{\mu}P_{L}v_{\nu}$, 5. V.

$\overline{v}_{\mu}P_{L}\cancel{p}_{1}\gamma^{\mu}\cancel{p}_{1}\gamma^{\nu}\cancel{p}_{1}v_{\nu}=0$, 6. VI.

$\overline{u}_{\mu}P_{L}\cancel{p}_{2}\gamma^{\mu}\cancel{p}_{2}\gamma^{\nu}\cancel{p}_{2}v_{\nu}=0$, 7. VII.

$\overline{u}_{\mu}P_{L}\cancel{p}_{1}\gamma^{\mu}\cancel{p}_{1}\gamma^{\nu}\cancel{p}_{2}v_{\nu}=0$, 8. VIII.

$\overline{u}_{\mu}P_{L}\cancel{p}_{1}\gamma^{\mu}\cancel{p}_{2}\gamma^{\nu}\cancel{p}_{2}v_{\nu}=0$, 9. IX.

$\overline{u}_{\mu}P_{L}\cancel{p}_{2}\gamma^{\mu}\cancel{p}_{1}\gamma^{\nu}\cancel{p}_{1}v_{\nu}=-4m_{3/2}k^{\mu}p^{\nu}\overline{u}_{\mu}P_{R}v_{\nu}$, 10. X.

$\overline{u}_{\mu}P_{L}\cancel{p}_{2}\gamma^{\mu}\cancel{p}_{2}\gamma^{\nu}\cancel{p}_{1}v_{\nu}=-4m_{3/2}k^{\mu}p^{\nu}\overline{u}_{\mu}P_{L}v_{\nu}$, 11. XI.

$\overline{u}_{\mu}P_{L}\cancel{p}_{2}\gamma^{\mu}\cancel{p}_{1}\gamma^{\nu}\cancel{p}_{2}v_{\nu}=0$, 12. XII.

$\overline{u}_{\mu}P_{L}\cancel{p}_{1}\gamma^{\mu}\cancel{p}_{2}\gamma^{\nu}\cancel{p}_{2}v_{\nu}=0$,

where the momentum for $\overline{u}_{\mu}$ and $v_{\nu}$ is assigned to be $\overline{u}_{\mu}(p)$ and $v_{\nu}(k)$.

For diagram C, we have

[TABLE]

where $c_{H}$ represents a constant factor coming from the mixing of Higgs bosons, and the vertex factor is defined as

[TABLE]

Using Eq. (44), we obtain

[TABLE]

where $B_{0}$ is the two-point scalar function defined as

[TABLE]

The coupling $c_{H}$ is found by recalling the relation between $H_{u,d}$ and their mass eigenstates101010We denote the lightest neutral Higgs $h^{0}$, heavier neutral Higgs $H^{0}$, CP-odd Higgs $A$, charged Higgs $H^{\pm}$, neutral Goldstone boson $G^{0}$, and charged Goldstone boson $G^{\pm}$. given as

[TABLE]

with $s_{x}\equiv\sin x$ and $c_{x}\equiv\cos x$. We find

[TABLE]

In the $m_{A}\gg m_{Z}$ limit and $\tan\beta\simeq 1$, we have $-s_{\alpha}c_{\alpha}\sim s_{\beta}c_{\beta}\sim 1/2$. Therefore, in the supersymmetric limit, all of the scalar masses are degenerate in mass at $\mu$, and the sum over all contributions identically vanishes. In the case that $m_{H^{0}}\sim m_{A}\sim m_{H^{\pm}}\sim\mu\gg m_{h^{0}}$ and $m_{G^{0}}=m_{G^{\pm}}=0$, we have

[TABLE]

so the contribution from diagram C turns out to be negligible.

For diagram D, the amplitude is given by

[TABLE]

where $p_{t}$ is the inflaton momentum, and

[TABLE]

After computing the trace in the amplitude, the loop integral is

[TABLE]

where $B_{1}$ is defined as

[TABLE]

The spinor piece in the amplitude then gives

[TABLE]

The first term may be written as

[TABLE]

where we have used the EOM. For the second term we obtain

[TABLE]

where we have used $\cancel{k}v_{\rho}(k)=-m_{3/2}v_{\rho}(k)$ and $\overline{u}_{\mu}(p)\cancel{p}=m_{3/2}\overline{u}_{\mu}(p)$. Therefore, the amplitude vanishes for diagram D.

Hence, the total amplitude, whose dominant contributions come from diagrams A and B, can be written as

[TABLE]

where the overall factor 2 comes from the exchange of final state gravitinos, and we define

[TABLE]

with the coefficients being

[TABLE]

By denoting $F_{i,j}\equiv{\cal M}_{i}^{*}{\cal M}_{j}$ with $i=1R,1L,2R,2L$, we have

[TABLE]

where $\tau\equiv 4m_{3/2}^{2}/M_{T}^{2}$.

Next, let us express all Higgs and Higgsinos as mass eigenstates.111111We denote the Higgsino-like chargino as $\widetilde{\chi}_{2}^{\pm}$, and $\widetilde{\chi}_{3,4}^{0}$ are the Higgsino-like neutralinos. The amplitude for A and B is given as

[TABLE]

for $I=A,B$, and thus we have,

[TABLE]

As a good approximation for high-scale supersymmetry, let us suppose that $m_{\widetilde{\chi}_{2}^{\pm}}\sim m_{\chi_{3,4}^{0}}\sim m_{H^{\pm}}\sim m_{H^{0}}\sim m_{A}\sim\mu\gg m_{h^{0}}$. Then, we may write the total amplitude as

[TABLE]

where ${\cal M}_{\rm tot}^{h}\equiv\hat{\cal M}_{A+B}(\widetilde{\chi}^{0}_{i},h^{0})=\hat{\cal M}_{A+B}(\widetilde{\chi}^{0}_{i},G^{0})=\hat{\cal M}_{A+B}(\widetilde{\chi}^{0}_{i},G^{\pm})$ and ${\cal M}_{\rm tot}^{H}\equiv\hat{\cal M}_{A+B}(\widetilde{\chi}_{2}^{\pm},H^{\pm})=\hat{\cal M}_{A+B}(\widetilde{\chi}^{0}_{i},H^{0})=\hat{\cal M}_{A+B}(\widetilde{\chi}^{0}_{i},A)$.

We can further simplify the amplitude by taking the limit $m_{3/2}\ll M_{T}\ll\mu$. In this limit, the leading amplitudes among $F_{i,j}$ are $F_{1R,1L},F_{2R,2L}$, and $F_{1R,2L}$ for $\tau\ll 1$, namely,

[TABLE]

Besides, as we can see from the dimensionality, $C_{00}$ and $C_{00i}$ in $c_{1R}$ and $c_{1L}$ remain constant (or logarithmic of the involved mass scales, i.e., $\mu$) as it is dimensionless, while $C_{12}^{B}$ in $c_{2R}$ and $c_{2L}$ scales as $1/\mu^{2}$. Indeed, we obtain analytic expressions of those functions in this limit:

[TABLE]

and $C_{001}^{I,h(H)}=C_{002}^{I,h(H)}$ for $I=A,B$. Therefore, the dominant contribution in the decay width is the term involving ${\rm Re}(c_{1R}^{*}c_{1L})F_{1R,1L}$ in the squared amplitude. Then, we obtain

[TABLE]

where $\bar{c}_{i}\equiv c_{i}^{h}+c_{i}^{H}$ for $i=1R,1L,2R,2L$, and we have used

[TABLE]

with the relation $C_{tSS}=2\mu C_{tFF}$.