Solution (source code)

= Solution

At large radius $V_c\to V_0$, so the spherical enclosed <mass> is $M(<r)\simeq V_0^2r/G$. A spherical <mass density> satisfies $dM/dr=4\pi r^2\rho(r)$, and hence
$$
\boxed{\rho(r)\simeq\frac{V_0^2}{4\pi Gr^2}\propto r^{-2}.}
$$
Equivalently, a local power law $\rho\propto r^{-p}$ with enclosed <mass> dominated by that extended profile has $M\propto r^{3-p}$ and $V_c^2\propto r^{2-p}$. A flat <galaxy rotation curve> requires $p=2$.

Differentiating the exact <mass> profile verifies the limit and the finite central <density>:
$$
\rho(r)=\frac{V_0^2}{4\pi Gr^2}\left[(1-e^{-x})^2+2x(1-e^{-x})e^{-x}\right],\qquad x=r/r_0.
$$
As $x\to0$, the bracket is $3x^2+O(x^3)$, so $\rho(0)=3V_0^2/(4\pi Gr_0^2)$. As $x\to\infty$, it tends to one, recovering the $r^{-2}$ tail. Extending the flat <galaxy rotation curve> to infinity makes total <mass> diverge linearly; a real finite halo must eventually depart from this idealization.