Paper deep dive
Discrete Wasserstein Flows for One-Step Generative Modeling
Alessandro Micheli, Andrea Zerio, Samir Bhatt
Intelligence
Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 85%
Last extracted: 10/3/2026, 2:55:30 AM
Summary
This paper introduces Discrete Drifting, a framework for one-step generative modeling on finite state spaces. It extends the concept of drifting models to discrete data by utilizing discrete Wasserstein geometry (Maas, 2011) to define a target-relative KL gradient flow over reversible Markov kernels. The method realizes this probability flow at the particle level through Markov jumps during training, allowing a latent-conditioned neural generator to amortize these transport updates. This enables one-step inference, where the iterative dynamics are bypassed after training. The authors provide theoretical validation of KL dissipation and consistency between particle dynamics and probability flow, along with empirical results showing that finite-capacity neural generators can track exact transport targets.
Entities (8)
Relation Signals (7)
Discrete Wasserstein Geometry → proposedby → Maas (2011)
confidence 99% · We build both from the discrete Wasserstein geometry of Maas (2011).
Discrete Drifting → enables → One-Step Generative Modeling
confidence 95% · We introduce a new framework for one-step generative modelling on finite state spaces.
Discrete Drifting → uses → Discrete Wasserstein Geometry
confidence 95% · we use discrete Wasserstein geometry to define a target-relative KL gradient flow
Discrete Drifting → realizesflowvia → Markov Jumps
confidence 92% · We realize this probability flow at the particle level through Markov jumps
Neural Generator → amortizes → Transport Updates
confidence 90% · amortize the resulting transport updates into a latent-conditioned generator
Discrete Drifting → extends → Drifting Models
confidence 90% · Extending this idea to discrete data is not straightforward... Drifting Models (Deng et al., 2026) amortize this evolution into training
Discrete Drifting → optimizes → KL Divergence
confidence 90% · define a target-relative KL gradient flow... verify KL dissipation
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:We introduce a new framework for one-step generative modelling on finite state spaces. To extend drifting beyond continuous domains, we use discrete Wasserstein geometry to define a target-relative KL gradient flow over the transitions of a reversible Markov kernel. We realize this probability flow at the particle level through Markov jumps and amortize the resulting transport updates into a latent-conditioned generator, so that the iterative dynamics are required only during training while inference remains one-step. In a controlled setting where the underlying distributions and transport dynamics can be computed exactly, we verify KL dissipation, consistency between the particle dynamics and the probability flow, and the predicted numerical scaling. We further show that a finite-capacity neural generator can track these exact transport targets while retaining one-step generation. These results validate the basic construction and provide a foundation for scaling Discrete Drifting to structured discrete data.
Tags
Links
- Source: https://arxiv.org/abs/2610.01355v1
- Canonical: https://arxiv.org/abs/2610.01355v1
Trouble viewing inline? Open PDF directly →
Full Text
88,713 characters extracted from source content.
Expand or collapse full text
Geometric Distributional Deep Learning (GDDL) Discrete Wasserstein Flows for One-Step Generative Modeling Alessandro Micheli Affiliation: Imperial College London Affiliation: London, United Kingdom Email: a.micheli19@imperial.ac.uk Andrea Zerio Affiliation: Department of Computer Science, Aalborg University, Copenhagen, Denmark; Affiliation: Centre for Frontier AI Research (CFAR), Affiliation: Institute of Advanced Intelligence and Computing (IAIC), Affiliation: A*STAR, Singapore anze@cs.aau.dk Samir Bhatt Affiliation: Imperial College London; University of Copenhagen Affiliation: London, United Kingdom; Copenhagen, Denmark Email: s.bhatt@imperial.ac.uk Abstract We introduce a new framework for one-step generative modelling on finite state spaces. To extend drifting beyond continuous domains, we use discrete Wasserstein geometry to define a target-relative KL gradient flow over the transitions of a reversible Markov kernel. We realize this probability flow at the particle level through Markov jumps and amortize the resulting transport updates into a latent-conditioned generator, so that the iterative dynamics are required only during training while inference remains one-step. In a controlled setting where the underlying distributions and transport dynamics can be computed exactly, we verify KL dissipation, consistency between the particle dynamics and the probability flow, and the predicted numerical scaling. We further show that a finite-capacity neural generator can track these exact transport targets while retaining one-step generation. These results validate the basic construction and provide a foundation for scaling Discrete Drifting to structured discrete data. 1 Introduction Diffusion and flow-based generators construct samples through iterative inference-time transport, whereas Drifting Models (Deng et al., 2026) amortize this evolution into training, enabling one-step generation. Extending this idea to discrete data is not straightforward. In Euclidean space, drifting relies on infinitesimal particle displacements x↦x+ηV(x)x x+η V(x); on a finite state space =x1,…,xNX=\x_1,…,x_N\, however, particles cannot move infinitesimally and probability must instead be transferred between distinct states. This obstruction is also geometric: when the underlying state space is finite, the ordinary Wasserstein-22 metric does not induce a useful infinitesimal transport geometry on the corresponding space of probability distributions. Indeed, transporting O(δt)O(δ t) mass across a positive distance produces W2=O(δt)W_2=O( δ t), yielding a divergent metric derivative as δt→0δ t→ 0 (Maas, 2011). A discrete analogue of drifting therefore requires both an appropriate transport geometry and a particle dynamics compatible with it. We build both from the discrete Wasserstein geometry of Maas (2011). An irreducible, reversible Markov kernel K defines the admissible local transitions on X, while the resulting geometry provides a notion of gradient flow for probability distributions on the finite state space. Choosing ℱp(q)=KL(q∥p),F_p(q)=KL(q\|p), with p the data distribution, yields a flow that redistributes probability along the transitions of K so as to decrease the target-relative KL. We interpret this redistribution at the particle level through a Markov jump process, providing the discrete counterpart of the continuous drift V(x)V(x). Crucially, these jumps are used only during training to construct updated targets: unlike discrete diffusion models (Campbell et al., 2022; Lou et al., 2024), generation does not require simulating a Markov process at inference time. Instead, a latent-conditioned generator is trained to absorb the successive transport updates and directly produce samples in a single forward pass. Our contributions are threefold: (1) we formulate drifting on finite state spaces using discrete Wasserstein geometry; (2) we derive a Markov-jump particle realization of the resulting probability flow that can be used to construct stop-gradient training targets; and (3) we provide a controlled empirical validation that separately tests the underlying discrete transport and its neural amortization. 2 Background We consider probability evolutions on a finite space =x1,…,xN,X=\x_1,…,x_N\, equipped with a Markov transition kernel K∈ℝN×NK ^N× N, that is, a non-negative, row-stochastic real matrix. We assume throughout that K is irreducible and reversible. In particular, irreducibility guarantees the existence of a unique invariant probability measure π∈ℝNπ ^N, with strictly positive entries, satisfying π⊤K=π⊤.π K=π . We denote by ()P(X) the set of probability densities on X relative to π. We further write ∗()P_*(X) for the subset of strictly positive densities. Throughout the paper, the symbol ρ is used for densities relative to π, whereas p denotes ordinary probability masses. The two representations are related by p(x)=π(x)ρ(x),x∈,p(x)=π(x)ρ(x), x , (1) or, equivalently, p=ρ⊙πp=ρ π, where ⊙ denotes componentwise multiplication. A canonical example of the probability evolutions considered in this work is the discrete heat equation dtρt=(K−I)ρt. ddt _t=(K-I) _t. (2) This dynamics is widely used as a forward noising process in discrete generative modelling (Campbell et al., 2022; Lou et al., 2024). Our objective is to interpret evolutions such as equation 2 as gradient flows of functionals defined on ()P(X), in analogy with variational constructions on continuous spaces (Jordan et al., 1998; Bunne et al., 2022; Terpin et al., 2024). In the Euclidean setting, probability evolutions can be characterized as gradient flows of an energy functional ℱF with respect to the 22-Wasserstein geometry. The Jordan–Kinderlehrer–Otto (JKO) scheme describes this evolution through the proximal recursion pk+1∈argminp∈(ℝd)ℱ(p)+12τW22(p,pk).p_k+1∈ *arg\,min_p (R^d) \F(p)+ 12τW_2^2(p,p_k) \. (3) A common choice is to consider ℱF of the free-energy form ℱ=+βℋ,F=V+ , (4) where V is a potential energy, ℋH is an entropy functional, and β≥0β≥ 0 controls the strength of the entropic contribution. The standard W2W_2 distance does not provide a suitable gradient-flow geometry on finite state spaces. This limitation appears even in the simplest non-trivial case of a two-point state space; see (Rancati et al., 2026, Lemma 2.1) and (Maas, 2011, Remark 2.1). To address this limitation, the transport geometry of Maas (2011) provides the basis for a discrete analogue of the JKO framework by replacing W2W_2 with a metric adapted to the transition structure of the Markov kernel K. The construction of this metric parallels the dynamic formulation of optimal transport due to Benamou and Brenier (2000). Its central ingredient is a mobility associated with each edge of the graph induced by K, which determines how the cost of transporting mass between x,y∈x,y depends on their current densities. In this work we use the logarithmic mean m(a,b):=∫01asb1−s,s,m(a,b):= _0^1a^sb^1-s,ds, and define the corresponding mobility by mρ(x,y):=m(ρ(x),ρ(y)).m_ρ(x,y):=m (ρ(x),ρ(y) ). This choice is particularly natural here because, with the corresponding transport metric, the discrete heat equation in equation 2 is recovered as the gradient flow of entropy. The logarithmic mean is one of several possible symmetric mobilities; alternatives include the geometric mean and suitable power-law mobilities (Maas, 2011; Erbar and Maas, 2014; Erbar et al., 2018). The mobility is therefore a modelling choice that determines which probability evolution is represented as a gradient flow under the resulting transport geometry. Using the logarithmic-mean mobility, Maas (2011) defines a dynamic transport distance WKW_K on ()P(X) by minimizing a discrete kinetic action over paths connecting two probability densities. Specifically, for ρ0,ρ1∈() _0, _1 (X), consider piecewise continuously differentiable curves ρ:[0,1]→()ρ:[0,1] (X) and measurable curves ψ:[0,1]→ℝNψ:[0,1] ^N satisfying, for almost every t∈[0,1]t∈[0,1], the discrete continuity equation ρ˙t(x)+∑y∈(ψt(y)−ψt(x))K(x,y)mρt(x,y)=0,x∈, ρ_t(x)+ _y ( _t(y)- _t(x) )K(x,y)m_ _t(x,y)=0, x , (5) with boundary conditions ρt=0=ρ0 _t=0= _0 and ρt=1=ρ1 _t=1= _1. The squared transport distance is then defined as the minimum action over all such admissible pairs: WK(ρ0,ρ1)2:=infρ,ψ12∫01∑x,y∈ W_K( _0, _1)^2:= _ρ,ψ \ 12 _0^1 _x,y (ψt(x)−ψt(y))2K(x,y)mρt(x,y)π(x)dt. ( _t(x)- _t(y) )^2K(x,y)m_ _t(x,y)π(x)\,dt \. (6) The quantity in equation 6 is the discrete counterpart of the kinetic action in the Benamou–Brenier formulation. Accordingly, WK(ρ0,ρ1)2W_K( _0, _1)^2 is the minimum kinetic energy required to transport ρ0 _0 into ρ1 _1, subject to conservation of mass through equation 5. The difference ψt(x)−ψt(y) _t(x)- _t(y) plays the role of a velocity along the edge connecting x and y, while K(x,y)mρt(x,y)π(x)K(x,y)m_ _t(x,y)π(x) determines the corresponding density-dependent transport weight. Maas (2011) shows that WKW_K defines a transport metric on ()P(X) and induces a Riemannian structure on its interior ∗()P_*(X). Crucially, when the mobility mρm_ρ is chosen to be the logarithmic mean, the discrete heat equation 2 is the gradient flow, with respect to WKW_K, of the relative entropy ℋ(ρ):=∑x∈ρ(x)logρ(x)π(x),H(ρ):= _x ρ(x) ρ(x)\,π(x), (7) defined relative to the invariant measure π, with the convention 0log0=00 0=0. a 3 Discrete Drifting on Finite Spaces Section 2 reviewed the discrete Wasserstein geometry of Maas (2011), under which the heat equation (equation 2) is the WKW_K-gradient flow of the relative entropy ℋ(ρq)=KL(q∥π)H( _q)=KL(q\,\|\,π). Consequently, the heat flow relaxes the model distribution q toward the invariant measure π. For generative modelling, however, the desired equilibrium is the data distribution p, which need not coincide with π. We therefore replace the entropy relative to π by the target-relative functional ℱp(ρq):=∑x∈ρq(x)logρq(x)ρp(x)π(x)=KL(q∥p),F_p( _q):= _x _q(x) _q(x) _p(x)\,π(x)=KL(q\,\|\,p), (8) where ρp:=p/π _p:=p/π and ρq:=q/π _q:=q/π. Throughout this section we assume ρp,ρq∈∗() _p, _q _*(X), so that the logarithmic density ratios are well defined. This is an entropy-plus-potential free energy with Vp=−logρpV_p=- _p, so the finite-state gradient calculus of Maas (2011) applies directly. The corresponding WKW_K-gradient flow is the discrete continuity equation with potential ψp,q=log(ρp/ρq) _p,q= ( _p/ _q). Equivalently, defining the edge current Jp,q(x,y):=π(x)K(x,y)mρq(x,y)(logρq(x)ρp(x)−logρq(y)ρp(y)),J_p,q(x,y):=π(x)K(x,y)m_ _q(x,y) ( _q(x) _p(x)- _q(y) _p(y) ), (9) the induced evolution of probability masses is q˙t(x)=−∑y∈Jp,qt(x,y). q_t(x)=- _y J_p,q_t(x,y). (10) By reversibility of K and symmetry of the logarithmic mean, Jp,q(x,y)=−Jp,q(y,x)J_p,q(x,y)=-J_p,q(y,x), so the flow conserves total probability. The standard gradient-flow dissipation identity specializes here to dtKL(qt∥p)=−ℐK(qt∥p), ddtKL(q_t\,\|\,p)=-I_K(q_t\,\|\,p), (11) where ℐK(q∥p):=12∑x,y∈π(x)K(x,y)mρq(x,y)(logρq(x)ρp(x)−logρq(y)ρp(y))2.I_K(q\,\|\,p):= 12 _x,y π(x)K(x,y)m_ _q(x,y) ( _q(x) _p(x)- _q(y) _p(y) )^2. (12) Since K is irreducible and ρp,ρq>0 _p, _q>0, ℐK(q∥p)=0I_K(q\,\|\,p)=0 only when log(ρq/ρp) ( _q/ _p) is constant on X; normalization then implies q=pq=p. Thus p is the unique stationary point of equation 10 in ∗()P_*(X). These are direct consequences of the Maas gradient-flow framework; we use equation 10 as the distribution-level evolution to be realized by discrete drifting. 3.1 A Normalization-Free Drift in the Logit Chart The gradient flow in equation 10 cannot be used as a training signal in the form given: it is written in terms of the normalized laws q and p, whereas a parametric model emits unnormalized logits and the target is typically known only up to a constant. We now show that the flow admits an exact representation involving only local differences of unnormalized log-densities, and is therefore free of both normalizing constants. Let g:→ℝg:X be a log-density parameter and set qg(x):=π(x)eg(x)Zπ(g),Zπ(g):=∑z∈π(z)eg(z),ρg(x):=qg(x)π(x)=eg(x)Zπ(g).q_g(x):= π(x)e^g(x)Z_π(g), Z_π(g):= _z π(z)e^g(z), _g(x):= q_g(x)π(x)= e^g(x)Z_π(g). (13) For an oriented edge (x,y)(x,y), define the model log-density difference and the target edge score δg(x,y):=g(y)−g(x),sp(x,y):=logρp(y)−logρp(x). _g(x,y):=g(y)-g(x), s_p(x,y):= _p(y)- _p(x). (14) Both depend only on local differences and are therefore independent of the corresponding normalization constants; in particular, the target need only be known up to a multiplicative constant. To express the logarithmic-mean mobility in terms of these differences, define L(r):=m(1,er)=er−1r,r≠0,1,r=0.L(r):=m(1,e^r)= cases e^r-1r,&r≠ 0,\\[8.00003pt] 1,&r=0. cases (15) Transporting the flow of equation 10 to this chart yields a drift that depends only on edge quantities. Define the normalization-free logit drift Δg(x):=∑y∈K(x,y)L(δg(x,y))(δg(x,y)−sp(x,y)), _g(x):= _y K(x,y)L ( _g(x,y) ) ( _g(x,y)-s_p(x,y) ), (16) a K-average over the neighbours of x of the model-minus-target edge discrepancy δg−sp _g-s_p, weighted by the logarithmic-mean factor L. Every quantity on the right-hand side is computable from logits and from the target up to normalization. Proposition 1 (Logit representation of the flow). Let I⊆ℝI be an interval, let t↦gt∈ℝt g_t ^X be continuously differentiable, and set qt:=qgt.q_t:=q_g_t. Then qtq_t solves the probability-coordinate gradient flow in equation 10 if and only if there exists a continuous scalar function c:I→ℝc:I such that g˙t=Δgt+ct, g_t= _g_t+c_t1, (17) where (x)=11(x)=1 for every x∈x . The proof is given in Appendix D. The drift equation 16 splits into two terms with distinct roles, Δg(x)=∑y∈K(x,y)L(δg(x,y))δg(x,y)⏟model repulsion−∑y∈K(x,y)L(δg(x,y))sp(x,y)⏟target attraction. _g(x)= _y K(x,y)L ( _g(x,y) )\, _g(x,y)_model repulsion\;-\; _y K(x,y)L ( _g(x,y) )\,s_p(x,y)_target attraction. (18) The first term depends only on the current model through δg _g; where the model already favours y (δg(x,y)>0 _g(x,y)>0), it raises g(x)g(x), opposing further concentration there. The second term is the target attraction: where the target favours a neighbour y (sp(x,y)>0s_p(x,y)>0), it lowers g(x)g(x), moving mass toward y. This is the finite-space analogue of the attraction–repulsion structure of Euclidean drifting. 3.1.1 Stop-gradient realization Proposition 1 determines the logit dynamics up to an additive gauge. Since qg+c=qgq_g+c1=q_g for every scalar c, all gauge choices induce the same probability law. We therefore choose the zero-gauge representative g˙t=Δgt. g_t= _g_t. Its explicit-Euler discretization with step size η>0η>0 is g+=g+ηΔg.g^+=g+η _g. (19) The next proposition shows that this chart update is first-order consistent with the probability-coordinate gradient flow and reproduces its instantaneous relative-entropy dissipation to first order. Proposition 2 (Local consistency and relative-entropy descent). For fixed g, as η→0η→ 0, qg+ηΔg=qg+ηqg⊙Δg+O(η2),q_g+η _g=q_g+η\,q_g _g+O(η^2), (20) and KL(qg+ηΔg∥p)=KL(qg∥p)−ηℐK(qg∥p)+O(η2).KL (q_g+η _g\,\|\,p )=KL(q_g\,\|\,p)-η\,I_K(q_g\,\|\,p)+O(η^2). (21) In particular, whenever qg≠pq_g≠ p, the update strictly decreases KL(qg∥p)KL(q_g\,\|\,p) for every sufficiently small positive η. The proof is given in Appendix D. The update equation 19 is realized by regression onto a frozen target. Writing g¯ g for the current iterate, define ℒsg(g,g¯):=12∑x∈(g(x)−sg[g¯(x)+ηΔg¯(x)])2.L_sg(g; g):= 12 _x (g(x)-sg [ g(x)+η _ g(x) ] )^2. (22) As in Euclidean drifting, the stop-gradient is what makes this a discretization of the flow rather than an ordinary loss. Propagating gradients through the target, with the same variable g appearing on both sides, collapses the objective to ℒno-sg(g):=12‖g−(g+ηΔg)‖22=η22‖Δg‖22.L_no -sg(g):= 12 \|g- (g+η _g ) \|_2^2= η^22 \| _g \|_2^2. (23) Gradient descent on ℒno-sgL_no -sg minimizes the squared magnitude of the drift through derivatives of Δg _g; it does not implement the vector field Δg _g itself. 3.2 A Minimal Particle Realization The flow in equation 10 is Eulerian: it specifies the net change of probability mass at each state, but does not prescribe how individual particles should move. On a finite state space, the natural particle analogue of a Euclidean displacement is a Markov jump. Our goal is therefore to construct nonnegative jump rates whose directional probability traffic reproduces the signed current Jp,qJ_p,q. For strictly positive probability mass functions p,qp,q on X, define δq(x,y):=logρq(y)−logρq(x),Δp,q(x):=−1q(x)∑y∈Jp,q(x,y). _q(x,y):= _q(y)- _q(x), _p,q(x):=- 1q(x) _y J_p,q(x,y). Using the positive 11-homogeneity of the logarithmic mean, the current in equation 9 admits the equivalent representation Jp,q(x,y)=q(x)K(x,y)L(δq(x,y))[sp(x,y)−δq(x,y)].J_p,q(x,y)=q(x)K(x,y)L ( _q(x,y) ) [s_p(x,y)- _q(x,y) ]. (24) Since jump rates must be nonnegative, we assign to each oriented edge the positive directional part of this current, ap,q(x,y):=[Jp,q(x,y)]+q(x)=K(x,y)L(δq(x,y))[sp(x,y)−δq(x,y)]+,x≠y.a_p,q(x,y):= [J_p,q(x,y)]_+q(x)=K(x,y)L ( _q(x,y) ) [s_p(x,y)- _q(x,y) ]_+, x≠ y. (25) These rates define the generator Ap,q(x,y)=ap,q(x,y),x≠y,−∑z≠xap,q(x,z),x=y.A_p,q(x,y)= casesa_p,q(x,y),&x≠ y,\\[3.00003pt] - _z≠ xa_p,q(x,z),&x=y. cases (26) Its off-diagonal entries are nonnegative and its rows sum to zero, so for fixed q it is a continuous-time Markov generator. Moreover, antisymmetry of Jp,qJ_p,q yields (qAp,q)(x)=−∑y∈Jp,q(x,y)=q(x)Δp,q(x).(qA_p,q)(x)=- _y J_p,q(x,y)=q(x) _p,q(x). (27) Thus, when the rates are recomputed from the evolving law, the nonlinear Kolmogorov equation q˙t=qtAp,qt q_t=q_tA_p,q_t coincides exactly with the relative-KL gradient flow in equation 10. Markov-jump realizations of graph-based probability flows are standard (Chow et al., 2011; Sun et al., 2023); the construction in equation 25 is the positive-part realization of the Maas current. Discrete drifting requires only a single step of this dynamics. We therefore freeze the rates at the current law q and take a forward-Euler step, Pη,p,q:=I+ηAp,q.P_η,p,q:=I+η A_p,q. (28) Under the step-size condition ηmax∑y≠x∈ap,q(x,y)≤1,η _x _y≠ xa_p,q(x,y)≤ 1, (29) Pη,p,qP_η,p,q is a Markov transition kernel, and qPη,p,q=q+ηq⊙Δp,q.qP_η,p,q=q+η\,q _p,q. (30) Hence one frozen stay-or-jump transition realizes exactly the forward-Euler step of the law-level gradient flow. These generator identities and the Markov-kernel property are verified in Lemma 5 of Appendix D. We now relate this particle update to the logit update. When q=qgq=q_g, write Pη,g:=Pη,p,qg.P_η,g:=P_η,p,q_g. We are now ready to state our next result. Proposition 3 (Particle realization of the logit update). Let q=qgq=q_g, and suppose that equation 29 holds with q=qgq=q_g. Then Pη,gP_η,g is a Markov transition kernel and, as η→0η→ 0, qg+ηΔg=qgPη,g+O(η2).q_g+η _g=q_gP_η,g+O(η^2). (31) Thus one frozen stay-or-jump transition realizes, to first order, the same evolution of the model law as the logit update g↦g+ηΔg g+η _g. The proof is given in Appendix D. Finally, the positive-part construction is edgewise minimal. Any other pair of nonnegative directional rates realizing the same signed current differs from equation 25 by equal counterflow in the two directions. This counterflow cancels from the marginal evolution while increasing the total particle traffic. Consequently, equation 25 is the unique realization minimizing directional traffic on each edge. Moreover, when q=pq=p, the current vanishes and all minimal jump rates are zero, so particles stop pathwise. These statements are formalized in Proposition 6 and Corollary 7 of Appendix D. 4 Controlled Validation Figure 1: Exact tabular flow on the 32×3232× 32 torus. First 1600 updates. Top: Relative density ρqt=||qt _q_t=|X|q_t (blue) at five stages; orange contours show the fixed target density. The model begins as a single bump separated from the target, transports mass across low-density regions, and converges to the three modes in the correct proportions. Bottom: KL(qt∥p)KL(q_t\,\|\,p) versus flow time on a logarithmic scale. Numbered markers associate each displayed distribution with its position along the trajectory. We design a two-stage synthetic experiment to separate the correctness of the discrete transport from the error introduced by neural amortization. Both stages use the same target distribution, transport geometry, initialization, and time discretization. In the first stage, we represent the model distribution exactly with one logit per state, eliminating estimation, representation, and optimization error and allowing us to study the proposed transport update in isolation. In the second stage, we replace this tabular representation with a shared neural generator and examine how accurately the resulting transport steps can be amortized. This controlled progression separates discrepancies in the underlying discrete flow from those introduced by neural approximation. The state space is the periodic grid =(ℤ/32ℤ)2X=(Z/32Z)^2, ||=1024|X|=1024, equipped with the nearest-neighbour random-walk kernel K and its uniform invariant measure π. The target distribution p is a mixture of three Gaussian bumps with unequal weights and a small uniform floor, ensuring the strict positivity assumed in Section 3. We initialize the model as a single bump separated from all three target modes, so reaching p requires transporting probability across low-density regions rather than merely performing local relaxation. Unless stated otherwise, both stages use the same step size η=0.8η=0.8. Full target parameters and remaining hyperparameters are given in Appendix B. Exact Tabular Validation. We first ask whether the proposed update reproduces the predicted discrete Wasserstein gradient flow when every distributional quantity is available exactly. At scale, this update is obscured by three additional sources of error: the relevant density ratios must be estimated from samples, the generator has limited representational capacity, and optimization is stochastic. We remove all three by representing the model with one logit per state, qg=softmax(g).q_g=softmax(g). Because π is uniform, these logits coincide with the chart in equation 13 and can represent any strictly positive distribution on X. Since both p and qgq_g are known exactly, all target scores, currents, mobilities, and jump rates can likewise be evaluated without estimation. Moreover, the update depends only on model-logit differences and target log-density ratios along edges, so no normalizing constant is required. Any discrepancy from the predictions can therefore be attributed to the transport update or its numerical discretization rather than to statistical estimation or neural amortization. We apply the exact tabular update for 64006400 steps, evaluating every quantity entering the update exactly. Each step is implemented through the stop-gradient objective (equation 22); with unit learning rate, one gradient step coincides with the explicit Euler update of equation 19. Figure 1 visualizes the first 16001600 steps of this trajectory. During this initial phase, the model transports mass across the torus, splits toward the three target regions, and recovers their unequal mixture structure, while KL(qt∥p)KL(q_t\,\|\,p) decreases sharply. Continuing the same exact trajectory to 64006400 steps drives the KL to a numerical floor near 1.3×10−101.3× 10^-10, while total variation decreases from 0.710.71 to 7.7×10−67.7× 10^-6. The complete trajectory and additional numerical controls are reported in Appendix C. Figure 2 shows the edgewise transport driving this evolution at initialization. The target-dependent component directs mass toward states favoured by p, while the model-dependent component spreads mass away from regions where q is concentrated. Their sum gives the net current of the target-relative Maas flow, carrying probability from the initial model surplus toward the target modes. Figure 2: Edgewise decomposition of the relative-entropy current at initialization. Left: Target-attraction current over the target density. Center: Model-repulsion current over the model density; by the logarithmic-mean identity, this component coincides with the linear diffusion current. Right: Net current over log(ρq/ρp) ( _q/ _p), transporting mass from model surplus toward target deficit. Particle Approximation Validation. The aforementioned tabular update of equation 19 acts on all |||X| states at once, which is possible only because the space is small enough to enumerate, whereas the particle realization of Section 3.2 instead moves individual samples; we therefore ask whether one frozen stay-or-jump transition of equation 28 reproduces it, as Proposition 3 claims. Freezing the minimal rates of equation 25 at states of the tabular trajectory, we sweep the step size across five decades, from well below the bound of equation 29 to past it, and evaluate both residuals in double precision (Figure 3 (left)). The law-level identity of equation 30 holds to machine precision throughout, so one transition is not an approximation of the flow of equation 10 but its exact explicit-Euler step. The residual of equation 31 instead decays with a fitted slope of two, so the mismatch between moving particles and moving logits is the curvature of the chart of equation 13 that Proposition 2 already ascribes to the logit update, and nothing beyond it. Second order rather than first is what makes the substitution safe, since the per-step error then accumulates to O(η)O(η) over a fixed horizon and vanishes with the step size, and the bound of equation 29 stays far from binding at the η used throughout. Neural amortization. We next replace the tabular representation with a finite-capacity neural generator and ask whether a single parameterized model can absorb the sequence of transport updates. In the tabular experiment, the logits were themselves the model parameters, so the Euler update equation 19 could be applied directly. A neural generator, by contrast, represents its marginal through a mixture of latent-conditional distributions and cannot be updated by directly perturbing its marginal logits. We therefore construct training targets using the equivalent frozen law-level transport of equation 28. This reintroduces representation and optimization error while keeping the model marginal, target scores, and transport kernel exactly computable, so no estimation error is present. The generator is a latent-conditioned MLP Gθ(⋅∣z)G_θ(· z), evaluated on a fixed bank of B=256B=256 latent codes. The grid, target, initial marginal, transport kernel, and step size are identical to the tabular experiment. At transport step k, we denote the conditional distributions by Rk,i:=Gθk(⋅∣zi),i=1,…,B,R_k,i:=G_ _k(· z_i), i=1,…,B, and form the generator marginal qk:=1B∑i=1BRk,i.q_k:= 1B _i=1^BR_k,i. The current marginal determines a single frozen transport kernel, which is then applied to every conditional: Tk,i:=Rk,iPη,p,qk,Pη,p,qk:=I+ηAp,qk.T_k,i:=R_k,iP_η,p,q_k, P_η,p,q_k:=I+η A_p,q_k. (32) Holding the teachers Tk,iT_k,i fixed, we take a small number of stochastic gradient steps on their conditional cross-entropy with Gθ(⋅∣zi)G_θ(· z_i), corresponding to the stop-gradient objective in equation 22, and then recompute the marginal and the next transport step. We repeat this procedure for 16001600 transport updates. This construction preserves the desired marginal dynamics. Indeed, linearity of the frozen kernel gives 1B∑i=1BTk,i=(1B∑i=1BRk,i)Pη,p,qk=qkPη,p,qk. 1B _i=1^BT_k,i= ( 1B _i=1^BR_k,i )P_η,p,q_k=q_kP_η,p,q_k. (33) Thus, if the conditional teachers were fitted exactly, the neural marginal would advance by exactly the frozen law-level Euler step equation 30, which agrees with the corresponding logit update to first order in η. In practice, the same network must successively approximate every teacher along the trajectory, so fitting errors may accumulate. We therefore measure its deviation from the exact transported trajectory throughout the rollout rather than only at the final iterate. Figure 3: Left: one frozen transition of equation 28 against the two updates it replaces: fitted slope 2.002.00 for the residual of equation 31, and machine precision for the Euler step of the law, which equation 30 makes an identity. The vertical rule marks ηCFL _CFL, the largest step size admitted by equation 29, past which no particle can make the move. Right: final gap to the conditional oracle against batch size b, at four step sizes with ηNη N fixed; medians over five seeds, band spanning seeds at η=0.8η=0.8. As a reference, we use a conditional oracle that starts from the same conditionals and applies equation 32 exactly at every transport step, without neural fitting. Figure 5 shows that the neural generator closely follows this oracle: the initial mass is transported across the low-density regions and distributed among the three target modes with the correct unequal weights, while the final marginal is visually indistinguishable from the target at the resolution shown. The fitted conditionals attain comparable error on training and held-out latent codes, indicating that the network learns the transport rule rather than merely memorizing the training bank. The right panel reports the total-variation gap between the neural rollout and the oracle together with the oracle’s distance to p. Increasing the per-step latent batch size b reduces the neural–oracle gap by more than an order of magnitude, with the smallest error obtained at b=Bb=B. This behavior indicates that the remaining discrepancy is primarily associated with fitting each transport teacher rather than with the underlying transport itself. Complete KL and total-variation results, together with the architecture, optimizer, and held-out evaluation protocol, are given in Appendices B and C. 5 Conclusion We introduced Discrete Drifting, which replaces Euclidean drift with a relative-KL gradient flow in Maas geometry and realizes each frozen flow step through Markov jumps. In early results, on an enumerable finite-state problem, the proposed transport follows the predicted KL-descending dynamics and can be amortized by a finite neural generator while preserving one-step inference. Future work will focus on scalable local ratio estimation and evaluation on real discrete data against strong discrete diffusion and one-step baselines. References Benamou and Brenier (2000) J. Benamou and Y. Brenier A computational fluid mechanics solution to the monge-kantorovich mass transfer problem. Numerische Mathematik 84 (3), p. 375–393. External Links: ISSN 0945-3245, Link, Document Cited by: §2. Bunne et al. (2022) C. Bunne, L. Papaxanthos, A. Krause, and M. Cuturi Proximal optimal transport modeling of population dynamics. In Proceedings of The 25th International Conference on Artificial Intelligence and Statistics, G. Camps-Valls, F. J. R. Ruiz, and I. Valera (Eds.), Proceedings of Machine Learning Research, Vol. 151, p. 6511–6528. External Links: Link Cited by: §2. Campbell et al. (2022) A. Campbell, J. Benton, V. De Bortoli, T. Rainforth, G. Deligiannidis, and A. Doucet A continuous time framework for discrete denoising models. In Advances in Neural Information Processing Systems, S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh (Eds.), Vol. 35, p. 28266–28279. External Links: Link Cited by: §1, §2. Chow et al. (2011) S. Chow, W. Huang, Y. Li, and H. Zhou Fokker–planck equations for a free energy functional or markov process on a graph. Archive for Rational Mechanics and Analysis 203 (3), p. 969–1008. External Links: ISSN 1432-0673, Link, Document Cited by: §3.2. Deng et al. (2026) M. Deng, H. Li, T. Li, Y. Du, and K. He Generative modeling via drifting. External Links: 2602.04770, Link Cited by: §1. Erbar et al. (2018) M. Erbar, J. Maas, and M. Wirth On the geometry of geodesics in discrete optimal transport. Calculus of Variations and Partial Differential Equations 58 (1). External Links: ISSN 1432-0835, Link, Document Cited by: §2. Erbar and Maas (2014) M. Erbar and J. Maas Gradient flow structures for discrete porous medium equations. Discrete and Continuous Dynamical Systems 34 (4), p. 1355–1374. External Links: ISSN 1078-0947, Document, Link Cited by: §2. Jordan et al. (1998) R. Jordan, D. Kinderlehrer, and F. Otto The variational formulation of the fokker–planck equation. SIAM Journal on Mathematical Analysis 29 (1), p. 1–17. External Links: ISSN 1095-7154, Link, Document Cited by: §2. Lou et al. (2024) A. Lou, C. Meng, and S. Ermon Discrete diffusion modeling by estimating the ratios of the data distribution. External Links: 2310.16834, Link Cited by: §1, §2. Maas (2011) J. Maas Gradient flows of the entropy for finite markov chains. Journal of Functional Analysis 261 (8), p. 2250–2292. Cited by: §1, §1, §2, §2, §2, §2, §3, §3. Rancati et al. (2026) D. Rancati, J. Maas, and F. Locatello Learning discrete diffusion on graphs via free-energy gradient flows. In Forty-third International Conference on Machine Learning, External Links: Link Cited by: §2. Sun et al. (2023) H. Sun, H. Dai, B. Dai, H. Zhou, and D. Schuurmans Discrete langevin samplers via wasserstein gradient flow. In Proceedings of The 26th International Conference on Artificial Intelligence and Statistics, F. Ruiz, J. Dy, and J. van de Meent (Eds.), Proceedings of Machine Learning Research, Vol. 206, p. 6290–6313. Cited by: §3.2. Terpin et al. (2024) A. Terpin, N. Lanzetti, M. Gadea, and F. Dorfler Learning diffusion at lightspeed. In The Thirty-eighth Annual Conference on Neural Information Processing Systems, External Links: Link Cited by: §2. Appendix A Neural Amortization of the Drifting Update Section 4 replaces the tabular logits with a shared generator and fits each frozen transition rather than applying it. Two facts make that substitution well posed, and we establish them here. First, the conditional teachers average to the exact law-level Euler step, so the generator is asked to follow the trajectory validated in Section 4 rather than a nearby one. Second, the distance between the fitted marginal and that step is controlled by the mean conditional fitting error, which separates the exact teacher construction from the amortization gap we measure. A.1 Conditional teachers At outer step k, let Rk,z R_k,z :=Gθk(⋅∣z), :=G_ _k(· z), (34) qk q_k :=z∼ν[Rk,z], :=E_z ν[R_k,z], Tk,z T_k,z :=Rk,zPη,p,qk. :=R_k,zP_η,p,q_k. Every conditional is advanced by the same kernel, built once from th marginal and then frozen. Linearity gives z[Tk,z]=z[Rk,z]Pη,p,qk=qkPη,p,qk,E_z[T_k,z]=E_z[R_k,z]P_η,p,q_k=q_kP_η,p,q_k, (35) which by equation 30 is exactly the law-l update. No approximation enters before fitting. A.2 Fitting and the amortization gap Holding the teachers fixed, we fit the generator by minimizing the mean conditional cross-entropy ℒk(θ):=−z∼ν∑x∈Tk,z(x)logGθ(x∣z),L_k(θ):=-E_z ν _x T_k,z(x) G_θ(x z), (36) the per-latent form of the stop-gradient objective equation 22. Because the teachers are c ℒk(θ)=z[H(Tk,z)+KL(Tk,z∥Gθ(⋅∣z))],L_k(θ)=E_z [H(T_k,z)+KL (T_k,z\,\|\,G_θ(· z) ) ], (37) so the only θ-dependent term is the mean divergence to the teachers. In the exact-fitting limit Gθk+1(⋅∣z)=Tk,zG_ _k+1(· z)=T_k,z learned marginal absorbs one exact Euler step. Finite fitting leaves a residual, whose effect on the marginal is bounded by convexity of total variation. Writing q¯k+1:=zGθk+1(⋅∣z) q_k+1:=E_zG_ _k+1(· z), TV(q¯k+1,qkPη,p,qk)≤zTV(Gθk+1(⋅∣z),Tk,z).TV\! ( q_k+1,q_kP_η,p,q_k ) _zTV\! (G_ _k+1(· z),T_k,z ). (38) The right-hand side is the local fitting error equation 51 and the left-hand side the induced marginal error equation 52, both reported in Appendix C.2. The teacher construction therefore contributes nothing to the amortization gap of Figure 3 (right), and every departure from the validated trajectory is attributable to the fitting budget. A.3 Training loop and cost For |E||E| directed edges and B latents, one outer step costs O(|E|)O(|E|) to build the marginal rates and O(B|E|)O(B|E|) to transport all conditionals, so teacher construction is linear in the bank size and independent of the inner budget. The Markov kernel is used only during training: sampling from the fitted generator is a single forward pass. Algorithm 1 Amortized discrete drifting. One outer step advances every conditional by the same frozen kernel; the inner loop fits the generator to the result. 1: Target scores sps_p, latent law ν, step size η, outer steps N, inner steps m 2: Initialize θ0 _0 3: for k=0,…,N−1k=0,…,N-1 do 4: Rk,z←Gθk(⋅∣z)R_k,z← G_ _k(· z) ⊳ Conditionals 5: qk←z∼ν[Rk,z]q_k _z ν[R_k,z] ⊳ Model marginal 6: ap,qk←a_p,q_k← rates of equation 25 from δqk _q_k and sps_p 7: Pη,p,qk←I+ηAp,qkP_η,p,q_k← I+η A_p,q_k ⊳ Markov under equation 29 8: Tk,z←sg[Rk,zPη,p,qk]T_k,z [R_k,zP_η,p,q_k ] ⊳ Frozen teachers 9: θk+1←m _k+1← m optimizer steps on ℒkL_k of equation 36, starting from θk _k 10: end for 11: θN _N Appendix B Experimental Setup and Reproducibility Details This section provides the complete setup for the controlled experiments in Section 4. Both experiments use the same enumerable state space, target distribution, initialization, transport graph, and time discretization. The tabular experiment evaluates the proposed transport without representation or optimization error. The neural experiment then introduces finite representation and fitting while retaining exact, enumerable transport teachers. B.1 Shared enumerable problem State space and transport graph. The state space is the periodic grid =ℤ32×ℤ32X=Z_32×Z_32, represented by x=(x1,x2)∈0,…,312x=(x_1,x_2)∈\0,…,31\^2, with ||=1024|X|=1024. Each state is connected to its four nearest neighbours, (x)=x±(1,0),x±(0,1)(mod32).N(x)=\x±(1,0),x±(0,1)\ 32. We use the nearest-neighbour random-walk kernel K(x,y)=1/4,y∈(x),0,otherwise,π(x)=11024.K(x,y)= cases1/4,&y (x),\\ 0,&otherwise, cases π(x)= 11024. (39) This kernel is reversible with respect to the uniform reference law π. Periodic Gaussian family. For a,b∈0,…,31a,b∈\0,…,31\, define the periodic distance d32(a,b) d_32(a,b) :=min|a−b|, 32−|a−b|, := \|a-b|,\,32-|a-b|\, d2(x,μ) d_T^2(x,μ) :=d32(x1,μ1)2+d32(x2,μ2)2. :=d_32(x_1, _1)^2+d_32(x_2, _2)^2. (40) The normalized periodic Gaussian with centre μ and scale σ is Zμ,σ Z_μ,σ :=∑y∈exp[−d2(y,μ)2σ2], := _y [- d_T^2(y,μ)2σ^2 ], φμ,σ(x) _μ,σ(x) :=Zμ,σ−1exp[−d2(x,μ)2σ2]. :=Z_μ,σ^-1 [- d_T^2(x,μ)2σ^2 ]. (41) Target and initialization. The target is a mixture of three periodic Gaussians with a uniform floor: p(x) p(x) =(1−ε)p~(x)+επ(x), =(1- ) p(x)+ π(x), (42) p~(x) p(x) :=∑j=13wjφμj,σj(x),∑j=13wj=1, := _j=1^3w_j _ _j, _j(x), _j=1^3w_j=1, q0(x) q_0(x) =φμ0,σ0(x). = _ _0, _0(x). The target modes have unequal weights and are separated by low-density regions. The initial mode is placed away from all three target modes, so reaching the target requires both splitting and transporting probability mass. The uniform floor ensures the strict positivity required by the flow. All centres and scales are given in lattice units, and take the values (μ1,σ1,w1) ( _1, _1,w_1) =((8.00,8.00), 2.88, 0.462), = ((8.00,8.00),\,2.88,\,0.462 ), (43) (μ2,σ2,w2) ( _2, _2,w_2) =((23.04,10.88), 2.56, 0.284), = ((23.04,10.88),\,2.56,\,0.284 ), (μ3,σ3,w3) ( _3, _3,w_3) =((16.00,24.96), 3.20, 0.254), = ((16.00,24.96),\,3.20,\,0.254 ), (μ0,σ0) ( _0, _0) =((26.24,24.96), 4.48),ε=0.168. = ((26.24,24.96),\,4.48 ), =0.168. The floor is set so that the least likely state still carries relative density ρp=0.168 _p=0.168, against 7.737.73 at the highest mode. The implementation builds the same law from unnormalized periodic Gaussians and a log-domain floor; equation 43 reports the equivalent parameters of the normalized form in equation 42. Because X is finite, every distribution, directed edge, current, transition probability, and conditional teacher can be evaluated by enumeration. Using the same problem throughout allows us to distinguish errors in the proposed flow from those introduced by particle sampling or neural amortization. B.2 Exact tabular validation The tabular model assigns one unconstrained logit g(x)g(x) to each state and represents qg=softmax(g).q_g=softmax(g). The initial logits represent q0q_0 from equation 42. At each outer step, we enumerate all states and directed nearest-neighbour edges to compute δg(x,y) _g(x,y), sp(x,y)s_p(x,y), L(δg(x,y))L( _g(x,y)), the current, the logit drift, and the frozen Markov transition. Numerically, the logarithmic-mean factor L(r)L(r) is evaluated as expm1(r)/rexpm1(r)/r away from zero and by its continuous limit at zero. Writing the frozen transition of equation 28 entrywise, Pη,p,q(x,y)=ηap,q(x,y),y∈(x),1−η∑z∈(x)ap,q(x,z),y=x,0,otherwise,P_η,p,q(x,y)= casesη\,a_p,q(x,y),&y (x),\\[1.99997pt] 1-η _z (x)a_p,q(x,z),&y=x,\\[1.99997pt] 0,&otherwise, cases (44) we verify that every entry is nonnegative before applying each update. This is the step-size condition of equation 29: the smallest stay probability over the tabular trajectory is 0.9500.950, attained at η=0.8η=0.8. Component Setting Model One unconstrained logit per state; qg=softmax(g)q_g=softmax(g) Initialization Exact logits representing q0q_0 Drift evaluation All 10241024 states and directed nearest-neighbour edges Outer updates 64006400 Flow step η=0.8η=0.8 Final flow time 51205120 Chart optimization One full-batch step on equation 22 with learning rate 11 Exact references Probability-space Euler update and frozen Markov law Particle counts M∈[2×103,2×105]M∈[2× 10^3,2× 10^5]; six repetitions per count Table 1: Configuration of the exact tabular validation. To quantitatively validate every step of the construction, from the law-level flow of equation 10 to the sampled stop-gradient update, we perform six numerical checks: • Chart versus density Euler. We compare the distribution represented by g+ηΔg+η _g with the explicit density update qg+ηqg⊙Δgq_g+η q_g _g. • Frozen Markov law versus density Euler. We compare qgPη,p,qgq_gP_η,p,q_g with qg+ηqg⊙Δgq_g+η q_g _g. • Particle histogram versus frozen Markov law. We compare the empirical post-jump histogram q^M+ q_M^+ with qgPη,p,qgq_gP_η,p,q_g. • Frozen transition versus logit update. We compare qgPη,p,qgq_gP_η,p,q_g with qg+ηΔgq_g+η _g, which is the residual of equation 31. • Stop-gradient control. We compare the intended frozen-target regression in equation 22 with the differentiate-through-target objective in equation 23. • Sampled state–edge control. We compare an unbiased sampled implementation with full enumeration. States are drawn uniformly with replacement, 40964096 per update, with two neighbours sampled per drawn state. The scaled loss is minimized by gradient descent at learning rate 0.50.5 for 1280012800 updates, so the sampled control covers the same flow time 51205120 as the enumerated run. Results for these checks are reported in Appendix C.1. Additionally the frozen-transition order test, also appears in Figure 3 (left) of the main text and in Table 2 below. Step-size and particle sweeps. The chart-versus-Euler and frozen-transition checks are order tests, so their sweeps are specified here. The production flow runs in single precision, whose probability floor near 10−710^-7 would hide a second-order signal, and both order tests are therefore recomputed on the host in double precision from the same frozen logits. The chart-versus-Euler control freezes the trajectory at outer step 12801280 and sweeps η∈1,2,5×10−4,10−3,10−2,10−1∪4×10−1η∈\1,2,5\×\10^-4,10^-3,10^-2,10^-1\∪\4× 10^-1\. The frozen-transition control freezes the trajectory at outer steps 00, 128128, 320320, and 720720, and extends the same grid to η∈0.8,1.6,3.2,6.4,12.8η∈\0.8,1.6,3.2,6.4,12.8\, which carries it past the step-size bound at the earliest state. The slope is fitted by least squares in logη η over the points that satisfy equation 29 and whose residual exceeds 10−1410^-14, so that neither an inadmissible step nor double roundoff enters the fit. Table 2 reports the four frozen states; Figure 3 (left) plots the first. The particle control freezes the trajectory at outer step 320320, uses η=0.8η=0.8, and draws M∈2,5×103,104∪105,2×105M∈\2,5\×\10^3,10^4\∪\10^5,2× 10^5\ particles, with six independent repetitions at each count. Flow time t KL(qt∥p)KL(q_t\,\|\,p) Step-size bound Fitted slope Law-identity residual 00 1.851.85 3.723.72 2.00142.0014 3.7×10−163.7× 10^-16 102.4102.4 3.08×10−13.08× 10^-1 7.197.19 2.00012.0001 1.2×10−161.2× 10^-16 256.0256.0 7.47×10−27.47× 10^-2 16.0616.06 2.00002.0000 9.3×10−179.3× 10^-17 576.0576.0 1.20×10−21.20× 10^-2 38.2038.20 2.00012.0001 6.9×10−176.9× 10^-17 Table 2: Frozen-transition order test at four states of the tabular trajectory. The step-size bound is the largest η admitted by equation 29, the fitted slope is that of ‖qg+ηΔg−qgPη,g‖1\|q_g+η _g-q_gP_η,g\|_1, and the last column is the largest observed violation of equation 30 over the whole sweep. The bound grows as the flow converges and the rates decay. B.3 Finite neural amortization Generator and initialization. The neural generator is a latent-conditioned MLP Gθ(⋅∣z)G_θ(· z), where z∼(0,I16)z (0,I_16). It has three hidden layers of width 128128, SiLU activations, and a 10241024-way softmax output. The headline experiment draws a fixed training bank of B=256B=256 latents and uses the uniform distribution over this bank as the finite latent prior. A separate fixed bank of 256256 held-out latents is never used for optimization. A statewise output bias is calibrated once so that the training-bank marginal matches the tabular initialization q0q_0, and is then frozen. The bias therefore supplies only the common initial marginal. All subsequent transport must be represented by the shared MLP weights. The bias is obtained by a damped fixed-point iteration: starting from logq0 q_0 centred to zero mean, each of 700700 iterations adds 0.75(logq0−logq¯)0.75 ( q_0- q ) to the current bias, where q¯ q is the bank marginal it induces, and recentres the result. The iteration count is fixed rather than tolerance-driven, so calibration is deterministic given the seed. Exact conditional teachers. For the fixed training bank, let Rk,i R_k,i =Gθk(⋅∣zi), =G_ _k(· z_i), (45) qk q_k =1B∑i=1BRk,i, = 1B _i=1^BR_k,i, Tk,i T_k,i =Rk,iPη,p,qk. =R_k,iP_η,p,q_k. We enumerate all 10241024 states and all B conditionals. Consequently, qkq_k, Pη,p,qkP_η,p,q_k, and every Tk,iT_k,i are exact for this finite latent prior. Moreover, 1B∑i=1BTk,i=qkPη,p,qk, 1B _i=1^BT_k,i=q_kP_η,p,q_k, (46) so the mean teacher is exactly the density-Euler update validated in the tabular experiment. At each outer step, all teachers are frozen and fitted jointly by minimizing their mean conditional cross-entropy. We take 1616 Adam steps and carry the optimizer moments across outer updates. Carrying the moments is useful near convergence, where the required teacher movements become small. Component Setting Latent distribution z∼(0,I16)z (0,I_16) Generator Three hidden layers of width 128128, with SiLU activations Output 10241024-way softmax Training latent bank B=256B=256, fixed for the complete run Held-out latent bank 256256, fixed and never optimized against Initialization Statewise output bias calibrated to q0q_0, then frozen Outer updates 16001600 Flow step η=0.8η=0.8 Final flow time 12801280 Inner optimization 1616 full-bank Adam steps per outer update Adam learning rate 2×10−32× 10^-3 Optimizer state First and second moments carried across outer updates Independent repetitions Independent deterministic seeds Table 3: Configuration of the finite neural amortization experiment. Conditional oracle. We compare the neural rollout with an oracle that applies the same frozen transition directly to each conditional: Rk+1,i∗=Rk,i∗Pη,p,qk∗,qk∗=1B∑i=1BRk,i∗.R^*_k+1,i=R^*_k,iP_η,p,q_k^*, q_k^*= 1B _i=1^BR^*_k,i. (47) This oracle contains no neural fitting error. As an implementation check, its marginal is also integrated independently using the probability-space Euler update. Across the rollout, the two implementations differ by at most 6.6×10−96.6× 10^-9 in any state’s probability. B.4 Finite-batch sweep The headline neural experiment enumerates its complete 256256-latent training bank. We separately examine how finite-batch fitting affects the accumulated amortization error. For this experiment, we fix a population of Bpop=4096B_pop=4096 latents and vary the finite-batch size over b∈16,32,64,128,256,512,1024,4096.b∈\16,32,64,128,256,512,1024,4096\. The b=4096b=4096 condition uses the complete latent population. Unless stated otherwise, the architecture, flow step, optimizer, number of inner updates, and number of outer updates match Table 3. Each condition is repeated over independent deterministic seeds. At each outer update we draw a fresh subset of b latents without replacement from the fixed population, and the same subset is reused for all 1616 inner steps of that update. The subset alone defines the marginal qkq_k, the frozen kernel Pη,p,qkP_η,p,q_k, and every conditional teacher, so a small b perturbs the update itself and not merely the gradient. Evaluation is unaffected by this choice: the reported neural marginal qk(b)q_k^(b) and the conditional-oracle marginal qk∗q_k^* are always averaged over the complete 40964096-latent population. Each subset is the length-b prefix of a permutation drawn from a seed-derived key, so conditions sharing a seed also share the latent population, the initial parameters, the optimizer state, and nested batches. Let qk(b)q_k^(b) denote the neural marginal produced with finite-batch size b, and let qk∗q_k^* denote the conditional-oracle marginal defined on the same latent population. For Nout=1600N_out=1600, the quantities reported in Figure 3 (right) are gfinal(b) g_final^(b) :=TV(qNout(b),qNout∗), :=TV (q_N_out^(b),q_N_out^* ), (48) gavg(b) g_avg^(b) :=1Nout+1∑k=0NoutTV(qk(b),qk∗). := 1N_out+1 _k=0^N_outTV (q_k^(b),q_k^* ). (49) The first measures endpoint agreement, while the second measures agreement over the complete rollout. Both marginals are recorded at every outer step, and the time average runs over all Nout+1N_out+1 of them, including the common initialization at k=0k=0. B.5 Metrics and statistical reporting For probability laws r and s on X, we report TV(r,s) (r,s) :=12∑x∈|r(x)−s(x)|, := 12 _x |r(x)-s(x)|, (50) KL(r∥s) (r\,\|\,s) :=∑x∈r(x)logr(x)s(x). := _x r(x) r(x)s(x). For the neural experiment, the local training error is etrain,k:=1B∑i=1BTV(Rk+1,i,Tk,i),e_train,k:= 1B _i=1^BTV(R_k+1,i,T_k,i), (51) and the induced marginal fitting error is emarg,k:=TV(qk+1,qkPη,p,qk).e_marg,k:=TV (q_k+1,q_kP_η,p,q_k ). (52) For the held-out bank, we construct teachers using the same frozen transition as the training bank: Tk,ih:=Rk,ihPη,p,qk.T^h_k,i:=R^h_k,iP_η,p,q_k. The held-out fitting error is eheld,k:=1Bh∑i=1BhTV(Rk+1,ih,Tk,ih).e_held,k:= 1B_h _i=1^B_hTV (R^h_k+1,i,T^h_k,i ). (53) Results for the headline neural experiment report means across seeds, with pointwise 95%95\% t intervals m¯±t0.975,n−1sn, m\ ±\ t_0.975,\,n-1 s n, (54) where n is the number of seeds and s the sample standard deviation across them. The finite-batch sweep has too few seeds per condition for a meaningful interval, so we report the median across seeds and the full seed range instead. To test whether the finite latent bank materially affects evaluation, each final generator is additionally reevaluated using fresh Gaussian banks of sizes 256256, 10241024, and 40964096. The 40964096-latent marginal serves as the reference. This diagnostic changes only the bank used to evaluate a fixed generator; it is distinct from the finite-batch sweep, which changes the fitting procedure during the rollout. Appendix C Additional Results C.1 Exact transport and particle realization The main text establishes that the tabular model transports mass toward the target and that one frozen particle transition is consistent with the corresponding logit update. Here, we test each link in this construction over the complete trajectory: whether the stop-gradient update follows the intended probability flow, whether its descent obeys the predicted dissipation identity, and whether the particle realization introduces any error beyond Monte Carlo sampling. The tabular update follows the intended flow. Figure 4 reports KL and total variation along the stop-gradient chart update, the probability-space Euler reference, the sampled state–edge update, and the differentiate-through-target control. The stop-gradient chart update closely tracks the Euler reference throughout the descent, with both approaching their numerical floors. Replacing full enumeration with unbiased state–edge sampling preserves the same convergence behaviour, showing that stochastic evaluation does not alter the underlying dynamics. By contrast, differentiating through the transport target produces only limited descent and remains far from the target. This confirms that the stop-gradient in equation 22 is essential: it implements the prescribed vector field, whereas the differentiated-target objective does not. The descent obeys the gradient-flow identity. Convergence toward p alone would not establish that the update follows the proposed Maas gradient flow: different dynamics can share the same equilibrium while following different transport paths and descent rates. dissipation prescribed by the Maas geometry. We measure the corresponding finite-difference rate as Dk:=−KL(qk+1∥p)−KL(qk∥p)ΔtD_k:=- KL(q_k+1\,\|\,p)-KL(q_k\,\|\,p) t Here, qkq_k and qk+1q_k+1 are the model laws before and after the update, p is the fixed target, and the negative sign makes DkD_k positive when KL decreases. Thus, DkD_k measures the observed KL decrease per unit flow time. The continuous-time Maas flow predicts the corresponding instantaneous rate to be ℐK(qk∥p)I_K(q_k\,\|\,p), the discrete relative Fisher information that aggregates variations of the log-density ratio across the edges of K, weighted by their mobility. As shown in Figure 4 (bottom left), the measured and predicted rates agree to within one percent until both become dominated by numerical precision. The remaining panels verify the algebra underlying this evolution throughout the rollout. The relative residuals for mass conservation, current antisymmetry, the attraction–repulsion decomposition, and the logarithmic-mean identity remain below 4.1×10−64.1× 10^-6. The absolute mass-conservation defect remains at approximately 10−1610^-16. Figure 4: End-to-end diagnostics for the exact tabular flow. Top left: relative KL for the stop-gradient chart update, probability-space Euler reference, sampled state–edge implementation, and differentiate-through-target control. Top right: the corresponding total-variation trajectories. Bottom left: predicted instantaneous dissipation and measured finite-difference KL decrease. Bottom right: numerical residuals for four structural identities. The intended chart update, probability-space reference, and sampled implementation agree throughout the resolved descent, whereas the differentiate-through-target control does not. The chart and particle updates exhibit the predicted scaling. We first test the two convergence predictions governing the numerical realization. Across five decades of step size, the discrepancy between the chart update and the explicit probability-space Euler step has fitted slope 2.002.00, verifying the O(η2)O(η^2) local consistency predicted by equation 20. Across the particle-count sweep, the empirical post-jump law converges to the exact transported law with fitted slope −0.51-0.51, matching the canonical O(M−1/2)O(M^-1/2) Monte Carlo rate. Thus, decreasing the step size recovers the intended probability flow, while increasing the particle population recovers its exact law-level transition. The complementary test of equation 31, reported in Figure 3 (left) and Table 2, verifies the same second-order agreement between the frozen transition and the logit update. Consequently, both the chart update and the frozen particle transition reproduce the probability-space Euler evolution to first order in η. The frozen transition realizes the Euler step. At the representative outer step used for the particle control, the production step size remains comfortably within the admissible range of equation 29. At the law level, the exact post-jump distribution agrees with the explicit probability-space Euler update to machine precision, directly verifying equation 30; the corresponding residuals are reported in Table 2. At the particle level, the observed root-mean-square error is within two percent of the scale predicted by multinomial sampling, the empirical two-standard-deviation coverage is within half a percentage point of its nominal value, and no state deviates by more than three predicted standard deviations. Together with the O(M−1/2)O(M^-1/2) convergence above, these diagnostics establish both parts of the particle construction: the frozen transition realizes the Euler update exactly at the law level, and finite particle populations recover that law at the expected Monte Carlo rate. C.2 Finite neural amortization The preceding experiments validate the transport update without representation or optimization error. We now isolate the approximation introduced when a single finite-capacity generator must absorb the successive conditional teachers. We examine this at three levels: the accuracy of each local fit, the accumulated deviation from the conditional oracle, and the sensitivity of the reported marginal to finite latent banks. The generator learns the local transport update. At the final outer step, the training, held-out, and induced marginal teacher errors are all approximately 2.1×10−42.1× 10^-4 in total variation. The training and held-out conditional errors are indistinguishable even though the held-out latents are never used for optimization, leaving no measurable generalization gap across latent codes. The equally small marginal error shows that this conditional accuracy transfers directly to the model law rather than being lost when the conditionals are averaged. The neural rollout preserves target-directed descent. Finite inner optimization does not guarantee that every fitted update decreases KL. Nevertheless, every seed retains the strong overall descent of the exact flow. Small one-step increases occur during the rollout, but their combined magnitude remains below one percent of the cumulative decrease. Finite fitting therefore introduces minor local fluctuations without disrupting the global movement toward the target. Consistently, Figure 5 shows that every final marginal recovers the three target modes and their unequal allocation of mass. The individual runs closely resemble their mean, ruling out an apparently accurate average produced by qualitatively different endpoints. Figure 5: Final neural marginals across five seeds. The panels show the target, the mean final marginal, and one representative seed. Orange contours mark the fixed target-density levels. Every seed recovers the three modes and their unequal mass allocation, with limited variation across runs (final total variation (5.12±1.09)×10−3(5.12± 1.09)× 10^-3). Larger fitting batches improve oracle tracking. The conditional oracle applies every transport teacher exactly, so its difference from the neural rollout isolates error accumulated through amortization rather than error in the transport itself. Figure 3 (right) shows that the final neural–oracle gap falls sharply as the fitting batch grows, across all tested step sizes, before levelling off at a small residual. The time-averaged gap follows the same pattern, showing that the improvement holds throughout the rollout rather than only at its endpoint. Because the full-population condition removes latent subsampling but retains a small gap, the remaining discrepancy can be attributed to finite representation and optimization rather than to the teacher construction. Endpoint estimates are stable across latent banks. We finally reevaluate each trained generator using progressively larger fresh latent banks, without changing its parameters. The resulting variation in the estimated marginal decreases with bank size and remains well below the final distance to the target even for the smallest evaluation bank. The reported endpoint error is therefore not an artefact of Monte Carlo variation in the evaluation marginal. This diagnostic is distinct from the preceding batch sweep: fitting-batch size changes the transport updates learned during the rollout, whereas evaluation-bank size changes only how the marginal of an already trained generator is estimated. Appendix D Proofs Throughout this appendix, the target probability mass function p is fixed and strictly positive. For any strictly positive probability mass function q on X, we write ℓp,q(x):=logρq(x)ρp(x),x∈, _p,q(x):= _q(x) _p(x), x , (55) and, for any h:→ℝh:X , ⟨h⟩q:=∑x∈q(x)h(x). h _q:= _x q(x)h(x). (56) The proofs rely repeatedly on three elementary identities: the differential of the logit chart, the representation of the Maas current in logit coordinates, and the fact that the resulting logit drift is centered under the current model law. We collect them first. D.1 Preliminary identities Let Φ:ℝ→∗() :R^X _*(X) denote the logit chart Φ(g)=qg (g)=q_g defined in equation 13. For g,h∈ℝg,h ^X, we write DΦ(g)[h]D (g)[h] for the directional derivative of Φ at g in the direction h, namely DΦ(g)[h]:=dεΦ(g+εh)|ε=0.D (g)[h]:= . dd (g+ h) |_ =0. Equivalently, for each x∈x , DΦ(g)[h](x)=dεqg+εh(x)|ε=0.D (g)[h](x)= . dd q_g+ h(x) |_ =0. Lemma 4 (Logit differential and current identities). For every g,h∈ℝg,h ^X, DΦ(g)[h](x)=qg(x)(h(x)−⟨h⟩qg),x∈.D (g)[h](x)=q_g(x) (h(x)- h _q_g ), x . (57) Moreover, the current admits the logit-coordinate representation Jp,qg(x,y)=qg(x)K(x,y)L(δg(x,y))(sp(x,y)−δg(x,y)),J_p,q_g(x,y)=q_g(x)K(x,y)L ( _g(x,y) ) (s_p(x,y)- _g(x,y) ), (58) and therefore −∑y∈Jp,qg(x,y)=qg(x)Δg(x).- _y J_p,q_g(x,y)=q_g(x) _g(x). (59) In particular, the logit drift is centered under the current model law: ⟨Δg⟩qg=∑x∈qg(x)Δg(x)=0. _g _q_g= _x q_g(x) _g(x)=0. (60) Proof. We first differentiate the logit parametrization. For a perturbation g+εhg+ h, qg+εh(x)=π(x)eg(x)+εh(x)Zπ(g+εh).q_g+ h(x)= π(x)e^g(x)+ h(x)Z_π(g+ h). Differentiating at ε=0 =0 gives DΦ(g)[h](x) D (g)[h](x) =π(x)eg(x)h(x)Zπ(g)−π(x)eg(x)Zπ(g)2∑z∈π(z)eg(z)h(z) = π(x)e^g(x)h(x)Z_π(g)- π(x)e^g(x)Z_π(g)^2 _z π(z)e^g(z)h(z) =qg(x)(h(x)−∑z∈qg(z)h(z)), =q_g(x) (h(x)- _z q_g(z)h(z) ), which is equation 57. The second term is the correction induced by normalization; in particular, DΦ(g)[]=0D (g)[1]=0, reflecting the gauge invariance qg+c=qgq_g+c1=q_g. We next express the current in the same coordinates. By positive 11-homogeneity of the logarithmic mean, mρg(x,y)=ρg(x)m(1,ρg(y)ρg(x)).m_ _g(x,y)= _g(x)m (1, _g(y) _g(x) ). Since ρg(y)ρg(x)=eg(y)−g(x)=eδg(x,y), _g(y) _g(x)=e^g(y)-g(x)=e _g(x,y), the definition of L yields mρg(x,y)=ρg(x)L(δg(x,y)).m_ _g(x,y)= _g(x)L ( _g(x,y) ). (61) At the same time, ℓp,qg(x)−ℓp,qg(y) _p,q_g(x)- _p,q_g(y) =(logρg(x)−logρg(y))−(logρp(x)−logρp(y)) = ( _g(x)- _g(y) )- ( _p(x)- _p(y) ) =sp(x,y)−δg(x,y). =s_p(x,y)- _g(x,y). Substituting these identities into equation 9, and using qg(x)=π(x)ρg(x)q_g(x)=π(x) _g(x), gives equation 58. Hence −∑y∈Jp,qg(x,y) - _y J_p,q_g(x,y) =qg(x)∑y∈K(x,y)L(δg(x,y))(δg(x,y)−sp(x,y)) =q_g(x) _y K(x,y)L ( _g(x,y) ) ( _g(x,y)-s_p(x,y) ) =qg(x)Δg(x), =q_g(x) _g(x), which proves equation 59. It remains to establish the centering identity. Reversibility of K and symmetry of the logarithmic mean give π(x)K(x,y)mρg(x,y)=π(y)K(y,x)mρg(y,x).π(x)K(x,y)m_ _g(x,y)=π(y)K(y,x)m_ _g(y,x). Since ℓp,qg(y)−ℓp,qg(x)=−(ℓp,qg(x)−ℓp,qg(y)), _p,q_g(y)- _p,q_g(x)=- ( _p,q_g(x)- _p,q_g(y) ), the current is antisymmetric: Jp,qg(y,x)=−Jp,qg(x,y).J_p,q_g(y,x)=-J_p,q_g(x,y). (62) Summing equation 59 over x therefore gives ∑x∈qg(x)Δg(x) _x q_g(x) _g(x) =−∑x,y∈Jp,qg(x,y) =- _x,y J_p,q_g(x,y) =−12∑x,y∈(Jp,qg(x,y)+Jp,qg(y,x)) =- 12 _x,y (J_p,q_g(x,y)+J_p,q_g(y,x) ) =0. =0. This is equation 60. ∎ D.2 Logit representation of the gradient flow Proof of Proposition 1. Let qt=qgtq_t=q_g_t. Applying the chain rule to the logit chart and using Lemma 4 gives q˙t(x)=qt(x)(g˙t(x)−⟨g˙t⟩qt). q_t(x)=q_t(x) ( g_t(x)- g_t _q_t ). (63) Suppose first that qtq_t satisfies the probability-coordinate gradient flow in equation 10. By equation 59, q˙t(x)=qt(x)Δgt(x). q_t(x)=q_t(x) _g_t(x). (64) Since qt(x)>0q_t(x)>0 for every x∈x , comparison with equation 63 yields g˙t(x)−⟨g˙t⟩qt=Δgt(x). g_t(x)- g_t _q_t= _g_t(x). Thus g˙t=Δgt+ct,ct:=⟨g˙t⟩qt. g_t= _g_t+c_t1, c_t:= g_t _q_t. The scalar ctc_t is continuous because t↦gtt g_t is continuously differentiable and g↦qgg q_g is smooth. Conversely, suppose that g˙t=Δgt+ct g_t= _g_t+c_t1 for some continuous scalar function ctc_t. The centering identity equation 60 implies ⟨g˙t⟩qt=⟨Δgt⟩qt+ct=ct. g_t _q_t= _g_t _q_t+c_t=c_t. Substituting this into equation 63 eliminates the gauge term and gives q˙t(x) q_t(x) =qt(x)Δgt(x) =q_t(x) _g_t(x) =−∑y∈Jp,qt(x,y), =- _y J_p,q_t(x,y), where the final identity again follows from equation 59. Hence qtq_t satisfies equation 10, completing the equivalence. ∎ D.3 Local consistency and entropy descent Proof of Proposition 2. Fix g∈ℝg ^X, and abbreviate q:=qg,Δ:=Δg,ℓ(x):=logρq(x)ρp(x).q:=q_g, := _g, (x):= _q(x) _p(x). We first establish the first-order evolution of the probability law. Recall that Φ(g)=qg (g)=q_g. By definition of the directional derivative, DΦ(g)[Δ]=dεqg+εΔ|ε=0.D (g)[ ]= . dd q_g+ |_ =0. Since Φ is smooth, Taylor’s theorem gives qg+ηΔ=qg+ηDΦ(g)[Δ]+O(η2).q_g+η =q_g+η D (g)[ ]+O(η^2). (65) By equation 57, DΦ(g)[Δ](x)=q(x)(Δ(x)−⟨Δ⟩q).D (g)[ ](x)=q(x) ( (x)- _q ). The centering identity equation 60 gives ⟨Δ⟩q=0 _q=0, and hence DΦ(g)[Δ](x)=q(x)Δ(x).D (g)[ ](x)=q(x) (x). Substituting this into equation 65 yields qg+ηΔg=qg+ηqg⊙Δg+O(η2),q_g+η _g=q_g+η q_g _g+O(η^2), which proves equation 20. We next study the corresponding change in relative entropy. Define ℋp(q):=KL(q∥p)=∑x∈q(x)logq(x)p(x),H_p(q):=KL(q\,\|\,p)= _x q(x) q(x)p(x), and (g):=ℋp(qg).G(g):=H_p(q_g). Thus =ℋp∘ΦG=H_p . For a perturbation v∈ℝv ^X, the directional derivative of ℋpH_p at q is Dℋp(q)[v] DH_p(q)[v] :=dεℋp(q+εv)|ε=0 := . dd H_p(q+ v) |_ =0 =dε∑x∈(q(x)+εv(x))logq(x)+εv(x)p(x)|ε=0 = . dd _x (q(x)+ v(x) ) q(x)+ v(x)p(x) |_ =0 =∑x∈v(x)(logq(x)p(x)+1). = _x v(x) ( q(x)p(x)+1 ). Now let h∈ℝh ^X be an arbitrary perturbation of the logits. Since qg+εh=q+εDΦ(g)[h]+o(ε),q_g+ h=q+ D (g)[h]+o( ), we obtain D(g)[h] DG(g)[h] :=dεℋp(qg+εh)|ε=0 := . dd H_p(q_g+ h) |_ =0 =Dℋp(q)[DΦ(g)[h]] =DH_p(q) [D (g)[h] ] =∑x∈DΦ(g)[h](x)(logq(x)p(x)+1). = _x D (g)[h](x) ( q(x)p(x)+1 ). (66) Because qgq_g is normalized for every g, ∑x∈qg(x)=1. _x q_g(x)=1. Taking the directional derivative of this identity in the direction h gives ∑x∈DΦ(g)[h](x)=0. _x D (g)[h](x)=0. Hence the constant +1+1 in equation 66 makes no contribution, and D(g)[h]=∑x∈DΦ(g)[h](x)logq(x)p(x).DG(g)[h]= _x D (g)[h](x) q(x)p(x). Since q(x)p(x)=π(x)ρq(x)π(x)ρp(x)=ρq(x)ρp(x), q(x)p(x)= π(x) _q(x)π(x) _p(x)= _q(x) _p(x), we have logq(x)p(x)=ℓ(x). q(x)p(x)= (x). Taking h=Δh= , and using DΦ(g)[Δ](x)=q(x)Δ(x)D (g)[ ](x)=q(x) (x), therefore gives D(g)[Δ] DG(g)[ ] =∑x∈q(x)Δ(x)ℓ(x) = _x q(x) (x) (x) =−∑x,y∈Jp,q(x,y)ℓ(x), =- _x,y J_p,q(x,y) (x), where the second equality follows from equation 59. We now use antisymmetry of the current. Relabelling x and y in the double sum gives ∑x,y∈Jp,q(x,y)ℓ(y) _x,y J_p,q(x,y) (y) =∑x,y∈Jp,q(y,x)ℓ(x) = _x,y J_p,q(y,x) (x) =−∑x,y∈Jp,q(x,y)ℓ(x). =- _x,y J_p,q(x,y) (x). Therefore, ∑x,y∈Jp,q(x,y)ℓ(x) _x,y J_p,q(x,y) (x) =12∑x,y∈Jp,q(x,y)(ℓ(x)−ℓ(y)). = 12 _x,y J_p,q(x,y) ( (x)- (y) ). Using the definition of the current equation 9, we obtain D(g)[Δ] DG(g)[ ] =−12∑x,y∈π(x)K(x,y)mρq(x,y)(ℓ(x)−ℓ(y))2 =- 12 _x,y π(x)K(x,y)m_ _q(x,y) ( (x)- (y) )^2 =−ℐK(q∥p). =-I_K(q\,\|\,p). Since G is smooth, Taylor’s theorem now gives KL(qg+ηΔg∥p) (q_g+η _g\,\|\,p ) =(g+ηΔ) =G(g+η ) =(g)+ηD(g)[Δ]+O(η2) =G(g)+η DG(g)[ ]+O(η^2) =KL(qg∥p)−ηℐK(qg∥p)+O(η2), =KL(q_g\,\|\,p)- _K(q_g\,\|\,p)+O(η^2), which proves equation 21. It remains to establish strict descent away from the target. Every summand in ℐK(q∥p)I_K(q\,\|\,p) is nonnegative. Moreover, since π(x)>0π(x)>0 for every x∈x and mρq(x,y)>0m_ _q(x,y)>0 whenever ρq(x),ρq(y)>0 _q(x), _q(y)>0, the identity ℐK(q∥p)=0I_K(q\,\|\,p)=0 implies that, for every pair x,y∈x,y such that K(x,y)>0K(x,y)>0, (ℓ(x)−ℓ(y))2=0. ( (x)- (y) )^2=0. Hence ℓ(x)=ℓ(y)whenever K(x,y)>0. (x)= (y) K(x,y)>0. By irreducibility, for any x,y∈x,y there exists a finite path x=x0,x1,…,xn=yx=x_0,x_1,…,x_n=y such that K(xi,xi+1)>0,i=0,…,n−1.K(x_i,x_i+1)>0, i=0,…,n-1. Hence ℓ(x0)=ℓ(x1)=⋯=ℓ(xn), (x_0)= (x_1)=·s= (x_n), and therefore ℓ(x)=ℓ(y) (x)= (y). Since x and y are arbitrary, ℓ is constant on X. Thus, for some a∈ℝa , logρq(x)ρp(x)=afor all x∈, _q(x) _p(x)=a all x , or equivalently q(x)=eap(x)for all x∈.q(x)=e^ap(x) all x . Since both q and p are probability distributions, 1=∑xq(x)=ea∑xp(x)=ea,1= _xq(x)=e^a _xp(x)=e^a, so a=0a=0 and q=pq=p. Hence ℐK(q∥p)=0⟺q=p.I_K(q\,\|\,p)=0 q=p. Consequently, if qg≠pq_g≠ p, then ℐK(qg∥p)>0I_K(q_g\,\|\,p)>0. For sufficiently small |η||η|, write the Taylor remainder as R(η)R(η), with |R(η)|≤Cη2|R(η)|≤ Cη^2 for some constant C>0C>0. Then KL(qg+ηΔg∥p)−KL(qg∥p) (q_g+η _g\,\|\,p )-KL(q_g\,\|\,p) =−ηℐK(qg∥p)+R(η) =- _K(q_g\,\|\,p)+R(η) ≤−ηℐK(qg∥p)+Cη2 ≤- _K(q_g\,\|\,p)+Cη^2 =−η(ℐK(qg∥p)−Cη). =-η (I_K(q_g\,\|\,p)-Cη ). For 0<η<ℐK(qg∥p)C,0<η< I_K(q_g\,\|\,p)C, the right-hand side is strictly negative. Thus the update strictly decreases KL(qg∥p)KL(q_g\,\|\,p) for every sufficiently small positive η. ∎ D.4 Particle realization and minimality The particle construction rests on the elementary observation that an antisymmetric signed current can be represented by nonnegative directional jump rates by assigning each edge only its positive current component. The following lemma records the resulting generator and its induced law-level evolution. Lemma 5 (Properties of the positive-part generator). Let p and q be strictly positive probability mass functions, and let Ap,qA_p,q be defined by equation 25–equation 26. Then Ap,qA_p,q is a continuous-time Markov generator and (qAp,q)(x)=−∑y∈Jp,q(x,y)=q(x)Δp,q(x).(qA_p,q)(x)=- _y J_p,q(x,y)=q(x) _p,q(x). (67) If equation 29 holds, then Pη,p,q=I+ηAp,qP_η,p,q=I+η A_p,q is a Markov transition kernel and qPη,p,q=q+ηq⊙Δp,q.qP_η,p,q=q+η q _p,q. (68) Proof. By construction, the off-diagonal entries of Ap,qA_p,q are nonnegative, while its diagonal entries are chosen so that every row sums to zero. Hence Ap,qA_p,q is a continuous-time Markov generator. We regard probability distributions as row vectors. For every x∈x , (qAp,q)(x) (qA_p,q)(x) =∑y∈q(y)Ap,q(y,x) = _y q(y)A_p,q(y,x) =∑y≠xq(y)ap,q(y,x)−q(x)∑y≠xap,q(x,y) = _y≠ xq(y)a_p,q(y,x)-q(x) _y≠ xa_p,q(x,y) =∑y≠x([Jp,q(y,x)]+−[Jp,q(x,y)]+). = _y≠ x ([J_p,q(y,x)]_+-[J_p,q(x,y)]_+ ). Using antisymmetry, Jp,q(y,x)=−Jp,q(x,y)J_p,q(y,x)=-J_p,q(x,y), together with the scalar identity [−r]+−[r]+=−r,[-r]_+-[r]_+=-r, we obtain (qAp,q)(x) (qA_p,q)(x) =−∑y≠xJp,q(x,y) =- _y≠ xJ_p,q(x,y) =−∑y∈Jp,q(x,y) =- _y J_p,q(x,y) =q(x)Δp,q(x), =q(x) _p,q(x), where Jp,q(x,x)=0J_p,q(x,x)=0. This proves equation 67. It remains to verify the discrete-time step. For x≠yx≠ y, Pη,p,q(x,y)=ηap,q(x,y)≥0,P_η,p,q(x,y)=η a_p,q(x,y)≥ 0, while Pη,p,q(x,x)=1−η∑y≠xap,q(x,y)≥0P_η,p,q(x,x)=1-η _y≠ xa_p,q(x,y)≥ 0 by equation 29. The rows of Pη,p,qP_η,p,q sum to one because the rows of Ap,qA_p,q sum to zero. Thus Pη,p,qP_η,p,q is a Markov transition kernel. Finally, qPη,p,q=q(I+ηAp,q)=q+ηqAp,q=q+ηq⊙Δp,q,qP_η,p,q=q(I+η A_p,q)=q+η qA_p,q=q+η q _p,q, which proves equation 68. ∎ Proof of Proposition 3. For q=qgq=q_g, define Ap,g:=Ap,qg,Pη,g:=Pη,p,qg=I+ηAp,qg.A_p,g:=A_p,q_g, P_η,g:=P_η,p,q_g=I+η A_p,q_g. (69) By Lemma 5 and the step-size condition equation 29, evaluated at q=qgq=q_g, Pη,gP_η,g is a Markov transition kernel. The law-level and logit-level drifts coincide at qgq_g. Indeed, ρg(x)=eg(x)Zπ(g) _g(x)= e^g(x)Z_π(g) implies δqg(x,y) _q_g(x,y) =logρg(y)−logρg(x) = _g(y)- _g(x) =g(y)−g(x) =g(y)-g(x) =δg(x,y). = _g(x,y). Using equation 24, we therefore obtain Δp,qg(x) _p,q_g(x) =−1qg(x)∑y∈Jp,qg(x,y) =- 1q_g(x) _y J_p,q_g(x,y) =∑y∈K(x,y)L(δg(x,y))(δg(x,y)−sp(x,y)) = _y K(x,y)L ( _g(x,y) ) ( _g(x,y)-s_p(x,y) ) =Δg(x). = _g(x). The exact law-level Euler identity from Lemma 5 consequently becomes qgPη,g=qg+ηqg⊙Δg.q_gP_η,g=q_g+η q_g _g. (70) On the other hand, Proposition 2 gives qg+ηΔg=qg+ηqg⊙Δg+O(η2).q_g+η _g=q_g+η q_g _g+O(η^2). (71) Comparing the two expressions yields qg+ηΔg=qgPη,g+O(η2),q_g+η _g=q_gP_η,g+O(η^2), which proves equation 31. ∎ The positive-part realization has an additional useful property: among all nonnegative pairs of directional rates producing the same signed current across an edge, it introduces the least total particle traffic. We make this precise next. Proposition 6 (Edgewise minimality of the positive-part realization). Fix strictly positive probability mass functions p,qp,q and distinct states x,y∈x,y . Suppose that b(x,y),b(y,x)≥0b(x,y),b(y,x)≥ 0 are alternative jump rates satisfying q(x)b(x,y)−q(y)b(y,x)=Jp,q(x,y).q(x)b(x,y)-q(y)b(y,x)=J_p,q(x,y). (72) Then there exists a unique rx,y≥0r_x,y≥ 0 such that q(x)b(x,y) q(x)b(x,y) =[Jp,q(x,y)]++rx,y, =[J_p,q(x,y)]_++r_x,y, (73) q(y)b(y,x) q(y)b(y,x) =[−Jp,q(x,y)]++rx,y. =[-J_p,q(x,y)]_++r_x,y. (74) Consequently, the total probability traffic across the edge x,y\x,y\ satisfies q(x)b(x,y)+q(y)b(y,x)=|Jp,q(x,y)|+2rx,y≥|Jp,q(x,y)|.q(x)b(x,y)+q(y)b(y,x)=|J_p,q(x,y)|+2r_x,y≥|J_p,q(x,y)|. (75) Equality holds if and only if b(x,y)=ap,q(x,y),b(y,x)=ap,q(y,x).b(x,y)=a_p,q(x,y), b(y,x)=a_p,q(y,x). Proof. Set u:=q(x)b(x,y),v:=q(y)b(y,x),J:=Jp,q(x,y).u:=q(x)b(x,y), v:=q(y)b(y,x), J:=J_p,q(x,y). Then u,v≥0u,v≥ 0, and the current constraint becomes u−v=J.u-v=J. We first characterize all nonnegative pairs (u,v)(u,v) satisfying this constraint. If J≥0J≥ 0, then [J]+=J,[−J]+=0,[J]_+=J, [-J]_+=0, and u−v=Ju-v=J implies u=J+vu=J+v. Hence, setting rx,y:=v≥0r_x,y:=v≥ 0, u=[J]++rx,y,v=[−J]++rx,y.u=[J]_++r_x,y, v=[-J]_++r_x,y. If J<0J<0, then [J]+=0,[−J]+=−J,[J]_+=0, [-J]_+=-J, and u−v=Ju-v=J is equivalently v=−J+uv=-J+u. Setting rx,y:=u≥0r_x,y:=u≥ 0 again gives u=[J]++rx,y,v=[−J]++rx,y.u=[J]_++r_x,y, v=[-J]_++r_x,y. Thus, in either case, rx,y=minu,v,r_x,y= \u,v\, which also shows that rx,yr_x,y is uniquely determined. Adding the two directional traffic terms yields u+v=[J]++[−J]++2rx,y=|J|+2rx,y.u+v=[J]_++[-J]_++2r_x,y=|J|+2r_x,y. Since rx,y≥0r_x,y≥ 0, the total traffic is bounded below by |J||J|, with equality if and only if rx,y=0r_x,y=0. When rx,y=0r_x,y=0, u=[J]+,v=[−J]+.u=[J]_+, v=[-J]_+. Recalling the definitions of u and v, this gives b(x,y)=[Jp,q(x,y)]+q(x)=ap,q(x,y),b(x,y)= [J_p,q(x,y)]_+q(x)=a_p,q(x,y), and, using antisymmetry of the current, b(y,x)=[−Jp,q(x,y)]+q(y)=[Jp,q(y,x)]+q(y)=ap,q(y,x).b(y,x)= [-J_p,q(x,y)]_+q(y)= [J_p,q(y,x)]_+q(y)=a_p,q(y,x). Therefore the positive-part construction is the unique edgewise realization attaining the minimum total directional traffic. ∎ The minimal realization also has the expected behavior at equilibrium: once the model law reaches the target, not only does the marginal distribution stop evolving, but every particle jump rate vanishes. Corollary 7 (Pathwise stationarity at the target). If q=pq=p, then Jp,p(x,y)=0,ap,p(x,y)=0J_p,p(x,y)=0, a_p,p(x,y)=0 for all x≠yx≠ y. Hence Ap,p=0,Pη,p,p=I,A_p,p=0, P_η,p,p=I, and particles remain at their current states almost surely. Proof. If q=pq=p, then ρq=ρp _q= _p, and hence logρq(x)ρp(x)=0for every x∈. _q(x) _p(x)=0 every x . The current equation 9 therefore vanishes identically. It follows immediately from equation 25 that every off-diagonal jump rate is zero, so Ap,p=0A_p,p=0 and Pη,p,p=IP_η,p,p=I. Thus the particle realization is pathwise stationary at the target. ∎