Flux-Space Flow Matching in 2D Compact U(1) with Spatial β-Conditioning

1. Introduction

Lattice Gauge Theory (LGT) provides the non-perturbative framework for studying gauge fields. The standard method for numerical estimation is Hybrid Monte Carlo (HMC). However, as the lattice spacing decreases towards the continuum limit, HMC algorithms exhibit critical slowing down and topological freezing, which increase the autocorrelation time of observables.

Generative neural sampling methods aim to mitigate this by producing independent configurations [1,2,3]. Normalizing Flows, specifically Gauge Equivariant Flows [4], have demonstrated the ability to sample without autocorrelation by building symmetries directly into the network architecture. Following proof-of-principle studies [4,5], recent work has focused on scalable architectures for high-dimensional manifolds [2,6]. However, the standard approach relies on equivariant coupling layers [7]. While these architectures rigorously preserve symmetry, they often incur prohibitive costs due to determinant computations or restricted expressivity.

In this work, we are motivated by classic duality and Villain-type analyses of the 1980s [8,9], which make the topological excitation content of compact Abelian models explicit and highlight the utility of flux/defect-based representations. Historically, such duality transformations have allowed physicists to map a system with complex constraints or strong couplings into a dual space where the degrees of freedom decouple or become weakly interacting, offering a significant computational advantage. We leverage this principle to circumvent the limitations of direct link-based sampling.

We propose Flux-Space Flow Matching (FFM), a framework for Compact $U\left(1\right)$ LGT that integrates three specific components.

First, we transform the problem usage from link variables to gauge-invariant dual flux variables. In this representation, the local Bianchi identity becomes vacuous, yet we explicitly acknowledge that the global constraint imposed by the toroidal geometry remains.

Second, to address critical slowing down, we implement Spatial $\beta$-Conditioning, replacing global scalar conditioning with a local field.

Third, we build on the Flow Matching framework [10], which regresses the neural vector field onto optimal transport paths, offering superior stability over likelihood-based training. This formulation relates to the Probability Flow ODE of Score-Based Models [11], identifying the generative process as a deterministic flow [12].

We demonstrate that FFM enables scalable sampling on $L=48$ lattices with stable Metropolis correction, and—by validating against thermodynamic, confinement, and topology-sensitive observables—show that it achieves high-fidelity ensembles, bridging geometric deep learning with lattice-sampling methodology grounded in statistical mechanics.

2. Background

2.1. Link to Flux Transformation

Transformation. Standard LGT is formulated in terms of link variables $U_{\mu }\left(x\right)=e^{i{A}_{\mu }\left(x\right)}\in U\left(1\right)$ defined on oriented edges $(x,\mu )$ of a periodic $L\times L$ square lattice $\mathsf{\Lambda }$ (a 2-torus). Here $x\in \mathsf{\Lambda }$ denotes a lattice site, $\mu ,\nu \in \{1,2\}$ are direction indices, and $\widehat{\mu }$ is the unit vector in direction $\mu$. The Wilson action reads

$$S\left[U\right]=-\beta \sum _{x\in \mathsf{\Lambda }}\sum _{\mu <\nu }\Re \left[U_{\mu }\left(x\right)U_{\nu }(x+\widehat{\mu })U_{\mu }^{†}(x+\widehat{\nu })U_{\nu }^{†}\left(x\right)\right].$$

(1)

This link formulation induces a dense interaction graph: each link participates in two adjacent plaquettes, producing strong local couplings that lead to a stiff energy landscape and critical slowing down.

Figure 1. Link to Flux Transformation. Visualization of the mapping from gauge-redundant link variables (edges) to gauge-invariant flux variables (plaquettes).

Figure 1. Link to Flux Transformation. Visualization of the mapping from gauge-redundant link variables (edges) to gauge-invariant flux variables (plaquettes).

To mitigate this, we change variables from links to plaquette degrees of freedom. Let $p=(x;\mu ,\nu )$ denote an oriented plaquette and define the Wilson plaquette variable

$$P_p\equiv P_{\mu \nu }\left(x\right)\in U\left(1\right),P_p=e^{i{\varphi }_p},{\varphi }_p\in (-\pi ,\pi ].$$

(2)

Throughout this work, we use the term flux to refer to the plaquette angle ${\varphi }_p$ (equivalently, to the plaquette phase $P_p$). In this basis the action becomes

$$S\left[\varphi \right]=-\beta \sum _{p}cos\left({\varphi }_p\right).$$

(3)

Let $V:=L^2$ denote the lattice volume; in $D=2$ on an $L\times L$ torus, V equals the total number of plaquettes (and also the number of sites). Ignoring global constraints and measure effects, the local statistical weight is well-approximated by the factorized form

$$\pi \left(\varphi \right)\propto \prod _{p=1}^V{e}^{\beta cos\left({\varphi }_p\right)},$$

(4)

which yields a locally factorized optimization and sampling landscape compared to the coupled link formulation. Throughout, our goal is computational efficacy and empirical correctness of lattice observables, not a one-to-one reproduction of a specific analytic dual construction; in particular, the link-to-plaquette map is many-to-one, and the induced plaquette measure includes gauge redundancy and (on the torus) global degrees of freedom.

