PaperlayerView exact revision on arXiv

Source-linked AI summary

Structure-Preserving Detailed-Balance Master-Equation Discretizations for Fokker--Planck Equations

Satish Chandran, Yiwei Wang

arXiv:2608.30121v1math.NAcs.NA

TL;DR

Structure-preserving Fokker–Planck discretizations must jointly satisfy mass, positivity, free-energy, and equilibrium requirements, while existing constructions often depend on particular PDE forms. The paper derives reversible master-equation schemes from discrete energy-dissipation and detailed-balance conditions, with numerical tests reporting equilibrium accuracy and expected spatial orders.

  • Problem

    Structure-preserving Fokker–Planck discretizations must jointly conserve mass, preserve nonnegativity, dissipate discrete free energy, and reproduce equilibrium, while existing derivations often depend on particular PDE forms.

    Source: §1 · p. 2

  • Method

    The paper represents density by a finite probability vector and determines edge jump rates using the discrete energy-dissipation flux and detailed balance.

    Source: §1 · p. 2, §3 · p. 12

  • Results

    Tests report better Gibbs-state accuracy than upwind and central differences in smooth potentials, and near-second-order transient and equilibrium orders in nonlinear saturation cases.

    Source: §6.1 · p. 26, §6.3 · p. 30

  • Takeaways & Limitations

    Backward Euler and semi-implicit variants provide linear, mass-conservative, positivity-preserving updates; nonlocal energy decay is unconditional for negative-semidefinite kernels, but otherwise requires a kernel-dependent time-step condition.

    Source: §7 · p. 33

    General convergence to equilibrium is not addressed and may remain open even for the continuous PDE, particularly with irregular interaction kernels.

    Source: §5.2 · p. 22

Abstract

from arXiv · show

We develop a variational--Markov construction of detailed-balance master-equation discretizations for Fokker--Planck equations directly from their energy--dissipation laws. Rather than discretizing the differential operator, we represent mass transfer between neighboring grid points by two directional jump rates. On each edge, the discrete energy--dissipation law determines the net flux, while detailed balance fixes the rate ratio; a local reconstruction of edge quantities from neighboring grid-point values then uniquely determines both rates. A logarithmic-mean reconstruction recovers the classical Scharfetter--Gummel/Wang--Peskin--Elston rates, while alternative reconstructions yield other reversible schemes. The resulting semi-discrete systems conserve mass, preserve nonnegativity, satisfy detailed balance, and dissipate a discrete free energy. The same edgewise construction extends to state-dependent and degenerate mobilities, nonlocal interaction energies and higher-dimensional problems. The corresponding implicit and semi-implicit fully discrete schemes are linear, mass-conservative, positivity-preserving, and energy-stable, with a time-step restriction required only for nonlocal interaction energies whose kernel is not negative semidefinite. Numerical experiments confirm the expected spatial convergence rates of all detailed-balance schemes and verify their equilibrium accuracy, positivity preservation, and free-energy decay with discontinuous potentials, saturation, nonlocal interactions, and two-dimensional problems.

1. Introduction

The introduction motivates discretizations that preserve mass, nonnegativity, free-energy dissipation, and correct equilibrium, then presents an edgewise variational–Markov construction based on energy dissipation and detailed balance.

1. Introduction

  • Existing schemes often start from PDE-specific forms, while mass conservation, nonnegativity, energy decay, and correct equilibrium remain distinct requirements.A method can dissipate energy yet converge to an inaccurate equilibrium, or reproduce equilibrium without comparable transient accuracy.Source: §1 · p. 2
  • In one-dimensional linear constant-diffusion settings with intervalwise constant drift, Il’in–Allen–Southwell, SG, CC, and WPE reduce to the same Bernoulli-function flux.The WPE scheme was also interpreted as a detailed-balance master equation [14].Source: §1 · p. 2
  • The method avoids discretizing the differential operator: energy dissipation fixes each edge’s net flux, detailed balance fixes the rate ratio, and local reconstruction uniquely determines both directional rates.The rate ratio determines the discrete equilibrium, while the common edge prefactor controls transient evolution.Source: §1 · pp. 2–3
  • A logarithmic-mean reconstruction recovers classical SG/WPE rates; other local means yield reversible schemes, and the construction extends to nonlinear mobilities, nonlocal energies, and multiple dimensions.The stated extensions also include degenerate mobilities and general f-divergence free energies.Source: §1 · p. 3
  • Implicit and semi-implicit updates are linear, mass-conservative, positivity-preserving, and energy-dissipative; general nonlocal interaction kernels require a time-step condition for energy stability.For fixed potentials the backward-Euler update has these properties for every time-step size; state-dependent mobilities retain them when edge mobilities are degenerate.Source: §1 · p. 3

