Convergence Rates for Discretized Monge–Ampère Equations and Quantitative Stability of Optimal Transport

Foundations of Computational Mathematics https://doi.org/10.1007/s10208-020-09480-x Convergence Rates for Discretized Monge–Ampère Equations and Quantitative Stability of Optimal Transport Robert J. Berman1 Received: 29 June 2018 / Revised: 11 September 2020 / Accepted: 11 September 2020 © The Author(s) 2020 Abstract In recent works—both experimental and theoretical—it has been shown how to use computational geometry to efficiently construct approximations to the optimal transport map between two given probability measures on Euclidean space, by discretizing one of the measures. Here we provide a quantitative convergence analysis for the solutions of the corresponding discretized Monge–Ampère equations. This yields H 1 -converge rates, in terms of the corresponding spatial resolution h, of the discrete approximations of the optimal transport map, when the source measure is discretized and the target measure has bounded convex support. Periodic variants of the results are also established. The proofs are based on new quantitative stability results for optimal transport maps, shown using complex geometry. Keywords Monge–Ampère equations · Optimal transport · Numerical analysis · Complex differential geometry Mathematics Subject Classification 35J60 · 90C08 · 65N99 · 53C56 1 Introduction The theory of optimal transport [47], which was originally motivated by applications to logistics and economics, has generated a multitude of applications ranging from meteorology and cosmology to image processing and computer graphics in more recent years [44,45]. This has led to a rapidly expanding literature on numerical methods to construct optimal transport maps, using an appropriate discretization scheme. From the PDE point of view, this amounts to studying discretizations of the second boundary Communicated by Hans Munthe-Kaas. B Robert J. Berman 1 Mathematical Sciences, Chalmers University of Technology and the University of Gothenburg, 412 96 Göteborg, Sweden 123 Foundations of Computational Mathematics value problem for the Monge–Ampère operator. The present paper is concerned with a particular discretization scheme, known as semi-discrete optimal transport in the optimal transport literature (see [6] and references therein for other discretization schemes, based on finite differences). This approach uses computational geometry to compute a solution to the corresponding discretized Monge–Ampère equation and exhibits remarkable numerical performance, using a damped Newton iteration [34,38]. The convergence of the iteration toward the discrete solution was recently settled in [31], and one of the main aims of the present paper is to establish quantitative convergence rates of the discrete solutions, as the spatial resolution h tends to zero. 1.1 Background Throughout the paper, we fix open bounded domains X and Y in Rn with Y assumed convex and a probability measure ν on Y with a density which is uniformly bounded from below: ν = 1Y g(y)dy, g ∈ L 1 (Rn ), δ := inf g > 0 Y (1.1) We recall that in the case when μ is a probability measure on X which is also absolutely continuous with respect to d x, i.e., μ ∈ Pac (X ), then a map T in L ∞ (X , Y ) is said to be a transport map (with source μ and target ν) if T∗ μ = ν and T is said to be an optimal transport map (with respect to the Euclidean cost function |x − y|2 ), denoted by Tμν , if it realizes the infimum defining the Wasserstein L 2 -distance W2 (μ, ν) [47]:  W2 (μ, ν) := inf T |x − T x|2 μ, X where the infimum ranges over all transport maps. By Brenier’s theorem [14], there exists a unique optimal transport map Tμν and it has the characteristic property of being a gradient map: Tμν = ∇φ (in the almost everywhere sense) for a convex function φ on X , called the potential of Tμν (see [47] for further background on optimal transport theory). The potential φ is the unique (modulo an additive constant) convex solution to the corresponding second boundary value problem for the Monge–Ampère operator: The sub-gradient image of φ is contained in the closure of Y , (∂φ)(X ) ⊆ Ȳ 123 (1.2) Foundations of Computational Mathematics and φ solves the equation M A g (φ) = μ, (1.3) where the Monge–Ampère measure M A g is the probability measure on X defined by M A g (φ) := g(∇φ) det(∇ 2 φ)d x (1.4) when φ is C 2 -smooth and the general definition, due to Alexandrov, is recalled in Sect. 2.1. We will say that (X , Y , μ, ν) is regular if the corresponding solution φ is in C 2 ( X¯ ). By Cafferelli’s regularity results [17–19], this is the case if X and Y are assumed strictly convex with C 2 -boundary and the densities of μ and ν are Hölder continuous and strictly positive on X̄ and Ȳ , respectively. Then, φ defines a classical solution of the corresponding PDE, and the corresponding optimal transport map ∇φ yields a diffeomorphism between the closures of X and Y . In fact, as is well known, for any probability measure μ, there exists a solution (in the weak sense of Alexandrov) of the corresponding second boundary value problem, which is uniquely determined up to normalization (Lemma 2.2). In the sequel, it will be convenient to use the normalization condition that the integral of a solution over (X , d x) vanishes. 1.1.1 Discretization Using Semi-discrete Transport A time-honored approach for discretizing Monge–Ampère equations, which goes back to the classical work of Alexandrov on Minkowski type problems for convex bodies and polyhedra, amounts to replacing the given probability measure μ with a sequence of discrete measures converging weakly toward μ (see [1, Thm 7.3.2 and Section 7.6.2] and [5, Section 17]). A standard way to obtain such a sequence is to first discretize X by fixing a sequence of “point clouds” (x1 , . . . , x N ) ∈ X N and a dual tessellation N . This means that the union of C cover X , x ∈ C and of X with N cells (Ci )i=1 i i i the intersection of different cells have zero Lebesgue measure. For example, given a point cloud the corresponding Voronoi tessellation of X provides a canonical dual tessellation of X . The “spatial resolution” of the discretization is quantified by h := max diam(Ci ), i≤N where diam(Ci ), denotes the diameter of the cell Ci . The corresponding discretization of the measure μ is then defined by setting μh := N  f i δxi , f i := μ(Ci ) (1.5) i=1 where we have used the subindex h to emphasize that we are focusing on the limit when h → 0 (see also Sect. 3.4 for other discretizations). This discretization scheme corresponds, from the point of view of optimal transport, to the notion of semi-discrete optimal transport (since it corresponds to optimal transport between the “continuous” 123 Foundations of Computational Mathematics measure ν and the discrete measure μh ; see [31,34] and Sect. 5.2). From the point of view of numerics, this kind of discretization scheme was first introduced in the different setting of the Dirichlet problem in [43]. 1.2 Convergence R (...truncated)