Discrete Spatial Diffusion: Intensity-Preserving Diffusion Modeling
Authors: Javier E. Santos
Agnese Marcato
Roman Colman
Nicholas Lubbers
Yen Ting Lin
Venue: NeurIPS 2025 spotlight ()
Keywords: Diffusion models
physics
spatial diffusion
continuous-time discrete-state Markov process
mass conservation
intensity conservation
scientific data
Primary area: machine_learning_for_sciences
TL;DR: This paper presents Discrete Spatial Diffusion, a generative model for discrete-state data that ensures mass conservation, enabling applications in scientific domains like materials science, also demonstrating results on popular image benchmarks.

Abstract: Generative diffusion models have achieved remarkable success in producing high-quality images. However, these models typically operate in continuous intensity spaces, diffusing independently across pixels and color channels. As a result, they are fundamentally ill-suited for applications involving inherently discrete quantities such as particle counts or material units, that are constrained by strict conservation laws like mass conservation, limiting their applicability in scientific workflows. To address this limitation, we propose Discrete Spatial Diffusion (DSD), a framework based on a continuous-time, discrete-state jump stochastic process that operates directly in discrete spatial domains while strictly preserving particle counts in both forward and reverse diffusion processes. By using spatial diffusion to achieve particle conservation, we introduce stochasticity naturally through a discrete formulation. We demonstrate the expressive flexibility of DSD by performing image synthesis, class conditioning, and image inpainting across standard image benchmarks, while exactly conditioning total image intensity. We validate DSD on two challenging scientific applications: porous rock microstructures and lithium-ion battery electrodes, demonstrating its ability to generate structurally realistic samples under strict mass conservation constraints, with quantitative evaluation using state-of-the-art metrics for transport and electrochemical performance.

Discrete Spatial Diffusion: Intensity-Preserving Diffusion Modeling
Javier E. Santos* Agnese Marcato' Roman Colman" Nicholas Lubbers'
Yen Ting Lin'
Abstract
Generative diffusion models have achieved remarkable success in producing high- quality images. However, these models typically operate in continuous intensity spaces, diffusing independently across pixels and color channels. As a result, they are fundamentally ill-suited for applications involving inherently discrete quantities such as particle counts or material units, that are constrained by strict conservation laws like mass conservation, limiting their applicability in scientific workflows. To address this limitation, we propose Discrete Spatial Diffusion (DSD), framework based on continuous-time discrete-state jump stochastic process that operates directly in discrete spatial domains while strictly preserving particle counts in both forward and reverse diffusion processes. By using spatial diffusion achieve particle conservation, we introduce stochasticity naturally through a discrete for- mulation. We demonstrate the expressive flexibility of DSD by performing image synthesis, class conditioning, and image inpainting across standard image bench- marks, while exactly conditioning total image intensity. We validate DSD challenging scientific applications: porous rock microstructures and lithium-ion battery electrodes, demonstrating its ability to generate structurally realistic sam- ples under strict mass conservation constraints, with quantitative evaluation using oftte-of-thhaar metries for transport and electrochemical performance.
Introduction
Diffusion-basse based generative models have emerged as powerful tools for high-quality image generation 36 66 69]
Typically, these models inject noise into the images, then learn to reverse this noise- adding process to recover meaningful structure. In most frameworks, this is based on an Itô Stochastic Differential Equation (SDE) with Gaussian noise. While effective for many vision tasks, these approaches inherently assume continuous pixel intensities, which can cause difficulty when dealing with the discrete nature of many datasets. Tasks beyond vision, such as those in the physical sciences have many applications which require discrete physical quantities, such as particle counts in a simulation, or phases in materials microstructure. Conservation of total quantities is critical for scientific applications, and so generative modeling which can operate under constrained, discrete pixel intensities would enable scientifically grounded and physically consistent synthesis. Such capability might also prove useful within vision tasks, such as inpainting and super-resolution.
Scientific and engineering studies of the natural world using computational techniques often involve discrete variables in space and/or time. On microscopic scales, everyday materials exhibit extremely complex structural patterns which encode the history of their formation, and play large role in how the material functions on a macroscopic level. An important and wide-reaching field of study is materials microstructure, which is used in materials design 355 forensic analysis, hydrology (7).
Los Alamos National Laboratory
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