Absence of Local Constraints. A potential concern is whether an arbitrary plaquette field $\left\{{\varphi }_p\right\}$ corresponds to a valid underlying link configuration. In two dimensions, the usual local Bianchi identity $d\varphi =0$ is vacuous because there are no 3-cells ($C^3=0$), so no local geometric consistency constraint couples neighboring plaquettes. Moreover, as detailed in Appendix A, on any simply connected patch $\mathsf{\Omega }$ with boundary, $H^2(\mathsf{\Omega };\mathbb{Z})=0$ implies local realizability: every plaquette configuration on $\mathsf{\Omega }$ can be reproduced by a (gauge-fixed) assignment of links. Thus, in $D=2$ there are no local constraints that re-introduce dynamical coupling beyond the factorized plaquette action.

Global Sector Constraint on the Torus. Consequently, the remaining stiffness is global and topological on $T^2$, manifesting as the quantization of the total flux (topological charge)

$$Q\equiv {\displaystyle \frac{1}{2\pi }}\sum _p{\varphi }_p\in \mathbb{Z},$$

(5)

with the precise lifting/branch structure discussed in Appendix A. We enforce this global sector constraint via a manifold-style projector (“Relax-and-Project”) and use an independent Metropolis–Hastings accept/reject step as a bias-control mechanism (Section 3.3).

2.2. Spatial $\beta$-Conditioning

A key architectural choice in our conditional generator is to avoid injecting a global conditioning vector into every layer and instead present the inverse coupling $\beta$ as a spatial channel. Concretely, we lift the scalar $\beta$ to a field $B\left(x\right)\in {\mathbb{R}}^{H\times W}$ (constant in the homogeneous case, but representable as spatially varying in general) and concatenate it with the flux representation at the network input:

$$\mathrm{Input}\left(x\right)=\left[cos\varphi \left(x\right),sin\varphi \left(x\right)\right]\oplus B\left(x\right).$$

(6)

This is a reasonable inductive bias for locality and translation-equivariance: $\beta$ is provided in the same tensor format as the physical degrees of freedom, and its influence on the generated configuration is mediated through the same local convolutional computations rather than through per-layer global modulation pathways (e.g., FiLM-style feature-wise affine transforms).

This design is particularly relevant near criticality. Markov-chain methods propagate information through a sequence of updates, so long-range correlations emerge only after many steps and suffer from critical slowing down. In contrast, our model produces an entire configuration in parallel. Although the conditioning enters locally via $B\left(x\right)$, the deep hierarchical network has an effective receptive field that can span the lattice (via depth and/or dilation), enabling the synthesis of the long-range correlations required at large correlation length within a single forward pass. This does not eliminate the need for an exactness mechanism when required (e.g., a global accept/reject correction), but it amortizes the construction of correlated configurations across $\beta$ values within one trained model.

3. Generative Architecture: FluxUNet

3.1. Network Architecture (FluxUNet)

To implement the FFM framework, we parameterize the time-dependent velocity field $v_t\left(\varphi \right)$ using FluxUNet, a specialized Residual U-Net designed for compact angular variables on a toroidal lattice. Several architectural choices are motivated by the geometry of the lattice and by numerical stability of continuous-time flows:

• Toroidal Topology ($T^2$). All convolutional layers employ circular padding, enforcing periodic boundary conditions and preventing boundary artifacts.
• Residual Parameterization. We use residual blocks, $h_{\mathrm{out}}=h_{\mathrm{in}}+f\left(h_{\mathrm{in}}\right)$, which bias the network toward incremental feature updates rather than abrupt re-mappings. This typically stabilizes optimization and reduces spurious high-frequency responses. In our setting, it yields smoother spatial predictions under periodic convolutions and improves numerical robustness.
• Smooth Nonlinearity and Normalization. We use SiLU activations and Group Normalization ($G=8$) instead of ReLU or BatchNorm. SiLU provides a smooth nonlinearity, and GroupNorm avoids dependence on batch statistics, which is beneficial when training vector fields with small or varying batch sizes.

3.2. Training Strategy: Ensemble Distillation

We generate the training corpus using massively parallel persistent local Metropolis chains with a checkerboard update scheme; all hyperparameters (warmup, decorrelation interval, proposal amplitude, and $\beta$ scheduling) are documented in Appendix B. While individual short chains may remain trapped in metastable regions of configuration space, aggregating many chains initialized from random hot starts ($N_{\mathrm{chains}}\gg 1$) yields broad empirical coverage of typical configurations.

We observe an implicit denoising effect consistent with the spectral bias of convolutional networks [13]: the model preferentially captures low-frequency, large-scale structure while being less prone to memorizing high-frequency artifacts present in imperfectly equilibrated training sets. We therefore treat FluxUNet as an ensemble distiller in the empirical sense that its generated configurations exhibit improved agreement of lattice observables with reference results. Any claim of asymptotic correctness is instead attributed to the Metropolis correction described below.

Prior flow-based samplers for LGT often employ stacks of discrete invertible coupling layers [4]. In contrast, we adopt a continuous-time formulation in flux space and train the dynamics with an ODE-style flow matching objective. This provides a smooth transport/denoising learning signal well suited to compact angular variables, and it prioritizes accurate reproduction of lattice observables rather than an explicit analytic dual construction.

3.3. Sector-Projected Metropolis-Corrected Sampling (IMH)

