Super heavy thermal dark matter

TL;DR

This paper introduces a novel mechanism for super heavy thermal dark matter with masses up to 10^{14} GeV, maintaining equilibrium through a chain of nearly degenerate particles interacting with the Standard Model, and explores its cosmological and experimental implications.

Contribution

It presents a new model for super heavy thermal dark matter that avoids unitarity bounds and involves a chain of nearly degenerate particles maintaining equilibrium.

Findings

Dark matter can have masses up to 10^{14} GeV within this framework.

The model predicts heavy relics decaying at different cosmological epochs.

Implications for CMB, structure formation, and cosmic ray experiments are discussed.

Abstract

We propose a mechanism of elementary thermal dark matter with mass up to $10^{14}$ GeV, within a standard cosmological history, whose relic abundance is determined solely by its interactions with the Standard Model, without violating the perturbative unitarity bound. The dark matter consists of many nearly degenerate particles which scatter with the Standard Model bath in a nearest-neighbor chain, and maintain chemical equilibrium with the Standard Model bath by in-equilibrium decays and inverse decays. The phenomenology includes super heavy elementary dark matter and heavy relics that decay at various epochs in the cosmological history, with implications for CMB, structure formation and cosmic ray experiments.

Equations

$${\chi }_i+\mathrm{sm}\leftrightarrow {\chi }_{i+1}+\mathrm{sm},$$

$${\chi }_N\to \mathrm{sm}+\mathrm{sm}.$$

$${\chi }_1+\mathrm{sm}\leftrightarrow {\chi }_2+\mathrm{sm},$$

$${\chi }_i+\mathrm{sm}\leftrightarrow {\chi }_{i+1}+\mathrm{sm},$$

$${\chi }_N\to \mathrm{sm}+\mathrm{sm}.$$