2. Preliminaries

This section connects Fokker–Planck equations and finite-state master equations through energy–dissipation laws and spatial discretization. It establishes detailed balance and derives transition rates using local stationary problems or weighted-diffusion transforms, including rates that depend on how the potential is averaged across neighboring grid points.

2.1. Variational Structure of Fokker–Planck-type equations.

  • For generalized diffusions, admissible states, free energy, dissipation, and transport kinematics determine the evolution, with flow-map kinematics encoding mass and nonnegativity constraints.The energetic variational approach combines least action and maximum dissipation to obtain a force balance and generalized Darcy law.Source: §2.1 · p. 4

2.2. Variational structure of Master Equations with detailed balance.

  • For master equations, detailed balance with a positive stationary distribution yields an energy-dissipation law for convex f-divergences, with the KL case giving a specific dissipation form.The discrete free energy is an f-divergence relative to the stationary distribution.Source: §2.2 · p. 6

2.3. Spatial discretization of a Fokker–Planck equation as a master equation.

  • A finite-difference discretization turns grid masses into a nearest-neighbor master equation, whose desired stationary state is discrete Gibbs and whose rates should satisfy detailed balance.The Fokker–Planck equation and master equation are also connected as continuum and discrete descriptions.Source: §2.3 · p. 8
  • Finite-volume schemes conserve mass by construction but require cell-averaged Gibbs equilibria or a coarse-grained potential, so the authors use pointwise finite differences.The pointwise choice uses ρ_i ≈ ρ(x_i, t) and a simple discrete equilibrium.Source: §2.3 · pp. 8–9
  • Central differences can lose reversibility for discontinuous potentials and positivity when |δV|>2; upwind rates stay positive but their detailed-balance mismatch persists across potential jumps.For smooth potentials, central and upwind rate-ratio discrepancies are asymptotically small, but need not vanish across discontinuities.Source: §2.3 · p. 9
  • In linear constant-diffusion settings, the Bernoulli-function flux also appears as the Il’in–Allen–Southwell [20, 21], Scharfetter–Gummel, and Chang–Cooper schemes.Source: §2.3 · p. 10
  • The Bernoulli-function transition rates satisfy detailed balance at every mesh size, even for arbitrarily large potential differences, and both rates are positive.The local stationary exponential profile motivates the term exponential fitting.Source: §2.3 · p. 10
  • The Slotboom transform converts the equation to weighted diffusion; discretizing it gives rates using an averaged neighboring potential, and different averages yield different schemes.The same reformulation has been used to construct positivity-preserving, free-energy-dissipating schemes for nonlinear and nonlocal diffusion equations [30, 31, 52].Source: §2.3 · p. 10

3. Variational Construction of Jump Rates

The section constructs conservative master-equation discretizations from the Fokker–Planck variational structure, using flux consistency, detailed balance, and local reconstructions to determine reversible edge rates. It establishes variational consistency for five reconstructed rates, with second-order accuracy in discrete free energy and dissipation.

3. Variational Construction of Jump Rates

  • Projected energy–dissipation fixes each edge’s net flux, while detailed balance fixes the forward-to-backward rate ratio; local reconstructions then determine the jump rates.The projection yields a conservative semi-discrete master equation for a Markov jump process.Source: §3 · pp. 10–12
  • Replacing the exact density-dependent edge coefficient with a potential-based approximation gives rates independent of the evolving density, keeping the linear master equation’s generator fixed.The approximation uses a local reconstruction of the potential at the edge.Source: §3 · p. 12

3.1. Representative reconstructions.

  • A logarithmic-mean reconstruction recovers classical SG/WPE rates, while midpoint and harmonic-mean choices yield ED and IWPE variants, respectively.The detailed-balance multiplicative symmetry ensures the sampled Gibbs state is an exact stationary state.Source: §3.1 · pp. 14–15