Our simulations impose periodic boundary conditions, so the lattice is topologically a torus and flux configurations decompose into integer total-flux sectors. We therefore include an explicit sector projection step that removes residual non-integer drift in the total flux by shifting only the global mean. We then apply an Independent Metropolis–Hastings (IMH) correction using the model’s estimated proposal density.

Let the target density be $\pi \left(\varphi \right)\propto e^{-S_{\mathrm{phys}}\left(\varphi \right)}$. One IMH step proceeds as:

1.

Base distribution. Sample $z\sim p_0\left(z\right)=\mathcal{U}{[-\pi ,\pi ]}^V$.

2.

Flow proposal (CNF). Integrate

$${\displaystyle \frac{d\varphi }{dt}}=v_t\left(\varphi \right),t\in [0,1],$$

(7)

using an adaptive Runge–Kutta solver with tolerance $(\mathrm{atol},\mathrm{rtol})={10}^{-5}$ to obtain ${\varphi }^{\prime }=G\left(z\right)$. In parallel, we track an estimated log proposal density $log\widehat{q}\left({\varphi }^{\prime }\right)$ by integrating the instantaneous change-of-variables equation for the flow (divergence of $v_t$), using a Hutchinson trace estimator.

3.

Toroidal wrapping. Map the angles to the principal branch

$${\varphi }_{\mathrm{wrap}}\left(x\right)=({\varphi }^{\prime }\left(x\right)+\pi mod2\pi )-\pi .$$

(8)

4.

Nearest-sector projection (global constraint). Compute the total wrapped flux

$${\mathsf{\Phi }}_{\mathrm{tot}}=\sum _{x\in \mathsf{\Lambda }}{\varphi }_{\mathrm{wrap}}\left(x\right),Q=\mathrm{round}\left({\displaystyle \frac{{\mathsf{\Phi }}_{\mathrm{tot}}}{2\pi }}\right),$$

(9)

and project by shifting only the global mean:

$$\tilde{\varphi }\left(x\right)={\varphi }_{\mathrm{wrap}}\left(x\right)-{\displaystyle \frac{1}{V}}\left({\mathsf{\Phi }}_{\mathrm{tot}}-2\pi Q\right).$$

(10)

This step preserves the nearest integer sector label Q and enforces ${\sum }_x\tilde{\varphi }\left(x\right)=2\pi Q$ exactly.

5.

IMH acceptance Accept $\tilde{\varphi }$ with probability $\alpha =min(1,r)$. The acceptance ratio is:

$$r={\displaystyle \frac{\pi \left(\tilde{\varphi }\right)\widehat{q}\left({\varphi }_{\mathrm{curr}}\right)}{\pi \left({\varphi }_{\mathrm{curr}}\right)\widehat{q}\left(\tilde{\varphi }\right)}}=exp\left(-S_{\mathrm{phys}}\left(\tilde{\varphi }\right)+S_{\mathrm{phys}}\left({\varphi }_{\mathrm{curr}}\right)-log\widehat{q}\left(\tilde{\varphi }\right)+log\widehat{q}\left({\varphi }_{\mathrm{curr}}\right)\right).$$

(11)

Note on proposal-density estimation and exactness. The sector projection (10) is a deterministic global constraint enforcement that maps unconstrained proposals onto a codimension-one manifold (and is piecewise-defined due to the nearest-sector rounding). As a result, the exact proposal density of the projected samples requires a measure-consistent treatment on the constrained manifold. In addition, the divergence entering $log\widehat{q}$ is estimated stochastically via Hutchinson, which is unbiased for the trace but does not yield an exact likelihood when exponentiated. For these reasons, the IMH step in (11) should be interpreted as a Metropolis-corrected, bias-limited correction rather than a strictly exact sampler. We provide a quantitative empirical validation of the method’s fidelity and efficiency in the following section.

4. Results: Empirical Validation

We focus our evaluation on Compact $U\left(1\right)$ LGT on $L=48$ lattices. As detailed in Section 3.2, we utilize a Distilled Ensemble strategy to train on $1.024\times {10}^5$ samples ($\beta \in [0.1,6]$) using a single NVIDIA A100 GPU for 20 epochs. All sampling benchmarks employ an adaptive Runge–Kutta integrator (DOPRI5/RK45). We analyze the results from five perspectives.

4.1. Thermodynamic Consistency Checks

We quantify precision by monitoring the Energy Density during training. As shown in Figure 2, the model converges rapidly to the physical manifold. The final ensemble achieves a Mean Absolute Error (MAE) of 0.00095 relative to the theoretical prediction. For a detailed validation against the Von Mises distribution, see Appendix C.

This allows for precise spectroscopy of the topological susceptibility ${\chi }_t=\langle Q^2\rangle /V$ (Figure 3). Notably, our measured susceptibility slightly exceeds the perturbative non-compact limit (${\chi }_t\propto 1/\beta$). This deviation is physically robust, reflecting the additional contribution from non-perturbative vortex excitations characteristic of the Compact $U\left(1\right)$ theory, confirming that the generative model captures full non-perturbative dynamics beyond the Gaussian approximation. The deviation at high $\beta$ reflects expected Finite Volume Effects as the correlation length $\xi$ approaches the lattice size L. Finite-$\theta$ reconstructions via phase reweighting and sector-aware reconstruction, reported in Appendix E, further confirm the pipeline’s internal consistency within a statistically controlled analysis window.

