Projected Gradient Method Combined with Homotopy Techniques for Volume-Measure-Preserving Optimal Mass Transportation Problems

Abstract

Optimal mass transportation has been widely applied in various fields, such as data compression, generative adversarial networks, and image processing. In this paper, we adopt the projected gradient method, combined with the homotopy technique, to find a minimal volume-measure-preserving solution for a 3-manifold optimal mass transportation problem. The proposed projected gradient method is shown to be sublinearly convergent at a rate of O(1/k). Several numerical experiments indicate that our algorithms can significantly reduce transportation costs. Some applications of the optimal mass transportation maps—to deformations and canonical normalizations between brains and solid balls—are demonstrated to show the robustness of our proposed algorithms.

Appendices

Area-Weighted Stretching Energy Minimization [35]

Let ${\mathcal {M}}$ be a simply connected tetrahedral mesh and $\rho$ be a density function on ${\mathcal {M}}$. Let $g:\partial {\mathcal {M}}\rightarrow {\mathbb {S}}^2$ be a piecewise affine map induced by ${\mathbf {g}}=[({\mathbf {g}}^1)^\top , ({\mathbf {g}}^2)^\top , ({\mathbf {g}}^3)^\top ]^\top$. We define the area-weighted stretching energy functional by

$$\begin{aligned} {\mathcal {E}}_{\mathbb {A}}(g) = \frac{1}{2} \sum _{s=1}^3 ({\mathbf {g}}^s)^\top L_{\mathbb {A}}(g) {\mathbf {g}}^s, \end{aligned}$$

where $L_{\mathbb {A}}(g)$ is the area-weighted Laplacian matrix with

