Source-linked AI summary

Gradient-Enhanced Proximal Algorithms for Mean Field Planning on Surfaces

Chengrun Jiang

arXiv:2609.10370v1math.NAcs.NA

TL;DR

Surface transport updates depend on derivatives used in a shared continuity-constraint projection, so potential accuracy alone does not determine the density and momentum corrections. The paper combines temporal and spatial gradient recovery with surface discretization and incorporates the construction into proximal algorithms. It distinguishes recovered updates from exact projections and gives a conditional finite-iteration estimate for ISTA and FISTA.

  • Problem

    The shared continuity-constraint projection requires a time–space Poisson solve, and potential accuracy alone does not determine the accuracy of its density and momentum updates.

    Source: §1 · p. 2
  • Method

    The paper combines finite differences in time, surface finite elements in space, and temporal PPR with spatial PPPR recovery for the constraint projection.

    Source: §1 · p. 3
  • Results

    The paper states a conditional finite-iteration estimate for basic ISTA and FISTA.

    Source: §5.3 · p. 19
  • Takeaways & Limitations

    The recovered update is distinct from an exact discrete projection, and feasibility of the discrete constraint set is a separate linear-algebra condition.

    Source: §4.3 · p. 13

    The projection comparison identity is not the system used by the recovered update, and feasibility of the discrete constraint set must be checked separately.

    Source: §4.3 · p. 13, §7 · p. 25

Abstract

from arXiv · show

Mean field planning on a surface prescribes initial and terminal densities and minimizes a transport energy subject to the continuity equation. Proximal algorithms for this problem repeatedly solve a time--space Poisson equation, whose temporal derivative and surface gradient determine the density and momentum corrections. We study a gradient-enhanced approximation of this constraint projection using finite differences in time, surface finite elements in space, temporal polynomial preserving recovery, and spatial parametric polynomial preserving recovery. The same construction is incorporated into ISTA, FISTA, and Douglas--Rachford splitting. We distinguish the recovered update from an exact discrete projection and derive residual identities and conditional finite-iteration perturbation bounds that retain data, boundary, and linear-solver errors. Existing derivative-recovery estimates identify a higher-order contribution under suitable regularity and mesh assumptions; they do not by themselves establish convergence of the outer optimization iteration. Available numerical illustrations on the sphere and a more complicated algebraic surface are discussed together with the limits of the recorded refinement data.

1 Introduction

Mean field planning connects prescribed endpoint densities through a continuity-constrained transport problem, and its shared constraint projection determines density and momentum corrections. The paper develops a gradient-enhanced approximation of that projection, uses it in several proximal algorithms, and analyzes its difference from an exact projection.

1 Introduction

  • Raw surface finite-element gradients can have one less power of h in L2 accuracy than the potential, even though the gradient helps determine updated density and momentum.The potential’s temporal derivative and surface gradient drive the corrections, so potential accuracy alone does not determine update accuracy.Source: §1 · p. 2
  • The approximation combines a time–space Poisson solver using finite differences in time and surface finite elements in space with temporal PPR and spatial PPPR recovery.Recovered operators enter both the projection equation and the subsequent density and momentum updates.Source: §1 · p. 3
  • The shared construction is used in ISTA, FISTA, and Douglas–Rachford splitting, with FISTA including backtracking, adaptive restarting, and monotone safeguards.Proximal-gradient iterations treat the energy explicitly, whereas Douglas–Rachford evaluates its proximal map.Source: §1 · pp. 2–3
  • The analysis distinguishes the recovered update from an exact projection and tracks derivative-recovery, data, boundary, and solver errors in updates and finite iteration sequences.These discretization estimates are distinct from convergence results for the outer optimization iteration.Source: §1 · p. 3
  • The numerical illustrations examine stepsize dependence, interaction costs, sphere and Enzensberger–Stern surface transport, and refinement behavior.The introduction describes these study areas but does not report comparative outcomes.Source: §1 · p. 3

2 Preliminaries

The preliminaries formulate mean field planning on a compact surface as minimization of a dynamic cost subject to continuity and endpoint constraints. They also identify dynamic optimal transport as a special case and distinguish the fixed-endpoint planning setting from a terminal-cost MFG formulation.

