GEORCE: A Fast New Control Algorithm for Computing Geodesics

Computing geodesics for Riemannian manifolds is a difficult task that often relies on numerical approximations. However, these approximations tend to be either numerically unstable, have slow convergence, or scale poorly with manifold dimension and number of grid points. We introduce a new algorithm called GEORCE that computes geodesics in a local chart via a transformation into a discrete control problem. We show that GEORCE has global convergence and quadratic local convergence. In addition, we show that it extends to Finsler manifolds. For both Finslerian and Riemannian manifolds, we thoroughly benchmark GEORCE against several alternative optimization algorithms and show empirically that it has a much faster and more accurate performance for a variety of manifolds, including key manifolds from information theory and manifolds that are learned using generative models.

💡 Research Summary

The paper introduces GEORCE (Geodesic Optimization using Control Techniques), a novel algorithm for computing geodesics on both Riemannian and Finsler manifolds. Traditional approaches—gradient‑based minimization of the continuous energy functional, Newton or quasi‑Newton methods, and boundary‑value ODE solvers—suffer from instability, slow convergence, or poor scalability with dimension and discretization size. GEORCE circumvents these issues by reformulating the discretized energy minimization problem as a discrete optimal control problem.

Starting from the continuous energy (E(\gamma)=\frac12\int_0^1 \dot\gamma^\top G(\gamma)\dot\gamma,dt), the authors discretize the curve into (T+1) points and approximate the velocity by forward differences, yielding the discrete functional (E(x_{0:T})=\sum_{t=0}^{T-1}(x_{t+1}-x_t)^\top G(x_t)(x_{t+1}-x_t)). By defining the control variables (u_t = x_{t+1}-x_t), the problem becomes: minimize (\sum_{t=0}^{T-1} u_t^\top G(x_t) u_t) subject to the linear dynamics (x_{t+1}=x_t+u_t) and fixed endpoints.

Applying Pontryagin’s maximum principle to this discrete‑time control problem yields necessary optimality conditions: (2G(x_t)u_t+\mu_t=0) and a backward recursion for the costate (\mu_t) involving the gradient (\nu_t=\nabla_y\big(u_t^\top G(y)u_t\big)|{y=x_t}). Although these equations are nonlinear in general, the authors propose an iterative scheme that freezes (G) and (\nu) at the current iterate, turning the system into a linear one that can be solved analytically. Proposition 2 gives closed‑form updates for (\mu{T-1}), each control (u_t), and the next state (x_{t+1}).

The algorithm (Algorithm 1) proceeds as follows: initialize the curve by linear interpolation, then repeat – compute the local metric matrices (G_t) and gradients (\nu_t); compute provisional controls (\tilde u_t) using the closed‑form formulas; perform an exact line search to find the optimal step size (\alpha^*) that minimizes the discrete energy along the direction (\tilde u_t); update controls and states; check convergence of the energy gradient. The line search guarantees a monotone decrease of the energy, establishing global convergence. Near a solution, the update equations coincide with a Newton step, yielding quadratic local convergence, which the authors prove rigorously.

Extensive experiments benchmark GEORCE against ADAM, BFGS, Newton, and several boundary‑value ODE solvers (RK45, RK23, DOP853, Radau, BDF, LSODA) on a variety of manifolds: Fisher‑Rao information manifolds, tori, spheres, and high‑dimensional spheres (S^n) with (n=2,3,5,10,20). Metrics include final geodesic length, number of iterations, runtime, and memory consumption. GEORCE consistently reaches comparable or lower energy values in far fewer iterations and dramatically lower wall‑clock time—often an order of magnitude faster. In high dimensions, ODE solvers either exceed the 24‑hour time limit or exhaust 10 GB of RAM, while GEORCE solves the problem in seconds to minutes.

The method’s generality stems from the fact that only the positive‑definite metric matrix (G(x)) is required; no Christoffel symbols or higher‑order derivatives are needed, making it directly applicable to Finsler manifolds where the metric may be direction‑dependent. The algorithm operates within a single chart; if iterates approach the chart boundary, a soft line search or chart transition can be employed, but the authors report that this is rarely needed in practice.

Limitations include the assumption that the entire geodesic lies within one coordinate chart and the computational cost of the line search, which scales with the number of discretization points. Future work suggested by the authors involves automatic chart switching, adaptive discretization, integration with learned metric tensors from data‑driven models, and extensions to non‑smooth or discrete manifolds.

In summary, GEORCE provides a theoretically sound, globally convergent, and quadratically fast algorithm for geodesic computation that outperforms existing optimization and ODE‑based methods, especially in high‑dimensional settings, and opens new avenues for efficient geometric processing in many scientific and engineering applications.