$$\begin{aligned} {[}L_{\mathbb {A}}(g)]_{i,j} = {\left\{ \begin{array}{ll} -w_{i,j}(g) \equiv -\frac{1}{2}\left( \frac{\cot \theta _{i,j}(g)}{\sigma _{g}([v_i,v_j,v_k])} + \frac{\cot \theta _{j,i}(g)}{\sigma _{g}([v_j,v_i,v_\ell ])}\right) &{} \text {if }[v_i,v_j]\in {\mathcal {E}}({\mathcal {M}}),\\ \sum _{\ell \ne i} w_{i,\ell }(g) &{} \text {if }j=i,\\ 0 &{} \text {otherwise} \end{array}\right. } \end{aligned}$$

(38)

in which $\sigma _{g}(\alpha ) = \frac{\rho (\alpha ) |\alpha |}{|g(\alpha )|}$ is the stretching factor of g on the triangular face $\alpha =[v_i,v_j,v_k]$ or $[v_j,v_i,v_\ell ]$, and $\theta _{i,j}(g)$ and $\theta _{j,i}(g)$ are two angles opposite to the edge $g([v_i,v_j])$. We summarize area-weighted stretching energy minimization (ASEM) for the computation of spherical area-measure-preserving parameterizations below.

Proof of Lemmas

1.1 Proof of Lemma 1

(i) Let ${\mathbf {g}}_{\tau } = {\mathbf {g}}^{(k)} + \tau ({\mathbf {g}}^* - {\mathbf {g}}^{(k)})$ be the line segment between ${\mathbf {g}}^{(k)}$ and ${\mathbf {g}}^*$ for $\tau \in [0,1]$. Then, we have

$$\begin{aligned}&\left\| {\mathbf {g}}^{(k+1)} - {\mathbf {g}}^{*} \right\| _2 \\&\quad = \left\| \mathrm {P}_{{\mathcal {C}}}(\overline{{\mathbf {g}}}^{(k)}) - \mathrm {P}_{{\mathcal {C}}}({\mathbf {g}}^{*}) \right\| _2 \\&\overset{(27)}{\le } \left\| \overline{{\mathbf {g}}}^{(k)} - {\mathbf {g}}^{*} \right\| _2 \overset{(24)}{=} \left\| {\mathbf {g}}^{(k)} - {\eta }_k \nabla c({\mathbf {g}}^{(k)}) - {\mathbf {g}}^* \right\| _2 \\&\quad = \left\| (I - {\eta }_k \int _0^1 \nabla ^2 c({\mathbf {g}}_\tau ) \,{\mathrm d}\tau ) ({\mathbf {g}}^{(k)} - {\mathbf {g}}^*) \right\| _2 ~\text {(Fundamental Theorem of Calculus)}\\&\quad \le \sup _{\tau \in [0,1]} \left\| I - {\eta }_k \int _0^1 \nabla ^2 c({\mathbf {g}}_\tau ) \,{\mathrm d}\tau \right\| _2 \left\| {\mathbf {g}}^{(k)} - {\mathbf {g}}^* \right\| _2 \\&\quad \le \left\| {\mathbf {g}}^{(k)} - {\mathbf {g}}^* \right\| _2 ~ \text {(from the convexity of }c(\cdot )\text { and }{\eta }_k L < 2). \end{aligned}$$

From property (26) and (24)–(25) with ${\mathbf {h}}={\mathbf {g}}^*$ or ${\mathbf {g}}^{(k)}$, we have

$$\begin{aligned}&\left( \mathrm {P}_{{\mathcal {C}}}(\overline{{\mathbf {g}}}^{(k)}) - \overline{{\mathbf {g}}}^{(k)} \right) ^\top \left( \mathrm {P}_{{\mathcal {C}}}(\overline{{\mathbf {g}}}^{(k)}) - {\mathbf {h}} \right) \\&\quad = \left( {\mathbf {g}}^{(k+1)} - ({\mathbf {g}}^{(k)} - {\eta }_k \nabla c({\mathbf {g}}^{(k)}))\right) ^\top \left( {\mathbf {g}}^{(k+1)}-{\mathbf {h}}\right) \\&\quad = {\eta }_k \left( \nabla c({\mathbf {g}}^{(k)}) - \frac{1}{{\eta }_k}({\mathbf {g}}^{(k)} - {\mathbf {g}}^{(k+1)}) \right) ^\top \left( {\mathbf {g}}^{(k+1)} - {\mathbf {h}} \right) \\&\quad \overset{(28)}{=} {\eta }_k \left( \nabla c({\mathbf {g}}^{(k)}) - d_{{\mathcal {C}}}({\mathbf {g}}^{(k)}) \right) ^\top \left( {\mathbf {g}}^{(k+1)} - {\mathbf {h}} \right) \le 0. \end{aligned}$$

$\square$

1.2 Proof of Lemma 2

For convenience, we now delete the superscript "k" and denote

$$\begin{aligned} {\mathbf {g}}:= & {} {\mathbf {g}}^{(k)}, ~~ {\mathbf {g}}^{+}:= {\mathbf {g}}^{(k+1)}, \end{aligned}$$

(39)

$$\begin{aligned} \overline{{\mathbf {g}}}:= & {} \overline{{\mathbf {g}}}^{(k)} = {\mathbf {g}}^{(k)} - {\eta }_k \nabla c({\mathbf {g}}^{(k)}) = {\mathbf {g}} - {\eta } \nabla c({\mathbf {g}}), \end{aligned}$$

(40)

and

$$\begin{aligned} {\mathbf {g}}^{+} = \mathrm {P}_{{\mathcal {C}}}(\overline{{\mathbf {g}}}), ~ d_{{\mathcal {C}}}({\mathbf {g}}) := \frac{1}{\eta }({\mathbf {g}} - {\mathbf {g}}^+). \end{aligned}$$

(41)

(i) From (39)–(41), it follows that

$$\begin{aligned} c({\mathbf {g}}^*) - c({\mathbf {g}}^+)&= c({\mathbf {g}}^*) - c({\mathbf {g}}) -(c({\mathbf {g}}^+)-c({\mathbf {g}}))\\&\ge \nabla c({\mathbf {g}})^\top ({\mathbf {g}}^* - {\mathbf {g}}) - \left( \nabla c({\mathbf {g}})^\top ({\mathbf {g}}^+-{\mathbf {g}})+\frac{L}{2}\Vert {\mathbf {g}}^+-{\mathbf {g}}\Vert _2^2\right) \\&\text {(From convexity and }L\text {-smoothness of }c(\cdot ))\\&= \nabla c({\mathbf {g}})^\top ({\mathbf {g}}^*-{\mathbf {g}}^+) - \frac{L}{2}\Vert {\mathbf {g}}^+-{\mathbf {g}}\Vert _2^2 \\&\ge d_{{\mathcal {C}}}({\mathbf {g}})^\top ({\mathbf {g}}^*-{\mathbf {g}}^+) - \frac{L}{2}\Vert {\mathbf {g}}^+-{\mathbf {g}}\Vert _2^2&\text {(Lemma } 1 \text { (ii))} \\&= d_{{\mathcal {C}}}({\mathbf {g}})^\top ({\mathbf {g}}^*-{\mathbf {g}}) + d_{{\mathcal {C}}}({\mathbf {g}})^\top ({\mathbf {g}}-{\mathbf {g}}^+) - \frac{L}{2}\Vert {\mathbf {g}}^+-{\mathbf {g}}\Vert _2^2 \\&\overset{(41)}{=} d_{{\mathcal {C}}}({\mathbf {g}})^\top ({\mathbf {g}}^*-{\mathbf {g}}) + {\eta } \Vert d_{{\mathcal {C}}}({\mathbf {g}})\Vert _2^2 - \frac{L}{2}{\eta }^2 \Vert d_{{\mathcal {C}}}({\mathbf {g}})\Vert _2^2 \\&= d_{{\mathcal {C}}}({\mathbf {g}})^\top ({\mathbf {g}}^*-{\mathbf {g}}) + {\eta } \left( 1 - \frac{{\eta }L}{2} \right) \Vert d_{{\mathcal {C}}}({\mathbf {g}})\Vert _2^2. \end{aligned}$$

This implies that

$$\begin{aligned} c({\mathbf {g}}^*)-c({\mathbf {g}}^+) \ge d_{{\mathcal {C}}}({\mathbf {g}})^\top ({\mathbf {g}}^*-{\mathbf {g}}). \end{aligned}$$

From the Cauchy–Schwarz theorem and Lemma 1 (i), we have

$$\begin{aligned} \Vert d_{{\mathcal {C}}}({\mathbf {g}})\Vert _2 \ge \frac{c({\mathbf {g}}^+)-c({\mathbf {g}}^*)}{\Vert {\mathbf {g}}-{\mathbf {g}}^*\Vert _2} \ge \frac{c({\mathbf {g}}^+)-c({\mathbf {g}}^*)}{\Vert {\mathbf {g}}^{(0)}-{\mathbf {g}}^*\Vert _2}. \end{aligned}$$

(ii) Similar to the proof of (i), by setting ${\mathbf {g}}^* = {\mathbf {g}}$, we obtain

$$\begin{aligned} c({\mathbf {g}})&\ge c({\mathbf {g}}^+) + {\eta } \left( 1-\frac{{\eta } L}{2} \right) \Vert d_{\mathcal {C}}({\mathbf {g}})\Vert _2^2 \\&\ge c({\mathbf {g}}^+) + {\eta }_0 \left( 1-\frac{{\eta }_* L}{2} \right) \Vert d_{\mathcal {C}}({\mathbf {g}})\Vert _2^2 \\&= c({\mathbf {g}}^+) + \frac{1}{2q} \Vert d_{\mathcal {C}}({\mathbf {g}})\Vert _2^2. \end{aligned}$$

$\square$