2 Preliminaries

  • The setting is a compact, connected, oriented, C3-smooth two-dimensional surface without boundary, embedded in R3 with its induced metric.Stronger regularity is specified when derivative recovery requires it.Source: §2 · p. 3
  • The dynamic MFP formulation minimizes a cost involving density and momentum subject to the continuity equation and prescribed endpoint densities.Here momentum is m = ρv, with v the velocity field; the cost includes a dynamic term and an interaction cost.Source: §2 · p. 3
  • Nonnegativity is imposed through the energy domain rather than the affine projection constraint.The energy is assumed proper, lower semicontinuous, and convex in the density–momentum pair.Source: §2 · pp. 3–4
  • For the cited choice F = 0, the MFP problem becomes the dynamical formulation of optimal transport; replacing its fixed terminal-density constraint with a convex terminal cost G gives a variational formulation of potential mean field games.A finite terminal cost generally does not enforce a specified terminal density.Source: §2 · p. 4

3 Algorithm

Section 3 develops optimization algorithms for mean field planning by representing the continuity constraint with a set indicator and projecting onto its constraint set. It describes Douglas–Rachford splitting’s local proximal update for the nonlinear transportation cost and global projection for the linear constraint, along with limits on intermediate densities, convergence assumptions, and termination checks.

3.1 Computation of the constraint proximal map

  • The constraint proximal map is the orthogonal projection onto C, obtained from a uniquely solvable weak potential problem; it is firmly nonexpansive, while general inputs need not satisfy homogeneous temporal Neumann conditions.For smooth fields the update is ρ⋆ = ρ + ∂tφ and m⋆ = m + ∇gφ.Source: §3.1 · pp. 5–6

3.2 Descent of Y(ρ, m)

  • Descent of Y uses either explicit gradient steps, as in ISTA and FISTA, or proximal splitting methods that apply ProxY alongside the constraint projection.The section identifies Douglas–Rachford, primal-dual, and ALG2/ADMM schemes as proximal-splitting approaches.Source: §3.2 · p. 6

3.3 ISTA

  • ISTA uses step sizes η(k) ≤ 1/L(u), but the transport kinetic energy is not Lipschitz smooth as ρ approaches zero, and clipping can disrupt mass and continuity constraints.A positive density lower bound prevents division by zero, but a Lipschitz gradient bound also requires control of momentum and the interaction Hessian; the safeguarded procedure rejects inadmissible trials.Source: §3.3 · p. 7

3.4 FISTA

  • FISTA extrapolates updated iterates and has an O(k^-2) objective-gap bound only for the standard convex problem with a Lipschitz gradient and exact proximal steps.That rate does not automatically extend to recovered projections or all adaptive variants.Source: §3.4 · p. 7

3.4.1 Backtracking and adaptive restarting

  • Adaptive FISTA backtracking reduces the trial step until a descent inequality holds, while restarting resets extrapolation and acceleration; these safeguards do not prove monotonicity of Y + χC.The safeguarded variant also tests monotonicity and discards failed trials.Source: §3.4.1 · p. 8

3.4.2 Monotone FISTA

  • Monotone FISTA compares the new proximal candidate's energy with that of the previously accepted iterate, but recovered updates require a separate constraint check.For the energy comparison to match standard exact MFISTA for Y + χC, the accepted point and both candidates must lie in C ∩ dom Y.Source: §3.4.2 · pp. 8–9

3.5 Proximal Splitting

  • Douglas–Rachford splitting treats Y and the continuity constraint implicitly by applying both proximal operators, with the unbarred sequence giving projected iterates.Its iteration uses a fixed proximal parameter γ > 0 and relaxation parameter α ∈ (0, 2).Source: §3.5 · p. 9
  • Douglas–Rachford alternates a local transportation-cost proximal update with a global constraint projection, reducing them to independent cubic equations and a time-space Poisson equation, respectively.The cubic reduction applies when the discrete energy and metric use matching diagonal weights.Source: §3.5 · p. 10
  • With an exact projection, standard Douglas–Rachford theory applies under its usual solvability assumptions; termination requires checking proximal-output proximity and the constraint residual.An intermediate projected density need not be nonnegative, and the kinetic proximal map has the appropriate extended-value domain.Source: §3.5 · p. 10
  • A general nonlocal interaction need not admit a pointwise proximal calculation.Source: §3.5 · p. 10

4 Discretization

