For squared Euclidean distance, cyclical monotonicity becomes ordinary convexity. This identifies the optimal map as a gradient and gives a multidimensional replacement for monotone rearrangement.
Convex subgradients
Definition. For a proper lower-semicontinuous convex function \(u:\mathbb R^d\to(-\infty,\infty]\), its subdifferential at \(x\) is \[ \partial u(x)=\{y:u(z)\ge u(x)+y\cdot(z-x)\ \text{for every }z\}. \] If \(u\) is differentiable at \(x\), then \(\partial u(x)=\{\nabla u(x)\}\).
Rockafellar theorem. A set \(\Gamma\subset\mathbb R^d\times\mathbb R^d\) is cyclically monotone for the pairing, meaning \[ \sum_i x_i\cdot y_i\ge\sum_i x_i\cdot y_{\sigma(i)}, \] if and only if it is contained in \(\partial u\) for some proper lower-semicontinuous convex \(u\).
For \(c(x,y)=\tfrac12|x-y|^2\), expanding the squares shows that \(c\)-cyclical monotonicity is exactly cyclical monotonicity of the pairing. The quadratic terms cancel under permutations.
Brenier's theorem
Brenier's theorem. Let \(\mu,\nu\in\mathcal P_2(\mathbb R^d)\), and suppose \(\mu\) is absolutely continuous with respect to Lebesgue measure. Then there exists a convex lower-semicontinuous \(u\) such that \[ T=\nabla u,\qquad T_\#\mu=\nu. \] The plan \((\operatorname{id},T)_\#\mu\) is the unique optimal plan for the cost \(\tfrac12|x-y|^2\). Consequently \(T\) is the unique optimal map \(\mu\)-almost everywhere.
Proof. Existence of an optimal plan \(\pi\) follows from page 3. Its support is \(c\)-cyclically monotone, so after expanding the quadratic cost it is cyclically monotone for \(x\cdot y\). Rockafellar's theorem gives a convex \(u\) with \(\operatorname{supp}\pi\subset\partial u\).
A convex function on \(\mathbb R^d\) is locally Lipschitz in the interior of its effective domain and therefore differentiable Lebesgue-almost everywhere by Rademacher's theorem. Since \(\mu\ll\mathcal L^d\), \(\partial u(x)=\{\nabla u(x)\}\) for \(\mu\)-almost every \(x\). Disintegrating \(\pi=\mu(\mathrm{d}x)\pi_x(\mathrm{d}y)\), the support inclusion forces \(\pi_x=\delta_{\nabla u(x)}\) almost surely. Hence \(\pi=(\operatorname{id},\nabla u)_\#\mu\) and \((\nabla u)_\#\mu=\nu\).
If \(\widetilde\pi\) is another optimal plan, the average \((\pi+\widetilde\pi)/2\) is optimal and the same argument makes it a graph plan. But the conditional law of this average is \(\tfrac12\delta_{\nabla u(x)}+\tfrac12\widetilde\pi_x\), which can be a Dirac mass only when \(\widetilde\pi_x=\delta_{\nabla u(x)}\) almost surely. Thus the optimal plan is unique.
Converse optimality
Proposition. If \(u\) is convex, \(\nabla u\) exists \(\mu\)-almost everywhere, and \(T=\nabla u\) pushes \(\mu\) to \(\nu\), then \(T\) is optimal for quadratic cost.
Proof. The graph of \(\partial u\) is cyclically monotone by summing the subgradient inequalities around a cycle. Therefore the graph of \(T\) is \(c\)-cyclically monotone for quadratic cost. The optimal support theorem from page 5, or equivalently the convex conjugate inequality \(u(x)+u^*(y)\ge x\cdot y\) with equality at \(y=\nabla u(x)\), supplies dual potentials touching the cost on the graph. Complementary slackness proves optimality.
Inverse transport and Legendre duality
If both \(\mu\) and \(\nu\) are absolutely continuous, the reverse optimal map is \(\nabla u^*\), where \[ u^*(y)=\sup_x\{x\cdot y-u(x)\}. \] At differentiability points, \[ y=\nabla u(x)\quad\Longleftrightarrow\quad x=\nabla u^*(y). \] Thus the forward and reverse maps are inverse almost everywhere with respect to the corresponding measures.
Gaussian transport
Proposition. Let \(\mu=N(m_0,\Sigma_0)\), \(\nu=N(m_1,\Sigma_1)\), with positive-definite covariances. The Brenier map is \[ T(x)=m_1+A(x-m_0), \quad A=\Sigma_0^{-1/2} (\Sigma_0^{1/2}\Sigma_1\Sigma_0^{1/2})^{1/2} \Sigma_0^{-1/2}, \] and \[ W_2^2(\mu,\nu) =|m_0-m_1|^2+ \operatorname{tr}\!\left( \Sigma_0+\Sigma_1 -2(\Sigma_0^{1/2}\Sigma_1\Sigma_0^{1/2})^{1/2} \right). \]
Proof. The matrix \(A\) is symmetric positive definite, hence \(T=\nabla u\) for the convex quadratic \[ u(x)=m_1\cdot x+\tfrac12(x-m_0)^{\mathsf T}A(x-m_0). \] Direct multiplication gives \(A\Sigma_0A^{\mathsf T}=\Sigma_1\), so \(T_\#\mu=\nu\). The converse proposition proves optimality. Expanding \(\mathbb E|X-T(X)|^2\) and using the covariance identity yields the distance formula.
Polar factorization
Polar factorization theorem. Let \(\mu\ll\mathcal L^d\) and let \(S\in L^2(\mu;\mathbb R^d)\). If \(\nu=S_\#\mu\), then \[ S=\nabla u\circ s \] \(\mu\)-almost everywhere, where \(\nabla u\) is the Brenier map from \(\mu\) to \(\nu\) and \(s_\#\mu=\mu\). Under the usual nondegeneracy assumptions the factorization is unique.
This is the nonlinear analogue of matrix polar decomposition: the measure-preserving map \(s\) plays the role of the orthogonal factor, and the gradient of a convex function plays the role of the positive factor.
Why absolute continuity matters
Example. If \(\mu=\delta_0\) and \(\nu=\tfrac12(\delta_{-1}+\delta_1)\), the optimal plan cannot be induced by a map. More subtly, a source supported on a lower-dimensional surface may meet a multivalued part of \(\partial u\) on a set of positive \(\mu\)-mass. Absolute continuity removes this ambiguity because convex functions are differentiable Lebesgue-almost everywhere.