Ordinary differential equations for regularized variational problems involving semi-discrete optimal transport

Abstract: We consider entropically regularized, semi-discrete versions of variational problems on the set of probability measures involving optimal transport as well as other terms. We prove that the solutions can be characterized by well-posed ordinary differential equations in the regularization parameter. The initial conditions for these equations, corresponding to solutions to completely regularized problems, are typically known explicitly. The ODE can then be solved to recover the solution for an arbitrary degree of regularization; we verify that the solution is continuous in the regularization parameter, implying that taking the limit of the trajectory yields the solution to the fully unregularized problem. We establish analogous results for a version of the problem when the non-optimal transport term is not scaled with the regularization parameter. We exploit our characterization to numerically solve several example problems using standard ODE methods; this strategy exhibits superior robustness to alternatives such as Newton's method, as arbitrary initializations are not required.

• The paper establishes an ODE-based framework to trace the evolution of dual potentials in entropic regularization of semi-discrete OT, ensuring continuous convergence from regularized to unregularized regimes.
• The paper demonstrates that the ODE solver achieves first- and second-order empirical convergence while being robust to initialization compared to traditional iterative methods.
• The paper’s numerical experiments validate the homotopic approach across varied OT scenarios, offering an efficient and scalable method for high-dimensional problems.

Ordinary Differential Equation Approach for Regularized Semi-Discrete Variational Optimal Transport

Problem Context and Formulation

This work develops a unified ODE characterization for a broad class of entropically regularized, semi-discrete variational problems involving optimal transport (OT) with additional convex functionals. The foundational setup considers a primal variational problem formulated as

$$\inf_{\substack{\nu\in\mathcal{P}(Y)\\gamma\in\Gamma(\mu,\nu)}} t \int_{X\times Y} c\; d\gamma + (1-t)\mathrm{Ent}(\gamma|\mu\otimes\sigma) + tF(\nu),$$

where $\mu$ is a given reference measure, $Y$ is a finite discrete set, and $F$ a convex functional. The entropic regularization parameter $t \in [0,1]$ interpolates between the fully regularized ($t=0$) and unregularized ($t=1$) regimes. This formulation encapsulates a range of applications including Wasserstein gradient flows, Cournot-Nash equilibria, hedonic pricing, and urban planning.

Practical computation motivates discretization of $Y$, leading to semi-discrete OT, and the addition of entropy regularization to obtain strictly convex programs. The dual problem is then expressed in finite dimensions, allowing variational and algorithmic methods from convex optimization.

ODE Characterization of the Solution Path

A central contribution is the derivation and analysis of an ordinary differential equation (ODE) governing the evolution of the dual potential $\psi(t)$ as a function of $t$. This is obtained by differentiating the first-order optimality conditions with respect to $\mu$0 and employing the implicit function theorem. The resulting ODE has the schematic form

$\mu$1

where $\mu$2 is the concave dual functional depending on the regularization parameter $\mu$3. The strong convexity of $\mu$4 ensures negative-definiteness and thus well-posedness except at $\mu$5, motivating a rescaled formulation.

Critically, the ODE remains well-posed up to $\mu$6, and under a uniform boundedness result for $\mu$7 (via coercivity and strong convexity), the solution path $\mu$8 is shown to be continuous at $\mu$9, recovering the unregularized dual optimum in the limit.

Figure 1: Uniform density visualization, highlighting the distribution $Y$0 on the continuous space $Y$1.

Parameter-Independent Functional Variant

A second variational scheme is analyzed, in which the convex potential $Y$2 is not scaled with $Y$3, leading to a distinct yet structurally similar dual ODE. Uniform boundedness and continuity at $Y$4 are established in this case as well. Both the scaled and unscaled ODEs support explicit computation of initial conditions due to closed-form solutions in the fully regularized regime ($Y$5). These properties collectively provide a robust homotopic framework for tracing the full solution curve from the entropic to the unregularized problem.

Numerical Methodology and Empirical Results

The ODE characterization allows the application of standard ODE integrators (e.g., third-order Runge-Kutta), bypassing the need for iterative fixed-point or line search algorithms such as Sinkhorn or Newton for each regularization parameter value. Numerical experiments encompass a suite of prototypical variational OT instances: regularized discrete transport with entropy, additional quadratic functionals, composite cost, and discrete Wasserstein barycenters.

Key findings from the extensive computational study include:

• First- and Second-Order Empirical Convergence: In multiple settings, the ODE-based solver achieves first- or second-order error decay with respect to $Y$6 in both one- and two-dimensional domains, conditional on quadrature accuracy for integrals over the continuous space.
• Robustness to Initial Guesses: Unlike Newton's method for the unregularized dual, the ODE scheme is insensitive to initialization, since the trajectory is uniquely determined by the well-posed Cauchy problem. This is particularly significant for high-dimensional or poorly conditioned cases, where Newton's method can fail or stagnate.
• Computation of the Entire Homotopy: The solver outputs the entire path $Y$7, providing fine-grained access to the evolution of induced Laguerre (power) cells and the measures assigned to each, supporting detailed geometric analysis and visualization.
• Numerical Benchmarks: In all regimes, residuals and CPU timings are reported for increasing $Y$8 (number of discrete target points) and step refinement. The ODE approach demonstrates robust and scalable performance, with error stagnation solely attributed to quadrature when $Y$9 and the integrand concentrates.
• Hybrid Strategies: The solution obtained via ODE at $F$0 can be used to robustly initialize Newton iterations for the final, unregularized problem, combining homotopic robustness with the fast asymptotic convergence of Newton in favorable cases.

Figure 2: Illustration of the evolution of Laguerre cell measures for increasing $F$1.

Figure 3: Visualization of the Laguerre cell configuration at $F$2 for benchmark geometry.

Figure 4: Visualization of the Laguerre cell configuration at $F$3 for an alternative density.

Theoretical Implications and Future Directions

The ODE interpretation for solution paths in regularized semi-discrete OT variational problems provides analytic insight into the stability and geometry of entropically regularized solutions. The continuity result at $F$4 ensures rigorous connection between regularized and unregularized settings, facilitating fine extrapolation and error certification in practical applications.

Potential future developments and research avenues include:

• Extension to General Cost and Functionals: The well-posedness analysis accommodates general cost functions $F$5 and convex potentials $F$6, suggesting straightforward generalization beyond the quadratic and entropy cases.
• Adaptivity in Solver Precision: Since the numerical bottleneck is in quadrature as $F$7, adaptive mesh refinement or Monte Carlo schemes could further improve empirical accuracy and runtime.
• Infinite-Dimensional Extensions: While the present framework is finite-dimensional, related ODE-based homotopy approaches may be adapted for infinite-dimensional convex variational problems, for example in function-space Wasserstein gradient flows.

Conclusion

This paper establishes an ODE-based framework for computing and analyzing the full regularization trajectory for a wide class of variational semi-discrete OT problems with entropy and convex penalties. The approach provides both rigorous theoretical guarantees and effective practical algorithms, combining robustness to initialization, strong convergence properties, and the ability to visualize the geometric evolution of partitions. The framework is broadly relevant for OT applications requiring scalable, regularization-continuation methods with rigorous guarantees on solution recovery and continuity (2604.03780).