The discretization combines a time–space finite-element framework with recovered derivatives to construct a constraint-projection update. It specifies the weighted metric, boundary treatment, residual checks, and safeguards needed to use that update in ISTA, FISTA, and Douglas–Rachford splitting.

4 Discretization

  • The discretization uses a shape-regular triangulated surface and a uniform piecewise-linear time grid, with trapezoidal temporal weights and a consistent spatial mass matrix.The optimization metric and local energy use matching lumped weights.Source: §4 · p. 10, §4.1 · p. 11

4.1 Inner products, energy, and temporal boundary conditions

  • The inner products define metric gradients, while temporal boundary conditions add ghost-point corrections to the Poisson right-hand side.Endpoint interpolants are normalized separately; that approximation error is retained in the analysis.Source: §4.1 · p. 11

4.2 Gradient-enhanced approximation of the projection

  • Temporal PPR and spatial PPPR provide recovered derivatives for the projection update, but the recovered divergence is not asserted to be the metric negative adjoint of the recovered gradient.The endpoint reset can differ from the prescribed Neumann datum, and the update includes neither density clipping nor time-slice mass renormalization.Source: §4.2 · pp. 12–13

4.3 Constraint residual and exact discrete projection

  • The exact W_h-orthogonal projection is defined by the weak constraint matrix as a comparison, not as the system solved by the recovered FDM–sFEM update.Centering ensures Poisson solvability but does not enforce continuity, so the discarded constant must be monitored as a compatibility defect.Source: §4.3 · p. 13, §4.2 · p. 13

4.4 Gradient-enhanced algorithms

  • The recovered update is incorporated into ISTA, FISTA, and Douglas–Rachford splitting, with domain, residual, and step safeguards.Explicit methods require a positive initial density; stopping diagnostics are not optimality certificates for the recovered iteration.Source: §4.4 · pp. 14–15

5 Discretization error analysis

Section 5 analyzes recovered Poisson-solve accuracy and its effect on outer iterations, deriving conditional finite-iteration perturbation bounds that retain data, boundary, solver, projection, and local-proximal defects. Its Douglas–Rachford bound is finite-step and does not establish defect decay, global convergence, or an unconditional accelerated rate.

5.1 Lifting, interpolation, and reference problems

  • For a smooth pair, interpolation bounds separate temporal Cτ²∥∂ttu∥L²(Q) and spatial Ch²∥u∥H²(T;H²(M)) terms, while the error also includes reconstructed tangent-plane error.An FDM–sFEM solution has an additional equation error and cannot be substituted for the interpolant without further analysis.Source: §5.1 · p. 16
  • The recovered-derivative estimate follows cited theory under its regularity and mesh hypotheses, not a new superconvergence theorem; uncovered mean-mode, boundary, or transfer errors remain in εrec,h.The cited bound requires the irregularity and chart conditions, not merely shape regularity; Proposition 5.1 is restricted to applicable mean-zero estimates.Source: §5.1 · pp. 16–17

5.2 Projection error, data perturbation, and residuals

  • For the affine recovered map, feasibility requires AhHh = 0 and Ahkh = bh, while idempotence additionally requires Hh^2 = Hh and Hhkh = 0; these conditions do not establish orthogonality.The residual identities provide matrix checks for feasibility and repeated application, not proof that the map is an exact projection.Source: §5.2 · p. 18

5.3 Finite-iteration bounds for basic ISTA and FISTA

  • For fixed K, conditional ISTA/FISTA bounds propagate projection and half-step errors; recovery contributes O(h^(1+min{1,σ}) + τ^2), while data, endpoint, consistency, and solver errors remain explicit.The estimate requires the stated regularity and trajectory assumptions; it does not establish a mesh-independent rate as K tends to infinity.Source: §5.3 · p. 19