4.2. Sampling Efficiency Benchmarks

Performance analysis.Table 1 demonstrates a clear separation between (i) mixing, quantified by ${\tau }_{\mathrm{int}}$, and (ii) absolute wall-clock throughput, quantified by ESS/s. At low coupling (e.g., $\beta =2.0$), HMC benefits from inexpensive leapfrog updates and achieves substantially higher ESS/s for Q (about $7\times$ higher than FFM in our benchmark). In this regime, our FFM sampler is not wall-clock competitive because the DOPRI5 proposal requires repeated model evaluations and divergence estimation. This is a compute-cost issue rather than a mixing issue, since ${\tau }_{\mathrm{int}}$ remains order unity.

At larger $\beta$, HMC exhibits pronounced topological freezing: at $\beta =6.0$, ${\tau }_{\mathrm{int}}\left(Q\right)$ becomes extremely large with substantial run-to-run variability ($782\pm 1055$), indicating practically non-ergodic behavior on the simulated timescale. In contrast, the conditional FFM sampler maintains ${\tau }_{\mathrm{int}}\left(Q\right)=\mathcal{O}\left(1\right)$ across the tested $\beta$ values, yielding a $>500\times$ reduction in ${\tau }_{\mathrm{int}}$ at $\beta =6.0$.

Crucially, this advantage is amortized across parameter space: the same trained conditional model is used for all $\beta$, whereas HMC typically requires re-thermalization (and often re-tuning of integrator hyperparameters) at each $\beta$. Consequently, despite a higher per-step cost, FFM crosses over in absolute wall-clock efficiency for the topological observable at $\beta \ge 5.0$ (ESS/s in Table 1), precisely where HMC begins to freeze.

4.3. Confinement Diagnostics: V(R) and Creutz Ratios

Using our FFM sampler, we successfully resolved the static potential $V\left(R\right)$ across two distinct physical regimes (Figure 4). The comparison reveals that FFM correctly captures the expected $\beta$-dependence of confinement observables:

• Strong Coupling ($\beta =2.2$): The system exhibits harsh confinement characteristics. The potential rises steeply ($\sigma \approx 0.32$), and Wilson loop signals decay rapidly, reflecting a "stiff" gauge field.
• Scaling Regime ($\beta =4.5$): As the system approaches the continuum limit, the string tension softens significantly to $\sigma \approx 0.12$. Remarkably, the model preserves the linear confining potential even in this weak coupling limit, consistent with confinement in 2D compact U(1) across couplings.

A striking feature of the static potentials is the Ground State Dominance. The flatness of effective potential plateaus suggests early plateau behavior / limited excited-state contamination within the accessible T-window for loop sizes with reliable signal.

4.4. Multiscale Sensitivity (ERF Analysis)

To probe how FluxUNet aggregates information across length scales, we computed the Effective Receptive Field (ERF) by measuring the sensitivity of representative internal activations to localized perturbations of the input flux field (Figure 5). The ERF maps reveal a clear separation between shallow, predominantly local processing and deep, globally aggregated processing:

• Local sensitivity (short-range features). In early and intermediate layers, the ERF is sharply localized: the response is concentrated in a compact neighborhood around the perturbed site. This indicates that these layers primarily encode short-range structure, such as local plaquette-angle fluctuations and short-wavelength textures. Such locality is consistent with the fact that many thermodynamic contributions are governed by local statistics, and it provides an inductive bias that avoids unnecessarily entangling distant regions when learning microscopic structure.
• Global sensitivity (long-range features). In the bottleneck layer, the ERF becomes broad and diffuse, indicating that deep features depend on information distributed across the lattice. This pattern is consistent with the network forming global summaries that cannot be inferred from a small patch alone. In particular, sector-dependent structure on a torus (e.g., correlations tied to the total-flux sector Q or other global modes) is inherently nonlocal; representing such effects requires access to long-range context. We emphasize that ERF demonstrates capacity for global dependence, not that the network explicitly computes Q.

4.5. Generalization to Spatial $\beta$-Conditioning

Standard algorithms typically require re-thermalization or careful annealing for every new potential landscape. In contrast, our model—trained with Spatial $\beta$-Conditioning—exhibits a qualitative zero-shot response to unseen conditions. To demonstrate this, we evaluate the model on a "Dipole Trap" configuration (Figure 6) without any retraining. We engineer a spatial coupling profile $\beta \left(x\right)$ with a vacuum background ($\beta =4.0$) and two localized "hot spots" ($\beta =0.5$). We then compute the local magnetic flux density directly from the generated plaquettes:

$$q\left(x\right)={\displaystyle \frac{1}{2\pi }}\varphi \left(x\right)\in (-0.5,0.5],$$

(12)

where $\varphi \left(x\right)$ is the principal branch of the plaquette angle.

