Geodesic Interpolation

A geodesic is the straightest possible curve on a manifold; the analogue of a straight line in Euclidean space. cartan exposes geodesic interpolation as a first-class trait for computing midpoints, slerp, and Fréchet means.

Definition

Definition null (Geodesic).

A smooth curve $\gamma : [0,1] \to M$ is a geodesic if its velocity field is parallel along itself:

Equivalently, $\gamma$ locally minimises arc length. The unit-speed geodesic from $p$ to $q$ is $\gamma(t) = \exp_p(t \cdot \log_p(q))$.

The `GeodesicInterpolation` Trait

View in rustdoc →

pub trait GeodesicInterpolation: Manifold {
  /// Point at fraction t ∈ [0,1] along the geodesic from p to q.
  /// γ(0) = p, γ(1) = q, constant speed.
  fn geodesic(
      &self,
      p: &Self::Point,
      q: &Self::Point,
      t: f64,
  ) -> Option<Self::Point>;

  /// Riemannian midpoint γ(0.5).
  fn midpoint(
      &self,
      p: &Self::Point,
      q: &Self::Point,
  ) -> Option<Self::Point> {
      self.geodesic(p, q, 0.5)
  }
}

Implemented by: Sphere<N>, Euclidean<N>, SO<N>, SE<N>, SPD<N>, Grassmann<N,K>

Returns None when p and q are antipodal (log undefined at cut locus).

Boundary Values and Constant Speed

Proposition null (Properties of γ(t)).

For $\gamma(t) = \exp_p(t \cdot \log_p(q))$:

$\gamma(0) = p$, $\gamma(1) = q$
is constant (constant-speed parametrisation)

Sphere Slerp

On , geodesic interpolation is slerp (Spherical Linear intERPolation):

For $\theta \approx 0$ (nearly coincident points), use the Taylor expansion $\sin(t\theta)/\sin\theta \approx 1; t(1-t)\theta^2/6$ to avoid division by a near-zero $\sin\theta$.

Example ?: Slerp on S²; interpolating two points

use cartan::prelude::*;
use cartan::manifolds::Sphere;
 
let s2 = Sphere::<3>;
let p  = [0.0_f64, 0.0, 1.0];   // north pole
let q  = [1.0_f64, 0.0, 0.0];   // equator point
 
// Midpoint along the great circle arc.
let mid = s2.midpoint(&p, &q).unwrap();
assert!((mid.norm() - 1.0).abs() < 1e-12);          // on S²
assert!((s2.dist(&p, &mid) - s2.dist(&mid, &q)).abs() < 1e-10);  // equidistant
 
// 10-point interpolation.
let path: Vec<_> = (0..=10)
    .map(|i| s2.geodesic(&p, &q, i as f64 / 10.0).unwrap())
    .collect();

Fréchet Mean

The Fréchet mean (Riemannian centre of mass) of points minimises:

On negatively curved manifolds (SPD, hyperbolic) the Fréchet mean is unique and the gradient descent

converges globally. On positively curved manifolds (, $SO(N)$) a unique mean exists when the points are contained in a convex geodesic ball of radius $< \pi/2$.

Cut Locus and Failure

geodesic returns None only when log is undefined, which occurs exactly at the cut locus of p. For all $t \in (0, 1)$ strictly interior the geodesic is unique; only the endpoints can trigger a None. In practice:

// Safe pattern: treat None as "points are antipodal, choose fallback".
let mid = s2.midpoint(&p, &q).unwrap_or_else(|| {
    // Arbitrarily break symmetry; pick a pole orthogonal to p.
    s2.project_point(&orthogonal_to(&p))
});