energy storage 665, and even medicine, such as in studies of bone structure [53]. For example, crystal grain shapes can give rise to complex stress patterns which affect the yield strength of a metal [10] Materials microstructures are frequently represented by small number of discrete phases that describe their underlying chemical structures. In sandstone, for example, the overall arrangement of nanocrystals is highly disordered, and gives rise to complex pore structures, through which subsurface water flows. This microstructure has an enormous influence on the rate of transport of fluids and contaminants [6]. Microstructure of electrodes is also known to have an immense impact on the characteristics of electrochemical devices [58]. Small changes in thermodynamic properties can cause drastic changes in microstructure,, such as in stainless steels [79), requiring the study of microstructure as a function of phase contents. Furthermore, gathering real-world data is often complex and expensive; decades of work have been applied to computational modeling of the generation and consequences of microstructure prior to the widespread popularization of machine learning [72].
While diffusion models have been constructed for microstructure 4 25 37 have not handled exact porosity (intensity) constraints.
46
47
52
80
they
In this work, we introduce Discrete Spatial Diffusion (DSD), a discrete-state Markov chain-based diffusion framework in which the forward process redistributes discrete units of intensity in space (Fig. [ Unlike previous diffusion models, DSD exactly preserves total intensity throughout both the forward and reverse phases, ensuring that global properties- -ucch as mass fractions- are exactly conserved by representing them in terms of conserved particles. We demonstrate that DSD enables scientific applications and can extend the capabilities of conventional image processing tasks, like image generation and inpainting in discrete domains. By directly modeling discrete transitions, DSD enables generative modeling under conservation laws, allowing models that specialize for constrained conditions in scientific applications and beyond.
2
Background
Among the body of literature on generative diffusion models, originating from the pioneering work
of Sohl-Dickstein et al. [66], the most relevant to our work fall into two broad categories: (1) those employing discrete-state Markov chains to introduce noise in the forward process 3
I
39
51 63
64 (7), and (2) those incorporating spatial dynamics into the forward diffusion process 5 39 61].
Generative diffusion modeling based discrete-state Markov chains has become active area of research in recent years. Early work, such as Austin et al. 33. Hoogeboom et al. 39, introduced discrete-state and discrete-time Markov chains as alternative to the Gaussian noise used in conventional diffusion models [36 66 (6). Campbell et al. generalized these formulations to continuous- framework, providing more rigorous theoretical foundation for discrete-state generative diffusion modeling. Santos et al. (63 employed operator algebraic analysis to formally establish the existence of the reverse-time dynamics and derived the stochastic generator for arbitrary discrete-state Markov processes. Similar formulations were independently developed by Sun et al. 770 and Lou et al. [51], with an emphasis on defining and estimating score functions for discrete-state systems. The Markov process operates in intensity space in all the aforementioned diffusion models, treating each pixel as an independent stochastic process (Fig. D a): Gaussian; Fig. Tb): Discrete).
This study focuses on a spatially correlated process for generative modeling for two reasons: (1) for structured images, it is more natural to incorporate spatial correlations into the generative process, and (2) spatially decorrelated noise makes it difficult to preserve total intensity. A spatially correlated approach has been explored for continuous systems. Cold Diffusion 5) introduced a deterministic blurring transformation, where image degradation follows a predefined forward process, and reconstruction is learned as an inverse mapping. However, lacking a probabilistic latent distribution (as in VAES (42), Cold Diffusion is not a true generative model. Inverse Heat Dissipation Model (IHDM, Rissanen al. ( uses the heat equation as corruption model. Since the heat equation is deterministic and reversible (except that the homogeneous solution at t - oc is singular), a naïve inversion would again result in deterministic reconstructions. Uncorrelated Gaussian noise was added to the heat equation to overcome this limitation, relaxing the deterministic process into a probabilistic Itō diffusion. Later, Blurring Diffusion Model (BDM, Hoogeboom and Salimans (38) recognized that IHDM could be recast as Gaussian diffusion model in the spectral domain. BDA extended IHDM and achieved state-of-the generative performance, validating the hypothesis that spatially structured diffusion processes can enhance image generation. Nevertheless, the probabilistic

---

formulation of IHDM and BDM only preserves intensity on average, not exactly per-sample, and their continuous-state nature makes it difficult apply to discrete datasets.
(b)
-11
p1, > s16)
(c)
(x,y) 4,, 82)
Figure 1: Schematic diagrams illustrating how intensity is modeled different diffusion frameworks. (a) Gaussian Diffusion relies on the OrnsteinUUhlenbeck process in continuous intensity space. (b) Prior discrete-state diffusion models apply a discrete-state Markov process to independent pixel intensities. (c) Discrete Spatial Diffusion (this work) relies on Markov jump process of the intensity units over discrete spatial lattice, redistributing particles while exactly conserving total intensity per color channel.
Our goal of generating samples with exactly conditioned total intensity aligns with conditional diffusion modeling. However, existing approaches all rely on some degree of approximation. Song and Ermon [68] proposed simple conditional sampling method by passing class labels into the neural network during training, but this does not guarantee exact enforcement of the condition in generated samples. A more structured approach was introduced by Chung and Ye (15), Chung et al.
[8), which interleaved projection steps with diffusion sampling to enforce linear constraints image generation. However, these projections disrupt the exactness of the forward corruption and reverse inference dynamics [1 63 leading to mismatch between the projected and true data manifolds. To address this, Chung al. [16] eliminated projection steps but instead relaxed deterministic constraints into probabilistic formulation via noisy measurement model. However, this method does not apply to deterministic constraints, as becomes singular in the limit of zero measurement noise. An alternative approach leverages Bayes' theorem for a posteriori conditional sampling, that is, p(S{C) x (CSSs,C)S) where "S" stands for samples and "C" for condition(s). Because pSS) is given by a trained unconditional diffusion model, conditioning can be performed one has pCC) S), which however intractable for arbitrary data distributions Existing methods approximate this term crudely by training separate classifier in Song et al. [69]. by Gaussian approximation with moment moment-matchine matching as in [3 27]. None of these methods guarantees that the generated samples exactly conditioned.
The principal contribution of our work that it provides a new capability for diffusion models to preserve intensity exactly in fully discrete-state context. The approach is based entirely on how the diffusion process is built, and how the model is trained; it is readily usable with existing diffusion model neural network (NN) architectures. To our knowledge, this is the first diffusion model to incorporate spatially correlated noise, achieved through stochastic jump process that allows units of intensity perform random walk. This fact also demonstrates that more complex noise processes can themselves be tractable. We furthermore demonstrate that such a model is effective for conventional image synthesis tasks. The relevance and power of the approach is then demonstrated through application to scientific data in the field of materials microstructure, where the ability to generate complex data-driven images with constrained total intensity is key to physically meaningful generative modeling and scientific reliability. To clarify DSD's advantage in scientific settings, we explicitly compare traditional frameworks in Appendix E
II is challenging because the constraint is imposed on the final samples at the end of the inference, but the conditioning "S" samples generated during the inference.

---

3 Methods
3.1 Corruption Process
In this manuscript, we adopt the language of image processing and consider 2-dimensional images, although the application context and spatial dimensionality of the data are not constrained by the mathematical framework provided here. We treat a digital image with discretized intensity
values L g.e E Zz0 at pixel (z, y) E {I,..., W} x {1,... (l...... ..... } in color channel c {l...., C}. Within the DSD framework, the intensity values in the image are treated as a spatially organized collection of particles, with one particle for each intensity unit. Below, we will interchangeably use "particles" and "intensity units" to denote these fundamentally discrete units. Specifically, L2..c n implies n particles of type c location (x.y), and the total number of particles of the system is Σd1 Σn= Σe L.s..- In the forward stochastic process with the time parameter t, each of the particles the system independently performs a continuous- time and discrete-state staree random walk:
(I, c)) 4 (a: + Mz + uyy c), === (Uz, My) E {(1, ( (-1,0),
(1a)
(lb)
where r is the transition rate of the particle jumping to one of their nearest neighbors, and v is a set of four directions the particles can hop to their nearest neighbors. A schematic diagram is shown in Fig. Dp. Note that the particles perform jumps in the (r, !) space at random times, but do not change their color coordinate c. Because of this, the forward process conserves the total number of particles Σ" Σyy La.y, in each color channel independently. We impose either no-flux boundaries, such that the transition rate for jumps out of the image domain are zero, or periodic boundaries, so that jump to x= = W + becomes a jump to x = 0, vice-versa, and analogously for y.
We refer to the spatial hopping process (Fq. I1) as the Discrete Spatial Diffusion (DSD), noting the "discreteness" refers to both the discretized intensity units and the discreteness of the spatial lattice W
{I....,If} _If where the particles allowed to reside. DSD, well as similar discrete-state random walks, have been extensively studied in non-equilibrium statistical physics and stochastic processes (Van Kampen 775), Gardiner 128), Giuggioli 33 and references therein). The evolution of the probability distribution of the single random walk in the continuum space limit, under the appropriate scaling of the transition rate (26), converges to the Fokker-Plann
Equation (FPE, Risken [60]. Van Kampen 175)), which is mathematically identical to the heat equation. Because of the duality between the probabilistic FPE and the deterministic heat equation [44], DSD can be considered as a microscopic description of the macroscopic heat dissipation that inspires IHDM and BDM. Notably, the correlated noise is built in DSD, in contrast to the heuristic addition of uncorrelated Gaussian noise in IHDM and BDM. Fig. Fig-Dillin iillustrates the application of DSD to a sample image. Due to the stochasticity of the random jumps, the limiting behavior (t + 20) of this process is a random configuration with no discernible structure or similarity to the original spatial organization aside from the conserved global particle counts in each color channel.
We use (X,, Yr, C₂) to denote the random process in (I, y, c) space, and (mo, yo, oo) are the initial condition of a specific particle. We use Ie to denote the randomly corrupted image at the time t, where [4el..g.c is the total number of particles at (z, y) in color channel c. The process can be represented in these two dual representations: with (X,, Yr, C2) the process is formulated in the frame of a moving particle (the Lagrangian frame), and with 44 the process formulated as histogram in space-time (the Eulerian frame). Below, we will use these two representations interchangeably.
The forward solution and the transition probabilities (z., y, c/ yo, cu) P(X, = x,Y,
y. C, = Xn = To. Yo = co. C. = c. can be computed by integrating the Master Equation 75 777 This corresponds to exponentiating the Markov transition matrix of the process defined (28
in Eq. I While the matrix exponential required numerically for no-flux boundaries is expensive, the solution can be stored and reused to corrupt images and to compute the reverse- -rrasse-tannitionn rates (see Sec. (.3) for learning. When periodic boundary conditions are imposed, the transition matrix is diagonal the discrete Fourier space, facilitating the efficient computation of pe(-I-) (see Sec. Cl.
3.2 Designing Noise Schedules by Structural Similarity Index Metric (SSIM)
Since the corruption process is time-homogeneous ( the noise applied to each partiele remains constant over time. However, it has been shown that inhomogeneous noise schedules can facilitate

---

learning [55]. We use the formulation of a recent study [62] identified the unique correspondence Caastein-UUlenbert process [73 and We follow the same philosophy as Santos and Lin [62 to construct a sequence of observation times 4g < l2 < ... < lT 1, at which we will generate random samples for learning. Here, T denotes the total number of discrete time steps used to generate corrupted sample images during training.
between non-uniform observation times in a homogeneous
noise schedule in conventional diffusion models [36, 67
We adopt a heuristic approach to construct the discrete times. The idea is to use a metric to quantify how much the "quality" of the images has been degraded up time f, and we aim to design tę such that the metric degrades from & = O to Similarity Index Metric (SSIM, Wang et al. (76]) between the corrupted image and the original image. We generalize a generic monotonic relation between to & proposed by Santos et al. 633
T as evenly as possible. We chose the Structural
(k: 1) 32) (T' 2 (e
(c Tate)
(2)
where ((p) := logp/(1 - p) is the logit function, TI and To are parameters used to construct the observation times. Note that Ir 1 in the above parametrization. Specifically, we tune 72 and the unit transition rate r process ( using subset of training samples, aiming to cover an degradation of the SSIM throughout observation times. We use Ti 7.5 and
2.5, and r 120-160, chosen to equalize SSIM degradation over timesteps. We remark that the choice of the functional form in Eq. (2) is picked empirically; treat Eq. (2) as versatile monotonic fitting function, whose corresponding SSIM degradation is more symmetric than polynomial and cosine
schedules [55] for the DSD process (see Appendix Fig. 6
(b)
Gauneian Noioe
73%
37%
00
0.10%
TEDA
0.00%
DSD
Figure 2: (a) The forward processes for Gaussian Diffusion [36 Inverse Heat Dissipation Model 61 and Discrete Spatial Diffusion (ours) applied on an image, sampled at discrete times. (b) Percentage change in intensity relative to the original image under the forward process.
3.3 Reverse-time process
Following the process formalism from 11 opposite time and whose joint probability distribution is identical to that of the forward process (T) Specifically, the reverse-time process corresponding process is:
63
there exists reverse-time process that evolves in
Pe(x- Peta.y.c M0.
(x, y, c)
(a + De, 4 1 c),
(3)
where the admissible reverse-time transitions & (Vas By) :.):: == ((-1,0) (1,0),
1)}
are the reversed direction of the forward jumps (D, = -vz: Dy = -2y). The framework ensures the same boundary condition to be imposed (no-flux or periodic, according to the forward process). We note that the reverse-timg process, and therefore the generated images, also conserve the total particle number per color channel.
note that the reverse transition rate depends on both the initial condition (xo- Voo cg) particle and the forward solution Pt (z, y, clæo, yo, co), V(I, y, c). This is analogous to conventional diffusion
5

---

models, where either the reverse-time drift [36 66 or the score function 66 formally depend on the initial sample and the solution of the forward process. However, during the inference, the initial particle configuration is not known, and as such, we train an NN to learn the reverse transition rates using samples I generated from the forward process ① at O. Additionally, the particles are indistinguishable, but the rate prescribed in Eq. is per-particle, raising the question: what is the appropriate per-pixel reverse transition rate that the NN ought to model? This question can be answered by performing the survival analysis of the many many-particle system in light of the independence of particle motion; see Appendix A for a derivation. The analysis shows that the reverse transition rate of the first jump of = [ill, particles is simply the sum of the instantancous transition rates:
pelar Pz. y B C
2
[1] [F] 11
MO
Pr(r, y, c c 90
b1
Eq. prescribes the rate for the first of the particles (which g..)) to jump to one of its neighboring pixels. Intuitively, this can also be derived by combining the first-reaction method [30] and inhomogeneous Poisson process (e.g., see Corbella et al. (20)). Eq. Jalso prescribes the rate that the NN will model. This rate is implicitly time-dependent through the dependence on the forward solution Pr, similar to standard continuous-time diffusion models.
3.4 Loss functions
Our goal is to provide the corrupted images 7, at a sampled time / >> 0 to a neural network (NN) and to train it to predict the reverse-time transition rates F E RA 4x of Eq. 4 The first dimension (of size four) represents rates for the four nearest-neighhor transitions. We denote the NN modeled rates as ,NN (4, t) with noisified image I at time t.
We use two approaches formulate the loss function. The first and more commonly used approach adopts a metric to heuristically match the NN prediction and the ground truth, as in DDPM (36), ccore-matching [69], and flow-matching [49]. We extend these schemes to rate-matching, where we minimize a chosen norm of the difference between the predicted and true rates FNN and F. For example, for using L1, a loss L:
Crate-matching Ei, [mean(] | NN - -FD]
(5)
Here, & E (1.... } is uniformly sampled, Le, is drawn from the random process ( at the sampled times, F = Fp.m.y.e (4, f/4)) is the true reverse-time transition rate ( NN = -NN NN-predicted reverse-time transition rate, and the mean is over all the indices (17, x, y, c). The second and more principled approach is through minimization of the negative likeliilddood L of the NN-induced process to predict the analytical reverse-ting process
f.x.y.cc () is the
63 66):
Likelibood log -EI, [/.0
Σ (FNN
F Flog NNN) dt
(6)
Because only observe the process at discrete times prescribed in Eq. approximate the continuous satinuous-time integration above by
log L = Eieg (Ix 12-1)
Σ (F - Flg ( NNN
(7)
where we again take expectation over randomly sampled ta and Ig. In this study, we experimented with both loss functions and observed no noticeable difference, giving evidence that the DSD forward process ( is not sensitive to the choice of the loss function. This stands in contrast to Gaussian diffusion models, where training loss choice significantly impacts performance Ho et al. [36). which used the heuristic approach improve over Sohl-Dickstein et [66]. which adopted the second approach. We focus on learning the transition rates of the reverse-time dynamics, which is distinct from the ratio-matching approach [51] 70 which focuses on learning the probability distribution pr()), although a similar formulation ("implicit score entropy") proposed by [5] can be regarded as the process likelihood 7 first proposed in Santos et Algorithm describes the DSD training pseudocode.

---

3.5 Sampling with an Adaptive Time Stepping using the Courant-Friedrichs-Leee condition
Once trained, the neural network will predict reverse-time rates 4 from the configuration of system Jt, at time 1. 2 0. Because the reverse rates are time-dependent, one could generate the exact sample paths of the inhomogeneous Poisson process by integrating the survival function of the first reaction on each pixel in each color channel (sce e.g., algorithms reported in Corbella et al. (20). However, this approach is not computationally efficient, so we resort to r-leaping [3I, an integrator that has been adopted by essentially all continuous-tiie time and discrete-state diffusson models I 59 63 78), analogous to the Euler's 's method for ordinary or partial differential equations and Itō SDEs. The central idea of r-leaping is to approximate the reverse-time transition rates f as fixed constant in a small enough window (t - T, t), assuming the time-dependent rates change slowly in the period, a condition often termed as the "leap condition"
for
[12] 31 With this assumption, the original r-leaping algorithm by Gillespie [31] generates Poisson random numbers to update the system's 's discrete states. However, this approach could sometimes lead
negative population of particles, which cannot happen in the process, due to violations of the leap condition. Mitigation strategies exist [12 13 32), however, some of them are limited to small reaction networks and not suitable for the DSD sampling task, which involves a very large number (4 x H x W x C) of transition rates to estimate.
As such, we propose more efficient (but arguably cruder) approach to select the step size T adaptively, (CFL) condition [21] to conservatively determine the adaptive step size T. Specifically, since the jump scale is fixed at the pixel length scale, the timescale fully determines the CFL condition. The idea is to choose a 7 such that the CFL number is fixed throughout the inference] To achieve this, we compute the reverse-time transition rates Fp for each pixel in each channel, noting that the probability of particle in that channel will jump to one of its neighboring locations is FGT. Then, we determine by fixing the largest probability across all the pixels and color channels a constant. Algorithm describes the DSD inference pseudocode.
Our idea is to combine the binomial r-leaping [14 7 and the
4 Computational Experiments
We employ the Noise Conditional Score Network (NCSN++) [1]69] with two modifications: the final convolutional layer outputs 4 times the number of input channels (c.g.. 3 for RGB) to represent four directions (up, down, left, right), and a SoftPlus activation function ensure non-negativity in the predicted rates. The hyperparameters are listed in Appendix H
4.1 Image synthesis benchmarks
While Discrete Spatial Diffusion (DSD) was developed to enable generative modeling under a strict intensity constraint, first evaluate its performance conventional image synthesis benchmarks. Specifically, we trained DSD models on MNIST [45), CIFAR- 10 [43], and CelebA [50), suggesting that DSD can approximate key statistical and structural properties of image distributions, despite operating under discrete-state and global conservation constraints.
We begin with MNIST, where intensity naturally corresponds to stroke thickness and digit area. In Fig. a), we show inpainting results with a fixed mask with no-flux boundary conditions. During inference, particles inside the masked region rearrange to complete the digit structure, while the surrounding region remains fixed. Holding the visible structure constant, we varied the number of particles allowed to move, revealing that the total intensity governs which digit is most likely to emerge. This highlights the ability of DSD to incorporate hard constraints in downstream tasks such as inpainting. Additionally, we trained a conditional DSD model that employed standard conditionnn [69]. In BigureB]) we illustrate the class-conditioned generated images with different total numbers of particles, varying from low, typical, to high total intensities. While these
EEvn though CFL condition is more commonly used in PDE integrators, the concept can be applied for our stochastic system. Suppose the reverse-time rate F. On average, the particle would move at timescale 1/r one of its neighbors, traveling Az. Then, the velocity e fåz. The CFL condition is then cAt/A Az where in our scheme At is the r; thus, the classical CFL convergence condition translates to the obvious bound of transition probability er < 1. This motivates us to ensure conservative estimation of r, but enforcing a small rr reduce the error.
7

---

do exhibit some artifacts, DSD notably learns the spatial structure of the digits and generates "Bolder" or "Lighter" digits without saturating the upper bound of the intensity (i.e. 255 for uint8). This would not have been precisely realizable using conventional diffusion models. Comparisons with other conditioning approaches can be found in Appendix
We also evaluated DSD on RGB datasets to explore scalability and generality. In Fig. c..d), we show unconditional generations from models trained CIFAR-10 and CelebA. Despite the discrete- state space and peeensttyprreervingg dynamics, the model captures complex semantic structures from lighting, and textures to animals, vehicles, and human facial features. Appendix D includes large grids of generated samples, Fréchet Inception Distance (FID), and spatial FID (sFID) metrics, sampling ablations, and post-process processing strategies for improving sample smoothness.
9
7
3
P
8
3 S
/ 2
3
4
5
9
8
9
/
2
Y
7
Revense-inee linee
Class
Figure 3: (a) Inpainting realizations on MNIST; 15% difference of conditioning intensity between consecutive rows. (b) Conditioned MNIST generations across different intensities and classes. (c) Unconditional CIFAR-10 generations. (d) Unconditional CelebA generations.
4.2 Subsurface rock microstructures
The microstructure of subsurface rocks governs wide range of physical processes, including fluid transport, electrical resistivity, and mechanical deformation [7]. This originates from connected pores on the nano- and micro-scale, which vary size, structure, and coordination degree across rock types. High-resolution 3D imaging via X-ray microtomography enables detailed pore-scale reconstructions, these scans are expensive and limited to sample sizes on the order millimeters to centimeters [99] While direct imaging of rock microstructure is costly, measuring porosity (defined here as average intensity over the image) across large formations inexpensive and can be performed without specialized equipment [48] 57). This enables large-scale field measurements of porosity, even when high-resolution microstructural data unavailable.
0.58
Berea Bandatione
Bavonnirrss Carbonatt
Massangia Limestane
Figure 4: Distributions of porosity and generated samples for three rock classes. Each plot shows the porosity distribution across the dataset (shaded curve) and corresponding DSD DSD-generated mi- crostructures selected quantiles. (Left) Berea Sandstone, (center) Savonnières Carbonate, and (right) Massangis Limestone. Samples illustrate the model's ability to generate structurally diverse images conditioned on exact global intensity. Quantitative metrics are presented Appendix F.2
To overcome this limitation, synthetic models are frequently used to generate representative pore structures for computational physics studies (56). However, conventional reconstruction techniques impose strong geometric assumptions that fail capture the heterogeneity observed in real rocks like Berea Sandstone, Savonnières and Massangis Carbonates. We trained DSD models using these rock samples, which represent broad spectrum of pore structures (including granular, fossiliferous, and

---

dissolution-dee driven features) across two lithologies. A description of the training datasets is provided in Appendix F.T Fig. presents representative outputs from our models trained on 256x256 binary images. The generated samples accurately reproduce the key morphological and statistical features of the original datasets, including two-point spatial correlations and pore size distributions which are primary microstructural descriptors that govern permeability, tortuosity, and diffusion behavior in porous media (Fig. [9]. Given that DSD allows for precise control over total porosity, one can generate synthetic microstructures that match the porosity measured in the field, enabling the reconstruction of representative pore-scale samples even in the absence of direct imaging.
4.3 Lithium-ion electrodes
Electrodes in lithium-ion batteries are porous materials with complex microstructure that governs key properties like ion transport and electrochemical performance. Nickel- manganese-cobalt cathodes, among the most common, are composed of three phases: the active material driving the electrochemical reaction, the carbon binder ensuring elec- trical conductivity and mechanical stability, and the pore space filled with electrolytes. The active material expensive, creating strong economic incentive to understand how its volume fraction and distribution influence electrode behavior. While tomographies are needed for studying microstructures and en- abling computational modeling, acquiring di- verse datasets is challenging [22]. To over- come this, researchers often rely on computa- tional methods to generate synthetic microstruc- tures [24] While GANs have been explored for this purpose, they do not control phase vol- ume ratio (29). We trained a DSD model on tomography data [74), where two color chan- nels were used to represent the carbon binder and active materials. The results, shown in Fig. 5, demonstrate DSD enables precise tun- ing of phase volume fractions. We then com- puted key morphological metrics characterizing the electrode structure [41] demonstrating the strong generative capability the DSD method, which generalizes across the training set more effectively than existing approaches in the liter- ature. For more details on datasets and recon- struction metrics applied to these samples, see Appendix G
Hlectritte
Actiow matirr
Carbor hinder
Active material volume fraction
57%
24%
8%
12%
16%
Carbon binder valume fraction
Figure 5: Generated cathode microstructures with varying phase volume fractions. The carbon binder domain appears in black, active material particles in gray, and electrolyte-file filled pore space in white.
5 Limitations
The computational cost of forward sampling during training and reverse-time sampling during inference in DSD scales linearly with the total intensity of the image. While this makes DSD highly efficient for low-bit-depth binary datasets, it may become less efficient than other techniques for higher-resolution images or datasets with higher intensity saturation, such as standard uint8 images. While the FLOP cost of the forward noising process is higher than in conventional Gaussian diffusion models, this is handled entirely by CPU workers during the data loading phase, and so the added cost effectively hidden completely behind the data pipeline and does not impact GPU throughput; DSD models were negligibly slower than their Gaussian counterparts, and models were trained in this work using one A100 GPU in < 750 hours per model (RGB images) and < 50 hours per model (microstructures. MNIST), Inference is slower (approximately 2x per), due to the need to sample individual particle transitions, but affordable across datasets of practical interest. Despite this cost scaling, many of the largest available microstructural datasets in subsurface modeling are on the order of 1000 x 1000 resolution. We trained DSD on such large-scale scientific data, demonstrating that the method remains computationally feasible in this regime, also suggesting the feasibility of extending

---

to 3D by carrying out implementation changes. Representative results from high-resolun regh-resolttion generation are presented in Appendix F3 Additionally, enforcing strict intensity conservation requires custom forward process code (Eq. (T) and novel sampling scheme, deviating from conventional Gaussian diffusion models. This introduces steeper learning curve for practitioners accustomed to standard diffusion approaches. However, we argue that these trade-offs are necessary to achieve exact intensity constraints, which is not possible with existing methods.
Conclusion
We introduced Discrete Spatial Diffusion (DSD)) fully discrete, antensity-preserviing generative model approach for images and scientific data. The foundation is the use of discrete-state, continuous time statistical processes incorporating jump dynamics, rather than SDEs, and in particular is the first discrete diffusion model to explore spatially correlated corruption. DSD demonstrates competitive quality on standard benchmarks while enabling exact global constraints on particle count (thereby conserving image intensity or mass fractions under various applications) that are critical in many scientific applications. By preserving these constraints in both forward and reverse processes, DSD enables exactly constrained data generation, which we explored on image synthesis and domain- specific datasets. also demonstrates that more complex statistical processes (in this case, random walks) can be used for diffusion modeling, opening the door for further models to exploit structure in their dynamics such as conservation laws and symmetries.
Code availability
All code supporting this study is Discretesspaiialiffffussoo
available
https://github. com/ lanl/
Acknowledgments and Disclosure of Funding
This work was performed Los Alamos National Laboratory (LANL), operated by Triad National Security, LLC, for the National Nuclear Security Administration of the U.S. Department of Energy (Contract No. 8923321 (9233218CNA000000 The authors are supported by the Laboratory Directed Re- search & Development project "Diffusion Modeling with Physical Constraints for Scientific Data (20240074ER)".
References
[1] Alexandre Adam. score_models: repository to store state of the score model architectures. https://githuut http////eeeaereeeda//eeeeeeeeeee
2023. Accessed: 2025-05-13.
[2] Brian D.O. Anderson. Reverse-time diffusion equation models. Stochastic Processes and their Applicaatoons, 12(3):31--326, 1982. ISSN 0304-4149. doi: https://doi.rg//0..100 0304-4149(82)900555-5
[3] Jacob Austin, Daniel Johnson, Jonathan Ho, Daniel Tarlow, and Rianne van den Berg. Structured denoising diffusion models in discrete state-spaces, In Mare' Aurelio Ranzato, Alina Beygelzimer, Yann Dauphin, Percy Liang, and Jennifer Wortman Vaughan, editors, Advances in Neural Information Processing Systems 34: Annual Conference on Neural Information Processing Systems 2021, NeurIPS 2021, December 6-14, 2021, virtual, pages 17981-17993, 2021.
[4] Erfan Azqadan, Hamid Jahed, and Arash Arami. Predictive microstructure image gener- ation using denoising diffusion probabilistic models. Acta Materialia, 261: 261:119406, 2023. ISSN 1359-6454. doi: oat.eccceee...ccccccooaarrrrrrea
19446
19406. URL https::tp///. //www
sciencedirect.comm
[5] Arpit Bansal, Eitan Borgnia, Hong-Min Chu, Jie S. Li, Hamid Kazemi, Furong Huang, Micah Goldblum, Jonas Geiping, and Tom Goldstein. Cold diffusion: Inverting arbitrary image transforms without noise, 2022.
10

---

[6] Jacob Bear. Dynamics of fluids in porous media. Courier Corporation, 2013.
[7] Martin J Blunt, Branko Bijeljie, Hu Dong, Oussama Gharbi, Stefan Iglauer, Peyman Mostaghimi, Adriana Paluszny, and Christopher Pentland. Pore-scale imaging and modelling, Advances in Water resources, 51:197-216, 2013.
[8] MA Boone. 3D mapping of water in oolithic limestone at atmospheric and vacuum saturation using x-ray micro-ct differential imaging. Materials Characterization, 97:150-166 2014.
[9] Tom Bultreys, Jeroen Van Stappen, Tim De Kock, Wesley De Boever, Marijn A Boone, Luc Van Hoorebeke, and Veerle Cnudde. Investigating the relative permeability behavior of microporosity- carbonates and tight sandstones with multiscale pore network models. Journal of Geophysical Research: Solid Earth, 121(11):7929-777 2016.
[10] Marion Calcagnotto, Yoshitaka Adachi, Dirk Ponge, and Dierk Raabe. Deformation and fracture mechanisms fine-and -ltrafieegraiied ferrite/martensite dual-phase steels and the effect of aging. Acta Materialia, 59(2):658-670, 2011.
[11] Andrew Camphell, Joe Benton, Valentin De Bortoli, Thomas Rainforth, George Deligiannidis, and Arnaud Doucet. A continuous time framework for discrete denoising models. S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh, editors, Advances in Neural Information Processing Systems, volume 35, pages 28266-28279. Curran Associates, Inc., 2022.
[12] Yang Cao, Daniel T. Gillespie, and Linda R. Petzold. Avoiding negative populations in explicit Poisson tau-leaping. The Journal of Chemical Physics, 123(5):054004, August 2005. ISSN 0021-9606. doi: 10.1063/1.1992473
[13] Yang Cao, Daniel T. Gillespie, and Linda R. Petzold, Efficient step size selection for the tau-leaping simulation method. The Journal of Chemical Physics, 124(4):044109, January 2006. ISSN 0021-9606. doi: 10. N..10631.2219468.
[14] Abhijit Chatterjee, Dionisios G. Vlachos, and Markos A. Katsoulakis. Binomial distribution based T-leap accelerated stochastic simulation. The Journal of Chemical Physics, 122(2), January 2005. ISSN 0021-9606. doi: 00.106//1.183335.
[15] Hyungjin Chung and Jong Chul Ye. Score-based diffusion models for accelerated MRI. Medical Image Analysis, 80:102479, August 2022. ISSN 1361-8415. doi: 10. (0166..media.2022.110
[16] Hyungjin Chung. Jeongsol Kim, Michael Thompson McCann, Marc Louis Klasky, and Jong Chul Ye. Diffusion Posterior Sampling for General Noisy Inverse Problems. In The Eleventh International Conference Learning Representations, September 2022.
[17] Hyungjin Chung. Byeongsu Sim, Dohoon Ryu, and Jong Chul Ye. Improving Diffusson Models for Inverse Problems using Manifold Constraints. In Advances in Neural Information Processing Systems, October 2022.
[18] Hyungjin Chung, Bycongsu Sim, and Jong Chul Ye. Come-Closer Sasser--Dffuee-Faae Accelerating Conditional Diffusion Models for Inverse Problems through Stochastic Contraction. In 2022 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), pages 12403- 12412, June 2022. doi: 10.1 10.110//CPRR226888 202.00.00....00009..
[19] Veerle Cnudde and Matthieu Nicolaas Boone. High-resolution x-ray computed tomography in geosciences: A review of the current technology and applications. Earth-Science Reviews, 123: 1-17, 2013.
[20] Alice Corbella, Simon E. F. Spencer, and Gareth O. Roberts. Automatic Zig-Zag sampling in practice. Statistics and Computing, 32(6):107, November 2022. ISSN 1573-1375. doi: 122077//112220022-1112222
[21] R. Courant, K. Friedrichs, and H. Lewy. Über die partiellen Differenzenllscchuneen der mathematischen Physik. Mathematische Annalen, 100(1):32-74, December 1928. ISSN 1432-1807 doi: 10.1007/BBO01448839,

---

[22] Zhe Deng, Xing Lin, Zhenyu Huang, Jintao Meng, Yun Zhong, Guangting Ma, Yu Zhou, Yue Shen, Han Ding, and Yunhui Huang. Recent progress on advanced imaging techniques for lithium-ion batteries. Advanced Energy Materials, 11(2):2000806, 2021.
[23] Pan Du, Meet Hemant Parikh, Xiantao Fan, Xin- Ying Liu, and Xun Wang. CoNFiLD: Conditional Neural Field Latent Diffusion Model Generating Spatiotemporal Turbulence, March 2024.
[24] Mare Duquesnoy, Chaoyue Liu, Diana Zapata Dominguez, Vishank Kumar, Elixabete Ayerbe, and Alejandro A Franco. Machine learning-assisted assisted multi-objective optimization of battery manufacturing from synthetic data generated by physics- -ased simulations. Energy Storage Materials, 56:50-61, 2023.
[25] Christian Düreth, Paul Seibert, Dennis Rücker, Stephanie Handford, Markus Kästner, and Maik Gude. Conditional iiffusiio-basee microstructure reconstruction. Materials Today Communica- tions, 35:105608, 2023. ISSN 2352-4928. doi: https://doi.org/11..00 http://1..1111...../.. /tt1//..mccommm222333333333 1056088 URL httph:tt://// www.sccdiieiccc..cc..mc.. htwpeeddrirrtttccceeeeoeooeeeeee com wtt.//occcccccea http:////...o.ooerrr..aaaa//////
[26] Albert Einstein. Über die von der molekularkinetischen Theorie der Wärme geforderte Bewe- gung von ruhenden Flüssigkeeten suspendierten Teilchen. Annalen der Physik, vol. 4, t. 17, 1905.
[27] Marc Anton Finzi, Anudhyan Boral, Andrew Gordon Wilson, Fei Sha, and Leonardo Zepeda- Nunez. User-defined Event Sampling and Uncertainty Quantification in Diffusion Models for Physical Dynamical Systems. In Proceedings of the 40th International Conference on Machine Learning, pages 10136-10152 PMLR, July 2023.
[28] Crispin W. Gardiner. Stochastic Methods: A Handbook for the Natural and Social Sciences. Number 13 in Springer Series in Synergetics. Springer, Berlin Heidelberg, 4th ed edition, 2009. ISBN 978-3-642-08962-6 978-3-540-70711-77
[29] Andrea Gayon-Lombardo, Lukas Mosser, Nigel P Brandon, and Samuel J Cooper. Pores for thought: generative adversarial networks for stochastic reconstruction 3d multi-phase electrode microstructures with periodic boundaries. npj Computational Materials, 6(1):82, 2020.
[30] Daniel T Gillespie. A general method for numerically simulating the stochastic time evolution of coupled chemical reactions. Journal of Computational Physics, 22(4):403-434, December 1976. ISSN 0021-9991. doi: 10. 10010166002-199911199990444
[31] Daniel T. Gillespie. Approximate accelerated stochastic simulation of chemically reacting systems. The Journal of Chemical Physics, 115(4):116--1733, 2001. doi: https://doi.org///
[32] Daniel T. Gillespic and Linda R. Petzold. Improved leap-size selection for accelerated stochastic simulation. The Journal of Chemical Physics, 1 119(1)):222988))) 8234, October 2003. ISSN 0021- 9606. doi: 10.1063/1.1 1...1666...11254.
[33] Luca Giuggioli. Exact Spatiotemporal Dynamics of Confined Lattice Random Walks in Arbitrary Dimensions: A Century after Smoluchowski and PVllya. Physical Review X, 10(2):021045, May 2020. doi: 10. 10..P 10.021045.
[34] Jeff Gostick, Zohaib A Khan, Thomas G Tranter, Matthew DR Kok, Mehrez Agnaou, Mohammadamin Sadeghi, and Rhodri Jervis. Porespy: A python toolkit for quantitative analysis of porous media images. Journal of Open Source Software, 4(37): 2019.
[35] Grace X Gu, Chun Thun-TTh&on Chen, Deon J Richmond, and Markus J Buehler. Bioinspired hierarchical composite design using machine learning: simulation, additive manufacturing, and experiment. Materials horizons., 5(5):939-945, 2018. ISSN 2051-6355 (clectronic).
[36] Jonathan Ho, Ajay Jain, and Pieter Abbeel. Denoising diffusion probabilistic models. In
H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 6840-6851. Curran Associates, Inc., 2020.
12

---

[37] Nathan Hoffmaan, Cashen Diniz, Dehao Liu, Theron Rodgers, Anh Tran, and Mark Fuge. Grainpaint: A multi-scale clffusion-based generative model for microstructure reconstruction of large-scale objects. Acta Materialia, 288:120784 2025. ISSN 1359-6454. doi: https://doi.
00g10000006/j..ctamtt2022..3 120784. URL att.cleepp11/13336444220000
http://www..ciinnneerrretttttt wtw.errtttttcceeeeeeoemmeeeeeeee cww/coo//c///o com/sciences
[38] Emiel Hoogeboom and Tim Salimans. Blurring Diffusion Models. In The Eleventh International Conference on Learning Representations, September 2022.
[39] Emiel Hoogeboom, Didrik Nielsen, Priyank Jaini, Patrick Forré, and Max Welling, Argmax flows and multinomial diffusion: Learning categorical distributions. M. Ranzato, A. Beygelz- imer, Y. Dauphin, P.S. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, volume 34, pages 12454-12465. Curran Associates, Inc., 2021.
[40] Michael Janner, Yilun Du, Joshua B. Tenenbaum, and Sergey Levine. Planning with diffusion for flexible behavior synthesis, 2022. lttpps:: org www.///.....0e.e 09991
[41] Steve Kench, Isaac Squires, and Samuel Cooper. Taufactor 2: A gpu accelerated python tool for microstructural analysis. Journal of Open Source Software, 8(88):5358, 2023.
[42] Diederik Kingma and Max Welling. Auto-Encoding Variational Bayes. In 2nd International Conference on Learning Representations, ICLR 2014, Banff, AB, Canada, April 14-16, 2014, Conference Track Proceedings, 2014.
[43] Alex Krizhevsky, Geoffrey Hinton, et al. Learning multiple layers features from tiny images. 2009.
[44] Gregory Lawler. Random Walk and the Heat Equation. Number volume 55 in Student Mathematical Library. American Mathematical Society, Providence, R.I, 2010. ISBN 978-0- 8218-4829-6.
[45] Yann LeCun, Corinna Cortes, and Christopher J.C. Burges. Mnist handwritten digit database, 2010. URL battp://yann lecun. www.co/cc///eee wtteeedbbmmmiieee
[46] Kang- Lee and Gun Jin Yun. Microstructure reconstruction using diffusion-based genera- tive models. Mechanics of Advanced Materials and Structures, 31(18):4443-446 2024.
[47] Kang-Hyun Lee, Hyoung Jun Lim, and Gun Jin Yun. A data-driven framework for designing microstructure of multifunctional composites with deep-Icarned diffusion- generative models. Engineering Applicattions Artificial Intelligence. 129; 129:107590, 2024. ISSN 0976-1976. doi: https:///oi.oooo org/ http://..0010/...appppia.0000000 107590. URL cattps:/// http://wwwwww.. www.sccinneeireetttt
com/science
[48] RW Leonard. Simplified porosity measurements. The Journal of the Acoustical Society of America, 20(1):39-41, 1948.
[49] Yaron Lipman, Ricky T. Q. Chen, Heli Matching for Generative Modeling. In The Eleventh International Conference on Learning Representations, September 2022.
Maximilian Nickel, and Matthew Le. Flow
[50] Ziwei Liu, Ping Luo, Xiaogang Wang, and Xiaoou Tang. Deep learning face attributes in the wild. Proceedings of the IEEE International Conference on Computer Vision (ICCV), pages 3730-3738, 2015.
[51] Aaron Lou, Chenlin Meng, and Stefano Ermon. Discrete diffusion modeling by estimating the ratios of the data distribution. In Ruslan Salakhutdinov, Zico Kolter, Katherine Heller, Adrian Weller, Nuria Oliver, Jonathan Scarlett, and Felix Berkenkamp, editors, Proceedings of the 41st International Conference Machine Learning, volume 235 of Proceedings of Machine Learning Research, pages 1289942 PMIR, 21-27 Jul 2024. URI. https //proceedings.co mlr. wwwpp/eeeevvae gpees/v235/lou244 html
[52] Xianrui Lyu and Xiaodan Ren. Microstructure reconstruction of 2d/3d random materials via diffusion-basse based deep generative models. Scientific Reports, 14(1):5041, 2024.
13

---

[53] Carolina Montoya, Yu Du, Anthony L. Gianforcaro, Santiago Orrego, Maobin Yang, and Peter I. Lelkes. On the road smart biomaterials for bone research: definitions, concepts. advances, and outlook. Bone Research, 9(1): February 2021. ISSN 2095-6231. doi:
URL
1038/ s414144411-1000 00010///0...11114411100000000000 020-0011-222
[54] Rodrigo F Neumann, Mariane BadseeAn Everton Lucas-Oliveira, Hugo Barbalho, Willian A Trevizan, Tito J Bonagamba, and Mathias B Steiner. High accuracy capillary network representation in digital rock reveals permeability scaling functions. Scientific reports, 11(1):11370, 2021.
[55] Alexander Quinn Nichol and Prafulla Dhariwal. Improved denoising diffusion probabilistic models. In Marina Meila and Tong Zhang, editors, Proceedings of the 38th International Conference on Machine Learning, volume 139 of Proceedings of Machine Learning Research, pages 8162-8171. PMLR, 18-24 Jul 2021.
[56] Pâl-Eric Øren and Stig Bakke. Process based reconstruction of sandstones and prediction of transport properties. Transport in porous media, (0(2)::11-3433 2002.
[57] QR Passey, S Creaney, JB Kulla, FJ Moretti, and JD Stroud. A practical model for organic richness from porosity and resistivity logs. AAPG bulletin, 74(12):1777-1794, 1990.
[58) Peeyush Phogat, Shreya Sharma, Ranjana Jha, and Sukhvir Singh. Microstructural influence on electrochemical devices. In Electrochemical Devices: Principles to Applications, pages 257-306. Springer, 2024.
[59] Yinuo Ren, Haoxuan Chen, Grant M. Rotskoff, and Lexing Ying. How Discrete and Continuous Diffusion Meet: Comprehensive Analysis of Discrete Diffusion Models via Stochastic Integral Framework, October 2024.
[60] Hermann Risken. The Fokker-Planck Equation. Springer-Verlag Berlin Heidelberg, 1984. ISBN 978-3-642-96809
[61] Severi Rissanen, Markus Heinonen, and Arno Solin. Generative modelling with inverse heat dissipation. In The Eleventh International Conference on Learning Representations, 2023. doi: 11.48550/mxxiv2226.1333
[62] Javier E. Santos and Yen Ting Lin. Understanding Denoising Diffusion Probabilistic Models and their Noise Schedules via the Ornstein-Uhlenbeck Process. In NeurlPS 2023 Workshop on Diffusion Models, October 2023.
[63] Javier E. Santos, Zachary Fox, Nicholas Lubbers, and Yen Ting Lin. Blackout diffusion: Generative diffusion models in discrete-state spaces. In Andreas Krause, Emma Brunskill, Kyunghyun Cho, Barbara Engelhardt, Sivan Sabato, and Jonathan Scarlett, editors, Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pages 9034-9059. PMLR, 23-29 Jul 2023. URL https:// proceedings
html
[64] Jiaxin Shi, Kehang Han, Zhe Wang, Arnaud Doucet, and Michalis Titsias Simpli- fied and generalized masked diffusion for discrete data. In A. Globerson, L. Mackey, D. Belgrave, A. Fan, U. Paquet, J. Tomczak, and C. Zhang, editors, Advances in Neu- ral Information Processing Systems, volume 37, pages 103131- Curran Associates, ptpppppeeeeeeeeerr.peppppppppppp
Inc., 2024. URL htt.lhtttss//// htttps://pproceedingss www.nnn.ui.... cc.pss.cccpppaaaa famper_ffleee www.s/lllpppeeee palls//aaarr//22//ieeeee
[65] Patrice Simon and Yury Gogotsi. Materials for electrochemical capacitors. Nature Materials, 7 (11):845-885 854, 2008. ISSN 1476-4660. doi: 10.1038/mma2297. URL http::t/t//ooii /tttps://doo.oor/100 1038/ nmmt2297
[66] Jascha Sohl-Dickstein, Eric Weiss, Niru Maheswaranathan, and Surya Ganguli. Deep unsu- pervised learning using nonequilibrium thermodynamics. In Francis Bach and David Blei, editors, Proceedings of the 32nd International Conference on Machine Learning, volume 37 of Proceedings of Machine Learning Research, pages 2256-2265 Lille, France, 07-09 Jul 2015. PMLR. URL https::/ mlr prw.s/ srress//3i//ohh-ddcttteeeee dww.stiinss.....mmmmmmmmmmm....m html
14

---

[67] Jiaming Song. Chenlin Meng, and Stefano Ermon. Denoising diffusion implicit models. In International Conference on Learning Representations, 2021.
[68] Yang Song and Stefano Ermon. Generative modeling by estimating gradients of the data distribution. In Advances in Neural Information Processing Systems, pages 11895-11907, 2019.
[69] Yang Song, Jascha Sohl-Dickstein, Diederik P. Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-Basee Generative Modeling through Stochastic Differential Equations, February 2021.
[70] Haoran Sun, Lijun Yu, Bo Dai. Dale Schuurmans, and Hanjun Dai. Score-based Continuous- time Discrete Diffusion Models. In The Eleventh International Conference on Learning Representations, September 2022.
[71] Tianhai Tian and Kevin Burrage. Binomial leap methods for simulating stochastic chemical kinetics. The Journal of Chemical Physics, 121(21): December 2004. ISSN 0021-9606. doi: 10. 10.1063/1.181044
[72] Salvatore Torquato, Random heterogeneous materials: microstructure and macroscopic proper- fies. Springer, 2002.
[73] G. E. Uhlenbeck and L. S. Ornstein. On the theory of the brownian motion. Phys. Rev. 36: 823-841, Sep 1930. doi: 10. (0.1103/hhysRev33.822. https:s::///iik aps. org/doi/10
1103/PhysRev.. 36. 8233
[74] Francois LE Usseglio- Viretta, Andrew Colclasure, Aashutosh N Mistry, Koffi Pierre Yao Claver, Fezzeh Pouraghajan, Donal P Finegan, Thomas MM Heenan, Daniel Abraham, Partha P Mukherjee, Dean Wheeler, al. Resolving the discrepancy in tortuosity factor estimation for li-ion battery electrodes through micro-macro modeling and experiment. Journal of The Electrochemical Society, (65(14):3343--A3466 2018
[75] G Van Kampen. Stochastic Processes in Physics and Chemistry. Elsevier Science B. B.V., Amsterdam, 2007.
[76] Zhou Wang, A.C. Bovik, H.R. Sheikh, and E.P. Simoncelli. Image quality assessment: from error visibility to structural similarity. IEEE Transactions on Image Processing, 14())600--12, 2004. doi: 10. 200.//T..20038119866
[77] Markus F. Weber and Erwin Frey. Master equations and the theory of stochastic path integrals. Reports on Progress in Physics, 80(4), 2017. ISSN 00344885. doi: 10.
[78] Ludwig Winkler, Lorenz Richter, and Manfred Opper. Bridging discrete and continuous state spaces: exploring the ehrenfest process in time-continuous diffusion models. In Proceedings of the s1 International Conference on Machine Learning, (CML'24. JMLR.org. 2024.
[79] Wei Xiong, Malin Selleby, Qing Chen, Joakim Odqvist, and Yong Du Phase equilibria and thermodynamic properties in the fe-cr system. Critical Reviews in Solid State and Materials Sciences, 32(2):125-152 2010.
[80] L. Zhu, B. Bijeljic, and M.J. Blunt. Diffusion model- generation of three-dimensional multiphase pore-scale images. Transport in Porous Media, 152:22, 2025. doi: 10.1007/ s 11242-025-02158-44 URL httpp:: h/tp://doi..oo.oo org 10. 1007/ 122222111111111 1221-122-2222-----222 02158-4
15

---

A Deriving reverse-time transition rates
Here, we derive the reverse-time transition rate. Because the particles are moving independently, it is sufficient to discuss n particles colocalized at (I, y) in channel and the conclusion applies to other locations and color channels. For brevity, we will drop the (ær, y. c) dependence in this section when the context is clear. Let us index the particles by = 1 1...0 of the n particles, the reverse-time transition rates moving to (ær + Zz y Fyy) where (Po, {(-1,0), (1,0) (0, 1) (0, 1)} is
[Tol..gee. For each
F (1) Pt x + Dz + By cl 11 c
(8)
8] [ [1]
Pt 2. y, c 2
yo
co
according the general theory of reverse-time dynamics for continuous-time Markov systems [1663, We now perform the survival analysis for the inhomogeneous process. Within time dt, the probability that particle leaves (æ, y, and moves (x Pay y Py, c) r (1) dl +0 (d()). e) at time - dt is - Σ FD (t)
As such, the probability of the particle remains at (z,, v) O (dt?). Thanks to the independence between the particle dynamics, the probability of all n particles remaining at (x, y, c) at time t - dt (recall that we are evolving the reverse-time dynamics) is
-E= EE To (t) + (dr²). Then, the probability of no particle leaving at a previous time t At, where At == Ndt is
N
11-E
1 (t
dt) de + O (dt)).
(9)
k=l
which by sending dt O leads to the continuous-time survival function:
Σ 2 (C)a 1.D
PTT >
exp
(10)
where T is the random time of the first particle moving out of (æ, M, c). the sum is over all possible directions and all particle index (1...n} n). Identifying the total rate Et, p (t")dt' and the reverse-time transition rate for each particle and in each direction, Eq. 8 arrived at Eq. 4
B Training and generation algorithms.
gives the training algorithm using standard gradient descent techniques, and Algorithm gives the inference algorithm used in this work.
Algorithm I
16

---

Algorithm DSD training
Given the full transition probabilities (a v/, e v. c) repeat
To < sample drawn from the training set
Draw an index from {l.... T uniformly for Each discrete intensity unit in [ola.,e do
Draw (x', y', d') ~ Pe (a', y', xz 1r,y,C)
Move the unit from (ar, y. c) to (o'.v' end for
It, the corrupted image Compute the reverse transition rate Eq. if Using L rate-matching then
Loss + Ex.yo.. J,NN -
else if Using likelihood loss then
Loss -oloog log defined in Eq.
end if'
Take gradient step on VgLoss
until Converged
Algorithm 2 DSD inference
and desired total intensities in the color channels, initiate image To with desired total intensities in the color channels
Given CFL condition number
for Each discrete intensity unit in [4o]a. do
Draw (x', y,, (') ~ Pi (x', sf c' e)
Move the unit from (I, y, c) to Vd, .d end for
I < the fully corrupted image, while > 0 do
Evaluate NN predicted reverse rates NN D.a.yy.
T + min f, ming,a.y. (F
for each (ar, y, c) do
Sample total moving particles:
NNN ny: ~ Binom (Tela.,,ec Ep P.t.n.e)
Sample direction for each moving particle:
NN Po....
ni ~ Multinomial
[1g Pe Ent P pavuo
Move no intensity units to (x + Ūz + Dyy)
end for
Advance time: t t
I, t the configuration after movements
end while
17

---

CelebA, Ours
CelebA, x
Celeby, x²
0.75
MNIST. Ours
0.75
MNIST, x²
0.75
MNIST, x²
S5IM
0.50
0.50
0.50
0.25
0.25
0.00
0.00
0.00
500
1000
1500
2000
500
1000
1500
2000
500
1000
1500
2000
CelebA, x²
CelebA,
CelebA, x
0.75
MNIST, x²
0.75
MNIST, x²
0.75
MNIST, x
0.50
0.50
0.50
0.00
0.00
0.00
500
1000
1500
2000
500
1000
1500
2000
500
1000
1500
2000
1.00
1.00
1.00
ColebA,
ColcbA,
CelebA, Cosine
0.75
MNIST, x²
0.75
0.75
MNIST, Casine
SSW
0.50
0.50
0.50
0.25
0.25
0.25
0.00
0.00
0.00
500
1000
1500
2000
500
1000
1500
2000
500
1000
1500
2000
k(T-2000)
KIT-2000)
(TT2200))
Figure 6: Structural Similarity Index Metric between the original and corrupted MNIST and CelebA images (1,000 samples evaluated at uniformly sampled random time & E {1... .0000))) with various noise scheduler. We fixed 120 for MNIST and 200 for CelebA in this analysis. For ours, we use Eq. (2) with the cosine schedule, we use the heuristic formula given in (55] We remark that although the cosine schedule has been shown to be superior in previous studies [55 62]. the conclusion is based on the Ornstein OrnsteinUUhlenbeee process, which is distinct from the spatial diffusion process
7.5 and 72 2.5. For the polynomials, & (k/T)" , n 1... 7. For
C Spectral representation of DSD with periodic boundary conditions
Let us consider the two-dimensional DSD with periodic boundary conditions, whose master equation
di p = + ((p..(())))..))))++++ rps. (t).,
(11)
Here, we use periodic boundary conditions, PktNN.E Pk,t, Pk. = Pk,e over the domain (k, E) [1.... Nr} f1....Ny Ny). We can formally treat Pk. as rank-2 tensor. We can write the one- dimensional discrete-space Laplace operator:
LLapp. 6ptL, 8p-1,y 20,
(12)
where Lg is N₂ x N₂ and Ly is Ny x Ny. 5 is the Kronecker delta. Then, we can write the two-dimensional master equation succinetly
d dt
p(t) = (L, D1+10L,))
(13)
We aim to identify the spectral representation of the full operator L. To achieve this, we first define the two-dimensional discrete Fourier and inverse Fourier transformations:
N
2ri
mik No
km.n (1)
pk,
(14a)
VN
N₂
Ny
(t)
Σ
2
ma v.)))
km.
(14b)
mm
18

---

To simplify the notation, we write Up and Uy as the one one-dimensional inverse discrete Fourier transformation
pM
U U...ppp
(15)
N
and again use the tensor notation. Then, Eqs.
14
can be succinetly represented as
k(1) (U, @ 1) (I U,) p(t). p(t) == (U2 @ 11 @ k(t)
(16a)
(16b)
The master equation in the spectral space can be derived now:
d k(t) == (U, 2 (I
dt p(t)
dt
=r(,,0))(8U,)) (L2 @ L,)(()
r(U,,N)((3,,) (L. 1I (U, (U, (I) LUL 00
(17)
because ULU,- U I.
Let us now explicitly compute U L Un:
N₂
U U2) Uela
U
Ukkkm Lelma Uuln
N₂
N
mk-nn mn
(6m+1n + 8m-1,n 26m,n)
N
n=1
Nr
c
2x
mik (m++i
+e2n 2nl mi mk- (mn
2c 2 mk-me
N
m=1
2ns m(k-() No
20it N:
c N -2
N
2nT
2nTl N
2 cos
2
5k.F
4sin²
(18)
This shows that the process in the spectral space is diagonal with cosine spectrum. The same analysis applies to Ly Uy. Thus, element-wise, in the discrete Fourier space, we have
d km (t) = 4r ssnn dt
Tm N
+ sin 2
T Ng
km,n (t),
(19)
which means
4 sin? ( wT tsin? (x%)
km (1) e
km.
(20)
In practice, we can always shift the partiele of interest at the origin, i.e., po.o(0) 1 and Po.j if (i.j) # (). Then, km,, (0) 1/ NN,,, and consequently
1
e sin (*m)) (x;))rr
km (1)
(21)
VN..N,,
so the analytical solution of Pa.e can be expressed as
N₂
N
(t)
Σ
2ni
Nr
isiu?
cm
+sin
[
(22)
N₂ y
==1
Note that the double sum can carried out in one two-dimensional Fast Fourier Transformation. The above equation allows us compute the transition probability efficiently.
19

---

D Human-centric datasets
While our primary focus on scientific data, we include results on standard vision benchmarks to contextualize DSD's performance to the computer vision community.
D.1
MNIST experiments
As an initial validation, trained an unconditional DSD model on MNIST and found that it reliably captures digit structure and diversity (Fig.
k=1000
k=900
k=800
k= 700
k=600
k=500
k=400
k=300
k=200
k=100
k=0
9
9
9
9
9
9
9
9
3
3
3
3
3 7
7
Figure 7: Unconditional MNIST generations.
Then, our additional MNIST experiments, we explored class-conditional and inpainting generation. These experiments are particularly notable due to their interactions with the property DSD. For class-conditioning, we introduced class embeddings into our model following the approach described in [69] Our model performed well at the task of class retrieval, consistently producing the desired class (Fig.
k=1000
ko900
km:800
ko=
km6600
ko:500
kc=600
lc=300
k==200
lc=100
lcs0
8
8
8
/
/
/
f
/
2
2
2
2
3
3
3
3
3
4
4
4
4
Figure 8: Class-conditional MNIST generations for digits 0 through 4. Each row corresponds to specified target class.
For our intensity-related experiment, we tested our model on its ability to generate all of the classes given different starting intensity. Because generative models struggle to extrapolate beyond training data, our model demonstrated poor performance for certain digits initial intensity values that were too high or too low. response this, we picked the '1" with the highest intensity for our high-intensity test, and the '0" with the lowest for our low intensity test, as had the lowest intensity
20

---

of any of the numbers, and O had the highest. Our model performed very well on this task, consistently producing the target class even with varying intensity. See Fig. (d) for results.
In training our model to perform inpainting, we shrunk the size of the transition matrix and held the rest of the image static. We obtained high quality generations very quickly, within 40K training steps. For our intensity-related experiment, we tested the model's reaction to increasing intensity within the inpainted region and were able different number generations from the same starting image E for comparisons with other conditioning approaches.
(Fig. D See section
k-0000
90-900
ko-000
I--000
km-600
I--500
lcm400
ko-300
Ie-200
I-100
kmO
S
5
5
5
5
3
3
3
3
3
3
3
3
3
3
3
3
3
3
9
9
Unconditional MNIST inpainting with progressively increasing total intensity in the masked region, Rows correspond to increasing intensity levels, illustrating how digit identity and structure evolve under fixed spatial context.
Figure 9:
D.2 CIFAR-10 dataset
We trained DSD on the CIFAR-10 dataset and evaluated generative performance under different sampling configurations. Specifically, tolerance parameter, E E 0.01, 0.05, 0.10, 0.15, which controls the reverse-time integration step size. As expected, lower values of E yield more accurate reverse-time integration and lead to improved Fréchet Inception Distance (FID).
tested four values of the Courant Friedrichs-LLevee (CFL)
While FID is not the primary focus of our model, the discrete-state pixel-hopping nature of DSD introduces local high-frequency fluctuations during sampling (Fig. [T] that reduce smoothness. These artifacts, though minor, negatively impact FID and perceptual metrics. To address this, we apply a lightweight bilateral filter as a post-processing step to each generated image. This filter is non- trainable, respects both spatial proximity and intensity similarity, and acts as smooth denoising operator that preserves semantic edges and attenuates visually incoherent fluctuations. Importantly, apply it with default hyperparameters across all images, without per-sample tuning. The filtering operation is computationally negligible and can implemented as fixed convolutional layer the end of the generation pipeline.
This simple step substantially improves both FID and spatial FID (sFID), without degrading the underlying structure of the samples (Fig [2). In Fig HT, we quantify the impact of filtering across CFL tolerances. As an additional refinement step, we also use pre-trained off-the-shelf CIFAR- 10 classifier to discard semantically incoherent samples. Filtering based on classifier confidence can be integrated early in the sampling loop and correlates strongly with out----disti or degenerate outputs (Fig. IT)).
21

---

A summary of results of these sampling experiments is shown in Table While an FID score of around 20 does not represent state-of-the-an performance, we hypothesize that incorporating larger architectures, advanced optimization strategies, multi-GPU trainings, and data augmentation could further reduce this gap. DSD may be a viable foundation for discrete, constraint-preserving generative modeling even in conventional RGB settings.
Figure 10: Unconditional CIFAR-I0 generations produced by DSD with CFL tolerance 0.01. Samples are generated without any post-processing. While structurally coherent, mild pixel-level irregularities are visible.
Table 1: Quantitative evaluation of CIFAR-10 generations under different sampling and filtering strategies (all scores computed on 50,000 images). The CFL tolerance & controls the reverse-time integration step size. Filtering denotes post-processing with bilateral filter. Classifier thresholding retains samples with 2 0.99.
Sampling Configuration
FID +
sFID J
CFI. E = 0.01
46.3
23.4
CFI. e = 0.01 filtering
28.2
16.4
Classifier thresholding (p 2 0.99)
18.1
15.7
22

---

CFL sampling
Bilateral filtering
Classifier thresholding
Tok TOW 67
Toleurce 14t
Thocllol
10%
35%
0.99
00-
7:000 Number of Samples
20.0002.5003 10032
5.000 20.000200011000022200 201.2002.0022.2002....
Number Samples
Number of Samples
Figure 11: FID scores computed on 20,000 CIFAR-10 10 generations under different sampling configu- rations. (Left:) Samples generated with different CFL tolerances (é). (Middle:) Same samples after applying a fixed bilateral filter. (Right:) Bilaterally filtered samples sampled with 1% further refined by discarding outputs with classifier confidence below certain threshold.
Figure 12: CIFAR-10 generations after bilateral filtering. These images correspond directly
the unfiltered samples shown in Fig. I0 post-processed using fixed, non-trainable bilateral filter. The operation smooths pixel-level irregularities introduced by the discrete sampling process while preserving structural and semantic content.
23

---

D.3
CelebA dataset
To demonstrate scalability to higher-resolution, human-centric data, we trained DSD on the Celeb 64x64 dataset. While the resulting FID of 29 is not competitive with state-of-the-art models, it remains consistent with our CIFAR-10 results and underscores the feasibility of applying DSD to RGB images with more complex structure. As illustrated in Figure 13 the model is capable of producing globally coherent facial samples, despite operating in a discrete, intensity-consert speculate that such discrete spatial models could be beneficial for downstream applications involving controlled image manipulation, such as recolorization, region editing, or consistency musissency-preserving transformations.
setting. We
Figure 13: Unconditional generations of the CelebA 64x64 dataset.
Table 2: Quantitative evaluation of CelebA 64x64 generations under different sampling strategies. The CFL tolerance controls the reverse-time integration step size. Filtering refers to post-processing with bilateral filter. These scores computed with 50,000 realizations.
Sampling Configuration
FID +
sFID
CFL 0.05
44.2
29.0
CFL 0.05 + filtering
28.9
25.7
E
Comparison with existing frameworks
E.1
Gaussian Diffusion with Conditioning
To benchmark against conventional generative models, we implemented a Gaussian diffusion baseline conditioned on the target total intensity (or mass). The architecture employed continuous-time variant of NCSN++ with explicit conditioning on total intensity per image via an auxiliary embedding channel. Fig. 14 shows the behavior of this model across a range of target intensities. Near the dataset mean, where conditioning remains within the empirical distribution, the model performs moderately well. However, closer examination reveals that the generated images often contain small negative values, even though training data was strictly bounded in [0, 1). While this behavior may be inconsequential in floating-port point representations, becomes problematic when the outputs are discretized and clipped for evaluation downstream use. In particular, thresholding negative values
24

---

and rounding to 8-bit precision systematically increases the intensity error. These issues become increasingly severe in the distribution tails, where the model frequently fails respect the specified conditioning.
In contrast, our proposed DSD framework guarantees exact intensity preservation by construction, without requiring post hoc corrections, which holds uniformly across the intensity distribution, including extreme quantiles, as shown in Fig. 14 Notably, when sampling at the lowest intensity levels, DSD generates extremely thin or skeletal "1"s, while very high intensities produce over- saturated, bold structures such as "fat" zeros that span much of the image domain. This behavior suggests that the model is not merely enforcing conservation, but learning meaningful structural adaptations accommodate different global constraints, consistent with underlying patterns in the data.
105
Mean Intensity of the Realizations (0-255)
75
45
Gaussian diffusion
15-
Gaussian diffusion (Discretized)
Discrete Spatial Diffusion (ours)
2.6
6.4
12.8
25.5
51.0
76.5
89.2
102.0
Target Intensity (0-255)
Figure 14: Sample-wise mean intensity across 1,000 realizations of mnist digits for each target intensity. We compare Gaussian diffusion (white), its discretized version (gray), and Discrete Spatial Diffusion (black horizontal lines). Each violin shows the distribution of mean intensities computed across 1,000 individual samples.
Our observations indicate that: (I) intensity can be statistically conditioned within the distributional range, but without guarantee of exactness; (2) rounding and clipping of generated samples introduce further bias, revealing that the model relies physically invalid negative intensities approxi- mate the constraint; and (3) extrapolation outside the data manifold results in complete failure of conditioning.
E.2 Gaussian Diffusion with Conditioning for Scientific Data
To directly compare DSD with the standard continuous Gaussian diffusion framework in scientific setting, we trained continuous-time NCSN+ ++ baseline with an explicit embedding of the target total intensity fraction (porosity) and evaluated it on the Savonniéres carbonate microstructure dataset (binary, 128x128). For each target we drew 100 samples, measured sample-wise porosity error, and inspected morphology. The Gaussian baseline tracks the target near the training range but degrades the tails (large relative errors with low variance across samples) and exhibits characteristic sall-anddddaae salt-end--eepeer speckling and "staircase" grain boundaries consistent with conditioning that is only approximate and not enforced per-sample. These results are presented in Figure [5] DSD realizations withing the training distribution and well beyond the training range are shown in Figure 16 and Figure
artifacts sataacts-sal--addd
E.3 Guidance Function Approach
In reinforcement learning settings, a guidance function can be naturally defined using task reward or cost structures 40 In generative modeling, however, such functions are generally unavailable, and often must be learned or heuristic in nature, as discussed in work on posterior sampling. Specifically, exact guidance is only possible when the conditional distribution exact
25

---

0.72%
Predic
2.32%
1.82%
2.84%
24.45%
32.34%
Target porosity
Figure 15: Gaussian diffusion baseline with porosity conditioning on Savonnières carbonate. (Binary, 128 128). The baseline approximately matches targets near the training regime but incurs large relative errors when extrapolating. Relative enrors >20% are omitted from the plot.
0.08%
0.02%
0.14%
0.1 9%
0.28%
0.48%
0.19%
O. 66%
drod borsatu
0.15%
0.09%
0.10%
0.14%
0.15%
0.1..9%
Figure 16: DSD realizations of the Savonnières carbonate. Samples generated with decreasing target porosity, extending well beyond the training range. Per-sample porosity closely follows the target (points lie on the identity line), indicating exact global intensity conservation by construction. Nevertheless, A negligible residual error arises from small fraction of overlapping pixels at the end of the reverse process. These appear as small white specks in the plots (visible in the online high- resolution version) and can be trivially corrected by post-processing step that relocates overlapping pixels to the nearest empty voxel.
0.70
0.25
B.10
0.A
I-12
0.53
8.57
0.10
0.57
0.71
0.76
0.ED
2.04
0.19
1-
# 3r
# 30
p + 7a
Intensity training images
Figure 17: Out distribution DSD realizations the Savonnières carbonate. DSD-generated samples conditioned total intensity from μ- to μ+70. Despite being trained near the data distribution center, the model produces coherent and physically plausible structures far beyond the training range, highlighting robust extrapolation enabled by exact intensity preservation.
posteriors (of generative samples given the constraint of interest) is known, which is not the case for total intensity constraints in image generation.
To explore whether heuristic guidance functions could serve surrogate, we conducted additional experiments applying mean-squared error (MSE) penalty between the generated and target total intensity. Our results showed that while guidance can reduce error somewhat for intensities near the
26

---

mean training data intensity, it fails to enforce conditioning reliably, and especially so at the tails of the training distribution. In these cases, samples frequently severely violate the intensity constraint.
F Subsurface microstructures
F1 Detailed description of X-ray scans of subsurface rocks
This sandstone sample from [54] provides a high-resolution image of the rock microstructures obtained through X-ray microtomography (X-ray pCT). In this process, the rock sample is rotated while being scanned by an X-ray beam, capturing series of 2D radiographs different angles. These projections are then computationally reconstructed into a 3D volume, where each voxel represents the X-ray attenuation of the material at that location. The X-ray microtomography scans were performed using a SkyScan 1272 system, operating at 50 kV and 200 μA, with CCD detector capturing projections at a resolution of 2.25 jm per voxel. The resulting dataset consists of grayscale images with voxel size of 2.25 jum, where variations in intensity distinguish between the solid matrix and the pore space. The solid matrix primarily consists of tightly packed mineral grains- quartz-while the pores are voids that can be occupied by fluids such as water or hydrocarbons. After preprocessing steps like contrast enhancement, noise reduction, and segmentation, the final dataset represents the pore network. The Berea sample has measured porosity of 18. 96% and permeability of 121 mD. This dataset is particularly useful for computational modeling, it enables direct comparison between numerical simulations and experimentally measured permeability, providing a rich testbed for learning based methods that seek to map complex microstructural information to macroscopic transport properties. This sedimentary rock is a well-characterized geological benchmark, widely used in fluid flow studies due to its homogeneous grain structure and consistent permeability properties, making it a good first benchmark for our study.
Berea Sandstone:
Savonnières Carbonate: This carbonate sample, described in 2), is a layered, oolithic grainstone with wide porosity and permeability varying from 115 to over 2000 mD, depending on local heterogeneities. The rock is characterized by a highly multimodal and interconnected pore structure, with distinct macropores and microporosity. X-ray microtomography (X-ray μCT) was used to image the sample at a resolution of 3.8 pm voxel size, revealing intricate pore geometries. The sample was scanned at the Ghent University Centre for X-ray Tomography (UGCT) using their HECTOR scanner, developed in collaboration with XRE, Belgium. The macropores include both intergranular voids and hollow ooids, while the microporosity is found within ooid shells and intergranular spaces. Micropores in the sample often serve as the primary pathways connecting poorly connected macropores, creating complex hierarchical network. After preprocessing steps, including noise reduction, anisotropic diffusion filtering, and watershed segmentation, multiscale pore network model was extracted. This dataset is particularly compelling due to its extreme heterogeneity, with pore sizes spanning orders of magnitude, and its ability to represent coupled serial and parallel flow pathways. Savonnières serves as a test case for studying the impact of complex samples in our workflow.
Massangis Limestone: This oolitic limestone sample from [ is a highly heterogeneous car- bonate rock with a complex, multimodal pore structure resulting from diagenetic alterations, including dolomitization and dedolomitization. The rock contains mix of intergranular and moldic macroporosity, along with microporosity concentrated in ooid rims and partially dissolved dolomite regions. Its porosity ranges from 9.5% to 13.8%, depending on local variations, and its permcability is highly anisotropic due to the interplay between connected macropores and poorly accessible microporosity. X-ray microtomography (X-ray PCT) was used to image the sample at voxel resolution of 4.54 μa, capturing the intricate connectivity of macro- and micropores. The sample was scanned at the Ghent University Centre for X-ray Tomography (UGCT) using a FeinFocus FXE160.51 transmission tube, in collaboration with Paul Scherrer Institute (PSI), Switzerland. Differential imaging was applied to enhance the detection of fluid-filled microporosity. revealing the rock's internal heterogeneities. Unlike more uniform carbonate samples, Massangis exhibits significant spatial variations in pore connectivity, leading zones of high permeability interspersed
27

---

with isolated pore networks. This dataset serves as another challenging benchmark for modeling porous media microstructure.
Berea Sandstone
Savonnières Carbonate
Massangis Carbonate
(Thiiingg
(Training)
(Training)
(eterration)
(Concration)
(Geseration)
(Geseration)
(Cencration)
Figure 18: Schematic representation of three rock types: Berea Sandstone, Savonnières Carbonate, and Massangis Carbonate. The first image (top left) for each rock type shows one training sample, while the second one (top right) displays the generated sample conditioned the mean intensity l of the training set. The third (bottom left) and fourth (bottom right) samples illustrate the generated samples conditioned on 2a and μ + 2o, respectively, where represents the standard deviation of the training set intensity distribution.
F22 Quantitative metrics of the generated images
In porous media analysis, characterizing the spatial arrangement and size distribution of pores is crucial for understanding transport properties, mechanical behavior, and overall seuuctureffucctinn relationships. To quantify these characteristics, size distribution (PSD) using PoreSpy 34, a Python-based toolkit for quantitative analysis of porous media images. The Pore Size Distribution (PSD) characterizes the variation of pore sizes within porous material, providing insights into connectivity, permeability, and flow dynamics. The most common method to determine PSD computationally the local thickness approach. Given binary image 7 (xx v). where pore space is represented as and solid space as 0, the pore size function f(r) is defined the probability density function (PDF) of the largest sphere that be inscribed at any point within the pore space. The PSD provides a statistical summary of pore connectivity and transport properties. Small pores dominate permeability, while large pores govern bulk flow. Both of these metrics for the training and generated samples are shown in Fig. 19 our realizations exhibit excellent agreement with the training data.
compute the spatial correlation function and pore
Normelized pore size distribution
TWo point correlation function
Training data Generated
Training data Generated
150-
0.6-
0.4
50
0.2
0.0
20
40
60
0
50
100
150
200
0.0
0.2
0.4
0.6
0.8
Distance [pixels]
Normalized pore size [pixels]
0.2
0.4 0.6
0.8
1.0
Figure 19: Quantitative comparison between training and generated rock samples. (Left) Two-point correlation function, showing excellent agreement of spatial features between training data (gray) and 100 randomly generated samples (blue). (Middle) Normalized pore size distribution schematic. with colors indicating relative pore sizes. (Right) Cumulative distribution function (CDF) of the normalized pore sizes, comparing the statistical distribution of training and generated samples.
28

---

Additionally, following the evaluation protocol in Lee and Yun 46 we computed Fréchet Inception Distance (FID) for our generated carbonate samples. obtaining score of 0.9, substantially lower than their reported value of 18.1 rocks of similar lithology. While this comparison is not strictly one-to-one due dataset and preprocessing differences, offers a compelling quantitative indication of the improved realism and fidelity of DSDD-generate microstructures. We note, however, that FID may have limited interpretability for binary porous media images, and should be interpreted alongside domain-specific metrics.
F.3
Large scale generations.
To evaluate the scalability of DSD in realistic scientific settings, we trained model on high-resolution binary microstructure data of Leopard Sandstone, an interesting sample from the geological standpoint. This dataset consists of 1000 1000 samples representing complex subsurface structures characterized by heterogeneity, anisotropy, and nontrivial pore connectivity. Generating such large samples is particularly relevant for studying representative elementary volumes.
Despite the increased resolution and number of particles, DSD remains computationally feasible in this regime. The forward process remains efficiently parallelizable across CPU workers, and reverse-time sampling scales linearly with intensity while maintaining tractable runtimes. A subset of generated samples is shown in Fig. 20) demonstrating the model's ability reproduce structurally diverse and geologically plausible microstructures at scale.
While extending DSD 3D requires modest modifications to the implementation (e.g., defining 6-connected lattice for particle jumps), we find the approach remains tractable for volumetric data, further supporting its suitability for scientific applications in porous media and materials modeling.
Figure 20: Realizations of Discrete Spatial Diffusion (DSD) on 1,000 1,000 images of Leopard Sandstone. These samples demonstrate that model scales high-resolution, structurally complex porous media. Leopard Sandstone exhibits heterogeneous, anisotropic pore structures representative of real-world geological variability.
G
X-ray scans of NMC cathodes
G.1 Dataset description
This dataset provides high-resolution 3D images of Li-ion battery cathode composed of active material
cobalt oxide, NMC), carbon black, and polymer binder [74). The cathode sample was imaged via X-ray microtomography (X-ray μCT) and tanoooommgryhhy (X-ray nano-CT) capture both the overall electrode architecture and fine-scale features of the carbon/binder domain (CBD). For micro-CT. a Zeiss Xradia Versa 520 system was operated 80 kV and μA, acquiring projections at effective isotropic voxel size of approximately 398 nm over a field of view of about 400 μm. The nano-CT CT scans were performed using Zeiss Xradia Ultra 810 system with chromium target (35 kV, 25 mA), yiclding isotropic voxel sizes the order of 126 nm across field of view of approximately 64 μm. In both cases, the 2D radiographs were reconstructed into 3D grayscale volumes using filtered back-projection algorithm, capturing the X-ray attenuation due to the dense NMC particles and the less attenuating pore/CBD regions.
These tomographic datasets reveal the hierarchical microstructure of the electrode, from tens-of- micrometers NMC active particles to nanometer-son scale pores within the percolated carbon network. After preprocessing such as non-local mean filtering, contrast enhancement, and slice-by-slice local thresholding serssooddnggxsegmentaiinn identifies three main phases: (1) the NMC active material, (2) the CBD, and (3) the pore space. Measured porosity values for these cathodes can exceed 30%, while the typical
29

---

0.10
Normaten mneettee laneth
0.030
0.05
0.025
0.94
Mhose
0.33
0.020
0.32
0.81
0.015
Cartion bincer
Elcttolyte
Blecrrlyte
Figure 21: Microstructural characterization metrics for 80 training samples and 80 generated samples. The boxes show the 25th-50th-75tt percentile, the whiskers the minimum and maximum values. Metrics computed using TauFactor [41].
volume fraction of active material is on the order of 40%. The overall areal loading of the active material is around 29.78 mg-cm²² corresponding to about 33 mAh-cm" in specific capacity. These 3D reconstructions enable computational modeling of transport properties (e.g., tortuosity factor) and electrochemical performance, facilitating direct comparisons with experimentally measured parameters. Because of the electrode's s well-defined spherical NMC particles and percolating carbon network, this dataset serves as robust benchmark for multi-scale modeling and data-driven methods that aim to link microstructural features macroscopic cell behavior.
G.2 Effective metrics
The analysis of NMC cathode tomography and the generated images was conducted using three metrics: interface length, triple-phase boundary, and relative diffusivity. These metrics essen- tial for quantifying the morphological and transport characteristics that influence the electrode's electrochemical performance. Below describe these metrics in detail.
Interface length refers to the total length of boundaries where two distinct phases, such as active material and pore electrolyte, intersect. A higher interface length indicates more active sites for electrochemical reactions and enhances ion transport pathways, thereby improving the electrode's overall performance. This metric is calculated by identifying and summing the perimeters of all phase boundaries in the segmented image.
Triple-Phase Boundary denotes the regions where three different phases-typically solid active material, electrolyte, and conductive phase or pore convere the microstructure. TPBs are crucial for facilitating efficient electrochemical reactions, they provide optimal sites where all necessary phases interact. The total TPB length is determined by locating points or lines where three phases meet and summing their lengths within the image.
Relative Diffusivity quantifies the reduction in ion transport within the porous cathode structure relative to an unobstructed medium. It is defined as the ratio of the effective diffusivity, Dar. through the porous medium to the intrinsic diffusivity, Du, the conductive phase: D Dat/ Dg. This reduction is primarily attributed to the geometric complexities of the microstructure, encapsulated by tortuosity factor, т, in fact Dre Den / Do Vy T,, where Vy is the volume fraction the phase under analysis.
We computed these metrics using the Python library TauFactor 41, and the comparisons between the real and generated images based on these metrics are illustrated in Fig. 22, while collection of the training data and generated images is in Fig
H Hyperparameters for experiments
In our experiments. we thoroughly tested our model on various hyperparameters using the MNIST dataset. The MNIST dataset was chosen as a baseline for hyperparameter testing due to its low com- putational training cost. We found that our model was very robust with respect to the hyperparameters used, consistently generating quality generations without hyperparameter tuning. Due to limited compute, only limited tests were performed on CelebA, but we hypothesize that our model would perform well with different hyperparameters than the ones used. For the choice of our 'r', chose
30

---

Training Data
Generated
Figure 22: (Top) Eight randomly picked samples from the NMC cathodes dataset. (Bottom) Random unconditional realizations of our model.
rate that was large enough to allow full degradation, enabling the model to learn to predict starting from full noise. See 3 for our hyperparameters used.
Table 3: Hyperparameters used in all our experiments. All models were run in single NVIDIA A100 (or similar).
Boundary Condition
CFL Tolerance
Channel Multiplier
Training Iterations
Dataset
Schedule
Loss
Notes
MNIST
120
Ours
Periodic
Eq.
(2.2.2)
100K
unconditional
MNIST
120
Ours
Periodie
Eq.
0.15
(2,2.2)
90K
unconditional
MNIST
120
Ours
No-flux
Eq. 6
0.15
(2.2.2)
80K
unconditional
MNIST
120
Ours
No-flux
14g. 5
0.11
(2,2.2)
70K
class- conditioned
MNIST
85
Ours
No-flux
Eq.
0.07
(2,2.2)
40K
inpainting (14x14)
CIFARIO
160
Ours
Periodie
Eq.
Fig.
II
(1.2.2.2)
10M
unconditional
CelebA
200
Ours
No-flux
Eq.
0.05
(1,2,2,2
3M
unconditional
Electrodes
200
Periodie
Eq.
0.01
(1,2,2.2
180k
Rocks
Periodic
Hq.0
0.1/0.2/0.05
(1,2,222)
lolerance avoids lapping particles
31

---

1. Claims
Question:
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: The abstract and introduction clearly state the main contributions of the paper, namely the introduction of a discrete-state diffusion framework (DSD) that strictly preserves global intensity (e.g.. mass or particle counts). These claims are substantiated throughout the paper with theoretical formulations and experimental validation on both synthetic and scientific datascts.
2. Limitations
Question:
Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: The limitations are explicitly discussed Section including computational
costs and the requirement for custom forward process and novel sampling scheme.
Guidelines:
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
The paper provides theoretical derivations of the forward and reverse processes, including transition rates and training objectives. Derivations are included in Section 3 and Appendix A, and all assumptions are clearly stated.
Justification:
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
he paper discloses key components necessary to reproduce results, including training and sampling algorithms (Algorithm I and 2), model architecture, noise schedules, and dataset details. Further implementation details are provided in appendices.
Justification:
5. Open access to data and code
Question: Does the paper provide open access the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification: The paper mentions training on public datasets like MNIST and CelebA, and domain-specific open-source datasets are described and linked in Appendices F and G. While the exact repository link is not included in the draft (for anonymity, but also, the code release is still pending institutional review), the code will made public in case of acceptance.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (c.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: The paper specifies dataset types, image sizes, model architectures (e.g.. NCSN++ with modified outputs), and noise schedule parameters (e.g. Tı, r). These sufficient to understand and replicate results. Table includes the detailed values of every experiment ran.
7. Experiment statistical significance
32

---

Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [No]
Justification: While experiments demonstrate qualitative fidelity and conditioning ability, error bars or statistical significance metrics (e.g., variance across runs or seeds) are not reported. The central claim of the paper- the conservation of particle count-is mathemati- cally grounded and does not require error bars.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification: The paper addresses resource concerns from both a theoretical and practical perspective in the limitations section
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics
Answer: [Yes]
Justification: The research does violate the NeurIPS Code of Ethies. Anonymity is preserved, and no ethical concerns arise from the datasets or methods described.
10. Broader impacts
Question: Does the paper discuss both potential positive socictal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification: This paper presents work whose goal is to advance the field of generative modeling for spatial data. There are many potential societal consequences for generative modeling research, however, these are largely unspecific to the research presented here. This paper advances the ability to generate image data under constraints, which could improve the capabilities of generative models, and in particular in settings which are treated mathematically, such as physical simulations; we believe this does not broaden the scope of ethical concerns associated with generative models.
11. Safeguards
Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk of misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Answer: [NA]
The paper presents new capabilities for controlled generation of microstructures which do have high risk of misuse. The capability for controlled generation of human images is generally similar to that found in widely available models and does not pose large specific risk relative to the landscape of publicly available models.
Justification:
12. Licenses for existing assets
Question: Are the creators original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: Public datasets such MNIST and CelebA were used, and their sources cited. An open-source implementation of the NCSN+ architecture was used and cited.
13. New assets
Question: Are new assets introduced in the paper well documented and the documentation provided alongside the assets?
33

---

Answer: [No]
Justification: We plan for our paper to be accompanied with an open-source code release for the method, however, this is as-of-yet still pending institutional review.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification: No human subject research or crowdsourcing.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: Our IRB, who we consulted, deemed this work not human subject research.
16. Declaration of LLM usage
Question: Does the paper describe the usage of LLMs if is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
Answer: [NA]
Justification: LLMs were not involved in this research.
34


Reviews:
[{"review_id": "JsgJ3CM9Fq", "reviewer": "Reviewer_ADQ6", "summary": "The authors proposes Discrete spatial diffusion (DSD), a continuous time, discrete-state diffusion framework that preserves total image intensity (i.e., particle counts) in both forward and reverse processes by modeling intensity units as spatially correlated random walks. Their proposed method achieves competitive performance compared to continuous diffusion models (in the natural image domain), and demonstrates strong domain specific performance by accurately reconstructing porous rock microstructures and Li-ion battery electrode morphologies under strict intensity conservation.", "questions": "The custom jump-process and adaptive sampler adds engineering complexity into well established diffusion frameworks. Are there plans to release modular code or wrappers that would make DSD easy to plug into existing open source codebases?", "limitations": "yes", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 3, "significance": 3, "originality": 4, "strengths_and_weaknesses": "Strengths\n1. The work has a clear problem motivation, addressing the need for generative models that strictly conserve total image intensity, which is needed for scientific and engineering applications where mass/energy preservation is non‐negotiable.\n2. The proposed method guarantees strict mass conservation per color channel, enabling scientifically meaningful generation under hard constraints. They also rigorously justify their work through theory. The derivations of forward and reverse markov jump processes, noise scheduling via SSIM, and principled loss functions (rate matching) are clear.\n\nWeaknesses\nTheir proposed method's forward and reverse sampling scales linearly with total intensity, leading to slower inference compared to vanilla diffusion models and potentially prohibits runtimes for high bit depth or very high resolution images. However, I acknowledge that these weaknesses are not within the scope of this work, and the authors also make an effor to acknowledge these weaknesses through engineering implementations and a section in their appendix (F.3) that directly addresses that their method remains computationally feasible.\n\nOverall, this work is technically solid with high impact in generative modeling under conservation laws, well-substantiated theory, and strong domain specific results. Minor efficiency and usability tradeoffs do not outweigh its contributions.", "was_revised": true, "final_justification": "My final justification is to keep my original score. On the matter of oral/spotlight, I believe a strong rating from all reviewers make it eligible for spotlight, but it might be focused on a very narrow field for an oral presentationl.", "rebuttal": "Thank you for the thoughtful and encouraging review, we're happy that the motivation and theoretical clarity came through.\n\n**[Q1]** We appreciate the comment. The codebase which we will release is relatively modular and designed using scipy, numpy, pytorch, and the scoremodels python package so as to best make use of existing tools and isolate the complexities associated with DSD in particular. We have routines for a) generating the Greens functions associated with the noising process, b) applying the noising process using dataloader workers, and c) sampling from a trained DSD network, which are where the additional engineering complexities lie within DSD. A user can easily swap in standard diffusion components, such as different architectures, and time (noise) schedules. Datasets in binarized or RGB formats can be incorporated with minimal effort, and it is possible to modify the integration parameters/strategies used during sampling. We will release the self-contained codebase on GitHub upon acceptance."}, {"review_id": "ekpPR62oaI", "reviewer": "Reviewer_n1in", "summary": "This paper proposes a novel variant of diffusion models that preserves particle counts throughout both the forward and reverse processes, which is a property desired for scientific applications where discrete quantity conservation is essential. The key idea is to reformulate the standard Gaussian corruption process as a continuous-time, discrete-state jump stochastic process that operates directly in the discrete spatial domain. The model’s effectiveness is demonstrated through experiments on two scientific applications, with additional validation conducted on standard image benchmarks for computer vision tasks.", "questions": "Please see the weaknesses mentioned above.", "limitations": "Yes.", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "**S1:** The problem addressed is niche but compelling in the context of scientific applications. The paper is well-structured and clearly presented, with both the motivation and core technical contributions effectively communicated.\n\n**S2:** The Markovian jumping process is an interesting and intuitive idea to achieve the preservation of intensity from the formulation level, thus presenting a rather elegant way to solve the challenge.\n\n**S3:** The experiments effectively demonstrate the benefits of the proposed intensity-preserving property, even though this characteristic may not appear immediately critical for general vision tasks.\n\n**W1:** As acknowledged in the limitations section, the proposed DSD framework introduces significant computational overhead compared to standard Gaussian diffusion models, both during training and inference. Given that Gaussian diffusion is already computationally intensive, this further limits the practical adoption of DSD in mainstream vision applications.\n\n**W2:** The experimental setup and evaluation protocols, particularly for the scientific use cases, remain somewhat unclear. For instance, is the goal to synthesize images of subsurface rock microstructures and lithium-ion electrodes from randomly initialized images, as outlined in Algorithm 2? If so, the inference process appears to require predefined condition numbers and intensity values to initiate generation. How feasible or practical is this requirement within real-world scientific workflows? Similarly, on the evaluation side, while the proposed DSD demonstrates improved intensity preservation, are the generated scientific samples rigorously valid within their respective domains? A more detailed discussion on the scientific fidelity and downstream utility of the generated data would strengthen the work.\n\nMinors: There are some formatting inconsistencies across subtitles, for example, some section headers begin with bolded first letters while others do not.", "was_revised": true, "final_justification": "Keep my positive rating after rebuttal.", "rebuttal": "Thank you for the detailed and constructive feedback, we're glad the motivation and elegance of the approach came through, and appreciate the opportunity to clarify our scientific evaluation protocols.\n\n\n**[W2]** Thank you for pointing this out. We apologize that aspects of the scientific dataset evaluation were unclear. For the human-centric datasets, we initialized inference from a random distribution to demonstrate our model’s capacity to generate human-centric RGB images. However, for the scientific domains (here, rock and battery microstructures), the inference allows for more practical use cases. After training, our framework can interpolate across intensity values not seen during training, enabling the exploration of intermediate or extrapolated physical states. For instance, in Appendix F.3 (MNIST), we show that the model produces coherent, semantically valid outputs even under out-of-distribution intensity conditions, suggesting learned representations that generalize and reflect meaningful physical structure.\n\nThe ultimate goal of scientific applications is to generate microstructures conditioned on specific parameters for downstream analysis. In earth sciences, the types of samples we present in this work are immensely valuable, not only because acquiring them can cost millions of dollars, but because they offer a window through which we can observe and measure the complexity of the geology below us. Extracting as much information as possible from these scarce datasets, and then being able to generate new realizations, opens up avenues like performing uncertainty quantification, or interpolate to physical regimes that were not directly observed during sampling.\n\nPractically, in the case of rock microstructures, the generated samples could inform decisions in applications such as aquifer remediation, CO₂ sequestration, or hydrogen storage. In battery science, the aim is to generate novel electrode morphologies that are optimized for energy storage performance, subject to physical and material constraints (largely driven by economics). With these in mind, in real-world workflows, summary descriptors like porosity (for rocks) can often be estimated without expensive 3D x-ray imaging for example, through bulk measurements, indirect methods at the field (e.g., well logs), rapid and inexpensive lab-based porosity tests, or conventional 2D imaging, making it feasible to condition generation on realistic parameters. This opens up avenues for generating physically plausible realizations consistent with known constraints.\n\nWe also conducted rigorous domain-specific evaluations. For the rock samples, we computed pore size distributions and two-point correlation functions using PoreSpy[1]. These showed excellent agreement between real and generated samples (Appendix F.2, Fig. 16). For battery electrodes, we evaluated interface length, triple-phase boundary, and relative diffusivity using TauFactor[2] (Appendix G.2, Fig. 19), all of which quantitatively confirm the structural and functional fidelity of our generated images. Additionally, we report a significantly improved FID score (0.9 vs. 18.1 in [3]), although we acknowledge FID’s limitations for binary images. These analyses support the scientific validity and potential utility of DSD samples in downstream modeling/design/ and uncertainty quantification tasks.\n\n [1] Gostick, J. T., Khan, Z. A., Tranter, T. G., Kok, M. D., Agnaou, M., Sadeghi, M., & Jervis, R. (2019). PoreSpy: A python toolkit for quantitative analysis of porous media images. Journal of Open Source Software, 4(37), 1296.\n\n[2] Cooper, S. J., Bertei, A., Shearing, P. R., Kilner, J. A., & Brandon, N. P. (2016). TauFactor: An open-source application for calculating tortuosity factors from tomographic data. SoftwareX, 5, 203-210.\n\n[3] Kang-Hyun Lee and Gun Jin Yun. Microstructure reconstruction using diffusion-based generative models. Mechanics of Advanced Materials and Structures, 31(18):4443–4461, 2024.\n\n**[Minors]** Thank you, we have identified some inconsistencies in the capitalization used within section titles and have resolved the matter; hopefully this is the inconsistency to which you refer. We will make every effort to hunt down any other formatting problems, and if you have further issues to note then just let us know."}, {"review_id": "Y8xkIjnQcg", "reviewer": "Reviewer_tBZm", "summary": "The paper proposes a discrete spatial diffusion (DSD) model, in which multiple particles are distributed in space and their total number is preserved throughout both the forward and reverse processes. This conservation property allows the model to control the number of generated particles by specifying the initial count in the reverse process. In both diffusion directions, each particle stochastically transitions to a neighboring region. A deep neural network (DNN) is employed to learn the transition probabilities at each time step of the reverse process. The training objective is conditioned on the observed data $I_0$ and is formulated analogously to standard diffusion models. The authors validate the effectiveness of the proposed approach through experiments on common image datasets (MNIST, CIFAR10, CelebA) and domain-specific data such as rock microstructures and lithium-ion battery electrodes.", "questions": "### Questions\n\n* Regarding Eq. (1a):\n    * This equation appears to describe the transition of a single particle in the forward process, with $r$ referred to as the transition rate. However, in Table 8, values such as $r = 120$ are used. What does it mean for the transition rate to be $r = 120$?\n    * If $r$ represents the number of particles transitioning from $I_{x,y,c}$, how is the transition handled when $I_{x,y,c} < r$?\n    * Does the model assume that transitions are made with equal probability to the four neighboring directions included in $\\nu$ (or $\\bar{\\nu}$ in Eq. (3))?\n\n\n### Comments\n* In Appendix H, the authors state that “we found that our model was very robust with respect to the hyperparameters used.” It would be helpful to include a comparison of the forward and reverse processes under different values of $r$ (e.g., on MNIST), to illustrate this robustness more concretely.\n\n* Throughout the paper, the abbreviation for the Inverse Heat Dissipation Model should be used consistently. It is currently inconsistent between “IHDM” and “IHD,” and should be unified.", "limitations": "See weaknesses above.", "rating": 5, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 4, "strengths_and_weaknesses": "### Strengths\n* A novel approach to modeling spatial distributions of particles using diffusion models while preserving the total number of particles.\n* Experimental validation on real-world applications, including rock microstructures and lithium-ion battery electrodes.\n\n### Weaknesses\n* Certain aspects of the forward and reverse processes remain unclear (see Questions below).\n* As noted by the authors, computational cost becomes a critical issue in high-intensity scenarios.", "was_revised": true, "final_justification": "The authors' rebuttal satisfactorily addressed my concern regarding the definition and role of the transition rate $r$. I appreciate the clear explanation, and this issue is now resolved.\n\nThere are no remaining concerns from my side. I maintain my original score, which reflects my overall positive assessment of the paper.", "rebuttal": "Thank you for carefully reading of our work. We appreciate your thoughtful questions about the transition dynamics and your attention to detail regarding notation and clarity. Please find our detailed responses to your comments and questions below.\n\n\n**[Questions]**\n\n- We follow the standard model specification for continuous-time stochastic processes, for example, in Van Kampen’s *Stochastic Processes in Physics and Chemistry*, Gardiner’s *Stochastic Methods*, and Durrett’s *Essentials of stochastic processes*.\n\n- By definition, a transition rate is the probability of a potential “reaction” to take place per unit time, in the infinitesimal-time limit. For our forward process, $r=1000$ would mean that within time $\\Delta t \\downarrow 0$, a specific particle will move to its right (as well as left, up, and down) with a probability $p=1000 \\Delta t$. The infinitesimal condition is imposed to eliminate multiple jump events, such as two consecutive jumps, which have negligible probability ($\\mathcal{O}(\\Delta t^2)$) in this infinitesimal time limit. Writing these operations with finite-time probability equations is far more complicated as it requires accounting for multiple-jump events even when the process is generated by simpler single-jump rates.  To improve the clarity of the manuscript we will explain this more carefully.\n\n- It is important to point out that in Eq. (1), $r$ stands for the transition for each of the particles located at $(x,y,c)$. Thus, even if there is only one particle, the transition rate is still $r=1000$. Since the transition rate is probability of a reaction to take place per unit time, we cannot compare its numerical value directly to the particle number: Even if we have only one particle at that location, it is fine to have $r=1000$; it just means that the particle has a higher propensity to make a move. In addition, we remark that when there are $I_{x,y,c}$ particles, the *total transition rate* would be $r \\times I_{x,y,c}$ due to independence between the particles; in this case, it states that within time $\\Delta t \\downarrow 0$, *one* of the $I_{x,y,c}$ particles will move to its right (as well as left, up, and down) with a probability $p=1000 I_{x,y,c} \\Delta t$. \n\n- In terms of the choice of $r$, we chose it so that at time $t=1$, the spatial diffusion process ensures sufficient mixing to a random configuration (so easy to sample during inference; see for example the SSIM metric in the Appendix). We remark that $r$ and terminal time (1 in our case) are co-design choices. We heuristically designed the noise scheduler using the SSIM principle. We observed that the choice of the noise scheduler does not seem to significantly affect the performance of the model. \n\n- Yes, our model uses transitions that are made with equal probability to the four neighboring directions. This is because the particle’s transition rates to each of the directions are equal. Another way to characterize our forward process is that it is a spatially symmetric random walk. That said, we also carried-out experiments with no-flux boundary conditions (reported on Table 3) where the particles cannot transition periodically across the edge of the images. We did not identify discernible differences between the two approaches, showing the robustness of the spatial diffusion formulation. \n\n**[Comments]**\n\n- Thank you for pointing out the inconsistency in the use of the acronym. We have reviewed and corrected all occurrences to ensure consistent and accurate usage throughout the manuscript."}, {"review_id": "xx3YpKEgGg", "reviewer": "Reviewer_RTj5", "summary": "The article is about generating artificial images (or other data) by reversing a random process with a neural network. Image generation with a Gaussian diffusion reversing neural network is a popular topic. In this process, an image is generated from Gaussian noise by iteratively reversing the diffusion. In this process, the pixel values of the images are supposed to be continuous. \n\nThe article proposes a new random process in which the image (or data) has discrete values, with the constraint that the total mass (sum of all pixel values) is preserved between each random step in each color channel. This process can be applied to RGB images because they are originally uint8, but it is particularly suited for generating other types of data, such as mineral tomography. This data is binary, and controlling the amount of porosity (number of zeros) and the shape of the cavities is important so that they look natural. Thus, porosity can be controlled at initialization by the number of zeros in the initial random image.\n\nThe article presents the forward discrete mass-preserving process in detail. This process can be represented in two ways: with the random positions of individual particles or with population histograms of each cell in the image. In the forward process, each particle has an equal probability of jumping into one of its four neighboring cells or of not moving. The backward process is similar but different. Because it is conditioned by the initial position of the particle, the directions are no longer equiprobable, as the particle will be pushed to return to its initial state.\n\nAs with Gaussian diffusion, the discrete process is reversed using a network. The network is trained to estimate the probability that a particle will reverse jump in each direction. The authors propose two loss functions: one that minimizes the L1 distance between the network guess and the true backward process and one that maximizes the likelihood of the network in the backward process.\n\nIn the \"Experiment\" section, the authors demonstrate that their process can generate RGB images similar to those produced by Gaussian diffusion. Then, they present results on tomography generation, demonstrating that their method enables them to control porosity while producing realistic data.", "questions": "(0) My main concern is that the results of classical Gaussian diffusion should be provided for comparison in experiments. Even if the results are similar, this new framework will still be interesting.\n\n(1) Regarding my criticism of the notation, I think it would be clearer if the forward equations (1a, b) and backward equation (4) used probabilistic notation and definitions. For example, I think it would be better to describe the  \"noise\" $\\nu$ as categorical random variables with:\n\\begin{equation}\n\\mathbb{P}(\\nu) = \\begin{cases}\nr & \\text{if } \\nu\\in \\{(0,1),(0,-1),(1,0),(-1,0)\\};\\\\\n1-4r & \\text{if } \\nu = (0,0); \\\\\n0     & \\text{otherwise}\n\\end{cases}\n\\end{equation}\n\nIt may help to understand the likelihood in equation (6). Personally, I have trouble seeing where it comes from.\n\n(2) Table 3 should be cited in the text of the experiments so that readers know the experimental setup is provided in the supplementary material.\n\n(3) Equations (5), (6), and (7) should state the network inputs and variables on which the sums and integrals are based to make them easier to understand. For example, in equations (6) and (7), it is unclear what the sums are.\n\n(4) line 657, I think there is a typo and that $p_{i,l}$ should be replaced with $p_{k,l}$.", "limitations": "yes", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths :\n\nThe article proposes an interesting alternative to Gaussian diffusion. The random forward and backward processes are explained quite well. The authors demonstrate that this process works on natural RGB images, even though that is not their main target application. They also acknowledge the limitations, explaining that this process is not efficient for RGB images because the number of random variables is higher than in Gaussian diffusion.\n\nThe article is well-written and covers all aspects of this type of solution (random process, forward-backward, loss, sheduler ..), with extensive citation of applicative and theoretical articles.\n\nThe article also provides extensive supplementary material, including all the necessary mathematical derivations and additional experiments. \n\nIn this sense, I find that the article provides a thorough description of the proposed method.\n\nWeaknesses :\n\nEven though the article proposes many interesting experiments, I think an important one is missing. The article claims that Gaussian diffusion is not suited for tomography generation and that the proposed process is better. Therefore, I think the experiment should highlight this fact by comparing the proposed solution with tomographies generated using the Gaussian process  ( with a conditioning on porosity and relaxation from discrete/binary to continuous). \n\nIn my opinion, since the article describes a random process, it should use more common probabilistic notations and definitions. This would make it easier to compare to other methods of this kind.", "was_revised": true, "final_justification": "The article proposes a new generative method based on an original stochastic process that allows for mass conservation. This is both original and useful in certain application scenarios. The article is detailed and well explained.\n\nI just thought that, in order to demonstrate the usefulness of the method, the article needed to draw a comparison with classical Gaussian diffusion for these specific applications. The authors carried out the experiments and made the small corrections to the text that I requested. I also raised the concern that the mathematical notations were not standard, but the authors showed me that I was mistaken.\n\nI no longer have any concerns about this article. I think it's good work that should be published, which is why I changed my rating to Accept.", "rebuttal": "Thank you for the constructive feedback and for recognizing the novelty and thoroughness of our method. We believe the new experimental results and clarifications further strengthen the case for this work’s significance, particularly in scientific domains where conservation laws are important.\n\n**[Q0]**\n\nThank you for encouraging us to include a comparison with a Gaussian diffusion framework. This was a valuable suggestion, and we have carried out the experiment accordingly. We conducted the experiment using our Savonnières Carbonate sample, which displays the widest range of heterogeneous geological attributes in our dataset, including  dissolution features, fossils, and grain inclusions. The results confirm key limitations of Gaussian diffusion for our target application. We will document this experiment thoroughly in the final version, although we cannot include sample visualizations here, we will do our best to describe the quantitative and qualitative results:\n\n\n- The generated samples are noticeably “grainy” and exhibit artifacts.\n\n\n- Grain boundaries are not perfectly smooth or round (compared to our samples in Figure 4 and 15) but appear as staircase-like structures (like those in Minecraft).\n\n\n- Salt-and-pepper noise is present both inside and outside the solid regions.\n\n\n- While the porosity (i.e., the mean pixel intensity over image size) is approximately correct (especially near the mean ~0.6), this is achieved through a bimodal distribution centered near 0 and 1, with values falling outside the [0,1] range (e.g., down to –0.2 and up to 1.2). As a result, binarization does not significantly worsen the porosity error, unlike with datasets like MNIST (Figure 14).\n\n- We trained the model using data whose porosity is between 53% and 72%. We deployed the model to a wider range of target porosity between 45.74% and 81.24%. The model confidently produces incorrect predictions (low standard deviation over 100 generated samples) outside of the training distribution, indicating poor generalization outside the center of the distribution. In contrast, DSD was able to generalize, since intensity is a hard constraint in DSD. \n\n\n\n\n\n\n\n\n\n\n| **Porosity \\[%]** | **Prediction \\[Mean ± Std]** | **Relative Error \\[%]** |\n| ----------------- | ---------------------------- | ----------------------- |\n| 45.74             | 62.26 ± 0.60                 | 36.11                   |\n| 49.29             | 61.94 ± 0.45                 | 25.67                   |\n| 52.84             | 54.85 ± 0.54                 | 3.80                    |\n| 56.39             | 57.08 ± 0.52                 | 1.22                    |\n| 59.94             | 60.33 ± 0.49                 | 0.66                    |\n| 63.49             | 63.88 ± 0.53                 | 0.61                    |\n| 67.04             | 67.61 ± 0.48                 | 0.84                    |\n| 70.59             | 71.30 ± 0.55                 | 1.00                    |\n| 74.14             | 62.87 ± 0.45                 | 15.20                   |\n| 77.69             | 64.13 ± 0.52                 | 17.45                   |\n| 81.24             | 62.20 ± 0.51                 | 23.44                   |\n\n\n\n\nThese findings highlight that while the Gaussian diffusion model may approximate average porosity reasonably well near the mean, it fails to preserve structure and accuracy in the mean intensity away from this. More importantly, the lack of spatial coherence and the introduction of noise compromises the generated samples. This is especially important because transport properties in porous media are highly sensitive to such artifacts, small variations or noise can clog pores (percolation threshold) and significantly alter effective properties.\n\n**[Q1]** We thank the reviewer for the suggestion. However, respectfully, we are using standard notation. Since  $r$ is the transition rate (not a probability, but probability per unit time), we could not equate the $\\mathbb{P}(\\nu)$ to $r$ in the reviewer’s suggested notation (more about that on the answer to Reviewer tBZm's questions). Secondly, in the context of chemical reactions, $\\nu$ is often to be interpreted as stoichiometric coefficients (which are not random variables) which prescribe the change of chemical species’ populations (in our context, we can think of discrete location $(x,y,c)$ as populations of three chemical species.) We would argue that the notation in the submitted version *is* the standard notation in existing technical literature, for example, see *“Solving the chemical master equation for monomolecular reaction systems and beyond: a Doi-Peliti path integral view”* by Vastola and *“Efficient analysis of stochastic gene dynamics in the non-adiabatic regime using piecewise deterministic Markov processes”* by Lin and Buchler. As such, we would like to retain our current model specification.  To improve the clarity of the manuscript we will explain the link between the transition rates and jump probabilities in the infinitesimal time limit.  \n\n**[Q2]** We thank the reviewer for their suggestion. We have modified the sentence “The hyperparameters used can be found in Appendix H.” to “The hyperparameters used can be found in Table 3 in Appendix H.”\n\n**[Q3]** We thank the reviewer for the suggestions. We have inserted the dependence of $\\bar{\\nu}, x, y, c$ and explicitly annotated the sum and mean over these variables.\n\n**[Q4]** Thank you. We have fixed the typo."}]

Author rebuttal:

Meta-review: The paper proposed discrete spatial diffusion (DSD), based on discrete-state jump process for image generation with total intensity preserved.  The proposed method performs comparably to the standard Gaussian diffusion on general image generation (although slower), and shows clear advantages in scientific applications including porous rock microstructures and lithium-ion battery electrodes.

All reviewers acknowledged impactful contributions, and some concerns, including clarity and lack of investigation, have been fully addressed in the rebuttal.   Reviewers recommended the paper to be a spotlight, because of high novelty on an important research area.

AC would suggest the following revision.  The original version first shows comparable result to the standard diffusion model in general image generation, and then show results on the scientific applications ONLY by the proposed method.  Consequently, readers cannot find the real advantage of the proposed method.   My suggestion is to put more weight on the scientific applications, and compare the proposed method with the standard diffusion models.  Then, reader will clearly see what are the problems when the standard diffusion is used for those scientific applications, and that the proposed methods solve those problems.  I believe that this would make the story of the paper clearer, and more convincing.
Decision: Accept (spotlight)