3.2. Variational Consistency.

  • All five reconstructed schemes satisfy detailed balance and the discrete energy–dissipation law; their coefficients approximate the same limiting coefficient to second order under the stated smoothness assumptions.The continuum-limit coefficient is γ(x) = e^-V(x).Source: §3.2 · pp. 16–17

3.2. Variational Consistency

  • The variational-consistency analysis covers the semi-discrete master equation with all five reconstructed rates.The section states that consistency is established for the five reconstructed rates and refers to any of the five reconstructed coefficients.Source: §3.2 · pp. 17–18
  • For the stated positive periodic ρ, Fh[Phρ] and Δh[Phρ] each approximate their continuous counterparts to O(h^2).The proposition gives Fh[Phρ] = F[ρ] + O(h^2) and Δh[Phρ] = Δ[ρ] + O(h^2).Source: §3.2 · p. 18

4. Temporal Discretizations

The temporal schemes trade an explicit CFL-type time-step restriction for implicit Euler's unconditional positivity preservation, mass conservation, and free-energy stability. Crank–Nicolson is second-order accurate in time and mass-preserving, but positivity and free-energy dissipation are not generally unconditional.

4. Temporal Discretizations

  • The semi-discrete system is a conservative linear system with a Metzler generator, and the forward Euler update can be interpreted as a Markov chain on grid states.The chain's probability mass at state i is represented by p_i^n.Source: §4 · pp. 18–19
  • Explicit forward Euler requires a CFL-type restriction, Δt ≤ 1/max_i(F_i+1/2 + B_i−1/2), to keep its update coefficients nonnegative.For constant potential, this reduces to Δt ≤ h^2/2.Source: §4 · p. 19
  • For any time step, implicit Euler preserves positivity and mass and makes the discrete free energy non-increasing, assuming a fixed Markov generator satisfying detailed balance.The generator's detailed balance is with respect to the discrete Gibbs state.Source: §4 · p. 19
  • Crank–Nicolson is second-order accurate in time and preserves mass, but positivity and free-energy dissipation may require a time-step restriction.Those properties are generally not unconditional for this scheme.Source: §4 · p. 20

5. Extensions to More General Settings

The edgewise framework extends the discrete energy–dissipation and detailed-balance construction to nonlinear mobilities, nonlocal interaction energies, higher-dimensional grids, unstructured meshes, and weighted graphs. Interpreting the resulting graph operator as a continuous Fokker–Planck discretization requires separate geometric and graph-to-continuum analysis.

5.1. Nonlinear Mobility.

  • For nonlinear mobility, a symmetric consistent mean controls its approximation while the reconstructed rates preserve detailed balance exactly.The frozen-coefficient semi-implicit update conserves mass, preserves nonnegativity, and dissipates free energy; if edge coefficients vanish, the chain may become reducible and the Gibbs state need not be unique.Source: §5.1 · p. 21

5.2. Nonlocal Interaction Energies.

  • For nonlocal interactions, freezing the self-consistent potential gives a linear semi-implicit update that preserves mass and nonnegativity for every time step.Its fixed points are self-consistent discrete equilibria, but general convergence to equilibrium is not addressed; the scheme is formally first-order in time and second-order in space, without a uniform-in-grid stability proof.Source: §5.2 · pp. 22–23
  • The nonlocal update dissipates the true discrete free energy for every time step when the interaction operator is negative semidefinite on the zero-mass subspace; otherwise, a time-step condition suffices.The frozen step dissipates a lagged free energy, and the difference from the true-energy increment is a quadratic remainder; the sufficient restriction in the general case depends on kernel strength.Source: §5.2 · p. 23

5.3. Higher-dimensional cases.

  • Higher-dimensional rates are constructed edge by edge in each coordinate direction, and directional splitting can replace a coupled multidimensional solve if each substep preserves the structural properties and discrete Gibbs state.Splitting introduces temporal error; an ordered composition is generally not reversible unless the directional operators commute or the composition is palindromic.Source: §5.3 · pp. 24–25
  • On each coordinate edge, the construction extends the one-dimensional structural properties to multidimensional schemes that preserve mass and nonnegativity and satisfy detailed balance.Under periodic or no-flux boundary conditions, detailed balance determines the discrete Gibbs state; dimensional splitting can avoid a coupled solve but introduces temporal splitting error and may lose reversibility unless operators commute or the composition is palindromic.Source: §5.3 · pp. 24–25
  • On weighted graphs, edgewise energy–dissipation and detailed balance are algebraically well defined, but continuum fidelity depends on geometric encoding; convergence is established for regular admissible finite-volume meshes.Node weights and edge transmissibilities encode domain geometry; identifying a continuum operator on a general graph remains a separate graph-to-continuum problem.Source: §5.3 · p. 25