$$\text{leavevmode}to6.67pt\text{vbox}to6.67pt\text{pgfpicture}\text{makeatletter}\text{lower}-3.33301ptto0.0pt\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{definecolor}pgfstrokecolorrgb\mathrm{0,0,0}\text{pgfsys@color@rgb@stroke}000\text{pgfsys@invoke}\text{pgfsys@color@rgb@fill}000\text{pgfsys@invoke}\text{pgfsys@setlinewidth}0.4pt\text{pgfsys@invoke}\text{nullfont}to0.0pt\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{feynman}[small,\; baseline=(a)]\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@transformcm}1.00.00.01.00.0pt0.0pt\text{pgfsys@invoke}\text{definecolor}pgfstrokecolorrgb\mathrm{0,0,0}\text{pgfsys@color@rgb@stroke}000\text{pgfsys@invoke}\text{pgfsys@color@rgb@fill}000\text{pgfsys@invoke}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{vertex}[above\; left\; =of\; a]\; (i1){\chi }_1;\text{vertex}[below\; left=of\; a]\; (i2)\mathrm{sm};\text{vertex}[above\; right\; =of\; a]\; (f1){\chi }_2;\text{vertex}[below\; right=of\; a]\; (f2)\mathrm{sm};\text{diagram}*(i1)\; –\; (a)–\; (f1),\; (i2)\; –\; (a)\; –\; (f2),;\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{hss}\text{pgfsys@discardpath}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{hss}\text{lxSVG@closescope}\text{endpgfpicture}\text{leavevmode}to6.67pt\text{vbox}to6.67pt\text{pgfpicture}\text{makeatletter}\text{lower}-3.33301ptto0.0pt\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{definecolor}pgfstrokecolorrgb\mathrm{0,0,0}\text{pgfsys@color@rgb@stroke}000\text{pgfsys@invoke}\text{pgfsys@color@rgb@fill}000\text{pgfsys@invoke}\text{pgfsys@setlinewidth}0.4pt\text{pgfsys@invoke}\text{nullfont}to0.0pt\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{feynman}[small,\; baseline=(a)]\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@transformcm}1.00.00.01.00.0pt0.0pt\text{pgfsys@invoke}\text{definecolor}pgfstrokecolorrgb\mathrm{0,0,0}\text{pgfsys@color@rgb@stroke}000\text{pgfsys@invoke}\text{pgfsys@color@rgb@fill}000\text{pgfsys@invoke}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{vertex}[above\; left\; =of\; a]\; (i1){\chi }_2;\text{vertex}[below\; left=of\; a]\; (i2)\mathrm{sm};\text{vertex}[above\; right\; =of\; a]\; (f1){\chi }_3;\text{vertex}[below\; right=of\; a]\; (f2)\mathrm{sm};\text{diagram}*(i1)\; –\; (a)–\; (f1),\; (i2)\; –\; (a)\; –\; (f2),;\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{hss}\text{pgfsys@discardpath}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{hss}\text{lxSVG@closescope}\text{endpgfpicture}\cdots \text{leavevmode}to6.67pt\text{vbox}to15.89pt\text{pgfpicture}\text{makeatletter}\text{lower}-12.55656ptto0.0pt\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{definecolor}pgfstrokecolorrgb\mathrm{0,0,0}\text{pgfsys@color@rgb@stroke}000\text{pgfsys@invoke}\text{pgfsys@color@rgb@fill}000\text{pgfsys@invoke}\text{pgfsys@setlinewidth}0.4pt\text{pgfsys@invoke}\text{nullfont}to0.0pt\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{feynman}[small,\; baseline=(a)]\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@transformcm}1.00.00.01.00.0pt0.0pt\text{pgfsys@invoke}\text{definecolor}pgfstrokecolorrgb\mathrm{0,0,0}\text{pgfsys@color@rgb@stroke}000\text{pgfsys@invoke}\text{pgfsys@color@rgb@fill}000\text{pgfsys@invoke}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{vertex}[above\; left\; =of\; a]\; (i1){\chi }_{{}_{N-1}};\text{vertex}[below\; left=of\; a]\; (i2)\mathrm{sm};\text{vertex}[above\; right\; =of\; a]\; (f1){\chi }_{{}_N};\text{vertex}[below\; right=of\; a]\; (f2)\mathrm{sm};\text{vertex}[right=.8cm\; of\; f1]\; (c)\; ;\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@beginscope}\text{pgfsys@invoke}\text{pgfsys@transformcm}1.00.00.01.00.0pt-9.22356pt\text{pgfsys@invoke}\text{definecolor}pgfstrokecolorrgb\mathrm{0,0,0}\text{pgfsys@color@rgb@stroke}000\text{pgfsys@invoke}\text{pgfsys@color@rgb@fill}000\text{pgfsys@invoke}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{vertex}[\; right=.75cm\; of\; b]\; (f3)\mathrm{sm};\text{vertex}[\; below=.55cm\; of\; b]\; (d)\; ;\text{vertex}[\; right=.2cm\; of\; d]\; (f4)\mathrm{sm};\text{diagram}*(i1)\; –\; (a)–\; (f1)\; –\; (b)\; –\; (f3),\; (i2)\; –\; (a)\; –\; (f2),\; (b)\; –\; (f4),;\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{hss}\text{pgfsys@discardpath}\text{pgfsys@invoke}\text{lxSVG@closescope}\text{pgfsys@endscope}\text{hss}\text{lxSVG@closescope}\text{endpgfpicture}$$

$${\chi }_1+\mathrm{sm}\leftrightarrow {\chi }_2+\mathrm{sm}.$$

$${\dot{n}}_1+3H{n}_1$$

$${\dot{n}}_2+3H{n}_2$$

$$Y_1^{\prime }=\frac{\lambda }{x^2}(-Y_1+Y_{\mathrm{eq}}),$$

$$Y_1(\infty )\approx \frac{45g_{\chi }}{2^{3/2}{\pi }^3{g}_{\star s}}\lambda e^{-2\sqrt{\lambda }}\equiv Y_{\infty }(\lambda ),$$

$$Y_1^{\prime }=-\frac{\lambda }{x^2}{Y}_1,$$

$$Y_1^{\prime }$$

$$Y_j^{\prime }$$

$$Y_N^{\prime }$$

