Solution (source code)

= Solution

For fixed endpoint values, independent variations of $u$ and $v$ give the <Euler-Lagrange equation>[Euler-Lagrange equations]
$$
\boxed{\frac d{dt}F_{\dot u}-F_u=0,\qquad
\frac d{dt}F_{\dot v}-F_v=0.}
$$
The <sphere> parametrization has metric $ds^2=a^2(d\theta^2+\sin^2\theta\,d\phi^2)$, obtained by differentiating its position vector and taking its Euclidean squared length. The displayed functional is therefore one half the integral of squared speed, the <geodesic energy>. Its fixed-endpoint extremals are affinely parametrized geodesics; energy stationarity chooses a constant-speed parametrization.

Its explicit equations, after dividing by $a^2$, are
$$
\boxed{\ddot\theta-\sin\theta\cos\theta\,\dot\phi^2=0,
\qquad\frac d{dt}(\sin^2\theta\,\dot\phi)=0.}
$$
The second gives $\sin^2\theta\,\dot\phi=h$. Since the Lagrangian has no explicit time dependence, its energy $\dot q\cdot L_{\dot q}-L=L$ is conserved. Thus
$$
\boxed{\sin^2\theta\,\dot\phi=h,\qquad
\dot\theta^2+\sin^2\theta\,\dot\phi^2=E.}
$$
The same second identity follows by differentiating it and substituting the two equations. These polar-coordinate formulas hold on a nonsingular coordinate patch; the <geodesic> curve itself can pass through a pole in another chart.