The geodesic equation on the sphere reduces to the coupled differential equations \begin{align} \dot\phi &= \frac{\ell}{r^2\sin^2\theta} \\ r^2 \dot\theta^2 &= 1 - \frac{\ell^2}{r^2\sin^2\theta} \end{align} where the constant $\ell$ is the angular momentum about the $z$-axis (per unit mass), and where we have used arclength (here labeled $\lambda$) as the parameter (so that the curve is traversed at unit speed). We now proceed to solve these equations.

To solve these equations, we start with \begin{equation} r \dot\theta = \pm\frac{\sqrt{r^2\sin^2\theta-\ell^2}}{r\sin\theta} \end{equation} and then write \begin{equation} \frac{d\phi}{d\theta} = \frac{\dot\phi}{\dot\theta} = \frac{\pm \ell}{\sin\theta\sqrt{r^2\sin^2\theta-\ell^2}} \end{equation} so that \begin{align} d\phi &= \frac{\pm \ell\,d\theta}{\sin\theta\sqrt{r^2\sin^2\theta-\ell^2}} \nonumber\\ &= \frac{\pm \ell\,\csc^2\theta\,d\theta} {\sqrt{r^2-\ell^2\frac{\cos^2\theta+\sin^2\theta}{\sin^2\theta}}} \nonumber\\ &= \frac{\pm \ell\,\csc^2\theta\,d\theta}{\sqrt{r^2-\ell^2-\ell^2\cot^2\theta}} \nonumber\\ &= \frac{\mp \ell\,d\cot\theta}{\sqrt{(r^2-\ell^2)-\ell^2\cot^2\theta}} \end{align} This expression can now be integrated to yield \begin{equation} \phi = \mp\arcsin\left( \frac{\ell}{\sqrt{r^2-\ell^2}} \right) + \hbox{constant} \end{equation} or equivalently \begin{equation} \sin(\phi-\alpha) = \mp\sqrt{\frac{\ell^2}{r^2-\ell^2}}\cot\theta \label{geosolsph} \end{equation}

Just as before, we can bring this expression to a more familiar form by expanding the left-hand side, resulting in \begin{equation} \sin\theta(\sin\phi\cos\alpha - \cos\phi\sin\alpha) \pm\sqrt{\frac{\ell^2}{r^2-\ell^2}}\cos\theta = 0 \end{equation} which has the form of the equation of a plane through the origin, since $x=r\sin\theta\cos\phi$, $y=r\sin\theta\sin\phi$, and $z=r\cos\theta$. Not surprisingly, geodesics on the sphere lie on such planes; they are great circles.

However, ($\ref{geosolsph}$) is not a very practical description of geodesics on a sphere. Working in three dimensions, there is an easier way.

Any point ($\theta$,$\phi$) on a sphere is of course also a point ($x$,$y$,$z$) in three dimensions, or equivalently a vector from the origin. The length of the circular arc connecting these points is just the radius of the sphere times the angle — which is easily obtained using the dot product.

Figure 1: Using 3-dimensional vectors to construct geodesics on the sphere.

This strategy can also be used to find an explicit expression for the geodesic connecting two given points on the sphere. We can easily find two orthogonal vectors whose endpoints are on the geodesic, either using projections (essentially Gram-Schmidt orthogonalization), or by using the cross product twice. (The first time to obtain a vector orthogonal to the plane of the geodesic, then again to obtain a second vector in that plane orthogonal to one of the given vectors.) But any circle can be parameterized in vector form as $\uu\cos\beta+\ww\sin\beta$ using any two orthogonal vectors $\uu$ and $\ww$ on it. This construction is illustrated in Figure 1, in which the geodesic through two points $\uu$ and $\vv$ on the unit sphere is constructed using $\nn=\frac{\uu\times\vv}{|\uu\times\vv|}$ and $\ww=\nn\times\uu$.

Although this construction of the geodesic equation is quite messy to use in calculations by hand, it lends itself well to programmed numerical computations.