$$({\partial }_{\tau }-D{\partial }_{\ell }^2)Y_{\ell }(\tau )=0,$$

$${\partial }_{\ell }Y{|}_{\ell =0}=0,Y_{\pi /2}(\tau )=Y^{\mathrm{eq}}(\tau )$$

$$Y_{i}(\infty)\simeq Y_{\infty}(D)\frac{4}{\pi}\cos\Big{(}\frac{\pi}{2}\frac{j}{N}\Big{)}\,,$$

$${\Gamma }_1=\frac{1}{2m}\int d{\Phi }_{2N}|\mathcal{M}{|}^2.$$

$$\int d{\Phi }_{2N}=\frac{1}{S}\frac{2\pi }{\Gamma (2N)\Gamma (2N-1)}\frac{m^{4(N-1)}}{(16{\pi }^2{)}^{2N-1}},$$

$$|\mathcal{M}{|}^2\simeq S^2\frac{1}{m^{4N-4}}|\mathcal{M}{|}_{\chi \mathrm{sm}\to \chi \mathrm{sm}}^{2(N-1)}|\mathcal{M}{|}_{{\chi }_N\to \mathrm{sm}+\mathrm{sm}}^2$$

$$\frac{{\Gamma }_1}{{\Gamma }_N}\simeq \frac{S}{(2N-1)!(2N-2)!}{\left(\frac{{\alpha }^2}{16{\pi }^3}\right)}^{N-1}.$$

$$(m_{\max }-m_{\min })/m_{\min }<x_{\mathrm{fin}}^{-1}=({\pi }^2\lambda /4N^2{)}^{-1},$$

$$\mathcal{L}\supset -\frac{1}{2}{m}_i{\chi }_i{\chi }_i-\delta m{\chi }_i{\chi }_{i+1}-yS{\chi }_i{\chi }_i-\mu S|H{|}^2-\tilde{y}{\chi }_{N}LH,$$

$$M_{ij}=m_i{\delta }_{ij}+\delta m({\delta }_{i,j+1}+{\delta }_{i,j-1}).$$

$${\xi }_{\mathrm{loc}}^{-1}\simeq \ln \frac{\Delta m}{2\delta m}-1.$$

$${\psi }_i\sim \sum _j{\chi }_j\exp \left[-\frac{|i-j|}{{\xi }_{\mathrm{loc}}}\right].$$

$$\mathcal{L}\supset y{e}^{-1/{\xi }_{\mathrm{loc}}}S{\psi }_i{\psi }_{i+1},,.$$

$${\left(\frac{{\Gamma }_1}{{\Gamma }_N}\right)}_{\mathrm{mixing}}\simeq e^{-2N/{\xi }_{\mathrm{loc}}}.$$

$${\xi }_{\mathrm{loc}}^{-1}\gtrsim \frac{1}{N}\left(62+\log \frac{m}{10^{10}\mathrm{GeV}}\right)$$

$${\mathbf{Y}}^{\prime }\approx \frac{\lambda }{x^2}\mathbf{M}\cdot \mathbf{Y}+\frac{\lambda }{x^2}{\mathbf{n}}_{N-1}{Y}_{\mathrm{eq}},$$

$$(\mathbf{M}{)}_{ij}={\alpha }_i({\delta }_{i,j-1}-{\delta }_{ij})+({\delta }_{i,j+1}-{\delta }_{ij}),$$

Full text

I Introduction

One of the biggest questions in fundamental physics is the nature of dark matter (DM). The possibility that DM is a thermal Weakly Interacting Massive Particle (WIMP), whose abundance is determined by $2\to 2$ annihilations into Standard Model (SM) bath particles, is exciting, but has alluded detection thus far. The WIMP is particularly intriguing because it is very predictive—its abundance is determined only by its interactions with the SM, which informs us how it may be detected.

