Part IV: Dark Matter | Chapter 9

WIMPs: Freeze-out and the Boltzmann Equation

Deriving the thermal relic abundance from the Boltzmann equation, the WIMP miracle, and strategies for direct detection

Overview

A weakly interacting massive particle (WIMP) initially in thermal equilibrium with the primordial plasma β€œfreezes out” when the expansion rate overtakes the annihilation rate. The resulting relic abundance depends on a single quantity β€” the thermally averaged annihilation cross section. Remarkably, a weak-scale cross section yields precisely the observed dark matter density, a coincidence known as the WIMP miracle.

1. Comoving Number Density and Yield

The number density $n$ of a species in an expanding universe is diluted by$a^{-3}$. We define the yield as the ratio of number density to entropy density:

$$Y \equiv \frac{n}{s}\,,\qquad s = \frac{2\pi^2}{45}g_{*S}\,T^3\,,$$

where $g_{*S}$ counts the effective entropic degrees of freedom. Since entropy is conserved in a comoving volume ($sa^3 = \text{const}$), the yield $Y$tracks the true comoving abundance. The equilibrium yield for a non-relativistic species is

$$Y_{\rm eq}(x) = \frac{45}{4\pi^4}\frac{g}{g_{*S}}\,x^{3/2}\,e^{-x}\,,\qquad x \equiv \frac{m}{T}\,.$$

2. The Boltzmann Equation

The evolution of the WIMP number density is governed by the Boltzmann equation, which in an FRW background reads

$$\frac{dn}{dt} + 3Hn = -\langle\sigma v\rangle\left(n^2 - n_{\rm eq}^2\right)\,.$$

Switching to the variables $Y$ and $x = m/T$, and using$dt/dx = 1/(Hx)$ during radiation domination, this becomes

$$\frac{dY}{dx} = -\frac{\langle\sigma v\rangle\,s}{H\,x}\left(Y^2 - Y_{\rm eq}^2\right)\,.$$

This is a Riccati equation that cannot be solved in closed form. At early times ($x \ll x_f$), annihilations are rapid and $Y \approx Y_{\rm eq}$. After freeze-out ($x > x_f$), the equilibrium yield is exponentially suppressed and the equation simplifies to

$$\frac{dY}{dx} \approx -\frac{\langle\sigma v\rangle\,s}{H\,x}\,Y^2\,.$$

3. Freeze-out Temperature

Freeze-out occurs when the annihilation rate drops below the expansion rate,$n_{\rm eq}\langle\sigma v\rangle \approx H$. Solving iteratively gives

$$x_f = \ln\left[\frac{0.038\,g\,M_{\rm Pl}\,m\,\langle\sigma v\rangle}{\sqrt{g_*\,x_f}}\right] \approx 20\text{--}25\,,$$

remarkably insensitive to the WIMP mass and cross section (logarithmic dependence). For a 100 GeV WIMP, freeze-out occurs at $T_f \sim 5$ GeV.

4. Relic Density and the WIMP Miracle

Integrating the late-time Boltzmann equation from $x_f$ to $\infty$and converting to the density parameter:

$$\Omega_{\rm DM}h^2 \approx \frac{1.07\times10^9\;\text{GeV}^{-1}}{M_{\rm Pl}}\frac{x_f}{\sqrt{g_*(x_f)}}\frac{1}{\langle\sigma v\rangle}\,.$$

Numerically, this evaluates to

$$\Omega_{\rm DM}h^2 \approx 0.12\times\frac{3\times10^{-26}\;\text{cm}^3\text{s}^{-1}}{\langle\sigma v\rangle}\,.$$

A weak-scale annihilation cross section ($\langle\sigma v\rangle \sim \alpha_W^2/m_W^2 \sim 3\times10^{-26}$ cm$^3$s$^{-1}$) reproduces the measured abundance β€” the β€œWIMP miracle.” This numerical coincidence motivates extensive experimental searches and connects dark matter to the electroweak scale.

5. Velocity Expansion

The thermally averaged cross section is expanded in partial waves:

$$\langle\sigma v\rangle = a + b\langle v^2\rangle + \mathcal{O}(v^4) = a + \frac{6b}{x} + \cdots$$

For s-wave annihilation ($a \neq 0$), the cross section is velocity-independent and indirect detection constraints from the present epoch ($v/c \sim 10^{-3}$) directly probe the freeze-out cross section. For p-wave dominated processes ($a = 0$), present-day annihilation is suppressed by $v^2 \sim 10^{-6}$, evading indirect limits.

6. Direct Detection

A WIMP scattering off a nucleus deposits recoil energy $E_R$. The differential event rate per unit detector mass is

$$\frac{dR}{dE_R} = \frac{\rho_0}{m_\chi\,m_N}\int_{v_{\min}}^{v_{\rm esc}} v\,f(\mathbf{v})\,\frac{d\sigma}{dE_R}\,d^3v\,,$$

where $\rho_0 \approx 0.3$ GeV/cm$^3$ is the local DM density,$f(\mathbf{v})$ is the WIMP velocity distribution, and the minimum velocity to produce recoil $E_R$ is

$$v_{\min} = \sqrt{\frac{m_N E_R}{2\mu^2}}\,,\qquad \mu = \frac{m_\chi m_N}{m_\chi + m_N}\,.$$

The spin-independent (SI) differential cross section includes nuclear coherence through the Helm form factor:

