Solution (source code)

= Solution

The radial and azimuthal Hamilton equations give
$$
\boxed{\dot r=\frac{p_r}{m},\qquad
\dot\phi=\frac{p_\phi}{mr^2\sin^2\theta},\qquad\dot p_\phi=0.}
$$
Also $\dot\theta=p_\theta/(mr^2)$ and $\dot p_\theta=p_\phi^2\cos\theta/(mr^2\sin^3\theta)$. Differentiating the angular-momentum square,
$$
\frac{dJ^2}{dt}=2p_\theta\dot p_\theta
-\frac{2p_\phi^2\cos\theta}{\sin^3\theta}\dot\theta=0.
$$
Therefore both $p_\phi$ and $\boxed{J^2=p_\theta^2+p_\phi^2/\sin^2\theta}$ are constants. <Rotational symmetry> conserves the full angular-momentum vector; these are its axial component and squared magnitude in spherical canonical variables.