Solution (source code)

= Solution

For constant density,
$$
m(r)=\frac{4\pi\rho_0r^3}{3}=M\frac{r^3}{R^3},
\qquad M=\frac{4\pi\rho_0R^3}{3}.
$$
Substitution in the pressure equation makes it separable. Integrating from $r$ to the surface, where $p(R)=0$, gives the <Interior Schwarzschild solution>
$$
\boxed{
p(r)=\rho_0
\frac{\sqrt{1-2Mr^2/R^3}-\sqrt{1-2M/R}}
{3\sqrt{1-2M/R}-\sqrt{1-2Mr^2/R^3}}}.
$$
Direct differentiation verifies the TOV equation and the surface boundary condition.