$$\frac{d\sigma_{\rm SI}}{dE_R} = \frac{m_N\,A^2\,\sigma_0^{\rm SI}}{2\mu^2 v^2}\,|F(E_R)|^2\,,$$

where $F(q) = 3j_1(qR_1)/(qR_1)\,e^{-q^2s^2/2}$ and $A$ is the atomic mass number providing the $A^2$ coherent enhancement.

7. Indirect Detection

WIMPs annihilating in regions of high dark matter density produce Standard Model final states. The differential photon flux from a solid angle $\Delta\Omega$ is

$$\frac{d\Phi_\gamma}{dE} = \frac{\langle\sigma v\rangle}{8\pi\,m_\chi^2}\sum_f B_f\frac{dN_\gamma^f}{dE}\times\underbrace{\int_{\Delta\Omega}\int_{\rm l.o.s.}\rho^2(\mathbf{r})\,d\ell\,d\Omega}_{J\text{-factor}}\,,$$

where $B_f$ is the branching ratio into channel $f$ and the$J$-factor encodes the astrophysical target. The Galactic Center, dwarf spheroidal galaxies, and galaxy clusters are prime targets. Fermi-LAT observations of dwarf spheroidals constrain $\langle\sigma v\rangle < 3\times10^{-26}$cm$^3$s$^{-1}$ for $m_\chi \lesssim 100$ GeV annihilating to $b\bar{b}$, directly probing the thermal relic cross section.

8. Collider Searches

At the LHC, WIMPs are pair-produced via EFT contact operators or simplified models with a mediator particle. The signature is missing transverse energy recoiling against an initial-state radiation jet (monojet):

$$pp \to \chi\bar\chi + j\,,\qquad \sigma \sim \frac{g_q^2\,g_\chi^2}{M_{\rm med}^4}\,s\,,$$

in the contact-operator limit ($M_{\rm med} \gg \sqrt{s}$). For lighter mediators, the cross section is resonantly enhanced when $\sqrt{\hat{s}} \approx M_{\rm med}$. LHC Run 2 monojet searches exclude mediator masses up to $\sim 2$ TeV for$g_q = g_\chi = 1$, providing complementary coverage to direct detection in the low-mass ($m_\chi < 10$ GeV) regime.

9. Spin-Dependent Interactions

While SI scattering benefits from $A^2$ coherent enhancement, spin-dependent (SD) scattering couples to the nuclear spin and lacks this enhancement:

$$\sigma_{\rm SD} = \frac{32}{\pi}G_F^2\,\mu^2\,\frac{J+1}{J}\left[a_p\langle S_p\rangle + a_n\langle S_n\rangle\right]^2\,,$$

where $\langle S_{p,n}\rangle$ are the expectation values of the proton/neutron spin content and $a_{p,n}$ are the effective WIMP-nucleon couplings. The best SD limits on neutron coupling come from xenon experiments, while proton coupling is most sensitively probed by IceCube via neutrinos from WIMP annihilation in the Sun.

10. Current Experimental Limits

The most stringent SI limits come from dual-phase xenon time projection chambers:

  • LZ (2022): $\sigma_{\rm SI} < 9.2\times10^{-48}$ cm$^2$ at $m_\chi = 36$ GeV
  • XENONnT (2023): $\sigma_{\rm SI} < 2.6\times10^{-47}$ cm$^2$ at $m_\chi = 28$ GeV
  • Neutrino fog: Below $\sigma \sim 10^{-48}$ cm$^2$, coherent elastic neutrino-nucleus scattering ($\text{CE}\nu\text{NS}$) from solar and atmospheric neutrinos becomes an irreducible background

Next-generation experiments (DARWIN, XLZD) aim to reach the neutrino fog across a broad WIMP mass range, definitively testing the thermal WIMP paradigm above $\sim 5$ GeV.

11. Coannihilation and Resonances

When a particle $\chi'$ is nearly degenerate with the WIMP ($\Delta m/m \lesssim T_f/m \sim 5\%$), coannihilation processes$\chi\chi' \to \text{SM}$ contribute to the effective cross section:

$$\sigma_{\rm eff} = \sum_{ij}\sigma_{ij}\frac{g_i g_j}{g_{\rm eff}^2}(1+\Delta_i)^{3/2}(1+\Delta_j)^{3/2}e^{-x(\Delta_i+\Delta_j)}\,,$$

where $\Delta_i = (m_i - m_\chi)/m_\chi$. In the MSSM, stau coannihilation ($\tilde\chi^0_1\tilde\tau \to \tau\gamma$) is a classic example. Near a resonance ($2m_\chi \approx m_A$ for a pseudoscalar Higgs), the cross section is enhanced by the Breit-Wigner factor, yielding the correct relic density for much heavier WIMPs than the naive thermal estimate suggests.

Key Equations

$$\frac{dY}{dx} = -\frac{\langle\sigma v\rangle s}{Hx}(Y^2-Y_{\rm eq}^2)\,,\quad \Omega_{\rm DM}h^2 \approx 0.12\times\frac{3\times10^{-26}\text{cm}^3\text{s}^{-1}}{\langle\sigma v\rangle}$$

$$\frac{dR}{dE_R} = \frac{\rho_0}{m_\chi m_N}\int v\,f(\mathbf{v})\frac{d\sigma}{dE_R}d^3v\,,\quad x_f \approx 20\text{--}25$$

Rate this chapter: