From Independent to Correlated Diffusion: Generalized Generative Modeling with Probabilistic Computers

Abstract: Diffusion models have emerged as a powerful framework for generative tasks in deep learning. They decompose generative modeling into two computational primitives: deterministic neural-network evaluation and stochastic sampling. Current implementations usually place most computation in the neural network, but diffusion as a framework allows a broader range of choices for the stochastic transition kernel. Here, we generalize the stochastic sampling component by replacing independent noise injection with Markov chain Monte Carlo (MCMC) dynamics that incorporate known interaction structure. Standard independent diffusion is recovered as a special case when couplings are set to zero. By explicitly incorporating Ising couplings into the diffusion dynamics, the noising and denoising processes exploit spatial correlations representative of the target system. The resulting framework maps naturally onto probabilistic computers (p-computers) built from probabilistic bits (p-bits), which provide orders-of-magnitude advantages in sampling throughput and energy efficiency over GPUs. We demonstrate the approach on equilibrium states of the 2D ferromagnetic Ising model and the 3D Edwards-Anderson spin glass, showing that correlated diffusion produces samples in closer agreement with MCMC reference distributions than independent diffusion. More broadly, the framework shows that p-computers can enable new classes of diffusion algorithms that exploit structured probabilistic sampling for generative modeling.

• The paper shows independent Bernoulli diffusion is exactly p-bit dynamics with zero couplings, while restoring interaction couplings produces a temperature-controlled correlated diffusion process.
• The method combines a per-site neural estimator with Gibbs-sampled candidate states and likelihood reweighting, reproducing energy, magnetization, and Parisi overlap statistics more accurately than independent diffusion on Ising benchmarks.
• The approach shifts computation toward repeated Gibbs sampling, with raw random-number benchmarks indicating roughly 100× FPGA and projected 10,000× sMTJ energy-efficiency gains over GPUs, although end-to-end inference speed remains unmeasured.

Overview

This paper reformulates discrete diffusion models by generalizing their stochastic transition kernel from independent, site-wise noise injection to Markov chain Monte Carlo (MCMC) dynamics that incorporate known interaction structure (2603.27996). The authors' central observation is that standard Bernoulli diffusion is exactly equivalent to probabilistic-bit (p-bit) dynamics with all couplings set to zero ($J_{ij}=0$), so that the entire noising schedule reduces to a temperature sweep. Restoring physical couplings $J_{ij}\neq 0$ yields a "correlated diffusion" framework whose forward and reverse processes respect the spatial correlations of the target system. Because the resulting algorithm shifts computational weight toward repeated Gibbs sampling, it maps naturally onto probabilistic computers (p-computers), which the paper benchmarks as offering roughly $10^2\times$ (FPGA) to $10^4\times$ (projected sMTJ) better sampling energy efficiency than GPUs.

The work is positioned against prior diffusion applications to Ising-type systems—thermodynamic-fidelity studies of generative models for Ising systems and diffusion reconstruction for diluted Ising models—in which correlations are captured implicitly by the neural network rather than encoded explicitly in the diffusion process. The present framework instead embeds the interaction structure directly into both noising and denoising kernels, with the network reduced to a per-site conditional estimator.

Probabilistic computing primitives

The framework rests on two stochastic units. The p-bit is a binary stochastic unit updated as $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$, where the local field $I_i = \sum_j J_{ij}s_j + h_i$ is determined by programmable couplings and biases; a network of p-bits performs Gibbs sampling of the Ising energy $$E(\mathbf{s}) = -\sum_{i<j}J_{ij}s_is_j - \sum_i h_i s_i$$. Physical realizations include stochastic magnetic tunnel junctions (sMTJs), Zener diodes, and single-photon avalanche diodes. For continuous-state models, the paper introduces g-bits—tunable Gaussian stochastic units implemented natively or as networks of p-bits representing $2^N$ states—with Gibbs update equations $g_i = \mathcal{N}(I_i, \sigma_i^2)$ and $I_i = b_i + \sigma_i\sum_j W_{ij}g_j/\sigma_j$. Since diffusion has Bernoulli and Gaussian variants, p-bits and g-bits map onto them directly.

The hardware benchmarking compares 32-bit random sample generation across platforms:

Platform	Throughput (Gsamples/s)	Power (W)	Efficiency (Gsamples/J)	Improvement
GPU (experimental)	1000	280	3.5	1×
FPGA (experimental)	15000	25	600	~10²×
sMTJ (projected)	31250	2	15625	~10⁴×

Two caveats are stated plainly: these figures measure raw stochastic sample generation, not end-to-end diffusion inference, and a "sample" is a single 32-bit random number rather than a full configuration or Gibbs sweep. The GPU baseline is deliberately favorable—a best-case fabric-local LFSR microbenchmark on an A100 sustaining ~283 W—which makes the FPGA/sMTJ improvement factors conservative estimates.

Diffusion as p-bit dynamics

For the independent case, each site flips with probability $J_{ij}\neq 0$0 per step. Defining $J_{ij}\neq 0$1 and $J_{ij}\neq 0$2, the cumulative forward kernel compresses into $J_{ij}\neq 0$3, and Bayes' rule gives a closed-form one-step reverse posterior. These closed forms exist only because independence factorizes the $J_{ij}\neq 0$4 transition matrix into $J_{ij}\neq 0$5 identical $J_{ij}\neq 0$6 kernels.