6. Numerical Experiments

The experiments test five detailed-balance schemes across smooth and discontinuous potentials, nonlinear saturation, nonlocal aggregation–diffusion, and a two-dimensional case.

6.1. Smooth potential case.

  • All five DB schemes converge at second order for smooth-potential transients and reach round-off Gibbs-state error, unlike upwind and central differences.The upwind and central-difference equilibrium errors are 1.18 × 10^-2 and 1.12 × 10^-4, respectively.Source: §6.1 · p. 26

6.2. Discontinuous potential. (

  • At the potential jump, IWPE remains second-order transiently while the other four DB schemes and upwind are about first-order; all DB schemes share accurate equilibrium.The reported finest-pair ℓ1 orders are 2.00 for IWPE, 1.06–1.10 for SG/WPE, ED, AM-e−V, and LM-e−V, and 0.99 for upwind.Source: §6.2 · p. 28

6.3. Nonlinear mobility and saturation.

  • Across subcritical and supercritical saturation, DB schemes achieve transient and equilibrium orders near two, preserve ρ≤α to round-off, and harmonic edge mobility suppresses transport into saturated nodes.The reported transient orders are 1.91–2.02 and equilibrium orders are 1.98–2.00 across both mass regimes.Source: §6.3 · p. 30

6.4. Nonlocal aggregation–diffusion.

  • For nonlocal aggregation–diffusion, all methods preserve positivity and mass and approach the fixed point; DB transient orders are 2.02–2.03 versus upwind 1.05.At T = 4, ℓ1 distances to the fixed point range from 5.52 × 10^-13 to 2.53 × 10^-12 across the six methods.Source: §6.4 · p. 31

6.5. Two-dimensional case.

  • In two dimensions, DB schemes are second-order transiently and at round-off equilibrium error; upwind is first-order transiently, while central differences are second-order but less accurate at equilibrium.At T = 1.5, the upwind and central-difference Gibbs errors are 7.1×10^-2 and 3.5×10^-3, respectively.Source: §6.5 · p. 32

7. Conclusion

The variational–Markov construction separates edgewise flux, detailed balance, and rate reconstruction, clarifying their distinct roles in transient dynamics and equilibrium. It yields reversible schemes with structural guarantees in semi-discrete and fully discrete settings, and extends to nonlinear and nonlocal models and multiple dimensions.

7. Conclusion

  • The construction separates each edge’s net flux, directional-rate ratio, and common symmetric coefficient among the projected energy–dissipation law, detailed balance, and local reconstruction.This separation distinguishes how transient dynamics and equilibrium are determined.Source: §7 · pp. 32–33
  • Five second-order reconstructions yield exact sampled Gibbs stationary states; SG/WPE additionally satisfies the classical B-scheme additive consistency relation.The five are ED, SG/WPE, IWPE, AM-e−V, and LM-e−V.Source: §7 · p. 33
  • Backward Euler for fixed potentials is linear and uniquely solvable, preserving mass, nonnegativity, and discrete free-energy dissipation for every time-step size.For nonlocal interactions, the lagged-potential scheme dissipates true energy unconditionally for negative-semidefinite kernels; general symmetric kernels require a kernel-dependent time-step condition.Source: §7 · p. 33
  • The edgewise construction extends to state-dependent and vanishing mobilities, nonlocal interactions, general f-divergence energies, and Cartesian grids in multiple dimensions.The experiments cover smooth and discontinuous potentials, nonlinear mobility, nonlocal interactions, and two-dimensional problems.Source: §7 · p. 33
  • Convergence theory, especially for discontinuous potentials, nonlinear mobilities, and nonlocal interactions, remains open; high-dimensional sampling is another proposed direction.The authors suggest energy–dissipation estimates for continuum-limit analysis and detailed-balance-preserving transition rates for sampling.Source: §7 · pp. 33–34
Loading 2608.30121v1…