Convergence of spectral discretization for the flow of diffeomorphisms

The Large Deformation Diffeomorphic Metric Mapping (LDDMM) or flow of diffeomorphism is a classical framework in the field of shape spaces and is widely applied in mathematical imaging and computational anatomy. Essentially, it equips a group of diffeomorphisms with a right-invariant Riemannian metric, which allows to compute (Riemannian) distances or interpolations between different deformations. The associated Euler–Lagrange equation of shortest interpolation paths is one of the standard examples of a partial differential equation that can be approached with Lie group theory (by interpreting it as a geodesic ordinary differential equation on the Lie group of diffeomorphisms). The particular group $\mathcal D^m$ of Sobolev diffeomorphisms is by now sufficiently understood to allow the analysis of geodesics and their numerical approximation. We prove convergence of a widely used Fourier-type space discretization of the geodesic equation. It is based on a regularity estimate, for which we also provide a new proof: Geodesics in $\mathcal D^m$ preserve any higher order Sobolev regularity of their initial velocity.

💡 Research Summary

**
The paper investigates the convergence of a widely used Fourier‑type spectral discretization for the geodesic equation governing the flow of diffeomorphisms, a central model in Large Deformation Diffeomorphic Metric Mapping (LDDMM). The authors work on the Sobolev diffeomorphism group (\mathcal D^{m}) on the flat torus (T^{d}), where (m>d/2+1). This group is equipped with a right‑invariant Riemannian metric induced by a positive‑definite Fourier multiplier (L) (for instance (L=(I-\Delta)^{m})). The associated Euler–Poincaré (EPDiff) equation describes the evolution of a velocity field (v(t,\cdot)) and a density (\rho(t,\cdot)) that minimise the kinetic energy (\int_{0}^{1}\langle Lv_t, v_t\rangle_{L^{2}},dt).

A central theoretical contribution is a new regularity result (Theorem 2) stating that if the initial velocity belongs to the Sobolev space (H^{m+k}) with (k\ge1), then the whole geodesic trajectory retains this higher regularity for all times. This “no‑loss‑no‑gain” property improves earlier results that required (k) to be at least (m-d/2-1). The proof exploits the right‑invariance of the metric: the Riemannian exponential map is smooth in the velocity direction, and its differential is bounded only along that direction. Consequently, the metric need not be as smooth as previously assumed, and the regularity of the solution follows directly from the smoothness of the exponential map.

Using this regularity theorem, the authors analyse the spectral discretization introduced by Zeng and Fletcher (2019). In that scheme the velocity and density fields are approximated by band‑limited Fourier series with a cut‑off frequency (R). Non‑linear terms are truncated to the same band‑limit, yielding a finite‑dimensional ODE system that can be integrated efficiently. The paper proves that, provided the initial data are sufficiently smooth (i.e., belong to (H^{m+k}) with (k\ge1)), the numerical solution converges to the exact solution as (R\to\infty). Moreover, for (k\ge2) the error in the (H^{m}) norm decays like (R^{1-k}). This rate matches the intuition that each additional derivative of regularity improves the convergence by one power of (R). Numerical experiments reported in the literature even suggest a slightly better rate, indicating that the theoretical bound may be conservative.

The convergence proof relies on Grönwall’s inequality combined with the regularity estimate: the difference between the exact and discrete solutions satisfies an ODE whose right‑hand side can be bounded by a constant times the current error plus a term of order (R^{1-k}). Integration yields the desired error bound. The authors also discuss why a fully structure‑preserving discretisation (i.e., one that respects the Lie‑algebra structure of the infinite‑dimensional group) cannot be achieved with a simple band‑limited subspace: the Lie bracket of two band‑limited vector fields generally produces frequencies beyond the cut‑off, breaking closure. Consequently, the discretisation is only an approximation of the true Lie algebra, not an exact subalgebra.

Finally, the paper remarks on extensions to other smoothing operators such as Gaussian convolution. For such operators the underlying diffeomorphism group is not yet known to be a Banach manifold, and the right‑invariant metric may fail to be strong enough to guarantee well‑posed geodesics. Hence the current analysis does not directly apply, and establishing the geometric foundations for these alternative models remains an open problem.

In summary, the work provides (i) a strengthened regularity theorem for Sobolev diffeomorphism geodesics, (ii) a rigorous convergence analysis of a practical Fourier‑spectral scheme, and (iii) a clear discussion of the limitations of structure‑preserving discretisations in this infinite‑dimensional setting. The results give a solid theoretical foundation for the widespread use of spectral methods in computational anatomy and related imaging applications.