Solution (source code)

= Solution

<Spatial homogeneity> makes $\rho$ independent of position and makes the <two-particle number density> depend only on the separation. <Isotropy> removes its directional dependence. Consequently
$$
n_2(\mathbf x,\mathbf x_\star,t)=\rho^2g_2(|\mathbf x_\star-\mathbf x|,t).
$$
Here $g_2$ is a <radial distribution function>, with the finite-$N$ normalization inherited from <1/a>; statistical independence is unnecessary. In particular, the translated pair in <1/c> has separation $\mathbf R$ for every $s$, so its $s$ integral is simply $\rho^2g_2(R,t)$.

<Isotropy> of the random <velocity> $\mathbf c=\mathbf v-\mathbf u$ makes its second <velocity moment> proportional to the identity. Taking its trace fixes the proportionality constant:
$$
\boxed{\theta=\frac1{3\rho}\int|\mathbf v-\mathbf u|^2f\,d^3v,\qquad P_{ij}=\rho\theta\delta_{ij}.}
$$
Thus $\theta$ is the one-component random <velocity> variance; for unit <mass> in <thermal equilibrium>, $\theta=k_BT$. A <Maxwell-Boltzmann velocity distribution> is not required for the <pressure tensor> to be isotropic.

For the interaction contribution, the angular integral is $\int d\Omega\,\widehat R_i\widehat R_j=(4\pi/3)\delta_{ij}$, by <isotropy> and its trace. With $d^3R=R^2\,dR\,d\Omega$, <1/c> becomes
$$
P^\phi_{ij}=-\frac{2\pi}{3}\rho^2\delta_{ij}\int_0^\infty R^3\phi'(R)g_2(R,t)\,dR.
$$
Combining the two <pressure tensors> gives the <virial pressure from pair correlations>:
$$
\boxed{\mathbf P+\mathbf P^\phi=p\mathbf I,\qquad p=\rho\theta-\frac{2\pi}{3}\rho^2\int_0^\infty R^3\phi'(R)g_2(R,t)\,dR.}
$$
For a repulsive <pair potential>, the outward radial <force> is $-\phi'(R)\geq0$. Since $g_2\geq0$, the second term in $p$ is \b[nonnegative], and positive when occupied separations feel nonzero repulsion. It expresses additional <momentum transfer> through interactions. The displayed integral assumes a differentiable <pair potential> and convergence; a hard core requires a contact formula or a limiting argument.