5.4 Douglas–Rachford perturbations

  • With proper, lower-semicontinuous convex Y and fixed parameters, the exact Douglas–Rachford auxiliary-variable map is nonexpansive.This statement uses the exact projection P; the supplied passage does not give a perturbation bound for the recovered map.Source: §5 · p. 15, §5.4 · p. 20
  • Resolving the projection defect requires its input to have the stated smoothness; this does not impose that smoothness on an arbitrary discrete auxiliary variable.Source: §5.4 · p. 20
  • The finite-step Douglas–Rachford bound makes accumulated projection and local-proximal defects explicit, including scalar-solver and local proximal-solver errors.These are treated as total lifted-map defects.Source: §5.4 · p. 20
  • The bound does not show that defects decrease or are summable on a fixed mesh; inexact proximal-gradient results impose specific error criteria and decay conditions [31].Source: §5.4 · p. 20
  • The finite-iteration estimates do not assert global convergence or an unconditional accelerated rate for the recovered algorithms.Source: §5.4 · p. 20

6 Numerical Experiments

The experiments illustrate gradient-enhanced ISTA, FISTA, mean field planning, and Douglas–Rachford splitting on spherical and algebraic surfaces. Recorded refinement values show faster error decrease for gradient-enhanced variants over several refinements, but do not establish a uniform advantage or reproducible convergence study.

6 Numerical Experiments

  • The retained numerical illustrations and refinement records do not constitute a reproducible convergence or timing study, or independently validate the corrected boundary treatment, safeguards, or perturbation bounds.Original solver configurations and raw error data are unavailable.Source: §6 · p. 21

6.1 Gradient enhanced ISTA

  • ISTA density patterns differ across three stepsizes, especially at η = 1.8 × 10−4, motivating careful stepsize selection and residual monitoring without quantifying accuracy loss.The displayed patterns do not identify the source of the discrepancy.Source: §6.1 · p. 21, Figure 1 · p. 21

6.2 FISTA

  • FISTA’s two displayed density sequences appear qualitatively similar, but differing endpoints and unavailable backtracking history prevent a controlled comparison of acceleration or stepsize independence.The record associates these calculations with backtracking but gives neither accepted stepsizes nor the backtracking history.Source: §6.2 · p. 22, Figure 2 · p. 22

6.3 Mean Field Planning

  • The spherical planning illustrations compare four interaction integrands—zero, ρ2/2, ρ log ρ, and 1/ρ—but do not isolate the effects of cost, weight, or numerical parameters.The entropy derivative is singular at zero, while the reciprocal integrand is infinite there, requiring different domain checks from zero-interaction OT.Source: §6.3 · p. 22

6.4 Proximal Splitting

  • Douglas–Rachford splitting illustrates density redistribution between indicator-type endpoint densities near high-curvature corners of the Enzensberger–Stern surface, but snapshots do not certify optimality or individual geodesics.The nonsmooth endpoints lie outside the smooth-data recovery assumptions, and reproduction requires support regions, normalization, and solver tolerances.Source: §6.4 · pp. 23–24, Figure 4 · p. 24

6.5 Recorded refinement values

  • Gradient-enhanced errors decrease faster over several displayed refinements, but all three variants have larger errors at N = 8 than their baselines; the records support neither uniform advantage nor a common first-order characterization of all baselines.The last ADMM pair has an observed order near two.Source: §6.5 · p. 24

7 Conclusions

The paper presents a shared gradient-enhanced approximation of the continuity-constraint projection for variational mean field planning on surfaces and incorporates it into three splitting algorithms. Its analysis distinguishes the recovered update from an exact projection and gives conditional finite-iteration bounds, while identifying controls needed before stronger convergence or discretization claims.

7 Conclusions

  • The construction retains time finite differences and surface finite elements, using temporal PPR and spatial PPPR to update density and momentum.The same recovered quantities are used in ISTA, FISTA, and Douglas–Rachford splitting.Source: §7 · p. 25
  • The analysis distinguishes the recovered update from an exact projection and gives a matched-data error decomposition, defect tests, and conditional finite-iteration bounds for ISTA/FISTA and Douglas–Rachford.Higher-order recovered derivatives contribute under the stated regularity and mesh hypotheses.Source: §7 · p. 25
  • Input differentiation, endpoint enforcement, geometric approximation, and solver errors must be controlled before inferring an overall discretization rate or a convergent outer iteration.Source: §7 · p. 25
  • The numerical illustrations show intended surface transport and interaction-cost settings, while recorded refinement values motivate further quantitative verification.Source: §7 · p. 25
  • A complete assessment requires error definitions, reference solutions, conservation residuals, and computational costs; uniform stability, nonsmooth or vacuum data, and adaptive branch changes remain open directions.Source: §7 · p. 25
Loading 2609.10370v1…