The results reveal a sharp statistical separation: the model consistently generates a concentration of positive flux ($q>0$) in one trap and negative flux ($q<0$) in the other. This emergence of a conjugate pair is not random; it allows the system to satisfy the global topological neutrality of the torus ($Q=0$) while entropically localizing the necessary fluctuations within the low-stiffness regions where the energetic cost is minimal ($\beta =0.5$). This behavior is consistent with the network having learned a local conditional dependence $P\left(\varphi \right|\beta )$. This represents a functional capability beyond standard algorithms: the model can simulate domain walls and flux defects zero-shot, whereas standard HMC would typically require running a new chain for every variation in the $\beta$ profile.

5. Discussion

Dual Flux Transformation

Our results suggest that, in 2D compact $U\left(1\right)$ at large $\beta$, topology-sensitive slowing down is severe for standard link-based HMC on accessible timescales, whereas operating in gauge-invariant flux (plaquette-angle) variables yields proposals with markedly improved mixing for topology-sensitive observables. This performance leverages the factorization of the action in the dual basis (Section 2), which effectively decouples the local dynamics.

Similarly, the “Relax-and-Project” mechanism (Section 3.3) ensures global topological validity without restricting the local generative freedom. We then wrap these proposals in an Independent Metropolis–Hastings (IMH) step using an estimated proposal density (via a stochastic trace/divergence estimator). While the stochastic estimator precludes strict detailed balance, our validation confirms that the model, operating in the transformed dual basis, achieves high fidelity, empirically reproducing sensitive thermodynamic, confinement, and topology-related observables to high accuracy. We note that the current implementation reflects a practical compromise between training throughput and numerical precision; further gains in observable fidelity may be achievable with higher-precision arithmetic at the cost of increased training time.

Sampling Efficiency.

A central empirical advantage of FFM is its robust topological mixing. As detailed in Section 4.2 (Table 1), the HMC baseline exhibits extremely large and unstable integrated autocorrelation times for Q, consistent with non-ergodic behavior. In contrast, the conditional FFM sampler maintains ${\tau }_{\mathrm{int}}\left(Q\right)=\mathcal{O}\left(1\right)$ across the tested $\beta$ values. This “near-i.i.d.” behavior suggests that the FFM proposal distribution has strong phase-space overlap across the relevant topological structure, allowing the chain to decorrelate in a small number of steps. Although the ESS gain varies by observable, we emphasize that our sampling method can be easily parallelized on modern GPUs to achieve superior performance.

Local Topological Control

The Statistical Dipole Trap experiment demonstrates the utility of Spatial $\beta$-Conditioning. The ability to spontaneously nucleate conjugate flux concentrations solely through spatially varying $\beta \left(x\right)$ profiles is consistent with the model having learned the local conditional statistics of the theory. Rather than treating the inverse coupling $\beta$ as a global scalar, the network utilizes it as a dynamic background field, enabling zero-shot qualitative control over the vacuum structure. As discussed in Section 2, a single trained model samples a continuous range of $\beta$, amortizing training cost across the tested interval and yielding practical efficiency gains significantly greater than single-point ESS comparisons imply.

Theoretical Boundaries: Dimensionality and the Bianchi Constraint

A rigorous assessment must acknowledge the dimensional boundaries of this approach. The efficacy of the dual flux representation described here relies on the "2D Privilege": in two dimensions, plaquettes tile the plane without enclosing 3-volumes. Consequently, the dual variables ${\varphi }_p$ are free from local geometric constraints.

In 3D and 4D, however, fluxes live on the faces of elementary cubes. They must satisfy a nontrivial local geometric constraint (discrete Bianchi identity) up to $2\pi$ integers. Naive independent sampling in 4D would violate this consistency, leading to non-physical states. Therefore, FFM is not claimed as a universal solution for 4D Gauge Theory in its current form, but rather as a competitive alternative for 2D Compact U(1). Within this domain, it provides a highly efficient mechanism for mitigating critical slowing down in the calculation of topological observables.

Scaling Behavior

We validate FFM on lattices up to $L=64$ (Appendix D). Thermodynamic fidelity (energy density MAE $\sim {10}^{-3}$) and topological susceptibility scaling remain robust at larger volumes. The acceptance rate decreases as expected for independent proposal methods in higher-dimensional phase spaces (7–$36\%$ at $L=64$ vs. 50–$70\%$ at $L=48$), while topological autocorrelation ${\tau }_{\mathrm{int}}\left(Q\right)$ remains $\mathcal{O}\left(1\right)$. Confinement observables, including string tension extraction via Creutz ratios, remain stable, confirming that the physics fidelity of the method scales with system size.

6. Conclusion

We have presented Flux-Space Flow Matching (FFM), a framework that resolves topological freezing in 2D Compact $U\left(1\right)$ LGT. By lifting the dynamics to the gauge-invariant dual flux representation, we diagonalize the interaction and render sampling tractable. Combined with a “Relax-and-Project” paradigm, FFM yields an IMH-corrected kernel achieving a $\sim 500\times$ improvement in topological tunneling rates over HMC. Additionally, Spatial $\beta$-Conditioning enables zero-shot synthesis of topological defects. Our results suggest that scalable neural sampling requires architectures that respect the native geometric structure of the theory.

Appendix A.1. Discrete Variables

Let $\mathsf{\Lambda }$ be a 2D square lattice. Link variables are $U_{\mu }\left(n\right)\in U\left(1\right)$ on oriented edges, parameterized as $U_{\mu }\left(n\right)=e^{i{\theta }_{\mu }\left(n\right)}$ with ${\theta }_{\mu }\left(n\right)\in (-\pi ,\pi ]$. The plaquette variable on an elementary face is