The paper's key equivalence is exact: setting each p-bit's local field to self-feedback $J_{ij}\neq 0$7 with $J_{ij}\neq 0$8, the flip probability $J_{ij}\neq 0$9 matches $10^2\times$0 when

$10^2\times$1

so increasing noise corresponds monotonically to increasing temperature. Independent diffusion is thus p-bit dynamics under a temperature sweep, with no approximation introduced by this identification.

Structure-aware correlated diffusion

With couplings restored, the forward process becomes one sequential Gibbs sweep over the interacting Ising Hamiltonian at inverse temperature $10^2\times$2 per diffusion step, using a fixed update order $10^2\times$3 throughout. The one-step reverse posterior retains its Bayesian form,

$10^2\times$4

but loses all closed forms: the prior $10^2\times$5 must be estimated numerically via MCMC. The paper first develops an exact $10^2\times$6 transition-matrix treatment for two coupled spins, then presents the scalable realization: a neural network $10^2\times$7 predicts only independent per-site probabilities $10^2\times$8 that clean spin $10^2\times$9 (trained with binary cross-entropy); a binary estimate $10^4\times$0 is sampled from these probabilities; an ensemble of $10^4\times$1 forward Gibbs chains initialized at $10^4\times$2 produces candidates for $10^4\times$3; and candidates are reweighted by the directly evaluable one-step likelihood $10^4\times$4 before one is sampled. Notably, $10^4\times$5 is not conditioned on the diffusion timestep—all temporal structure enters through the Gibbs dynamics—and is a simple two-hidden-layer MLP, choices made explicitly to isolate the effect of the correlated kernel.

The cost profile is substantial: one generated sample requires approximately $10^4\times$6 forward Gibbs sweeps, i.e., 49,500 sweeps for $10^4\times$7 and $10^4\times$8. This is precisely what motivates dedicated sampling hardware: the correlated kernel's dominant operation is native p-computer behavior.

Results on benchmark systems

2D ferromagnetic Ising model. On a $10^4\times$9 lattice at criticality (10,000 equilibrium MCMC configurations, 80/20 split, $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$0, 100 reverse trajectories with $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$1), the independent kernel produces spatially fragmented intermediate states inconsistent with the correlation length at the corresponding temperature, while the correlated kernel coarsens domains gradually and tracks the Gibbs reference. Quantitatively, the supplementary analysis shows the independent kernel develops a systematic magnetization bias at small $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$2 even when energy appears comparable—the correlated kernel avoids this bias in both $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$3 and $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$4 across the full trajectory. This magnetization discrepancy is a notable finding: energy alone can mask failures of independent diffusion to preserve long-range order.

3D Edwards–Anderson spin glass. On a $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$5 lattice ($$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$6 spins, random $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$7 couplings, single disorder instance at criticality, 20,000 training configurations from chains equilibrated over $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$8 sweeps), the reverse trajectory tracks the forward energy profile, and the energy-per-spin distribution of generated samples agrees with Monte Carlo reference. The most stringent test is the Parisi overlap distribution $$s_i = \operatorname{sgn}[\tanh(\beta I_i) + \operatorname{rand}(-1,1)]$$9 across independently generated pairs, which probes ensemble-level structure rather than single-configuration observables; agreement with the reference indicates the framework reproduces the coexistence of many equilibrium states characteristic of the spin-glass phase near the transition.

Together, these results support the paper's claim that encoding physical knowledge into the diffusion dynamics outperforms relying on the network to learn all correlations implicitly—at least for equilibrium Ising systems where the couplings are fully known.

Limitations and open questions

Several limitations are conceded explicitly. First, the efficiency comparison concerns raw random-number generation, not full inference, so end-to-end speedups remain unquantified. Second, the reverse process relies on importance-weighted selection from a finite chain ensemble ($I_i = \sum_j J_{ij}s_j + h_i$0); more efficient strategies such as parallel tempering or learned proposals could reduce this cost but are left unexplored. Third, the conditional estimator is a plain MLP without timestep conditioning; whether architectures exploiting spatial structure (e.g., graph neural networks) improve accuracy is deferred to future work. Fourth, the demonstrations are restricted to equilibrium Ising systems with fully known couplings at criticality, small lattices, and a single spin-glass disorder instance; extension to partially known or inferred couplings, non-equilibrium dynamics, and continuous states via g-bits remains open. Finally, the exact transition-matrix formulation scales only to very small systems, so all large-$I_i = \sum_j J_{ij}s_j + h_i$1 results depend on the approximate importance-weighted posterior rather than an exact Bayesian inversion.

Conclusion

The paper establishes that independent diffusion is the $I_i = \sum_j J_{ij}s_j + h_i$2 special case of a generalized framework in which the stochastic kernel incorporates known interaction structure, derives the exact temperature-schedule equivalence between flip probability and inverse temperature, and demonstrates that correlated Gibbs kernels yield samples in closer agreement with MCMC references than independent kernels on both ordered and frustrated spin systems—including reproduction of the Parisi overlap distribution. Its broader contribution is methodological: the division of labor between neural-network estimation and structured stochastic sampling is presented as an underexplored design axis in diffusion modeling, one whose computationally intensive points become accessible given hardware whose native operation is fast Gibbs sampling.