The WIMP paradigm has been a guide towards the properties of DM, such as its mass and interactions. In particular, within the WIMP freezeout mechanism, there is a an upper bound on the DM mass from perturbative unitarity, of ${\cal O}(100)\,{\rm TeV}$ Griest:1989wd . The reason is that the WIMP annihilation rate is proportional to an exponentially decreasing DM density, and so the amount of dark matter that can be annihilated away before freezeout is limited by the theoretical size of the cross-section.

Of course, DM may be much heavier than this bound, if it is not a WIMP. Such models include non-thermal dynamics, decoupled dark sectors, inflationary and gravitational production, nonstandard cosmological histories, and large entropy production Hui:1998dc ; Kolb:1998ki ; Chung:1998rq ; Chung:2001cb ; Feng:2008mu ; Harigaya:2014waa ; Davoudiasl:2015vba ; Randall:2015xza ; Dev:2016xcp ; Harigaya:2016vda ; Berlin:2016vnh ; Berlin:2016gtr ; Bramante:2017obj ; Berlin:2017ife ; Hamdan:2017psw ; Cirelli:2018iax ; Babichev:2018mtd ; Hashiba:2018tbu ; Hooper:2019gtx . In none of these cases, however, is the DM abundance solely determined by its interactions with the SM. The common lore is that an elementary DM candidate that is thermally coupled to the SM, within a standard cosmological history, cannot have mass well above the WIMP perturbative unitarity bound. Exceptions include composite DM, but still do not go much beyond the above bound Harigaya:2016nlg ; Geller:2018biy ; Smirnov:2019ngs .

In this Letter we present a new freezeout mechanism within a standard cosmological history. The DM is an elementary particle that is thermalized with the SM at high temperatures, its relic abundance is determined via its freezeout from the SM bath, and the DM mass can be as high as $10^{14}$ GeV for s-wave processes, without violating the perturbative unitarity limit. In future work, we show how Planck-scale DM can be reached for velocity-dependent processes future1 .

The general idea is as follows. The dark matter consists of $N$ approximately degenerate states, $\chi_{i}$ ($i=1..N$). These states co-scatter DAgnolo:2017dbv off of the SM bath, but only in a chain of nearest-neighbor interactions

[TABLE]

while the $N^{\rm th}$ state co-decays Bandyopadhyay:2011qm ; Dror:2016rxc ; Kopp:2016yji ; Okawa:2016wrr ; Dror:2017gjq ; Dery:2019jwf in equilibrium with the SM,

[TABLE]

Here $\chi_{1}$ is the DM candidate. This setup is summarized in Fig. 1. The processes in Eq. (1) can maintain chemical equilibrium much longer than annihilations can, because the interaction rate for the scattering $\Gamma_{\rm sctr}=n_{\rm SM}\left<\sigma v\right>_{\rm sctr}$ never becomes exponentially suppressed. The in-equilibrium decay allows for the whole system to have vanishing chemical potential for a long time—if the $N^{\rm th}$ particle was annihilating with the bath, the system would inherit the unitarity bound from co-annihilations. Finally, the chain, which will typically require $N>5-20$, depending on the DM mass, ensures the stability of the DM $\chi_{1}$.

II General Idea

Consider a DM particle $\chi_{1}$, whose density changes via scattering with a light SM bath particle,

[TABLE]

where $\chi_{2}$ has similar mass to $\chi_{1}$. This process can maintain chemical equilibrium much longer than annihilations with the same interaction strength, because the interaction rate for the scattering does not depend on the DM density, and therefore the rate does not become exponentially suppressed. The DM is able to maintain equilibrium to smaller temperatures, becoming more Boltzmann suppressed than a WIMP for the same size cross-section.

However, $\chi_{1}$ can only reduce exponentially for as long as $\chi_{2}$ reduces exponentially (for instance, by maintaining chemical equilibrium with the SM bath). Thus in order to go beyond the unitarity bound on annihilations, some other process is still needed to reduce the $\chi_{2}$ density. If this process is annihilations (such as $\chi_{2}+\chi_{2}\leftrightarrow~{}{\rm sm}+{\rm sm}$), then it will freezeout too early, and the unitarity bound on annihilations will apply again. If instead, $\chi_{2}$ decays equilibrium with the SM bath, chemical equilibrium can be maintained for much longer. Thus, our proposed mechanism is a combination of co-decay and co-scattering dark matter.