$$P\left(n\right)\equiv U_1\left(n\right)U_2(n+\widehat{1})U_1{(n+\widehat{2})}^{-1}{U}_2{\left(n\right)}^{-1}\in U\left(1\right).$$

(A1)

Equivalently, in terms of angles,

$$P\left(n\right)=e^{i\varphi \left(n\right)},\varphi \left(n\right)\equiv {\mathsf{\Delta }}_1{\theta }_2\left(n\right)-{\mathsf{\Delta }}_2{\theta }_1\left(n\right)(mod2\pi ),$$

(A2)

where ${\mathsf{\Delta }}_{\mu }$ is the forward difference operator.

Appendix A.2. Absence of Local Bianchi Constraints in D=2

On a lattice, the Bianchi identity is the vanishing of the exterior derivative of the 2-cochain $\varphi$: $d\varphi =0$, where $d\varphi$ is a 3-cochain evaluated on elementary 3-cells (cubes). In $D=2$, there are no 3-cells; hence the condition $d\varphi =0$ is vacuous. Therefore, unlike $D\ge 3$, there is no local geometric consistency constraint coupling neighboring plaquettes.

Appendix A.3. Local Surjectivity (Constructive Validity on a Patch)

Consider any simply connected local patch $\mathsf{\Omega }\subset \mathsf{\Lambda }$ (e.g., a receptive-field-sized region) with open boundary. Because $\mathsf{\Omega }$ has no nontrivial 2-cycles, $H^2(\mathsf{\Omega };\mathbb{Z})=0$, so every 2-cochain is exact. Concretely, one may explicitly construct links realizing an arbitrary plaquette field $P\left(n\right)$ on $\mathsf{\Omega }$ by an axial-gauge choice. For example, set

$$U_2\left(n\right)\equiv 1\forall n\in \mathsf{\Omega },$$

(A3)

then the plaquette constraint reduces to $P\left(n\right)=U_1\left(n\right)U_1{(n+\widehat{2})}^{-1}$. Fix $U_1$ on the bottom boundary (e.g., $U_1(x,0)=1$) and define recursively for $y\ge 0$,

$$U_1(x,y+1)\equiv P{(x,y)}^{-1}{U}_1(x,y).$$

(A4)

This yields a valid link configuration on $\mathsf{\Omega }$ reproducing exactly the prescribed plaquettes. Hence, locally the map $\left\{U_{\mu }\right\}\mapsto \left\{P\right\}$ is surjective: arbitrary local plaquette fields correspond to some underlying links.

Appendix A.4. Global Topology on the Torus T²

Our simulations use periodic boundary conditions, so the lattice is topologically a torus. In this case $H^2(T^2;\mathbb{Z})\cong \mathbb{Z}$, and plaquette configurations decompose into topological sectors characterized by an integer total flux (first Chern number). To state this without ambiguity from branch cuts, it is convenient to work with a real-valued lift: choose integers $m\left(n\right)\in \mathbb{Z}$ such that

$$\mathsf{\Phi }\left(n\right)=\varphi \left(n\right)+2\pi m\left(n\right)\in \mathbb{R},Q\equiv {\displaystyle \frac{1}{2\pi }}\sum _{n\in \mathsf{\Lambda }}\mathsf{\Phi }\left(n\right)\in \mathbb{Z}.$$

(A5)

Thus, the only constraint on plaquettes imposed by the torus topology is the quantization of the total flux Q. In contrast, there are also two global holonomies (Wilson/Polyakov loops along the two non-contractible cycles) which are not determined by plaquettes, i.e. they are global degrees of freedom required to reconstruct links but do not constrain the local plaquette field.

Appendix A.5. Conclusion

In $D=2$, the absence of local Bianchi constraints implies that the plaquette field behaves locally as an unconstrained scalar field on faces, justifying direct generative modeling of $P\left(n\right)$. Global topology only induces an integer-valued total-flux sector label on $T^2$; no local geometric consistency constraints analogous to $D\ge 3$ are present.

Appendix D.1. Thermodynamic Fidelity

Figure A2 shows the energy density across $\beta \in [0.1,6.0]$. The model achieves MAE $=0.0009$, comparable to the $L=48$ result (MAE $=0.00095$).

Figure A2. Energy density check for $L=64$. MAE $=0.0009$ confirms thermodynamic fidelity scales with lattice size.

Figure A2. Energy density check for $L=64$. MAE $=0.0009$ confirms thermodynamic fidelity scales with lattice size.

Figure A3 further validates the marginal flux distributions at representative $\beta$ values. At $\beta =2.0$, the model achieves an energy error of just $0.008\%$, demonstrating exceptional precision in the intermediate coupling regime.

Figure A3. Von Mises distribution comparison for $L=64$ at $\beta =2.0$ (left, $0.008\%$ error) and $\beta =4.0$ (right, $0.203\%$ error). The tight alignment confirms the model learns the correct Boltzmann statistics.

Figure A3. Von Mises distribution comparison for $L=64$ at $\beta =2.0$ (left, $0.008\%$ error) and $\beta =4.0$ (right, $0.203\%$ error). The tight alignment confirms the model learns the correct Boltzmann statistics.

