Structure preserving schemes for a class of Wasserstein gradient flows

TL;DR

This work develops two time-discretization strategies for Wasserstein gradient flows that preserve core physical structures: mass, positivity, and energy dissipation. The first approach (S1) generalizes existing schemes and yields unconditional solvability with energy dissipation in specific settings; the second (S2) uses energy splitting and a scalar auxiliary variable (SAV) to achieve robust energy dissipation with a modified energy, applicable to a broader class of energies and including Onsager-type flows. The authors establish theoretical properties (mass conservation, positivity, unique solvability) and provide second-order variants, together with a finite-difference spatial discretization that preserves the discrete analogs of integration-by-parts. Extensive numerical experiments on Barenblatt PME, Fokker–Planck, and Fisher–KPP-type problems demonstrate accuracy, efficiency, and energy behavior, including positivity maintenance even for challenging regimes (e.g., large $m$). The schemes offer practical, structure-preserving tools for simulating a wide range of Wasserstein gradient flows with potential applications in physics and engineering.

Abstract

We introduce in this paper two time discretization schemes tailored for a range of Wasserstein gradient flows. These schemes are designed to preserve mass, positivity and to be uniquely solvable. In addition, they also ensure energy dissipation in many typical scenarios. Through extensive numerical experiments, we demonstrate the schemes' robustness, accuracy and efficiency.

Figures (7)

• Figure 1: Visualizations of the PME evolution with S1 and S2 at a cross section $y=0$ and energy dissipation comparison for $m=3$.
• Figure 2: Comparative analysis of the average number of Newton iterations required for schemes S1 and S2 every 50 computational steps over the simulation time from $t_0=0$ to $T=1$, demonstrating the relative computational efficiency of each scheme.
• Figure 3: Solution profiles for the Fokker-Planck equation. The simulation begins with the initial condition $\rho(x,y,1)$ as described in the equation \ref{['heat kernel']}. The reference steady state is given by $\rho_{\infty}$. Additionally, cross-sectional views of the solution at $y=0$ are presented in Figure \ref{['fkcrosss1']} and \ref{['fkcrosss2']}.
• Figure 4: The average number of Newton iterations for schemes S1 and S2, measured every 200 computational steps from $t_0=0$ to $T=4$.
• Figure 5: On the left, the initial state is set to a random configuration using the seed(1). On the right, we depict the potential function given by $V(x,y) = 1 - \sin(5\pi x)\sin(3\pi y)$.

Theorems & Definitions (12)

• Theorem 2.1
• Proof 1
• Theorem 2.2
• Proof 2
• Theorem 2.3
• Proof 3
• Theorem 2.4
• Proof 4
• Theorem 2.5
• Proof 5