= Solution
In a local orthonormal frame, a massless <neutrino> has $p^{\hat a}=E(1,\mathbf e)$ and $E=\epsilon/a$. With degeneracy $g_\nu$, kinetic theory gives the <stress-energy tensor>
$$
T^{\hat a\hat b}_\nu=
g_\nu\int\frac{d^3p}{(2\pi)^3}\,
\frac{p^{\hat a}p^{\hat b}}{E}\,f,
\qquad
\bar\rho_\nu=\frac{g_\nu}{a^4}\int
\frac{\epsilon^3d\epsilon\,d\Omega_e}{(2\pi)^3}\bar f(\epsilon).
$$
Coordinate components follow by inserting the orthonormal-frame vectors on the indices.
For the perturbation $\delta f=-\Theta\epsilon\bar f'(\epsilon)$, integration by parts gives
$$
-\int_0^\infty\epsilon^4\bar f'(\epsilon)\,d\epsilon
=4\int_0^\infty\epsilon^3\bar f(\epsilon)\,d\epsilon,
$$
assuming the boundary term vanishes. Consequently the physical trace-free spatial stress is
$$
\pi_\nu^{\hat i\hat j}
=4\bar\rho_\nu\int\frac{d\Omega_e}{4\pi}
\Theta\left(e^{\hat i}e^{\hat j}-\frac13\delta^{\hat i\hat j}\right).
$$
The paper uses the opposite signed source, $\Pi^{\hat i\hat j}=-\pi_\nu^{\hat i\hat j}$, consistently with its tensor equation having $-16\pi Ga^2\Pi_{ij}$ on the right. In that convention the <neutrino tensor anisotropic stress> is
$$
\boxed{\Pi^{\hat i\hat j}
=-4\bar\rho_\nu\int\frac{d\Omega_e}{4\pi}
\Theta\left(e^{\hat i}e^{\hat j}-\frac13\delta^{\hat i\hat j}\right)}.
$$
Back to article page