Appendix D.2. Topological Susceptibility

Figure A4 confirms ${\chi }_t\propto {\beta }^{-1}$ scaling, with values consistently above the non-compact perturbative limit due to vortex contributions.

Figure A4. Topological susceptibility scaling for $L=64$. The power-law ${\chi }_t\propto {\beta }^{-1}$ is preserved.

Figure A4. Topological susceptibility scaling for $L=64$. The power-law ${\chi }_t\propto {\beta }^{-1}$ is preserved.

Appendix D.3. Confinement Diagnostics

Figure A5 shows the static potential and Creutz ratios at $\beta =4.5$. The string tension $\sigma \approx 0.128$ is consistent with the $L=48$ result at the same coupling ($\sigma \approx 0.12$, see Figure 4), and plateaus remain stable for $R\le 6$.

Figure A5. Static potential $V\left(R\right)$ (left) and Creutz ratios (right) for $L=64$ at $\beta =4.5$. Average string tension $\sigma \approx 0.128$, consistent with the $L=48$ benchmark.

Figure A5. Static potential $V\left(R\right)$ (left) and Creutz ratios (right) for $L=64$ at $\beta =4.5$. Average string tension $\sigma \approx 0.128$, consistent with the $L=48$ benchmark.

Appendix D.4. Sampling Efficiency

Table A2 reports FFM sampling efficiency at $L=64$. While acceptance rates are lower than $L=48$ (as expected for independent proposals in larger phase spaces), the topological autocorrelation ${\tau }_{\mathrm{int}}\left(Q\right)$ remains $\mathcal{O}\left(1\right)$.

Table A2. FFM sampling efficiency on $L=64$ 2D $U\left(1\right)$ lattices. Values are mean ± std over three random seeds.

$\beta$	Acc.	${\tau }_{\mathrm{int}}\left(Q\right)$	ESS/s
2.0	$0.36$	$2.9\pm 2.3$	$0.22$
4.0	$0.23$	$4.8\pm 4.5$	$0.15$
5.0	$0.14$	$4.8\pm 2.0$	$0.10$
6.0	$0.07$	$8.4\pm 1.7$	$0.05$

Table A2. FFM sampling efficiency on $L=64$ 2D $U\left(1\right)$ lattices. Values are mean ± std over three random seeds.

$\beta$	Acc.	${\tau }_{\mathrm{int}}\left(Q\right)$	ESS/s
2.0	$0.36$	$2.9\pm 2.3$	$0.22$
4.0	$0.23$	$4.8\pm 4.5$	$0.15$
5.0	$0.14$	$4.8\pm 2.0$	$0.10$
6.0	$0.07$	$8.4\pm 1.7$	$0.05$

Appendix D.5. Conclusion

These results demonstrate that FFM scales successfully from $L=48$ to $L=64$, with physics fidelity preserved across all diagnostics: thermodynamic observables (MAE $\sim {10}^{-3}$), marginal distributions ($<0.3\%$ error), topological susceptibility scaling, and confinement behavior (string tension consistent with $L=48$). While acceptance rates decrease with increasing volume—an expected limitation of independent proposal methods—the topological autocorrelation remains $\mathcal{O}\left(1\right)$, consistent with the $L=48$ benchmarks. The successful scaling to $L=64$ suggests that FFM can extend to even larger lattices with appropriate computational resources.

Appendix E.1. Reweighting Protocol and CP-Even/CP-Odd Extraction

All finite-$\theta$ observables are reconstructed from a $\theta =0$ ensemble using phase reweighting. We define the complex ratio estimator

$$\begin{array}{cc}\hfill {\langle O\rangle }_{\theta }^{\left(\mathrm{raw}\right)}& \equiv {\displaystyle \frac{{〈O{e}^{i\theta Q}〉}_{\theta =0}}{{〈e^{i\theta Q}〉}_{\theta =0}}},\hfill \end{array}$$

(A6)

$$\begin{array}{cc}\hfill Z\left(\theta \right)& \equiv {〈e^{i\theta Q}〉}_{\theta =0}.\hfill \end{array}$$

(A7)

Here $Q\in \mathbb{Z}$ is the (projected) global topological charge and $P\left(Q\right)$ is its empirical distribution in the $\theta =0$ ensemble (sector decomposition).

In principle, CP-even observables yield real expectation values, whereas CP-odd quantities are purely imaginary in the CP-symmetric limit. Accordingly, we extract the reported values as:

$$\begin{array}{cc}\hfill {\langle O_{\mathrm{even}}\rangle }_{\theta }& \equiv \Re \left[{\langle O_{\mathrm{even}}\rangle }_{\theta }^{\left(\mathrm{raw}\right)}\right],\hfill \end{array}$$

(A8)

$$\begin{array}{cc}\hfill {\mathcal{O}}_{\mathrm{odd}}\left(\theta \right)& \equiv \Im \left[{\langle O_{\mathrm{odd}}\rangle }_{\theta }^{\left(\mathrm{raw}\right)}\right],\hfill \end{array}$$

(A9)

$$\begin{array}{cc}\hfill {\chi }_t\left(\theta \right)& \equiv \Re \left[{\displaystyle \frac{1}{V}}\left({\langle Q^2\rangle }_{\theta }^{\left(\mathrm{raw}\right)}-{\left({\langle Q\rangle }_{\theta }^{\left(\mathrm{raw}\right)}\right)}^2\right)\right].\hfill \end{array}$$

