= Solution
The energy density from a <phase-space distribution function> with $g=2$ internal states is
$$
\rho_\nu=\frac{g}{(2\pi)^3}
\int_{\mathbb R^3}\frac{\sqrt{p^2+m_\nu^2}}{e^{p/T_\nu}+1}\,d^3p
=\frac{g}{2\pi^2}\int_0^\infty
\frac{p^2\sqrt{p^2+m_\nu^2}}{e^{p/T_\nu}+1}\,dp.
$$
In the ultrarelativistic limit $T_\nu\gg m_\nu$, set $x=p/T_\nu$. The <Fermi-Dirac distribution> then gives
$$
\rho_\nu=\frac{T_\nu^4}{\pi^2}
\int_0^\infty\frac{x^3}{e^x+1}\,dx
=\frac{T_\nu^4}{\pi^2}
3!\left(1-2^{-3}\right)\zeta(4).
$$
Since $\zeta(4)=\pi^4/90$,
$$
\boxed{\rho_\nu=\frac{7\pi^2}{120}T_\nu^4}.
$$
Back to article page