However, one can easily see that the combination of the scattering process Eq. (3), and the in-equilibrium decay process $\chi_{2}\to{\rm sm}+{\rm sm}$, will necessarily lead to the fast decay of $\chi_{1}$ via an off-shell $\chi_{2}$. In the presence of a co-scattering chain, such that scatters take place only for nearest neighbors

[TABLE]

which codecay in equilibrium with $\chi_{N}$,

[TABLE]

$\chi_{1}$ can instead be long lived; here $\chi_{1}$ is still the DM, but its decay width is suppressed due to the large phase space needed to decay to the SM ($\chi_{1}$ decays to $2N$ SM particles). This is summarized in Fig. 1.

III Two particle case

As a toy example, we first work through the simplest case of $N=2$ dark matter particles with mass $m$. While this example is not cosmologically viable, it demonstrates the basic idea of the mechanism. Consider two degenerate states $\chi_{1}$ and $\chi_{2}$, which co-scatter DAgnolo:2017dbv off of SM bath particles via the process

[TABLE]

The dark matter candidate is $\chi_{1}$, while $\chi_{2}$ decays in equilibrium with the SM bath.

The Boltzmann equations for the system are

[TABLE]

where $\langle\sigma v\rangle$ is a thermally averaged cross section for the scattering process, $n_{\rm eq}$ is the equilibrium number density of $\chi$, and $\Gamma_{2}$ is the decay rate of $\chi_{2}$ in the thermal bath. This system is similar to the co-scattering scenario DAgnolo:2017dbv , but here the DM bath is kept in equilibrium with the SM bath via decays and inverse-decays, rather than annihilations. Ultimately, since decays and inverse-decays can stay in equilibrium much longer, this will allow for much heavier DM.

Unlike the case for freezeout via annihilations, the instantaneous freezeout approximation will not give a good estimate of the relic abundance. This is because the rate for $\chi_{1}$ scattering, $\Gamma\sim n_{\rm sm}\langle\sigma v\rangle$, is not dropping off exponentially fast with the expansion, and therefore freezeout takes a long time. However, an approximate analytic solution to the relic abundance can still be determined from the Boltzmann equations, as we now detail.

If the decay rate $\Gamma_{2}$ is larger than the Hubble expansion parameter when the temperature $T$ of the universe is equal to the DM mass, the number density of $\chi_{2}$ closely follows its equilibrium value. An approximate solution to the $\chi_{1}$ density can then be found by considering the single equation

[TABLE]

where $Y_{i}=n_{i}/s$ with $s$ the entropy density, $x=m/T$ and $\lambda=(n_{\rm sm}\langle\sigma v\rangle/H)|_{x=1}$. We have also assumed that the thermally averaged cross-section is velocity independent and took the number of relativistic degrees of freedom $g_{\star}=g_{\star s}$ to be constant. The asymptotic value of the relic abundance is

[TABLE]

where $g_{\chi}$ is the number of internal degrees of freedom of a $\chi$ particle.

From the Boltzmann equation (8), $\chi_{1}$ departs equilibrium when ${\lambda}/{x_{\rm fo}^{2}}\sim 1$. The fact that the relic abundance then scales as $Y_{1}(\infty)\sim e^{-2x_{\rm fo}}$ and not as $e^{-x_{\rm fo}}$ (as one finds in the instantaneous freezeout approximation), is the result of the slow freezeout. At this point $\chi_{1}$ creation stops, but $\chi_{1}$ can continue to scatter away. Neglecting the inverse process, solving

[TABLE]

yields that $Y_{1}(x)=Y_{1}(x_{\rm fo})e^{\frac{\lambda}{x}-\frac{\lambda}{x_{\rm fo}}}{\to}\,e^{-2x_{\rm fo}}$. One also sees from this solution that the abundance stops changing significantly when $x_{\rm fin}\sim\lambda=x_{\rm fo}^{2}$.

From the above estimation, one can find $\langle\sigma v\rangle$ that reproduces the observed relic abundance of dark matter, while satisfying the unitarity bound. For instance, parameterizing the cross-section as $\langle\sigma v\rangle\equiv\alpha^{2}/m_{\chi}^{2}$, one finds that for $\alpha\simeq 1$, ${m_{\chi}=6\times 10^{14}}$ GeV reproduces the correct relic abundance (assuming the DM scatters off of 4 SM degrees of freedom and that $g_{\star s}=106.75$).

IV $N>2$ Degenerate Case

Although the simplest $N=2$ example does not work, it indicates a clear path forward to a viable setup: suppressing the decay rate of $\chi_{1}$. We thus consider $N>2$ particles with the same type of nearest-neighbor interactions as before. Similarly, we assume that only $\chi_{N}$ is able to decay into SM particles, that the masses are degenerate, and that the cross section for each interaction is the same. The equations of the system are

[TABLE]

where the second equation is valid for $j=2,\cdots,N-1$, $\lambda=(n_{\rm sm}\langle\sigma v\rangle/H)|_{x=1}$, and $\lambda_{d}=(\Gamma_{N}/H)_{x=1}$. Here prime denotes a derivative with respect to $x$. We also ignore the time variation of $g_{\star}$ and $g_{\star s}$, and assume that the cross section is constant in the non-relativistic limit. Note also that we are interested in parameter space where $\lambda,\,\lambda_{d}\gg 1$.

This system can be solved in similar fashion to the $N=2$ case, by assuming $Y_{N}=Y_{\rm eq}$ and diagonalizing the differential equations. However, we can more simply obtain a solution with the new $N$ scaling by taking the large $N$ limit, and taking flavor space to be continuous. For $N\gg 1$, the system behaves as a random walk between $\chi_{j}$ states, with a reflecting wall at $j=1$ and an absorbing cliff at $j=N$ (corresponding to the decay). Taking the indices of $j$ as a continuous variable $\ell=\pi j/[2(N-1)]$, the system is described by the diffusion equation

[TABLE]

where $\tau=-1/x$ and $D=\pi^{2}\lambda/[4(N-1)^{2}]$ is a diffusion coefficient, and the boundary conditions are

[TABLE]

From the form of the diffusion equation, we see that the system depends on the cross-section only through the combination $D\simeq\pi^{2}\lambda/4N^{2}$. From solving the boundary conditions, the asymptotic profile ${Y_{\ell}(\infty)\propto\cos\ell}$. The solution of the diffusion equations gives the relic abundance

[TABLE]

where $Y_{\infty}$ is defined in Eq. (9).

In Fig. 2 we plot the yield of the $\chi_{i}$ particles for $N=10$, solving the full Boltzmann equations, Eqs. (11). We find that the relic abundances match the analytical estimate of Eq. (14).

Now, one must check that $\chi_{1}$ is sufficiently long lived to be a viable DM candidate. A model-independent bound on the DM lifetime is $\tau>5\times 10^{18}$ sec Audren:2014bca (for earlier studies, see Refs. Ichiki:2004vi ; PalomaresRuiz:2007ry ; Lattanzi:2007ux ; Lattanzi:2008ds ; Gong:2008gi ; DeLopeAmigo:2009dc ; Peter:2010au ; Peter:2010sz ; Wang:2010ma ; Peter:2010jy ; Huo:2011nz ; Aoyama:2011ba ). However, if the DM decays to SM particles, constraints on the lifetime may be much stronger, $\tau>10^{27}$ sec (see e.g. Refs. Cirelli:2012ut ; Essig:2013goa ; Blanco:2018esa ).

The decay of $\chi_{1}$ to SM particles takes place only through $N-1$ off-shell $\chi_{i}$ particles, into $2N$ SM particles, with decay width

[TABLE]

Treating the SM particles as massless, the total $2N$-body phase space is Kleiss:1985gy

[TABLE]