(A10)

The explicit $\Re [·]$ in Eq. (A10) removes a small non-physical imaginary component that may arise from finite-statistics asymmetries in the empirical $P\left(Q\right)$.

Figure A6. Reweighting signal strength $\left|Z\right(\theta \left)\right|$ and noise floor, overlaid with effective string tension ${\sigma }_{\mathrm{eff}}\left(\theta \right)$. The signal crosses the noise floor near $\theta \approx 1.0$.

Figure A6. Reweighting signal strength $\left|Z\right(\theta \left)\right|$ and noise floor, overlaid with effective string tension ${\sigma }_{\mathrm{eff}}\left(\theta \right)$. The signal crosses the noise floor near $\theta \approx 1.0$.

Appendix E.2. Analysis Window, Noise Floor, and Interpretation Scope

The reweighting denominator magnitude $\left|Z\right(\theta \left)\right|$ decays rapidly with $\theta$ due to phase cancellations across topological sectors. Once $\left|Z\right(\theta \left)\right|$ becomes comparable to the statistical noise scale $\mathcal{O}(1/\sqrt{N})$, ratio estimators may become noise-dominated and develop spurious structure. For completeness, the measured topological-charge statistics for this run are:

$$\begin{array}{cc}\hfill \mathrm{Mean}\left(Q\right)& =-0.20,\hfill \\ \hfill \mathrm{Std}\left(Q\right)& =3.65,\hfill \\ \hfill \mathrm{Range}\left(Q\right)& \in [-15,16].\hfill \end{array}$$

(A11)

For clarity, we therefore refer to the reported range as an analysis window rather than a “trust window”:

$$\theta \in [0,1.6].$$

(A12)

Within this window, the sub-interval $\theta \lesssim 1.0$ corresponds to the region where $\left|Z\right(\theta \left)\right|$ remains at or above the noise floor for the present statistics (Figure A6), and we treat it as the regime in which reweighted estimates are comparatively well-conditioned. By contrast, the region $\theta \gtrsim 1.0$ should be interpreted primarily as a numerical stability test of the pipeline under severe phase cancellations, not as a domain for precision finite-$\theta$ physics. In particular, near the upper end ($\theta \gtrsim 1.5$), the reconstruction is strongly noise-dominated, consistent with the observed instability in susceptibility-like diagnostics and in larger-loop observables.

Appendix E.3. Observables Within the Analysis Window

Within the analysis window, we examine four representative quantities:

1.

Energy density $E\left(\theta \right)=1-{\langle cos\varphi \rangle }_{\theta }$: The action density, serving as a bulk thermodynamic check. The reconstructed values remain smooth and vary modestly over $[0,1.6]$ (Figure A7).

Figure A7. Energy density $E\left(\theta \right)$ via phase reweighting at $\beta =5.0$ ($N=100,000$).

Figure A7. Energy density $E\left(\theta \right)$ via phase reweighting at $\beta =5.0$ ($N=100,000$).

2.

Wilson loops ${\langle W(R\times R)\rangle }_{\theta }$ for $R\in \{2,4,6\}$: Probes of confinement at multiple length scales. Small and intermediate loops are broadly stable; larger loops exhibit increased sensitivity near $\theta \gtrsim 1.5$ (Figure A8).

Figure A8. Wilson loops ${\langle W(R\times R)\rangle }_{\theta }$ for $R\in \{2,4,6\}$ at $\beta =5.0$.

Figure A8. Wilson loops ${\langle W(R\times R)\rangle }_{\theta }$ for $R\in \{2,4,6\}$ at $\beta =5.0$.

3.

Effective confinement proxy ${\sigma }_{\mathrm{eff}}\left(\theta \right)$: Extracted from a fixed $4\times 4$ loop ($R=T=4$), providing a finite-volume proxy for string tension. Over the majority of the trust window, ${\sigma }_{\mathrm{eff}}\left(\theta \right)$ evolves smoothly (see Figure A6).

4.

CP-odd / topological response: Including ${\langle Q\rangle }_{\theta }$ and susceptibility-like proxy ${\chi }_t\left(\theta \right)$. These quantities are particularly sensitive to phase cancellations; the visible breakdown near $\theta \gtrsim 1.5$ is consistent with the analysis-window boundary (Figure A9).

Figure A9. CP-odd diagnostics (${\langle Q\rangle }_{\theta }$, ${\chi }_t\left(\theta \right)$) via reweighting. Instability near $\theta \gtrsim 1.5$ marks the trust-window boundary.

Figure A9. CP-odd diagnostics (${\langle Q\rangle }_{\theta }$, ${\chi }_t\left(\theta \right)$) via reweighting. Instability near $\theta \gtrsim 1.5$ marks the trust-window boundary.

Appendix E.4. Takeaway

These checks show that, despite the rapid decay of the reweighting signal, our learned sampler combined with sector-aware reconstruction yields qualitatively stable finite-$\theta$ observables across $\theta \in [0,1.6]$. We use these experiments strictly as fidelity tests of the sampling-and-reconstruction pipeline, and we avoid physics-facing interpretation in the noise-dominated regime.