Solution (source code)

= Solution

Use $E=mc^2+p^2/(2m)+\cdots$ and replace the Fermi denominator by its Boltzmann form. The Gaussian <integral> gives
$$
\boxed{n_X=\left(\frac{2\pi m_XkT}{h^2}\right)^{3/2}e^{-m_Xc^2/kT}.}
$$