where $S$ is a symmetry factor accounting for identical particles in the final state. We may approximate the squared matrix element of the decay as 111We approximate internal propagators as $1/m^{4}$. One might wonder if the largest contribution arises from internal propagators being almost on-shell. Such possibility arises only with a small phase space density, so the parametric dependence of the decay rate on $m$ does not change. Nevertheless, internal propagators generally lead to larger $\Gamma_{1}$, but this changes our estimate on $N$ only by $\Delta N\sim 1\textrm{--}2$ for the range of $N$ of our interest.

[TABLE]

where the factor $S^{2}$ accounts for expected number of diagrams. The decay rate of $\chi_{1}$ is then estimated as

[TABLE]

Assuming the $\chi_{N}$ decay is in equilibrium $\Gamma_{N}>H|_{x=1}$, requiring the stability of the DM, $\tau>10^{27}$ sec, places a lower bound on the number of dark matter particles $N$. For instance, assuming that the DM decays to $2N$ identical SM particles corresponding to $S=(2N)!$, we find for $m=10^{7}$ GeV that $N\geq 5$, and for $m=10^{13}$ GeV that $N\geq 14$. In Fig. 3, we plot the minimum value of $N$ as a function of $m_{\chi}$, satisfying the lifetime requirement.

V Non-degenerate masses and couplings

Up until now, we have considered the case of exactly degenerate masses, but it is natural to consider small mass splittings between the particles. This will have two important effects: creating forbidden channels and allowing the heavier states to decay directly into the lighter states.

In standard two-particle forbidden annihilations, the mass splitting is negligible if it is smaller than the temperature when the abundance is still evolving, i.e., $\delta m/m<x_{\rm fo}^{-1}$. For the chain interactions, a similar condition is found:

[TABLE]

where $m_{\rm max}$ and $m_{\rm min}$ is the maximum and minimum particle mass in the chain.

If the mass splitting is larger than this, then the evolution will enter a forbidden regime. The abundance of the lightest particle (the DM candidate) decouples when the Boltzmann suppression factor due to the mass splitting becomes significant and the chain reaction departs from chemical equilibrium. After this point, the abundance of the heavier particles continues its Boltzmann suppression, $Y_{\rm heavy}\sim e^{-\Delta m/T}$. Consequently, the relic abundances of the heavier states are much smaller than that of the lightest state when in the forbidden regime. We leave a detailed study of this forbidden regime to future work.

The second important effect of the mass splitting is that the heavier states can quickly cascade decay to the lightest state (which is still very long lived). Depending on the size of the mass-splitting, the decays may be in- or out-of-equilibrium. If the decays are in-equilibrium, the calculation of the relic-abundance will mostly remain unchanged. If the decays are out of equilibrium, the heavier states will simply transfer their abundances to the lightest state, increasing the DM abundance by a factor $\mathcal{O}(N)$.

It is also reasonable to consider that the interaction strengths may vary across the chain. In this case, the Boltzmann equations are best solved by diagonalizing the system of differential equations. The relic abundance is approximately given by the degenerate case, but with the cross-section replaced by the smallest eigenvalue cross-section in the chain.

VI Phenomenology

One of the main phenomenological signatures of the mechanism are decays of long-lived relics. Late time decaying DM ($\tau>10^{10}$ yr) can potentially be probed for all masses. DM decays can change the CMB anisotropies and spectrum Slatyer:2016qyl , but a detailed study for very heavy dark matter has not been performed. Decays of super heavy DM can source ultra-high-energy cosmic rays (UHECR), with energies $\gtrsim 10^{9}$ GeV. UHECR can be observed in diffuse gamma ray satellites, such as FERMI-LAT Ackermann:2014usa , high energy neutrino experiments, such as IceCube Abbasi:2011ji , and dedicated UHECR observatories, such as Auger Abreu:2011zze . For studies of the indirect detection of heavy decaying DM, see Refs. Berezinsky:1997hy ; Kuzmin:1998uv ; Barbot:2002gt ; Esmaili:2012us ; Murase:2012xs ; Feldstein:2013kka ; Rott:2014kfa ; Aloisio:2015lva ; Kalashev:2016cre ; Kalashev:2017ijd ; Blanco:2018esa ; Alcantara:2019sco . Prospects for discovering UHECR are expected to improve with the proposed POEMMA mission Olinto:2017xbi .

