Two Stage Algorithm for Semi-Discrete Optimal Transport on Disconnected Domains

1. Introduction

We begin by recalling the semi-discrete optimal transport problem. Let $X⊂\Rn$ , $n \geq 2$ be compact and $Y:={yi}Ni=1⊂\Rn$ a fixed collection of finite points, along with a cost function $c:X×Y→\R$

. We also fix Borel probability measures

μ, ν

with

\sptμ⊂X

and

\sptν⊂Y

. Furthermore

μ

is absolutely continuous with respect to Lebesgue measure and we denote its density by

ρ

. We set

λ = (ν (y_{1}), \dots, ν (y_{N}))

We want to find a measurable mapping $T : X \to Y$ such that $T_{#} μ (E) := μ (T^{- 1} (E)) = ν (E)$ for any measurable $E \subset Y$ , and $T$ satisfies

(1)

∫Xc(x,T(x))dμ=min\tiT#μ=ν∫Xc(x,\tiT(x))dμ.

In this paper our goal is to propose and show convergence of a two stage algorithm. The first stage is a regularized gradient descent which achieves global linear convergence with respect to the regularized problem. The second stage is a damped Newton algorithm similar to that of [KMT19] which has superlinear local convergence. The main difference between our algorithm and that of [KMT19] is ours doesn’t require $\sptμ$ to be connected. Without a connected support the algorithm in [KMT19] may not converge at all.

Furthermore, we will show that the Laguerre cells, $T^{- 1} (y_{i})$ , for the approximate optimizers constructed by our algorithm converge both in the $L^{1} (μ)$ sense and in the sense of Hausdorff distance. This result is analogous to that of [BK19]. While the $L^{1} (μ)$ convergence follows directly from the results there, the Hausdorff convergence result in [BK19] heavily relies on the connectedness of $\sptμ$ .

2. Setup

2.1. Notations and Conventions

Here we collect some notations and conventions for the remainder of the paper. We fix positive integers $N, M, n$ with $N, n \geq 2$ and a collection $Y:={yi}Ni=1⊂\Rn$ . We will write $\onevect$

to denote the vector in

\RN

whose components are all

1

. For any vector

V∈\Rk

, we will write its components as superscripts so

V^{i}

is the

i

-th component of

V

. We reserve the notation

\normV

for the

l^{2}

-Euclidean norm, i.e.

\normV=√∑ki=1\absVi2

. We use

\normV1,\normV∞

for the

l^{1}

and

l^{\infty}

Euclidean norms, i.e.

\normV1=∑ki=1\absVi

and

\normV∞=maxi∈{1,…,k}\absVi

. We use

L

to denote

n

-dimensional Lebesgue measure and

H^{k}

to denote

k

-dimensional Hausdorff measure. Also we use the notation

Λ={λ∈\RN:λi∈[0,1],N∑i=1λi=1}

for the set of all admissible weight vectors.

We will assume that the cost function $c$ satisfies the following standard conditions:

(Reg)		$c (\cdot, y_{i})$	$\in C^{2} (X), \forall i \in {1, \dots, N},$
(Twist)		$\nabla_{x} c (x, y_{i})$	$\neq \nabla_{x} c (x, y_{k}), \forall x \in X, i \neq k .$

We also assume the following condition, originally studied by Loeper in [Loe09]. The cost, $c$ , is said to satisfy Loeper’s condition if for each $i \in {1, \dots, N}$ there exists a convex set $Yi⊂\Rn$ and a $C^{2}$ diffeomorphism $\cExpi⋅:Yi→X$ such that

(QC)

∀ t∈\R, 1≤k,i≤N, {p∈Yi∣−c(\cExpip,yk)+c(\cExpip,yi)≤t} is convex.

See Remark 2.1 below for further discussion of these conditions.

We also say that a set $X⊂\Rn$ is $c$ -convex with respect to $Y$ if $\invcExpiX$ is a convex set for every $i \in {1, \dots, N}$ .

For any $ψ∈\RN$ and $i \in {1, \dots, N}$ , we define the $i$ th Laguerre cell associated to $ψ$ as the set

\Lagi(ψ):={x∈X∣c(x,yi)+ψi=minic(x,yi)+ψi}.

We also define the function $G:\Rn→\weightvectset$ by

G(ψ):=(G1(ψ),…,GN(ψ))=(μ(\Lag1(ψ)),…,μ(\LagN(ψ))),

and denote for any $ϵ \geq 0$ ,

Kϵ:={ψ∈\RN∣Gi(ψ)>ϵ, ∀i∈{1,…,N}}.

The above conditions (Reg), (Twist), (QC) are the same ones assumed in [KMT19] and [BK19]. As is also mentioned in those papers, (Reg) and (Twist) are standard conditions in optimal transport. Furthermore, (QC) holds if $Y$ is a finite set sampled from from a continuous space, and $c$ is a $C^{4}$ cost function satisfying what is known as the Ma-Trudinger-Wang condition (first introduced in a strong form in [MTW05], and in [TW09] in a weaker form).

If $μ$ is absolutely continuous with respect to Lebesgue measure (or even if $μ$ doesn’t give mass to small sets), then (Twist) implies that the Laguerre cells are pairwise $μ$ -almost disjoint. In this case the generalized Brenier’s theorem [Vil09, Theorem 10.28], tells us that for any vector $ψ∈\RN$ , the map $T_{ψ} : X \to Y$ defined by $T_{ψ} (x) = y_{i}$ whenever $x∈\Lagi(ψ)$ is a minimizer in the optimal transport problem (1), where the source measure is $μ$ and the target measure is defined by $ν = ν_{G (ψ)}$ .

We define $C∇:=supx∈X,y∈Y\norm∇xc(x,y)$ .

Throughout the rest of the paper we will we assume that $X$ is compact and $c$ -convex with respect to $Y$ and that $ρ$ , the density of $μ$ , is $α$ -Hölder continuous.

In order to obtain our convergence results we will need $\sptμ$ to be at least locally connected in a quantitative way. First we recall a definition.

A non-zero measure $μ$ on a compact set $Z⊂\Rn$ satisfies a $(q, 1)$ -Poincaré-Wirtinger inequality for some $1 \leq q \leq \infty$ if there exists a constant $\Cpw>0$ such that for any $f \in C^{1} (Z)$ ,

\normf−∫ZfdμLq(μ)≤\Cpw\norm∇fL1(μ).

For brevity, we will write this as “ $μ$ satisfies a $(q, 1)$ -PW inequality on $Z$ ”. If $μ$ is defined on a larger set then $Z$ then the phrase “ $μ$ satisfies a $(q, 1)$ -PW inequality on $Z$ ” is understood to mean that $μ (Z) > 0$ and the restriction of $μ$ to $Z$ satisfies a $(q, 1)$ -PW inequality.

We remark that some version of all of our results hold with only $(1, 1)$ -PW inequalities but $q > 1$ will give better constants.

We now define a local version of the PW inequality that will serve as our quantitative measure of the local connectivity of $\sptμ$ .

A probability measure $μ$ on a compact set $X⊂\Rn$ satisfies a local $(q, 1)$ -Poincaré-Wirtinger inequality if there are compact $X_{i} \subset X$ which are pairwise disjoint and $μ$ -almost cover $X$ , i.e. $μ (X ∖ (X_{1} \cup \dots \cup X_{M})) = 0$ such that $μ$ satisfies a $(q, 1)$ -PW inequality on each $X_{i}$ with constant $κ_{i}$ .

We define the local PW constant $κ$ to be the $q$ -th power mean of the $κ_{i}$ , i.e. $κ = (\frac{1}{M} \sum_{i = 1}^{M} κ_{i}^{q})^{1 / q}$ . Note that $κ \leq {max}_{i} κ_{i}$ .

When we have $M = 1$ we recover the global connectivity assumption of [KMT19, Definition 1.3]. Throughout the paper we use the notation $χ^{i} = μ (X_{i})$ and the vector $χ=(χ1,…,χM)∈\RM$ .

In order to obtain convergence of our algorithm we require that optimal map $T$ , “splits up” the $X_{i}$ . In order to assure this we assume that the subsets sums of $χ$ and $λ$ are separated. For any vector $V∈\Rk$ we introduce the notation $S_{V}$ for the set of subset sums of $V$ , i.e

S_{V} = {\sum a \in A a : A \subset {V^{1}, \dots, V^{k}}}

and the notation $S_{V}^{p r o p}$ for the set of proper subset sums of $V$ , i.e.

S_{V}^{p r o p} = {\sum a \in A a : A ⊊ {V^{1}, \dots, V^{k}}, A \neq \emptyset} .

With this notation, we say that the subset sums of $χ$ and $λ$ are separated if $S_{λ}^{p r o p} \cap S_{χ} = \emptyset$ .

Note that for any fixed $X$ with fixed decomposition $X_{1}, \dots, X_{M}$ , the collection of all $λ \in Λ$ such that our condition is not satisfied (i.e. $S_{λ}^{p r o p} \cap S_{χ} \neq \emptyset$ ) is a “small set” in the sense that it has Hausdorff dimension $n - 2$ whereas $Λ$ has Hausdorff dimension $n - 1$ . In particular if $λ$ is randomly chosen from $Λ$ (with respect to any probability measure that is absolutely continuous with respect to $H^{n - 1}$ ), then with probability 1, we have $S_{λ}^{p r o p} \cap S_{χ} = \emptyset$ . We defer a formal proof of this till Proposition 8.

For the remainder of the paper we use the notation $d_{λ} := dist (S_{λ}^{p r o p}, S_{χ})$ for any $λ \in Λ$ . Note that since $N > 1$ , $S_{λ}^{p r o p}$ is never empty and so $d_{λ}$ is always defined. We defer a discussion of methods to compute and/or bound $d_{λ}$ to section 8.

It is well-known that the optimal transport problem has a dual problem with strong duality (we refer the reader to [Vil03, Chapter 1] for background),

(2)

min\tiT#μ=ν∫Xc(x,\tiT(x))dμ=supψ∈\RN∫X(minjc(x,yj)+ψj)\dmu−\innerψλ.

Furthermore given any dual optimizer $ψ$ the map $T_{ψ}$ from Remark 2.1 is a optimal transport map for the primal problem. Because of this we are able to optimize the finite dimensional dual problem instead of the infinite dimensional primal problem. We define Kantorovich’s functional as in [KMT19, Theorem 1.1]

Φ(ψ):=∫X(minic(x,yi)+ψi)\dmu−\innerψλ=∑i∫\Lagi(ψ)(c(x,yi)+ψi)\dmu−\innerψλ.

Our algorithm will find approximate maximizers of $Φ$ . We also recall that $\nabla Φ (ψ) = G (ψ) - λ$ .

2.2. Previous Results

There are many papers that apply a Newton-type algorithm to solve semi-discrete optimal transport problems and Monge-Ampère equations. In [OP88] the authors prove local convergence of a Newton Method for solving a semi-discrete Monge-Ampère equation. Global convergence is proved in [Mir15].

For the semi-discrete optimal transport problem the authors of [KMT19] give a Newton method similar to our Algorithm 1, in the case when the source measure is given by a continuous probability density and has connected support. An analogous result for the case when $μ$ is supported on a union of simplexes is given in [MMT18], however this paper also requires a connectedness condition on the support of $μ$ .

In [BK19], Kitagawa and the author present a Newton method to solve a generalization of the semi-discrete optimal transport problem in which there is a storage fee. Our results do not require a connectedness condition on the support of $μ$ , however we require a strict convexity assumption on the storage fee function that doesn’t hold for the classical semi-discrete optimal transport problem. Furthermore we show that the Laguerre cells associated to the approximate transport maps converge quantitatively both in the sense of $μ$ -symmetric distance and in the sense of Hausdorff distance.

2.3. Outline of the Paper

In section 3 we obtain quantitative bounds on the strong concavity of $Φ$ . In section 4 we exploit these bounds to obtain local convergence of a damped Newton algorithm. In section 5 we obtain global convergence of a gradient descent algorithm applied to a regularized version of $Φ$ . In section 6 we discuss convergence of our two stage algorithm. In section 7 we obtain quantitative convergence of the associated Laguerre cells. In section 8

we describe some methods to compute and/or estimate the parameter

d_{λ}

from the initial data.

3. Spectral Estimates

In this section we work toward quantitative bounds on the strong concavity of $Φ$ .

To start off we recall some notation in order to remain consistent with [KMT19]. We will write $\interior(Xi)$ to denote the interior of the set $X_{i}$ . Given an absolutely continuous measure $μ = ρ d x$ , where $ρ$ is continuous, and a set $A \subset X_{i}$ with Lipschitz boundary, we will write

\abs∂Aρ,Xi:

=∫∂A∩\interior(Xi)ρdHn−1,\absAρ:=μ(A).

The next lemma is virtually identical to [BK19, Lemma 7.4]. We omit the proof.

Suppose that $μ = ρ d x$ satisfies a $(q, 1)$ -PW inequality on $X_{i}$ . Then

infA⊂Xi\abs∂Aρ,Ximin(\absAρ,\absXi∖Aρ)1/q≥121/qκi,

where the infimum is over $A⊂\interior(Xi)$ whose boundary is Lipschitz with finite $H^{n - 1}$ -measure, and $min(\absAρ,\absXi∖Aρ)>0$ .

Recall from [KMT19] that $D G$ is negative semidefinite (as it is the Hessian of a concave function). Furthermore, recall that $\onevect$ is in the kernel of $D G$

, so the largest eigenvalue of

D G

0

. We now work toward the following estimate on the second largest eigenvalue.

Suppose that $μ$ satisfies a local $(q, 1)$ -PW inequality on $X$ with constant $κ$ . Let $ψ∈\RN$ be fixed. Assume that the subset sums of $G (ψ)$ and $χ$ are separated i.e., $S_{G (ψ)}^{p r o p} \cap S_{χ} = \emptyset$ . Then the second largest eigenvalue of $D G (ψ)$ is bounded away from zero. In particular it is bounded by $- \frac{2^{3 - 1 / q}}{C_{\nabla} M^{1 / q} N^{4} κ} d_{G (ψ)}^{1 / q} < 0$ where $d_{G (ψ)} = dist (S_{G (ψ)}^{p r o p}, S_{χ}) > 0$ .

At this point, given some $ψ∈\RN$ we recall the construction of $W$ from [KMT19, Section 5.3]. $W$ is the simple weighted graph whose vertex set is the collection $Y$ , and for any $y_{i}$ and $y_{j}$ , $i \neq j$ there is an edge between $y_{i}, y_{j}$ of weight $w_{i j}$ given by

wij:=DiGj(ψ)=DjGi(ψ)=∫\Lagi,j(ψ)ρ(x)\norm∇xc(x,yi)−∇xc(x,yj)dHn−1(x)

where we have used the notation

\Lagi,j(ψ):=\Lagi(ψ)∩\Lagj(ψ)

for $i$ , $j \in {1, \dots, N}$ . Under the conditions of Theorem 3, $W$ is connected by edges of weight at least $\frac{2^{1 - 1 / q}}{C_{\nabla} M^{1 / q} N^{2} κ} d_{G (ψ)}^{1 / q}$ .

Proof.

Suppose by contradiction that the proposition is false. This implies that removing all edges with weight strictly less than $\frac{2^{1 - 1 / q}}{C_{\nabla} M^{1 / q} N^{2} κ} d_{G (ψ)}^{1 / q}$ yields a disconnected graph. In other words, we can write $W = W_{1} \cup W_{2}$ where $W_{1}$ , $W_{2} \neq \emptyset$ and are disjoint, such that every edge connecting a vertex in $W_{1}$ to a vertex in $W_{2}$ has weight strictly less than $\frac{2^{1 - 1 / q}}{C_{\nabla} M^{1 / q} N^{2} κ} d_{G (ψ)}^{1 / q}$ . Letting $A=∪yi∈W1\Lagi(ψ)$ we see that

\abs∂Aρ,X≤∑{(i,j)∣yi∈W1, yj∈W2}2C∇wij<22−1/qM1/qN2κd1/qG(ψ)\absW1\absW2≤22−1/qM1/qN2κd1/qG(ψ)N24=2−1/qM1/qκd1/qG(ψ).

Now for each $i$ we claim that $min(\absA∩Xiρ,\absXi∖Aρ)<κqiMκqdG(ψ)$ . If $\absA∩Xiρ=0$ or $\absXi∖Aρ=0$ then the claim is trivial. Otherwise the claim follows from Lemma 3, after noting that $\abs∂Aρ,Xi≤\abs∂Aρ,X$ .

Now note that by construction

(3)

μ(A)=∑yi∈W1μ(\Lagi(ψ))=∑yi∈W1G(ψ)i∈SG(ψ).

On the other hand we can partition the $X_{i}$ into two types. Let $I$ be the index set so that $i \in I$ if and only if $\absXi∖Aρ<κqiMκqdG(ψ)$ . We see that

μ (A) = M \sum i = 1 μ (A \cap X_{i}) = \sum i \in I μ (A \cap X_{i}) + \sum i \notin I μ (A \cap X_{i}) = \sum i \in I μ (X_{i}) - μ (X_{i} ∖ A) + \sum i \notin I μ (A \cap X_{i}) .

Hence

	$\absμ(A)−∑i∈Iμ(Xi)$	$\leq \sum i \in I μ (X_{i} ∖ A) + \sum i \notin I μ (A \cap X_{i})$
		$< \sum i \in I \frac{κ_{i}^{q}}{M κ^{q}} d_{G (ψ)} + \sum i \notin I \frac{κ_{i}^{q}}{M κ^{q}} d_{G (ψ)} = \frac{1}{M κ^{q}} M \sum i = 1 κ_{i}^{q} d_{G (ψ)} = d_{G (ψ)} .$

By definition $\sum_{i \in I} μ (X_{i}) \in S_{χ}$ . Furthermore since both $W_{1}$ and $W_{2}$ are nonempty we have $μ (A) \in S_{G (ψ)}^{p r o p}$ , recalling (3). Hence we see that

d_{G (ψ)} = dist (S_{G (ψ)}^{p r o p}, S_{χ}) < d_{G (ψ)}

which is a clear contradiction. ∎

At this point Theorem 3 follows via the exact same method of proof as [BK19, Theorem 8.1]. We omit the proof.

4. Convergence of Newton Algorithm

First we recall the standard damped Newton Algorithm (similar to the one proposed in [KMT19]).

Input: A tolerance

ζ > 0

and an initial

ψ0∈\RN

while $\norm∇Φ(ψk)≥ζ$ do

Step 1:

Compute ${\to d}_{k} = - [D^{2} Φ (ψ_{k})]^{- 1} (\nabla Φ (ψ_{k}))$

Step 2:

Determine the minimum $ℓ∈\N$ such that $ψ_{k + 1, ℓ} := ψ_{k} + 2^{- ℓ} {\to d}_{k}$ satisfies

\norm∇Φ(ψk+1,ℓ)≤(1−2−(ℓ+1))\norm∇Φ(ψk).

Step 3:

Set $ψ_{k + 1} = ψ_{k} + 2^{- ℓ} {\to d}_{k}$ and $k \leftarrow k + 1$ .

end while

Algorithm 1 Damped Newton’s algorithm

We remark that unlike the algorithm proposed in [KMT19] we do not impose a condition on the cells not collapsing. This is because we will prove local convergence for the above algorithm and within the zone of convergence the error reduction requirement will automatically imply that the cells don’t collapse.

Our main result for this section is that Algorithm 1 converges locally with superlinear rate.

Suppose that $μ$ satisfies a local $(q, 1)$ -PW inequality on $X$ . If $\norm∇Φ(ψ0)≤dλ2√N$ then Algorithm 1 converges with linear rate and locally with rate $1 + α$ .

We start with a simple lemma.

For $λ_{1}, λ_{2} \in Λ$ , we have $\absdλ1−dλ2≤\normλ1−λ21$ .

Proof.

Let $x \in S_{λ_{1}}^{p r o p}$ . We claim that $dist(x,Spropλ2)≤\normλ1−λ21$ . Indeed if $x = \sum_{i \in I} λ_{1}^{i}$ then

\absx−∑i∈Iλi2=\abs∑i∈I(λi1−λi2)≤∑i∈I\absλi1−λi2≤\normλ1−λ21.

If we apply this to $x=argmin\tix∈Spropλ1dist(\tix,Sχ)$ then we get that there is $y \in S_{λ_{2}}^{p r o p}$ so that $\absx−y≤\normλ1−λ21$ . Hence

	$dist (S_{λ_{2}}^{p r o p}, S_{χ})$	$\leq dist (y, S_{χ})$
		$≤\absx−y+dist(x,Sχ)$
		$=dist(Spropλ1,Sχ)+\absx−y$
		$≤dist(Spropλ1,Sχ)+\normλ1−λ21$

and so we get $dλ2−dλ1=dist(Spropλ2,Sχ)−dist(Spropλ1,Sχ)≤\normλ1−λ21$ . A symmetric argument proves the opposite inequality.

∎

With the above lemma, Theorem 3 gives us that $Φ$ is locally well-conditioned near its maximum.

Suppose that $μ$ satisfies a local $(q, 1)$ -PW inequality on $X$ and $η < d_{λ}$ . On the set $Jη={ψ∈\RN:\norm∇Φ(ψ)1<η}$ we have that $Φ$ is uniformly $C^{2, α}$ and strongly concave, except in the direction $1$ . In particular on $J^{η}$

\normΦC2,α<C(dλ−η)−2

and

D^{2} Φ \leq - \frac{2^{3 - 1 / q}}{C_{\nabla} M^{1 / q} N^{4} κ} (d_{λ} - η)^{1 / q} P

where $P$

is the orthogonal projection onto the hyperplane perpendicular to

1

and

C

is a constant depending only on

\normρ∞,\diam(X)

, and

c

(the cost).

Proof.

Fix $ψ \in J^{η}$ . Since $0 \in S_{χ}$ and each $λ^{i} \in S_{λ}^{p r o p}$ , we have that each $λ^{i} \geq d_{λ}$ . Furthermore since

\absGi(ψ)−λi≤\norm∇Φ(ψ)∞≤\norm∇Φ(ψ)1<η

we obtain that $G^{i} (ψ) - λ^{i} > - η$ and so $G^{i} (ψ) > d_{λ} - η$ . In other words

(4)

J^{η} \subset K^{d_{λ} - η} .

Hence the first result follows by careful tracing through the proof of [KMT19, Theorem 4.1].

Now for the second claim Lemma 4 gives that

\absdG(ψ)−dλ≤\normG(ψ)−λ1=\norm∇Φ(ψ)1≤η.

In particular $d_{G (ψ)} \geq d_{λ} - η$ . The second result now follows from Theorem 3.

∎

We briefly note that this proposition gives us uniqueness of the maximizer of $Φ$ up to a multiple of $\onevect$ .

If $ψ_{1}, ψ_{2}$ are two maximizers of $Φ$ then there is $r∈\R$ so that $ψ1−ψ2=r\onevect$ .

Proof.

Since $Φ$ is concave and $ψ_{1}, ψ_{2}$ are two maximizers, $Φ$ must be constant along the line segment joining $ψ_{1}$ and $ψ_{2}$ . In particular if $t \in [0, 1]$ then $Φ (t ψ_{1} + (1 - t) ψ_{2}) = Φ (ψ_{1})$ and so $t ψ_{1} + (1 - t) ψ_{2}$ is a maximizer of $Φ$ . Hence $\nabla Φ (t ψ_{1} + (1 - t) ψ_{2}) = 0$ for all $t \in [0, 1]$ . In particular we can conclude $[D^{2} Φ (ψ_{1})] (ψ_{2} - ψ_{1}) = 0$ .

Now since $ψ_{1}$ is a maximizer we have $\nabla Φ (ψ_{1}) = 0$ so the conditions of Proposition 4 are satisfied. Hence we have $D^{2} Φ (ψ_{1}) \leq - \frac{2^{3 - 1 / q}}{C_{\nabla} M^{1 / q} N^{4} κ} (d_{λ} - η)^{1 / q} P$ where $P$ is the orthogonal projection onto the hyperplane perpendicular to $1$ . Hence we get $P (ψ_{2} - ψ_{1}) = 0$ which implies the claim.

∎

Given any maximizer, $ψ$ , of $Φ$ we see that $ψ+r\onevect$ is also a maximizer. In particular given any maximizer, $ψ$ , we can construct the maximizer $ψmax:=ψ−(∑Ni=1ψi)\onevect$ which is a maximizer of $Φ$ that satisfies $\sum_{i = 1}^{N} ψ_{m a x}^{i} = 0$ .

Conversely 4 tells us that there is a unique $ψ_{m a x}$ that maximizes $Φ$ and satisfies $\sum_{i = 1}^{N} ψ_{m a x}^{i} = 0$ . For the remainder of the paper we shall refer to this unique maximizer as $ψ_{m a x}$ .

It is well-known that Proposition 4 gives local $1 + α$ convergence of the standard damped Newton Algorithm. However for convenience of the reader we will give a self contained proof.

If $\norm∇Φ(ψ0)<dλ2√N$ then the iterates of Algorithm 1 satisfy

\norm∇Φ(ψk+1)≤(1−\conjτk2)\norm∇Φ(ψk)

where

\conjτk=min(β1+1αδL1α\norm∇Φ(ψ)21α,1)

and

	$δ$	$:= \frac{1}{2} min λ^{i}$
	$L$	$\leq C δ^{- 2}$
	$β$	$\geq \frac{2^{3 - 2 / q}}{C_{\nabla} M^{1 / q} N^{4} κ} d_{λ}^{1 / q}$

where $C$ is a constant depending only on $\normρ∞,\diam(X)$ , and $c$ (the cost).

Furthermore once we have $\conjτk=1$ we get

\norm∇Φ(ψk+1)≤L\norm∇Φ(ψ)1+αβ1+α.

Proof.

The proof is very similar to that of [KMT19, Proposition 6.1]. However we will give all the details.

Note that $\norm∇Φ(ψ0)<dλ2√N$ implies $\norm∇Φ(ψ0)1<dλ2$ .

For this proof define $\eps:=dλ2,δ:=12minλi$ . Note $\eps≤δ$ and so, by (4), $J\eps⊂K\eps⊂Kδ$ (Note that it is NOT true in general that $J^{η} \subset K^{η}$ ; our choice of $\eps$ is the largest choice for which this is guaranteed to work). Furthermore $L$ is the $C^{1, α}$ norm of $G$ on $K^{δ / 2}$ and $- β$ is the bound of the second largest eigenvalue of $D G$ on $J\eps$ . By Proposition 4 we have

β \geq \frac{2^{3 - 1 / q}}{C_{\nabla} M^{1 / q} N^{4} κ} (\frac{d_{λ}}{2})^{1 / q} = \frac{2^{3 - 2 / q}}{C_{\nabla} M^{1 / q} N^{4} κ} d_{λ}^{1 / q} .

Also by carefully tracing through the proofs in [KMT19] we have $L \leq C δ^{- 2}$ .

We analyze a single iteration of the Algorithm. Define $ψ:=ψk∈J\eps⊂K\eps⊂Kδ$ .

Let $v := [D^{2} Φ (ψ)]^{+} (\nabla Φ (ψ))$ . We see that

\normv≤\norm∇Φ(ψ)β,

as $Φ$ is $β$ -concave in the direction orthogonal to $\onevect$ .

Also define $ψ_{τ} = ψ - τ v$ . Let $τ_{1}$ be the first exit time from $K^{δ / 2}$ . We have that

δ2≤\normG(ψτ1)−G(ψ)≤Lτ1\normv≤Lβ\norm∇Φ(ψ)τ1

and so

τ1≥βδ2L\norm∇Φ(ψ).

Applying Taylor’s formula to $\nabla Φ$ we get

(5)

\nabla Φ (ψ_{τ}) = \nabla Φ (ψ - τ v) = \nabla Φ (ψ) - τ (D \nabla Φ (ψ)) v + R (τ) = \nabla Φ (ψ) - τ \nabla Φ (ψ) + R (τ)

where

	$\normR(τ)$	$=\norm∫τ0(D∇Φ(ψσ)−D∇Φ(ψ))vdσ$
		$=\norm∫τ0DG(ψσ)v−DG(ψ)vdσ$
		$≤∫τ0L\normψσ−ψα\normvdσ$
		$=∫τ0L\normσvα\normvdσ$
		$=L\normvα+1∫τ0σαdσ$
		$=L\normvα+1τα+1α+1$
(6)			$≤L\norm∇Φ(ψ)1+αβ1+ατ1+α$

for $τ \leq 1$ .

Finally we establish the error reduction estimates. (5) gives

\nabla Φ (ψ_{τ}) = (1 - τ) \nabla Φ (ψ) + R (τ)

so we have

\norm∇Φ(ψτ)≤(1−τ2)\norm∇Φ(ψ)

provided that

\normR(τ)≤τ2\norm∇Φ(ψ).

Again using (6) this will be true provided

τ≤min(τ1,β1+1αL1α\norm∇Φ(ψ)21α)=:τ2.

Hence we see that if we set $\conjτk:=τ2$ , then the claim is true. Furthermore as the error goes to zero, eventually we must have $τ_{k} = 1$ . When this happens (5) gives

\nabla Φ (ψ_{1}) = R (1)

and so we get the super-linear convergence.

∎

5. Convergence of Gradient Descent Algorithm

In order to remedy the fact that the above algorithm only has local convergence we propose a regularized version in order to get within a close enough error.

We define the regularized Kantorovich’s functional

\tiΦ(ψ)=Φ(ψ)−γ2\normψ2

where $γ > 0$ is some parameter.

We see that

∇\tiΦ(ψ)=∇Φ(ψ)−γψ=G(ψ)−λ−γψ

and

D2\tiΦ(ψ)=D2Φ(ψ)−γI=DG(ψ)−γI.

In particular we note that $\tiΦ(ψ)$ is strongly $γ$ -concave. If we let $C_{L}$ be the Lipschitz constant of $G$ ( $C_{L} = C N$ for some universal constant $C$ ) then $\tiΦ$ has Lipschitz gradient with constant $C_{L} + γ$ . Hence $\tiΦ$ is well-conditioned.

The projection of $K^{0}$ in the direction $\onevect$ is bounded by

M0:=√N\(maxx∈X(maxy∈Yc(x,y)−miny∈Yc(x,y))\).

In other words if $ψ \in K^{0}$ and $\sum_{i} ψ^{i} = 0$ then $\normψ≤M0$ .

Proof.

Since $\sum_{i} ψ^{i} = 0$ there are some indices $j_{1}, j_{2} \in {1, \dots, N}$ so that $ψ^{j_{1}} \leq 0$ and $ψ^{j_{2}} \geq 0$ . Suppose for sake of contradiction that $\normψ>M0$ . Then we must have $\normψ∞>M0√N=maxx∈X(maxy∈Yc(x,y)−miny∈Yc(x,y))$ . In particular there is some index $j$ so that $\absψj>maxx∈X(maxy∈Yc(x,y)−miny∈Yc(x,y))$ .

We first look at the case where $ψ^{j} > {max}_{x \in X} ({max}_{y \in Y} c (x, y) - {min}_{y \in Y} c (x, y))$ . In this case

	$\Lagj(ψ)$	$= {x \in X ∣ c (x, y_{j}) + ψ^{j} = min i c (x, y_{i}) + ψ^{i}}$
		$\subset {x \in X ∣ c (x, y_{j}) + ψ^{j} \leq c (x, y_{j_{1}}) + ψ^{j_{1}}}$
		$\subset {x \in X ∣ ψ^{j} \leq c (x, y_{j_{1}}) - c (x, y_{j})} = \emptyset$

and so we get $G (ψ)^{j} = 0$ which contradicts $ψ \in K^{0}$ . A similar argument handles the case where $ψ^{j} < - {max}_{x \in X} ({max}_{y \in Y} c (x, y) - {min}_{y \in Y} c (x, y))$ .

∎

We remark that since $X$ is compact and $c$ is continuous, $M_{0} < + \infty$ and so the bound isn’t vacuous.

Since $\tiΦ$ is strongly concave it has a unique maximizer. Although it may not be true that this maximizer is in $K^{0}$ we show that it is still bounded.

Let $ψ^{*}$ be the unique maximizer of $\tiΦ$ . Then $\normψ∗≤M0$ .

Proof.

Let $ψ_{m a x}$ be the maximizer of $Φ$ constructed in Remark 4 so that $\sum ψ_{m a x}^{j} = 0$ . Then $G (ψ_{m a x}) = λ$ and so $ψ_{m a x} \in K^{0}$ . In particular $\normψmax≤M0$ . But since $ψ^{*}$ is the maximizer of $\tiΦ$ we have $\tiΦ($

Two Stage Algorithm for Semi-Discrete Optimal Transport on Disconnected Domains

Authors

A Newton algorithm for semi-discrete optimal transport with storage fees and quantitative convergence of cells

Optimal transport: discretization and algorithms

Quantitative Stability and Error Estimates for Optimal Transport Plans

An algorithm for optimal transport between a simplex soup and a point cloud

Symmetrized semi-discrete optimal transport

Optimal partitions and Semi-discrete optimal transport

Computation of Optimal Transport with Finite Volumes

1. Introduction

2. Setup

2.1. Notations and Conventions

2.2. Previous Results

2.3. Outline of the Paper

3. Spectral Estimates

Proof.

4. Convergence of Newton Algorithm

Proof.

Proof.

Proof.

Proof.

5. Convergence of Gradient Descent Algorithm

Proof.

Proof.

Two Stage Algorithm for Semi-Discrete Optimal Transport on Disconnected Domains

Authors

Related Research

A Newton algorithm for semi-discrete optimal transport with storage fees and quantitative convergence of cells

Optimal transport: discretization and algorithms

Quantitative Stability and Error Estimates for Optimal Transport Plans

An algorithm for optimal transport between a simplex soup and a point cloud

Symmetrized semi-discrete optimal transport

Optimal partitions and Semi-discrete optimal transport

Computation of Optimal Transport with Finite Volumes

This week in AI

1. Introduction

2. Setup

2.1. Notations and Conventions

2.2. Previous Results

2.3. Outline of the Paper

3. Spectral Estimates

Proof.

4. Convergence of Newton Algorithm

Proof.

Proof.

Proof.

Proof.

5. Convergence of Gradient Descent Algorithm

Proof.

Proof.