Solution (source code)

= Solution

Fix the two conserved momenta $p_\psi=\lambda$ and $p_\phi=\mu$. The remaining canonical variables $(\theta,p_\theta)$ form a two-dimensional phase space with Hamiltonian
$$
H_{\lambda,\mu}(\theta,p_\theta)
=\frac{p_\theta^2}{2I}+V_{\lambda,\mu}(\theta),
$$
where
$$
V_{\lambda,\mu}(\theta)
=\frac{(\mu-\lambda\cos\theta)^2}{2I\sin^2\theta}
+\frac{\lambda^2}{2J}+Mga\cos\theta.
$$

An upright spinning solution has $\theta=0$, $p_\theta=0$, and necessarily $\mu=\lambda$, since the vertical and body-axis angular momenta then coincide. On that momentum level,
$$
V_{\lambda,\lambda}(\theta)
=\frac{\lambda^2}{2I}\tan^2\frac\theta2
+\frac{\lambda^2}{2J}+Mga\cos\theta.
$$
Near zero,
$$
V_{\lambda,\lambda}(\theta)
=\text{constant}
+\left(\frac{\lambda^2}{8I}-\frac{Mga}{2}\right)\theta^2
+O(\theta^4).
$$
Thus the upright configuration is a strict local minimum of the reduced Hamiltonian, and hence stable, when
$$
\boxed{\lambda^2>4IMga=3IMgl}.
$$
Here $\lambda$ is the angular momentum about the symmetry axis. This proves stability for sufficiently large axial angular momentum.