Depending on the parameters of the theory (especially the size of the mass splittings), some of the $\chi_{i}$ may decay at various cosmological epochs. For instance, a component of the DM could dissociate light elements during Big Bang Nucleosynthesis (see Ref. Kawasaki:2017bqm and references therein). Additionally, DM decay can lead to spectral distortions in the CMB Hu:1993gc ; Chluba:2011hw ; Poulin:2016anj , which would be probed by the proposed PiXiE experiment Kogut:2011xw . Partially decaying DM has also been shown to alleviate several cosmological tensions, such as the Hubble tension Vattis:2019efj ; Pandey:2019plg and small scale structure puzzles (see e.g. Refs. Wang:2014ina ; Enqvist:2015ara ). For further studies of multicomponent dark matter with varying lifetimes, see Refs. Dienes:2011ja ; Dienes:2018yoq . A detailed exploration of the phenomenology of our framework will be presented in upcoming work future .

VII A toy model

Having described the general mechanism, we now present a simple model realization. Consider the Lagrangian

[TABLE]

where $\chi_{i}$ ($i=1,\cdots,N$) are left-handed Weyl fermions, $S$ is a SM singlet scalar field, $L$ is the lepton doublet, and $H$ is the SM Higgs. The model is technically natural, since it respects a $Z_{2}^{N}$ symmetry that is broken by $\delta m$. Any correction to off-diagonal masses must then depend on some power of $\delta m$. For simplicity, we take $\delta m$ to be flavor independent; relaxing the assumption does not qualitatively alter the results.

The mass matrix is given by

[TABLE]

The above structure is similar to a tight-binding model along a one-dimensional wire in a quantum mechanical system. It is well-known that the wavefunctions in such system are localized at each site PhysRev.109.1492 . Assuming ${\Delta m\equiv(m_{\rm max}-m_{\rm min})\gg\delta m}$, the localization length is IZRAILEV2012125 ; Craig:2017ppp

[TABLE]

Due to the localization, the mass eigenstate $\psi_{i}$ can be approximated as

[TABLE]

In the mass basis, the generated nearest neighbor interactions

[TABLE]

and non-nearest neighbor interactions are exponentially suppressed.

The mixing can generate the direct decay ${\psi_{1}\to H+L}$, with width

[TABLE]

Taking $\Gamma_{N}=H(m)$, we find

[TABLE]

to ensure that $\tau\gtrsim 10^{27}$ sec. The ratio between $\Delta m$ and $\Delta m$ controls both the relic abundance and the lifetime of DM in this toy model. For the longevity of DM, a small localization length is preferred, while, for the correct relic abundance, the localization should be not so small as it could suppress the nearest neighbor interaction. In Fig. 3, we plot the minimum number of heavy particles $\chi_{i}$ needed for the stability and the observed abundance of dark matter.

VIII Summary

In this Letter, we presented a new freezeout mechanism for super heavy DM that freezes out with the SM within a standard cosmological history. The relic abundance is determined solely via its interactions with the SM. For a velocity-independent cross-section, we showed the DM mass could be as large as $m\sim 10^{14}$ GeV within the perturbative unitary limit. In an upcoming paper we show how velocity-dependent interactions, such as if the scattering was mediated by a light mediator, allows for DM to be as heavy as the Planck scale future1 .

Acknowledgements.

Acknowledgments — We are grateful to Tim Cohen, Raffaele D’Agnolo, Jared Evans, Yuval Grossman, Kenny C.Y. Ng, Jinhong Park, Josh Ruderman, Juri Smirnov, and the Weizmann Institute HEP lunch group for useful discussions. We especially thank Yonit Hochberg for useful discussions and comments on the manuscript. The work of EK is supported by the Israel Science Foundation (grant No.1111/17), by the Binational Science Foundation (grant No. 2016153) and by the I-CORE Program of the Planning Budgeting Committee (grant No. 1937/12).