  We propose an efficient numerical algorithm for the solution of diffeomorphic
image registration problems. We use a variational formulation constrained by a
partial differential equation (PDE), where the constraints are a scalar
transport equation.
  We use a pseudospectral discretization in space and second-order accurate
semi-Lagrangian time stepping scheme for the transport equations. We solve for
a stationary velocity field using a preconditioned, globalized, matrix-free
Newton-Krylov scheme. We propose and test a two-level Hessian preconditioner.
We consider two strategies for inverting the preconditioner on the coarse grid:
a nested preconditioned conjugate gradient method (exact solve) and a nested
Chebyshev iterative method (inexact solve) with a fixed number of iterations.
  We test the performance of our solver in different synthetic and real-world
two-dimensional application scenarios. We study grid convergence and
computational efficiency of our new scheme. We compare the performance of our
solver against our initial implementation that uses the same spatial
discretization but a standard, explicit, second-order Runge-Kutta scheme for
the numerical time integration of the transport equations and a single-level
preconditioner. Our improved scheme delivers significant speedups over our
original implementation. As a highlight, we observe a 20$\times$ speedup for a
two dimensional, real world multi-subject medical image registration problem.
