Entering the Era of Discrete Diffusion Models: A Benchmark for Schrödinger Bridges and Entropic Optimal Transport

## Abstract

Abstract The Entropic Optimal Transport (EOT) problem and its dynamic counterpart, the Schrödinger bridge (SB) problem, play an important role in modern machine learning, linking generative modeling with optimal transport theory. While recent advances in discrete diffusion and flow models have sparked growing interest in applying SB methods to discrete domains, there is still no reliable way to evaluate how well these methods actually solve the underlying problem. We address this challenge by introducing a benchmark for SB on discrete spaces. Our construction yields pairs of probability distributions with analytically known SB solutions, enabling rigorous evaluation. As a byproduct of building this benchmark, we obtain two new SB algorithms, DLightSB and DLightSB-M, and additionally extend prior related work to construct the α \alpha -CSBM algorithm. We demonstrate the utility of our benchmark by evaluating both existing and new solvers in high-dimensional discrete settings. This work provides the first step toward proper evaluation of SB methods on discrete spaces, paving the way for more reproducible future studies.

## 1 Introduction

The Entropic Optimal Transport (Cuturi, 2013, EOT) problem and its dynamic counterpart, the Schrödinger bridge (Schrödinger, 1931, SB), have recently attracted significant attention in the machine learning community due to their relevance for generative modeling and unpaired learning. A variety of methods have been developed to solve these problems in *continuous data spaces* such as (Daniels et al., 2021; Gushchin et al., 2023a; 2024b; Mokrov et al., 2024; Vargas et al., 2021; Chen et al., 2021; Shi et al., 2023; De Bortoli et al., 2024; Korotin et al., 2024; Gushchin et al., 2024a).

At the same time, much real world data are *discrete by nature*, including text (Austin et al., 2021; Gat et al., 2024), molecular graphs (Vignac et al., 2022; Qin et al., 2024; Luo et al., 2024), and protein sequences (Campbell et al., 2024). Others are *discrete by construction*, such as vector-quantized representations of images and audio (Van Den Oord et al., 2017; Esser et al., 2021).

Given the prevalence of such discrete data and the rapid progress in discrete diffusion/flow models (Hoogeboom et al., 2021; Austin et al., 2021; Campbell et al., 2022; Lou et al., 2023; Sahoo et al., 2024; Campbell et al., 2024; Gat et al., 2024), research on SBs has attracted growing attention in recent years. For instance, several recent works have already taken first steps in this direction (Kim et al., 2024, DDSBM;Ksenofontov & Korotin, 2025, CSBM), adapting diffusion methodologies from (Austin et al., 2021, D3PM;Vignac et al., 2022, DiGress), respectively.

The development of rigorous benchmarks is crucial for proper evaluation. These benchmarks enable us to determine whether SB methods actually solve the intended mathematical problem, separating true algorithmic performance from artifacts of specific parameterizations, regularization schemes, and other implementation choices. While this has recently become possible in the continuous setting of Schrodinger Bridges (Gushchin et al., 2023b), no such approach exists for discrete data, leaving it unclear how closely SB solvers approximate the true solution of the SB problem on discrete domains. To address this gap, we make the following contributions:

- •
Theory & Methodology. We introduce a method for constructing benchmark pairs of probability distributions on discrete spaces with analytically known SB solution by design (\wasyparagraph[3](https://arxiv.org/html/2509.23348v1#S3)).
- •
Algorithms. The introduction of our benchmark, combined with insights from the continuous-space LightSB methods (Korotin et al., 2024; Gushchin et al., 2024a), allows us to instantly derive two new discrete space SB solvers: DLightSB and DLightSB-M (\wasyparagraph[4.3](https://arxiv.org/html/2509.23348v1#S4.SS3) and \wasyparagraph[4.4](https://arxiv.org/html/2509.23348v1#S4.SS4)). Additionally, we introduce a new solver, $\alpha$-CSBM (\wasyparagraph[4.2](https://arxiv.org/html/2509.23348v1#S4.SS2)), by extending the recent CSBM method (Ksenofontov & Korotin, 2025) with the $\alpha$-DSBM principle used in continuous settings (De Bortoli et al., 2024).
- •
Practice. We use these benchmark pairs to evaluate both existing and newly introduced SB solvers in high-dimensional settings

#### Notation.

We consider a discrete state space $\mathcal{X}=\mathbb{S}^{D}$, where $\mathbb{S}=\{0,1,\dots,S-1\}$ is the set of $S$ categories and $D$ is the dimensionality. Each $x\in\mathcal{X}$ is a $D$-dimensional vector $x=(x^{1},\dots,x^{D})$. Time is discretized as $\{t_{n}\}_{n=0}^{N+1}$ with $0=t_{0}<t_{1}<\dots<t_{N}<t_{N+1}=1$. This gives $N+2$ time points and defines the *path space* $\mathcal{X}^{N+2}$ with the tuple $x_{\text{in}}\stackrel{{\scriptstyle\text{def}}}{{=}}(x_{t_{1}},\dots,x_{t_{N}})\!\in\!\mathcal{X}^{N}$ collecting the intermediate states. The set $\mathcal{P}(\mathcal{X}^{N+2})$ comprises all discrete time stochastic processes on the path space, with $\mathcal{M}(\mathcal{X}^{N+2})\!\subset\!\mathcal{P}(\mathcal{X}^{N+2})$ denoting the subset of *Markov processes*. Any $q\!\in\!\mathcal{M}(\mathcal{X}^{N+2})$ admits forward and backward representations: $q(x_{0},x_{\text{in}},x_{1})=q(x_{0})\prod_{n=1}^{N+1}q(x_{t_{n}}|x_{t_{n-1}})=q(x_{1})\prod_{n=1}^{N+1}q(x_{t_{n-1}}|x_{t_{n}})$, where $q(\cdot|\cdot)$ denotes conditional probabilities.

## 2 Background: Problem Statement

We establish the theoretical foundation for our benchmark by recalling discrete dynamic and static SB and their relationship. First, we present the dynamic SB and its reduction to a static problem (\wasyparagraph[2.1](https://arxiv.org/html/2509.23348v1#S2.SS1)). Next, we analyze diffusion-type reference processes (\wasyparagraph[2.2](https://arxiv.org/html/2509.23348v1#S2.SS2)) that yield practical cost functions, linking to the EOT framework in \wasyparagraph[2.3](https://arxiv.org/html/2509.23348v1#S2.SS3). Finally, we introduce our problem setting (\wasyparagraph[2.4](https://arxiv.org/html/2509.23348v1#S2.SS4)).

### 2.1 Dynamic and Static Schrödinger Bridges on Discrete Spaces

#### Dynamic Schrödinger Bridge.

The original SB problem (Schrödinger, 1931; 1932; Léonard, 2013) seeks to find a process $q^{*}\in\mathcal{P}(\mathcal{X}^{N+2})$ interpolating between an initial distribution $p_{0}$ at $t_{0}=0$ and a final distribution $p_{1}$ at $t_{N+1}=1$. This distribution is found by minimizing the Kullback-Leibler (KL) divergence with respect to a given Markov reference process $q^{\text{ref}}\in\mathcal{M}(\mathcal{X}^{N+2})$ under the following marginal constraints $p_{0}(x_{0})=q(x_{0})$ and $p_{1}(x_{1})=q(x_{1})$. One finds

$$ $q^{*}=\operatorname*{arg\,min}_{q\in\Pi_{N}(p_{0},p_{1})}\text{KL}\left(q(x_{0},x_{\text{in}},x_{1})\|q^{\text{ref}}(x_{0},x_{\text{in}},x_{1})\right),$ (1) $$

where $\Pi_{N}(p_{0},p_{1})\subset\mathcal{P}(\mathcal{X}^{N+2})$ denotes the subset of $\mathcal{X}$-valued stochastic processes which have $p_{0}$ and $p_{1}$ as marginals at times $t_{0}=0$ and $t_{N+1}=1$, respectively. In other words, the dynamic SB problem seeks the stochastic process $q^{*}$ that minimally deviates from a reference process $q^{\text{ref}}$ while respecting the boundary distributions $p_{0}$ and $p_{1}$.

#### Static Schrödinger Bridge.

We now introduce the static formulation of the SB. This begins with observing that ([1](https://arxiv.org/html/2509.23348v1#S2.E1)) admits the following decomposition:

$$ $\min_{q\in\Pi_{N}(p_{0},p_{1})}\Bigr[\text{KL}\left(q(x_{0},x_{1})\|q^{\text{ref}}(x_{0},x_{1})\right)+\mathbb{E}_{q(x_{0},x_{1})}\text{KL}\left(q(x_{\text{in}}|x_{0},x_{1})\|q^{\text{ref}}(x_{\text{in}}|x_{0},x_{1})\right)\Bigr].$ (2) $$

We further note that the optimal SB $q^{*}$ lies in the *reciprocal class* of the reference process $q^{\text{ref}}$, denoted $\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2})\subset\mathcal{P}(\mathcal{X}^{N+2})$. A reciprocal process is any $q\in\mathcal{P}(\mathcal{X}^{N+2})$ whose intermediate dynamics are given by the reference process $q^{\text{ref}}$, i.e., $q(x_{\text{in}}|x_{0},x_{1})=q^{\text{ref}}(x_{\text{in}}|x_{0},x_{1})$. Thus, restricting $q$ to $\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2})$ makes the conditional KL term vanish, leaving only the first term in the optimization. This yields the static SB problem, defined as

$$ $q^{*}(x_{0},x_{1})=\operatorname*{arg\,min}_{q\in\Pi(p_{0},p_{1})}\text{KL}\left(q(x_{0},x_{1})\|q^{\text{ref}}(x_{0},x_{1})\right),$ (3) $$

where $\Pi(p_{0},p_{1})\in\mathcal{P}(\mathcal{X}^{2})$ is the set of joint distributions $q(x_{0},x_{1})$ whose marginals are $p_{0}$ and $p_{1}$.

### 2.2 Formulations of the Reference Process

The key ingredient in both SB formulations is the reference process $q^{\text{ref}}$, which can be defined in either continuous time ($N=\infty$) or discrete time ($N<\infty$). In continuous time, $q^{\text{ref}}$ is specified through transition rates. In practice, modifying these rates is less direct than probabilities, and although any continuous-time process can be discretized, the converse does not hold, meaning discrete time processes form a strictly richer class. For these reasons, we focus on the discrete setting, which is well-suited for our benchmark design.

The reference process $\textstyle q^{\text{ref}}\in\mathcal{M}(\mathcal{X}^{N+2})$ is modeled as a *discrete-state diffusion process*, i.e., a discrete-time Markov chain defined by transition matrices $\textstyle Q_{n}^{\text{ref}}\in[0,1]^{|\mathcal{X}|\times|\mathcal{X}|}$, where $\textstyle q^{\text{ref}}(x_{t_{n}}^{d}|x_{t_{n-1}}^{d})=[Q_{n}^{\text{ref}}]_{x_{t_{n-1}}^{d},x_{t_{n}}^{d}}$. Assuming time-homogeneity ($\textstyle Q_{n}^{\text{ref}}=Q^{\text{ref}}$ for all $\textstyle n$), the $\textstyle n$-step transition probabilities are given by the matrix power $\textstyle\overline{Q}_{n}^{\text{ref}}=[Q^{\text{ref}}]^{n}$. To define $\textstyle Q$, we further restrict to $\textstyle D=1$ for clarity, noting that for $\textstyle D>1$ the transition probabilities are obtained as a product over dimensions. We now introduce two diffusion-like transitions: uniform (Hoogeboom et al., 2021; Campbell et al., 2022) and Gaussian-like (Austin et al., 2021).

#### Uniform Reference Process ( q unif q^{\text{unif}} ).

For unordered data, where no relation exists between categories, a natural choice is a so-called uniform transition matrix. For each dimension $d$, the elements of the transition matrix $Q^{\text{ref}}$ are defined by

$$ $\displaystyle[Q^{\text{ref}}]_{x^{d}_{t_{n-1}},x^{d}_{t_{n}}}=\begin{cases}1-\gamma,&\text{if }x^{d}_{t_{n}}=x^{d}_{t_{n-1}},\\ \frac{\gamma}{S-1},&\text{if }x_{t_{n}}^{d}\neq x_{t_{n-1}}^{d},\end{cases}$ (4) $$

where $\gamma\in[0,1]$ is an stochasticity parameter. This reference process introduces randomness independently of the distance between categories, assigning equal probability to transitions into any category and disregarding any inherent ordering or relationships.

#### Gaussian Reference Process ( q gauss q^{\text{gauss}} ).

For ordered data, where categories are expected to exhibit meaningful relations, a Gaussian-like transition matrix is more appropriate. With the stochasticity parameter $\gamma>0$ and the maximum category distance $\Delta=S-1$, the transition probabilities are

$$ $[Q^{\text{ref}}]_{x^{d}_{t_{n-1}},x^{d}_{t_{n}}}=\frac{\exp\!\left(-\frac{4(x^{d}_{t_{n}}-x^{d}_{t_{n-1}})^{2}}{(\gamma\Delta)^{2}}\right)}{\sum\limits_{\delta=-\Delta}^{\Delta}\exp\!\left(-\frac{4\delta^{2}}{(\gamma\Delta)^{2}}\right)},\qquad x^{d}_{t_{n}}\neq x^{d}_{t_{n-1}}.$ (5) $$

The diagonal entries take the remaining probability so that each row sums to 1.

### 2.3 Entropic Optimal Transport on Discrete Spaces

The proper construction of a Markov reference process (§[2.2](https://arxiv.org/html/2509.23348v1#S2.SS2)) allows us to formulate the static SB problem (§[3](https://arxiv.org/html/2509.23348v1#S2.E3)) in a manner that is equivalent to the entropic optimal transport (EOT) problem (Cuturi, 2013). To be precise, if we consider $q^{\text{ref}}(x_{0},x_{1})=q^{\text{ref}}(x_{0})q^{\text{ref}}(x_{1}|x_{0})$ and put $q^{\text{ref}}(x_{0})=p_{0}(x_{0})$,
the minimization problem in equation ([3](https://arxiv.org/html/2509.23348v1#S2.E3)) can be rewritten as

$$ $\displaystyle\min_{q\in\Pi(p_{0},p_{1})}$ $\displaystyle\text{KL}\left(q(x_{0},x_{1})\|q^{\text{ref}}(x_{0},x_{1})\right)=$ $\displaystyle=\min_{q\in\Pi(p_{0},p_{1})}\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log\frac{q(x_{0},x_{1})}{q^{\text{ref}}(x_{0})q^{\text{ref}}(x_{1}|x_{0})}$ $\displaystyle=\min_{q\in\Pi(p_{0},p_{1})}-H(q)-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})-\!\!\!\!\!\!\underbrace{\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{ref}}(x_{0})}_{=\sum_{x_{0}}p_{0}(x_{0})\log p_{0}(x_{0})=-H(p_{0})}$ (6) $\displaystyle=\min_{q\in\Pi(p_{0},p_{1})}\mathbb{E}_{q(x_{0},x_{1})}\!\bigl[-\log q^{\text{ref}}(x_{1}|x_{0})\bigr]-H(q)-\text{const}$ $\displaystyle=\min_{q\in\Pi(p_{0},p_{1})}\mathbb{E}_{(x_{0},x_{1})\sim q}\!\bigl[c(x_{0},x_{1})\bigr]-H(q)-\text{const},$ $$

where $H(q)$ is the entropy of $q$, while $H(p_{0})$ remains constant when minimizing over $q$. Thus, the static SB formulation becomes equivalent to the entropy-regularized optimal transport problem with cost $c(x_{0},x_{1})=-\log q^{\text{ref}}(x_{1}|x_{0})$. This perspective establishes a direct correspondence between SB and EOT, which we use in the design of our benchmark and methodological framework in \wasyparagraph[3](https://arxiv.org/html/2509.23348v1#S3).

Since the conditional distribution $q^{\text{ref}}(x_{1}|x_{0})$ is obtained by taking the $(N\!+\!1)$-th power of $Q^{\text{ref}}$, it admits the following closed-form expression in the uniform case:

$$ $\overline{Q}^{\text{ref}}_{N+1}=\left(1-\gamma\frac{S}{S-1}\right)^{N+1}\mathbb{I}+\frac{1-\left(1-\gamma\frac{S}{S-1}\right)^{N+1}}{S}\mathbf{1}\mathbf{1}^{\top},$ (7) $$

where $\mathbf{1}=[1,\dots,1]^{\top}\in\mathbb{R}^{S}$ is a vector full of ones. From here it can be seen that $\overline{Q}^{\text{ref}}_{N+1}$ converges to $(1/S)\mathbf{1}\mathbf{1}^{\top}$ when $(N\!+\!1)\rightarrow\infty$, that is a uniform distribution over the number of categories $S$, the derivation of ([7](https://arxiv.org/html/2509.23348v1#S2.E7)) can be found in Appendix [A](https://arxiv.org/html/2509.23348v1#A1). In the case of the Gaussian reference process, the closed-form expression can also be obtained, but it is much more complex.

### 2.4 Problem Setup for Discrete Schrödinger Bridges

In this section, we recall the *generative SB task on discrete spaces*, a well-established problem in the SB and OT literature (Kim et al., 2024; Ksenofontov & Korotin, 2025). In short, the goal is to learn an SB process or coupling that performs transport between probability distributions on discrete spaces using available empirical data samples. Formally, we consider the following learning setup:

Despite recent progress in the development of SB methods that solve this task, there remains no standard methodology for performance evaluation, mainly due to the absence of ground-truth distribution pairs $(p_{0},p_{1})$. In this work, we propose a benchmark construction, inspired by (Gushchin et al., 2023b), that enables standard evaluation of such methods on datasets built from SB pairs $(x_{0},x_{1})$ with known $q^{*}(x_{1}|x_{0})$. Such datasets provide more informative metrics and offer a consistent framework for assessing the performance of SB methods on discrete spaces.

Remark. Our paper is not related to the discrete EOT, which includes solvers such as the Sinkhorn algorithm (Cuturi, 2013) or gradient-based methods (Dvurechensky et al., 2018). These approaches are designed for a non-generative problem setting, see (Ksenofontov & Korotin, 2025, \wasyparagraph2.3). They treat samples as empirical distributions $p_{0}(x_{0})\!=\!\tfrac{1}{|I_{0}|}\!\sum_{i\in I_{0}}\delta_{x_{0}^{(i)}}$, $p_{1}(x_{1})\!=\!\tfrac{1}{|I_{1}|}\!\sum_{j\in I_{1}}\delta_{x_{1}^{(j)}}$. The resulting coupling is then a bi-stochastic $|I_{0}|\times|I_{1}|$ matrix, which does not support out-of-sample generation. While some extensions attempt to provide inference for unseen data (Hütter & Rigollet, 2021; Pooladian & Niles-Weed, 2021; Manole et al., 2024; Deb et al., 2021), they are designed for continuous spaces ($\mathcal{X}=\mathbb{R}^{D}$) rather than the discrete spaces ($\mathcal{X}=\mathbb{S}^{D}$) considered in our work.

## 3 Benchmark

This section outlines our theoretical and practical foundations necessary for constructing the benchmark for the SB. We introduce our benchmark construction in \wasyparagraph[3.1](https://arxiv.org/html/2509.23348v1#S3.SS1). This theoretical result requires a specific parameterization, which we explore in \wasyparagraph[3.2](https://arxiv.org/html/2509.23348v1#S3.SS2). This construction and parameterization is lately used to construct our high-dimensional Gaussian mixture benchmark \wasyparagraph[3.3](https://arxiv.org/html/2509.23348v1#S3.SS3).

### 3.1 Main Theorem for Benchmark Construction

Given an initial distribution $p_{0}\in\mathcal{P}(\mathcal{X})$, we aim to construct a target distribution $p_{1}\in\mathcal{P}(\mathcal{X})$ such that the static SB $q^{*}(x_{0},x_{1})$ between them is known by our construction. The resulting pair $(p_{0},p_{1})$ together with $q^{*}$ can then be used as benchmark data for evaluating SB methods. Our following theorem plays the key role in the construction of benchmark pairs.

###### Theorem 3.1 (Benchmark pair construction for SB on discrete Spaces) .

Let $p_{0}\in\mathcal{P}(\mathcal{X})$ be a given source distribution on a discrete space $\mathcal{X}$ and $v^{*}:\mathcal{X}\rightarrow\mathbb{R}$ be a given scalar-valued function. Let $q^{*}\in\mathcal{P}(\mathcal{X}^{2})$ be a joint distribution for which for all $x_{0}\in\mathcal{X}$ it holds that $q^{*}(x_{0})=p_{0}(x_{0})$ and

$$ $q^{*}(x_{1}|x_{0})\propto v^{*}(x_{1})q^{\textup{ref}}(x_{1}|x_{0}),$ (8) $$

Let $p_{1}\!\in\!\mathcal{P}(\mathcal{X})$ be the second marginal of $q^{*}$, i.e., $q^{*}(x_{1})\stackrel{{\scriptstyle\text{def}}}{{=}}p_{1}(x_{1})$. Then $q^{*}(x_{0},x_{1})$ is the static SB ([3](https://arxiv.org/html/2509.23348v1#S2.E3)) between $p_{0}$ and $p_{1}$. In turn, $q^{*}(x_{0},x_{\textup{in}},x_{1})\stackrel{{\scriptstyle\text{def}}}{{=}}q^{*}(x_{0},x_{1})q^{\textup{ref}}(x_{\textup{in}}|x_{0},x_{1})$ is the dynamic SB ([1](https://arxiv.org/html/2509.23348v1#S2.E1)).

Theorem [3.1](https://arxiv.org/html/2509.23348v1#S3.Thmtheorem1) establishes that any pair $(p_{0},v^{*})$ can be used to construct $(p_{0},p_{1})$ for the SB problem, thereby yielding a known solution $q^{*}$. The construction considers conditional distributions $q^{*}(x_{1}|x_{0})$ in an unnormalized form, so we further write

$$ $q^{*}(x_{1}|x_{0})=\frac{1}{c^{*}(x_{0})}v^{*}(x_{1})q^{\text{ref}}(x_{1}|x_{0}),$ (9) $$

where $c^{*}(x_{0})\stackrel{{\scriptstyle\text{def}}}{{=}}\sum_{x_{1}\in\mathcal{X}}v^{*}(x_{1})q^{\text{ref}}(x_{1}|x_{0})$ is the normalization constant.

Our benchmark construction idea may be non-trivial to implement in practice. Specifically, working in the high-dimensional space $\mathcal{X}=\mathbb{S}^{D}$ makes computing the normalization constant and sampling from $q^{*}$ computationally expensive. To address this, we introduce a parameterization that enables efficient computation and sampling, as detailed in the next section.

###### Theorem 3.1 (Benchmark pair construction for SB on discrete Spaces) .

### 3.2 Practical Parameterization of the Scalar-valued Function V ∗ V^{*}

We parameterize the scalar-valued function $v^{*}$ using a rank-1 Canonical Polyadic (CP) decomposition, which captures interactions across dimensions and provides a compact yet expressive representation. Such decompositions act as universal approximators, capable of modeling complex functions when the rank is sufficiently large (Cohen et al., 2016; Basharin et al., 2025). Thus, $v^{*}$ is given by

$$ $v^{*}(x_{1})\;=\;\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}r_{k}^{d}[x_{1}^{d}].$ (10) $$

Expression ([10](https://arxiv.org/html/2509.23348v1#S3.E10)) defines a mixture of $K$ *factorizable distributions*, each with weight $\beta_{k}\geq 0$. For each mixture component $k$ and dimension $d$, probabilities are given by non-negative vectors $r_{k}^{d}\in\mathbb{R}^{S}_{+}$, referred to as *CP cores*, where $r_{k}^{d}[x_{1}^{d}]$ denotes the probability of state $x_{1}^{d}$. The key advantage of this parameterization is that the factorization across dimensions makes both the normalizing constant $c(x_{0})$ and the conditional distribution $q^{*}(x_{1}|x_{0})$ computationally tractable. Specifically, the product structure allows efficient ancestral sampling by drawing each dimension independently.

###### Proposition 3.1 (Tractable Parameterization of Conditional Distributions) .

Given the CP decomposition of the scalar-valued function $v(x_{1})\!=\!\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}r_{k}^{d}[x_{1}^{d}]$ and a factorizable reference process $q^{\textup{ref}}(x_{1}|x_{0})\!=\!\prod_{d=1}^{D}q^{\textup{ref}}(x_{1}^{d}|x_{0})$, the optimal conditional distribution satisfies:

$$ $q^{*}(x_{1}|x_{0})\!=\!\tfrac{1}{c(x_{0})}\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}\Big[r_{k}^{d}[x_{1}^{d}]q^{\textup{ref}}(x_{1}^{d}|x_{0})\Big];\quad$ (11) $$

$$ $c(x_{0})\!=\!\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}\!\left(\!\sum_{x_{1}^{d}=0}^{S-1}r_{k}^{d}[x_{1}^{d}]q^{\textup{ref}}(x_{1}^{d}|x_{0})\!\right)$ (12) $$

where $c(x_{0})$ is the normalization constant. This formulation expresses $q^{*}(x_{1}|x_{0})$ as a mixture of $K$ factorizable distributions, each weighted by a scalar coefficient $\beta_{k}$.

###### Proposition 3.1 (Tractable Parameterization of Conditional Distributions) .

### 3.3 High-dimensional Gaussian Mixtures Benchmark Construction

We set $p_{0}$ as a noise distribution (uniform or discretized Gaussian) on $D\in\{2,16,64\}$ dimensions with $S=50$ categories. For $v^{*}$, we use $K=4$ components with uniformly initialized weights $\beta\in\mathbb{R}^{K}$, and the CP cores are initialized by setting their logarithms to the log-density of discretized Gaussians with varying means and fixed variance.
Given $p_{0}$ and $v^{*}$, we then construct $p_{1}$ (Theorem [3.1](https://arxiv.org/html/2509.23348v1#S3.Thmtheorem1)). This initialization produces a target $p_{1}$ resembling a discretized Gaussian mixture with a clear visual structure. Moreover, our benchmark formulation further allows the generation of an unlimited number of samples for training.

We construct pairs under different reference processes $q^{\text{ref}}$: Gaussian $q^{\text{gauss}}$ with $\gamma\in\{0.02,0.05\}$ and uniform $q^{\text{unif}}$ with $\gamma\in\{0.005,0.01\}$, using $N+1=128$ for both, see Figure [1(b)](https://arxiv.org/html/2509.23348v1#S5.F1.sf2) to visualize ground truth benchmark pairs.

## 4 Solvers for Evaluation

The field of discrete SB solvers remains in early development, with limited methods available for evaluation. We assess four approaches: the *Categorical Schrödinger Bridge Matching (CSBM)* method (Ksenofontov & Korotin, 2025), designed specifically for categorical distributions; our *$\alpha$-CSBM* extension, which applies the online methodology of (De Bortoli et al., 2024) to CSBM; new *Discrete Light Schrödinger Bridge (DLightSB)* solver, constructed using our benchmark framework (§[3](https://arxiv.org/html/2509.23348v1#S3)) and adapting ideas from (Korotin et al., 2024) to discrete settings; and finally new *DLightSB-M*, which extends DLightSB to dynamic setups following (Gushchin et al., 2024a). Further details about methods can be found in Appendix [B](https://arxiv.org/html/2509.23348v1#A2).

### 4.1 Categorical Schrödinger Bridge Matching (CSBM)

In (Ksenofontov & Korotin, 2025, Theorem 3.1), the discrete space SB problem is addressed by extending the discrete-time existence theorem of (Gushchin et al., 2024b, Theorem 3.6) to the discrete space/time setting, thereby establishing convergence of the *discrete time Iterative Markovian Fitting (D-IMF) procedure*. This constructive method uses the fact that the dynamic SB $q^{*}$ is both reciprocal ($q^{*}\in\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2})$) and Markov ($q^{*}\in\mathcal{M}(\mathcal{X}^{N+2})$). The D-IMF algorithm alternates between projections onto these two sets, starting from an initial process $q^{0}(x_{0},x_{1})q^{\text{ref}}(x_{\text{in}}|x_{0},x_{1})$, where $q^{0}(x_{0},x_{1})\in\Pi(p_{0},p_{1})$, e.g., $p_{0}(x_{0})p_{1}(x_{1})$, and converges to the SB $q^{*}$ in KL. Namely,

$$ $q^{2l}$ $q^{2l+2}$ $\operatorname{proj}_{\mathcal{M}}$ $\operatorname{proj}_{\mathcal{R}^{\textup{ref}}}$ $l=0,1,\dots$ $$

where

$$ $\displaystyle[\operatorname{proj}_{\mathcal{R}^{\textup{ref}}}(q)](x_{0},x_{\text{in}},x_{1})=\!\!\!\!\operatorname*{arg\,min}_{r\in\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2})}\!\!\!\!\text{KL}\left(q(x_{0},x_{\text{in}},x_{1})\|r(x_{0},x_{\text{in}},x_{1})\right),\quad\forall q\in\mathcal{P}(\mathcal{X}^{N+2}),$ (13) $\displaystyle\!\!\![\operatorname{proj}_{\mathcal{M}}(q)](x_{0},x_{\text{in}},x_{1})=\!\!\!\!\operatorname*{arg\,min}_{m\in\mathcal{M}(\mathcal{X}^{N+2})}\!\!\!\!\text{KL}\left(q(x_{0},x_{\text{in}},x_{1})\|m(x_{0},x_{\text{in}},x_{1})\right),\quad\!\!\forall q\in\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2}).$ (14) $$

#### Loss.

Because ancestral sampling makes the reciprocal part straightforward, the challenge lies in the Markov projection, for which the authors propose minimizing an alternative objective function.

$$ $\text{KL}\left(q(x_{0},x_{\text{in}},x_{1})\|m(x_{0},x_{\text{in}},x_{1})\right)=\mathbb{E}_{q(x_{0},x_{1})}\Bigg[\sum_{n=1}^{N}\mathbb{E}_{q^{\textup{ref}}(x_{t_{n-1}}|x_{0},x_{1})}\[-5.69054pt] \textup{KL}\!\left(q^{\textup{ref}}(x_{t_{n}}|x_{t_{n-1}},x_{1})\,\|\,m(x_{t_{n}}|x_{t_{n-1}})\right)-\mathbb{E}_{q^{\textup{ref}}(x_{t_{N}}|x_{0},x_{1})}\!\left[\log m(x_{1}|x_{t_{N}})\right]\Bigg].$ (15) $$

In practice, the D-IMF procedure is implemented in a bidirectional manner. That is, it first applies the Markovian projection using both forward and backward representations.

#### Remark.

A continuous-time IMF was introduced in the Discrete Diffusion Schrödinger Bridge Matching (Kim et al., 2024, DDSBM) paper, which performs the Markovian projection ([14](https://arxiv.org/html/2509.23348v1#S4.E14)) by matching the generator matrices of continuous-time Markov chains. As it reduces to the same loss and inference process due to the neccesity to discretize time, we report results only for CSBM.

### 4.2 α \alpha -Categorical Schrödinger Bridge Matching ( α \alpha -CSBM)

We adopt *the $\alpha$-IMF* procedure (De Bortoli et al., 2024, Equation 4), defined in discrete time as

$$ $q^{2l+3}=(1-\alpha)q^{2l+1}+\alpha\operatorname{proj}_{\mathcal{R}^{\textup{ref}}}\left(\operatorname{proj}_{\mathcal{M}}(q^{2l+1})\right).$ (16) $$

This procedure can be viewed as a relaxed formulation, where the exact projections [13](https://arxiv.org/html/2509.23348v1#S4.E13) and [14](https://arxiv.org/html/2509.23348v1#S4.E14) are replaced by partial updates (see (De Bortoli et al., 2024, Equation 9)). This avoids iterating each projection to convergence while still moving the distribution toward the double projection $\operatorname{proj}_{\mathcal{R}^{\textup{ref}}}(\operatorname{proj}_{\mathcal{M}}(\cdot))$, whose unique fixed point is the SB (De Bortoli et al., 2024, Theorem 3.1). Notably, it is closely related to an independently developed formulation in (Peluchetti, 2024). We therefore treat the discrete-time case in ([16](https://arxiv.org/html/2509.23348v1#S4.E16)) as a heuristic counterpart to the original procedure.

Loss. Practically, this extends the CSBM bidirectional setup (\wasyparagraph[4.1](https://arxiv.org/html/2509.23348v1#S4.SS1)) by updating forward and backward models jointly, with the loss computed for both parameterizations as:

$$ $L(\overrightarrow{m},\overleftarrow{m})=\tfrac{1}{2}\Big(\text{KL}\left(\overrightarrow{r_{\mathrm{sg}}}(x_{0},x_{\text{in}},x_{1})\|\overleftarrow{m}(x_{0},x_{\text{in}},x_{1})\right)\[-8.53581pt] +\text{KL}\left(\overleftarrow{r_{\mathrm{sg}}}(x_{0},x_{\text{in}},x_{1})\|\overrightarrow{m}(x_{0},x_{\text{in}},x_{1})\right)\Big),$ (17) $$

where $\rightarrow$ and $\leftarrow$ denote the direction of representations (forward and backward, respectively), and $r_{\mathrm{sg}}$ denotes $\operatorname{proj}_{\mathcal{R}^{\textup{ref}}}(m)$ evaluated with the stop-gradient operation.

### 4.3 Discrete Light Schrödinger Bridge (DLightSB)

We introduce a new solver for discrete spaces, derived as a byproduct of our benchmark design. Building on the CP parametrization ([10](https://arxiv.org/html/2509.23348v1#S3.E10)), our approach adapts the continuous solver of (Korotin et al., 2024) to discrete settings by optimizing a surrogate loss that directly learns the parameters $\{\beta_{k},r_{k}^{d}\}$ of the parameterized scalar-valued function $v_{\theta}$. Specifically, we minimize the KL divergence $\text{KL}\left(q^{*}\|q_{\theta}\right)$ between the true SB $q^{*}$ and our parametric model $q_{\theta}$, such that $q_{\theta}(x_{1}|x_{0})\stackrel{{\scriptstyle\text{def}}}{{=}}\frac{1}{c_{\theta}(x_{0})}v_{\theta}(x_{1})q^{\textup{ref}}(x_{1}|x_{0})$.

#### Loss.

Following (Korotin et al., 2024), we derive a discrete surrogate objective $\text{KL}\left(q^{*}\|q_{\theta}\right)$.

###### Proposition 4.1 (Feasible Discrete Reformulation of the KL Minimization.) .

For the characterization ([9](https://arxiv.org/html/2509.23348v1#S3.E9)) of $q(x_{1}|x_{0})$, it holds that the alternative KL objective $\text{KL}\left(q^{*}\|q\right)$ admits the representation $\text{KL}\left(q^{*}\|q_{\theta}\right)=\mathcal{L}(\theta)-\mathcal{L}^{*}$ where

$$ $\mathcal{L}(\theta)=\sum_{x_{0}\in\mathcal{X}}\log c_{\theta}(x_{0})p_{0}(x_{0})-\sum_{x_{1}\in\mathcal{X}}\log v_{\theta}(x_{1})p_{1}(x_{1}),$ (18) $$

and $\mathcal{L^{*}}\in\mathbb{R}$ is a constant value not depending on $\theta$, therefore, it can be omitted.

###### Proposition 4.1 (Feasible Discrete Reformulation of the KL Minimization.) .

### 4.4 Discrete Light Schrödinger Bridge Matching (DLightSB-M)

Inspired by (Gushchin et al., 2024a), we propose a matching approach for solving the SB problem under our parametrization in \wasyparagraph[3.2](https://arxiv.org/html/2509.23348v1#S3.SS2). This approach enables obtaining the SB in a single step, which is referred to as the *optimal projection*. Specifically, we restate the Markovian projection ([14](https://arxiv.org/html/2509.23348v1#S4.E14)) as the projection of a reciprocal process $q\in\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2})$ onto the set of all SBs:

$$ $\displaystyle{\mathcal{S}}(\mathcal{X}^{N+2})\stackrel{{\scriptstyle\text{def}}}{{=}}\Big\{q^{\text{SB}}\in\mathcal{P}(\mathcal{X}^{N+2}):\exists q^{\text{SB}}_{0},q^{\text{SB}}_{1}\in\mathcal{P}(\mathcal{X})$ (19) $\displaystyle\text{for which}\ q^{\text{SB}}=\operatorname*{arg\,min}_{q^{\text{SB}}\in{\mathcal{S}}(\mathcal{X}^{N+2})}\text{KL}\left(q\|q^{\text{SB}}\right)\ \text{and}\ q\in\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2})\Big\}.$ (20) $$

We show that (Gushchin et al., 2024a, Theorem 3.1) can be generalized to an arbitrary reference process $q^{\text{ref}}$, thereby enabling the application of the optimal projection in discrete space settings.

###### Proposition 4.2 (Optimal Projection with an Arbitrary Reference Process) .

Let $q\in\mathcal{R}^{\textup{ref}}(\mathcal{X}^{N+2})$ be a reciprocal process defined with a reference process $q^{\textup{ref}}\in\mathcal{M}(\mathcal{X}^{N+2})$ and a joint distribution $q(x_{0},x_{1})\in\Pi(p_{0},p_{1})$. Then, the optimal projection of $q$ onto the set of all SBs ${\mathcal{S}}(\mathcal{X}^{N+2})$ is the SB $q^{*}$ between the desired marginals $p_{0}$ and $p_{1}$, i.e.,

$$ $q^{*}=\operatorname*{arg\,min}_{q^{\textup{SB}}\in{\mathcal{S}}(\mathcal{X}^{N+2})}\text{KL}\left(q\|q^{\textup{SB}}\right).$ (21) $$

The main requirement is to define $q^{\text{SB}}$ such that the minimization is restricted to $q^{\text{SB}}\in\mathcal{S}(\mathcal{X}^{N+2})$. The following proposition establishes this characterization of SB transitions and, through its CP cores $r^{d}_{k}$, directly connects this approach to DLightSB (\wasyparagraph[4.3](https://arxiv.org/html/2509.23348v1#S4.SS3)).

###### Proposition 4.3 (The SB’s Transition Distributions with CP Decomposition) .

Let $q^{\textup{ref}}$ be a reference Markov process on a discrete space $\mathcal{X}$ with transition matrix $Q^{\textup{ref}}$. Using the CP decomposition of the scalar-valued function $v^{*}$ ([10](https://arxiv.org/html/2509.23348v1#S3.E10)), the marginal transition distributions of the SB are given by

$$ $q^{\textup{SB}}(x^{d}_{t_{n}}|x_{t_{n-1}})=q^{\textup{ref}}(x^{d}_{t_{n}}|x_{t_{n-1}})\sum_{k=1}^{K}\beta_{k}u^{d}_{k,t_{n}}\left[x^{d}_{t_{n}}\right]\prod_{\begin{subarray}{c}j=1\\ j\neq d\end{subarray}}^{D}u^{j}_{k,t_{n-1}}\left[x^{j}_{t_{n-1}}\right],$ (22) $$

where $u^{d}_{k,t_{n}}\left[x^{d}_{t_{n}}\right]=\sum_{x^{d}_{1}}[\overline{Q}^{\textup{ref}}_{\,N+1-n}\bigr]_{x^{d}_{t_{n}},\,x_{1}^{d}}r^{d}_{k}\left[x^{d}_{1}\right]$. Sampling is done via ancestral sampling.

Loss. The loss ([15](https://arxiv.org/html/2509.23348v1#S4.E15)) could be applied directly to train the SB transitions $q^{\text{SB}}$.

###### Proposition 4.2 (Optimal Projection with an Arbitrary Reference Process) .

###### Proposition 4.3 (The SB’s Transition Distributions with CP Decomposition) .

## 5 Evaluation

We first present our evaluation metrics (\wasyparagraph[5.1](https://arxiv.org/html/2509.23348v1#S5.SS1)), given the analogous problem structure, we adopt metrics from tabular data analysis (Zhang et al., 2024). Then we use them to assess the experimental setups from \wasyparagraph[3.3](https://arxiv.org/html/2509.23348v1#S3.SS3), and report the results in \wasyparagraph[5.2](https://arxiv.org/html/2509.23348v1#S5.SS2).

### 5.1 Metrics for Evaluation

Evaluating generative models on discrete data is challenging since common metrics (e.g., perplexity for text, FID for images (Heusel et al., 2017)) are domain-specific. Following work on tabular data evaluation (Zhang et al., 2024; Shi et al., 2025), we adopt the Shape Score and Trend Score metrics.

#### Shape Score.

This metric measures how well the predicted data preserves the marginal (per-dimension) distributions of the real data. We consider a dataset with $|I_{R}|$ real samples $x$ and corresponding predicted samples $\tilde{x}$. We compute a per-dimension score for the empirical distributions (expressed in $\delta$-delta notation) and report the average across all dimensions:

$$ $\text{SSM}_{d}=1-\frac{1}{2}\sum_{s=0}^{S-1}\left|\frac{1}{|I_{R}|}\sum_{i=1}^{|I_{R}|}\delta(s-x^{(i)}_{d})-\frac{1}{|I_{R}|}\sum_{j=1}^{|I_{R}|}\delta(s-\tilde{x}^{(j)}_{d})\right|,\quad\text{SSM}=\frac{1}{D}\sum_{d=1}^{D}\text{SSM}_{d}.$ $$

#### Trend Score.

This metric evaluates whether pairwise dimension dependencies in the real data are preserved in the predictions. For a dataset with ${|I_{R}|}$ real samples $x^{(k)}$ and corresponding predicted samples $\tilde{x}^{(k)}$. We compute a trend score and report the average across all dimension pairs:

$$ $\text{TSM}_{d_{i},d_{j}}=1-\frac{1}{2}\sum_{s_{i}=0}^{S-1}\sum_{s_{j}=0}^{S-1}\left|\frac{1}{{}^{|I_{R}|}}\sum_{k=1}^{{}^{|I_{R}|}}\delta(s_{i}-x^{(k)}_{d_{i}})\delta(s_{j}-x^{(k)}_{d_{j}})-\frac{1}{{}^{|I_{R}|}}\sum_{k=1}^{{}^{|I_{R}|}}\delta(s_{i}-\tilde{x}^{(k)}_{d_{i}})\delta(s_{j}-\tilde{x}^{(k)}_{d_{j}})\right|,$ $$

where ${x}^{(k)}_{d_{i}}$ represents the $d_{i}$-th dimension of the $k$-th sample in this case.

All metrics range from 0 to 1, with higher values indicating better performance. Conditional metrics are computed by repeatedly sampling several $x_{1}$ from each $x_{0}\!\sim\!p_{0}$ using both the benchmark and the method, and comparing the resulting $x_{1}$’s.

### 5.2 Results

We use our benchmark pair constructor differently for training and testing. For *training*, we randomly sample $x_{0}^{\textup{train}}\sim p_{0}$ and generate $x_{1}^{\textup{train}}\sim p_{1}$ via our benchmark theorem, allowing infinite sample generation. Training is performed in an unpaired manner. For *testing*, we use a fixed set of $20\ 000$ precomputed samples, which we provide to facilitate benchmarking new discrete SB solvers. We also use different training setups, first by varying $N$ across CSBM, $\alpha$-CSBM, and DLightSB-M. Since these methods share the same primary loss, we experiment with two loss functions: KL, defined in ([15](https://arxiv.org/html/2509.23348v1#S4.E15)), and a variant based on the mean-squared error (MSE); see (Ksenofontov & Korotin, 2025, Appendix C.1) for equivalence details. Further details are provided in Appendix [C](https://arxiv.org/html/2509.23348v1#A3).

#### High-Dimensional Gaussian Mixtures.

In this section, we report results on the high-dimensional Gaussian mixture benchmark constructed as in \wasyparagraph[3.3](https://arxiv.org/html/2509.23348v1#S3.SS3) using the methods from \wasyparagraph[4](https://arxiv.org/html/2509.23348v1#S4). Visual results are shown in Figure [1](https://arxiv.org/html/2509.23348v1#S5.F1) for $q^{\textup{gauss}}$ ($\gamma\!=\!0.02$) and $q^{\textup{unif}}$ ($\gamma\!=\!0.005$). See Appendix [D](https://arxiv.org/html/2509.23348v1#A4) for additional plots.

Figure: (a) Input/Target
Refer to caption: images/benchmark/pairs_d2_g002.png

Tables [1](https://arxiv.org/html/2509.23348v1#S5.T1) and [2](https://arxiv.org/html/2509.23348v1#S5.T2) show that DLightSB consistently achieves the best performance on Conditional Shape Score and Trend Score metrics, respectively. We attribute this to the benchmark pairs being built on the same principle used by the DLightSB solver. DLightSB-M, which incorporates this inductive bias as well, achieves similar results with a slight drop in metrics, likely due to error accumulation in the iterative sampling. Interestingly, our results resemble those on continuous data (Korotin et al., 2024, Table 2; Gushchin et al., 2024a, Table 1), showing comparable performance with a slight drop for the DLightSB-M. Unconditional metrics are reported in Tables [3](https://arxiv.org/html/2509.23348v1#A4.T3) and [4](https://arxiv.org/html/2509.23348v1#A4.T4).

**Table 1: Conditional Shape Score metric ($\uparrow$) on the high-dimensional Gaussian mixture benchmark. The best-performing method is highlighted in bold, and the second is underlined. Color code threshold: red for $<0.7$, yellow for $[0.7,0.9)$, and green for $\geq 0.9$.**
|  |  |  | $D\!=\!2$ | $D\!=\!16$ | $D\!=\!64$ |  |  |  |  |  |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
|  |  |  | gaussian | uniform | gaussian | uniform | gaussian | uniform |  |  |  |  |  |  |
| Method | Loss | $N\!+\!1$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ |
| DLightSB | – | – | 0.979 | 0.976 | 0.974 | 0.983 | 0.972 | 0.980 | 0.970 | 0.981 | 0.966 | 0.980 | 0.980 | 0.973 |
| CSBM | KL | 16 | $\pagecolor{myyellow}{0.849}$ | $\pagecolor{myyellow}{0.733}$ | $\pagecolor{mygreen}{0.919}$ | $\pagecolor{myyellow}{0.892}$ | $\pagecolor{myyellow}{0.884}$ | $\pagecolor{myyellow}{0.806}$ | $\pagecolor{myyellow}{0.841}$ | $\pagecolor{myyellow}{0.810}$ | $\pagecolor{mygreen}{0.929}$ | $\pagecolor{mygreen}{0.938}$ | $\pagecolor{mygreen}{0.918}$ | $\pagecolor{mygreen}{0.922}$ |
| 64 | $\underline{\pagecolor{mygreen}{0.934}}$ | $\pagecolor{myyellow}{0.888}$ | $\pagecolor{mygreen}{0.958}$ | $\pagecolor{mygreen}{0.958}$ | $\underline{\pagecolor{mygreen}{0.944}}$ | $\pagecolor{mygreen}{0.933}$ | $\pagecolor{mygreen}{0.933}$ | $\pagecolor{mygreen}{0.927}$ | $\underline{\pagecolor{mygreen}{0.934}}$ | $\pagecolor{mygreen}{0.963}$ | $\pagecolor{mygreen}{0.926}$ | $\underline{\pagecolor{mygreen}{0.949}}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.721}$ | $\pagecolor{myyellow}{0.700}$ | $\pagecolor{myyellow}{0.824}$ | $\pagecolor{myyellow}{0.846}$ | $\pagecolor{myyellow}{0.854}$ | $\pagecolor{myyellow}{0.783}$ | $\pagecolor{myyellow}{0.839}$ | $\pagecolor{myyellow}{0.745}$ | $\pagecolor{mygreen}{0.915}$ | $\pagecolor{mygreen}{0.932}$ | $\pagecolor{myyellow}{0.893}$ | $\pagecolor{myyellow}{0.896}$ |  |
| 64 | $\pagecolor{myred}{0.444}$ | $\pagecolor{myyellow}{0.841}$ | $\pagecolor{myyellow}{0.818}$ | $\pagecolor{myyellow}{0.780}$ | $\pagecolor{myyellow}{0.885}$ | $\pagecolor{mygreen}{0.902}$ | $\pagecolor{myyellow}{0.890}$ | $\pagecolor{myyellow}{0.894}$ | $\pagecolor{myyellow}{0.854}$ | $\pagecolor{mygreen}{0.942}$ | $\pagecolor{myyellow}{0.867}$ | $\pagecolor{mygreen}{0.928}$ |  |  |
| $\alpha$-CSBM | KL | 16 | $\pagecolor{myyellow}{0.829}$ | $\pagecolor{myyellow}{0.738}$ | $\pagecolor{mygreen}{0.927}$ | $\pagecolor{mygreen}{0.918}$ | $\pagecolor{myyellow}{0.881}$ | $\pagecolor{myyellow}{0.836}$ | $\pagecolor{myyellow}{0.873}$ | $\pagecolor{myyellow}{0.825}$ | $\pagecolor{mygreen}{0.930}$ | $\underline{\pagecolor{mygreen}{0.972}}$ | $\pagecolor{mygreen}{0.929}$ | $\pagecolor{mygreen}{0.943}$ |
| 64 | $\pagecolor{mygreen}{0.902}$ | $\pagecolor{myyellow}{0.896}$ | $\pagecolor{mygreen}{0.952}$ | $\pagecolor{mygreen}{0.958}$ | $\pagecolor{mygreen}{0.936}$ | $\underline{\pagecolor{mygreen}{0.963}}$ | $\pagecolor{mygreen}{0.932}$ | $\pagecolor{mygreen}{0.941}$ | $\pagecolor{mygreen}{0.927}$ | $\pagecolor{mygreen}{0.959}$ | $\pagecolor{mygreen}{0.924}$ | $\pagecolor{mygreen}{0.942}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.803}$ | $\pagecolor{myred}{0.695}$ | $\pagecolor{myyellow}{0.841}$ | $\pagecolor{myyellow}{0.890}$ | $\pagecolor{myyellow}{0.865}$ | $\pagecolor{myyellow}{0.820}$ | $\pagecolor{myyellow}{0.861}$ | $\pagecolor{myyellow}{0.815}$ | $\pagecolor{mygreen}{0.908}$ | $\pagecolor{mygreen}{0.943}$ | $\pagecolor{myyellow}{0.884}$ | $\pagecolor{mygreen}{0.910}$ |  |
| 64 | $\pagecolor{mygreen}{0.908}$ | $\pagecolor{myyellow}{0.896}$ | $\pagecolor{myyellow}{0.858}$ | $\pagecolor{myyellow}{0.875}$ | $\pagecolor{mygreen}{0.908}$ | $\pagecolor{mygreen}{0.924}$ | $\pagecolor{myyellow}{0.881}$ | $\pagecolor{mygreen}{0.911}$ | $\pagecolor{myyellow}{0.883}$ | $\pagecolor{mygreen}{0.925}$ | $\pagecolor{myyellow}{0.859}$ | $\pagecolor{mygreen}{0.913}$ |  |  |
| DLightSB-M | KL | 16 | $\pagecolor{mygreen}{0.926}$ | $\underline{\pagecolor{mygreen}{0.956}}$ | $\underline{\pagecolor{mygreen}{0.969}}$ | $\underline{\pagecolor{mygreen}{0.970}}$ | $\pagecolor{myyellow}{0.894}$ | $\pagecolor{mygreen}{0.930}$ | $\pagecolor{mygreen}{0.961}$ | $\pagecolor{mygreen}{0.952}$ | $\pagecolor{mygreen}{0.931}$ | $\pagecolor{mygreen}{0.929}$ | $\underline{\pagecolor{mygreen}{0.954}}$ | $\pagecolor{mygreen}{0.905}$ |
| 64 | $\pagecolor{mygreen}{0.907}$ | $\pagecolor{mygreen}{0.954}$ | $\pagecolor{mygreen}{0.967}$ | $\pagecolor{mygreen}{0.968}$ | $\pagecolor{myyellow}{0.878}$ | $\pagecolor{mygreen}{0.953}$ | $\underline{\pagecolor{mygreen}{0.962}}$ | $\underline{\pagecolor{mygreen}{0.967}}$ | $\pagecolor{mygreen}{0.910}$ | $\pagecolor{mygreen}{0.942}$ | $\pagecolor{mygreen}{0.950}$ | $\pagecolor{mygreen}{0.942}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.782}$ | $\pagecolor{mygreen}{0.951}$ | $\pagecolor{myyellow}{0.881}$ | $\pagecolor{mygreen}{0.926}$ | $\pagecolor{myyellow}{0.726}$ | $\pagecolor{mygreen}{0.921}$ | $\pagecolor{mygreen}{0.942}$ | $\pagecolor{mygreen}{0.951}$ | $\pagecolor{myyellow}{0.718}$ | $\pagecolor{mygreen}{0.918}$ | $\pagecolor{myyellow}{0.891}$ | $\pagecolor{myyellow}{0.850}$ |  |
| 64 | $\pagecolor{myyellow}{0.717}$ | $\pagecolor{mygreen}{0.942}$ | $\pagecolor{myyellow}{0.892}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{myred}{0.685}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{mygreen}{0.953}$ | $\pagecolor{mygreen}{0.943}$ | $\pagecolor{myred}{0.632}$ | $\pagecolor{mygreen}{0.906}$ | $\pagecolor{myyellow}{0.730}$ | $\pagecolor{myyellow}{0.879}$ |  |  |

**Table 2: Conditional Trend Score ($\uparrow$) on the high-dimensional Gaussian mixture benchmark. The best-performing method is highlighted in bold, and the second is underlined. Color code threshold: red for $<0.7$, yellow for $[0.7,0.9)$, and green for $\geq 0.9$.**
|  |  |  | $D\!=\!2$ | $D\!=\!16$ | $D\!=\!64$ |  |  |  |  |  |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
|  |  |  | gaussian | uniform | gaussian | uniform | gaussian | uniform |  |  |  |  |  |  |
| Method | Loss | $N\!+\!1$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ |
| DLightSB | – | – | 0.968 | 0.970 | 0.967 | 0.975 | 0.943 | 0.967 | 0.956 | 0.967 | 0.919 | 0.956 | 0.955 | 0.950 |
| CSBM | KL | 16 | 0.793 | 0.654 | 0.884 | 0.856 | 0.803 | 0.694 | 0.732 | 0.676 | 0.853 | 0.895 | 0.830 | 0.861 |
| 64 | 0.911 | 0.854 | 0.932 | 0.923 | 0.886 | 0.890 | 0.874 | 0.874 | 0.859 | 0.936 | 0.848 | 0.901 |  |  |
| MSE | 16 | 0.611 | 0.631 | 0.752 | 0.781 | 0.739 | 0.653 | 0.725 | 0.612 | 0.835 | 0.883 | 0.799 | 0.823 |  |
| 64 | 0.331 | 0.775 | 0.735 | 0.729 | 0.808 | 0.831 | 0.812 | 0.821 | 0.767 | 0.891 | 0.777 | 0.863 |  |  |
| $\alpha$-CSBM | KL | 16 | 0.773 | 0.651 | 0.898 | 0.876 | 0.810 | 0.744 | 0.783 | 0.724 | 0.854 | 0.945 | 0.847 | 0.891 |
| 64 | 0.874 | 0.855 | 0.921 | 0.913 | 0.878 | 0.934 | 0.877 | 0.903 | 0.852 | 0.929 | 0.845 | 0.896 |  |  |
| MSE | 16 | 0.728 | 0.603 | 0.756 | 0.829 | 0.771 | 0.716 | 0.769 | 0.710 | 0.818 | 0.883 | 0.781 | 0.821 |  |
| 64 | 0.861 | 0.855 | 0.797 | 0.807 | 0.829 | 0.863 | 0.795 | 0.846 | 0.798 | 0.848 | 0.747 | 0.817 |  |  |
| DLightSB-M | KL | 16 | 0.878 | 0.943 | 0.952 | 0.956 | 0.738 | 0.914 | 0.932 | 0.930 | 0.862 | 0.900 | 0.920 | 0.674 |
| 64 | 0.856 | 0.940 | 0.951 | 0.953 | 0.716 | 0.923 | 0.928 | 0.936 | 0.833 | 0.901 | 0.648 | 0.820 |  |  |
| MSE | 16 | 0.701 | 0.933 | 0.838 | 0.904 | 0.551 | 0.877 | 0.897 | 0.917 | 0.575 | 0.853 | 0.773 | 0.555 |  |
| 64 | 0.640 | 0.922 | 0.852 | 0.889 | 0.503 | 0.856 | 0.903 | 0.910 | 0.464 | 0.818 | 0.498 | 0.700 |  |  |

On the other hand, CSBM and $\alpha$-CSBM perform noticeably worse than DLight methods. Notably, $\alpha$-CSBM achieves similar quality to CSBM while halving computational cost, making it a more efficient alternative. Regarding $N$ and the loss function, increasing $N$ mostly improves metrics. For the loss function, KL consistently outperforms MSE, likely because MSE minimizes pointwise squared error and produces over-smoothed solutions that blur modes (see Figure [2](https://arxiv.org/html/2509.23348v1#A4.F2)).

## 6 Discussion

Our work fills a key gap in discrete SB research by introducing the first standardized benchmark for these methods. This contribution provides the community with ground truth data and standard evaluation metrics, enabling meaningful comparisons across methods. The benchmark reveals fundamental limitations of current approaches: CP-based solvers (DLightSB, DLightSB-M) face severe memory constraints in high dimensions, while matching-based methods (CSBM, $\alpha$-CSBM) struggle with parameter sensitivity and long training times. These findings highlight the need for more scalable architectures as well as more stable training procedures. Looking forward, our benchmark facilitates rigorous evaluation of discrete SB methods, guiding future research.

#### Reproducibility.

We provide the experimental details in Appendix [C](https://arxiv.org/html/2509.23348v1#A3) and the code to reproduce the conducted experiments in the supplementary materials (see readme.md).

#### LLM Usage.

Large Language Models (LLMs) were used only to assist with rephrasing sentences and improving the clarity of the text. All scientific content, results, and interpretations in this paper were developed solely by the authors.

## References

- Austin et al. (2021)
Jacob Austin, Daniel D Johnson, Jonathan Ho, Daniel Tarlow, and Rianne Van Den Berg.
Structured denoising diffusion models in discrete state-spaces.
*Advances in Neural Information Processing Systems*, 34:17981–17993, 2021.
- Basharin et al. (2025)
Artem Basharin, Andrei Chertkov, and Ivan Oseledets.
Faster language models with better multi-token prediction using tensor decomposition, 2025.
URL [https://arxiv.org/abs/2410.17765](https://arxiv.org/abs/2410.17765).
- Campbell et al. (2022)
Andrew Campbell, Joe Benton, Valentin De Bortoli, Thomas Rainforth, George Deligiannidis, and Arnaud Doucet.
A continuous time framework for discrete denoising models.
*Advances in Neural Information Processing Systems*, 35:28266–28279, 2022.
- Campbell et al. (2024)
Andrew Campbell, Jason Yim, Regina Barzilay, Tom Rainforth, and Tommi Jaakkola.
Generative flows on discrete state-spaces: Enabling multimodal flows with applications to protein co-design.
In Ruslan Salakhutdinov, Zico Kolter, Katherine Heller, Adrian Weller, Nuria Oliver, Jonathan Scarlett, and Felix Berkenkamp (eds.), *Proceedings of the 41st International Conference on Machine Learning*, volume 235 of *Proceedings of Machine Learning Research*, pp. 5453–5512. PMLR, 21–27 Jul 2024.
URL [https://proceedings.mlr.press/v235/campbell24a.html](https://proceedings.mlr.press/v235/campbell24a.html).
- Chen et al. (2021)
Tianrong Chen, Guan-Horng Liu, and Evangelos Theodorou.
Likelihood training of Schrödinger bridge using forward-backward SDEs theory.
In *International Conference on Learning Representations*, 2021.
- Cohen et al. (2016)
Nadav Cohen, Or Sharir, and Amnon Shashua.
On the expressive power of deep learning: A tensor analysis, 2016.
URL [https://arxiv.org/abs/1509.05009](https://arxiv.org/abs/1509.05009).
- Cuturi (2013)
Marco Cuturi.
Sinkhorn distances: Lightspeed computation of optimal transport.
*Advances in neural information processing systems*, 26, 2013.
- Daniels et al. (2021)
Max Daniels, Tyler Maunu, and Paul Hand.
Score-based generative neural networks for large-scale optimal transport.
*Advances in neural information processing systems*, 34:12955–12965, 2021.
- De Bortoli et al. (2024)
Valentin De Bortoli, Iryna Korshunova, Andriy Mnih, and Arnaud Doucet.
Schrödinger bridge flow for unpaired data translation.
In *The Thirty-eighth Annual Conference on Neural Information Processing Systems*, 2024.
URL [https://openreview.net/forum?id=1F32iCJFfa](https://openreview.net/forum?id=1F32iCJFfa).
- Deb et al. (2021)
Nabarun Deb, Promit Ghosal, and Bodhisattva Sen.
Rates of estimation of optimal transport maps using plug-in estimators via barycentric projections.
*Advances in Neural Information Processing Systems*, 34:29736–29753, 2021.
- Dvurechensky et al. (2018)
Pavel Dvurechensky, Alexander Gasnikov, and Alexey Kroshnin.
Computational optimal transport: Complexity by accelerated gradient descent is better than by Sinkhorn’s algorithm.
In *International conference on machine learning*, pp. 1367–1376. PMLR, 2018.
- Esser et al. (2021)
Patrick Esser, Robin Rombach, and Bjorn Ommer.
Taming transformers for high-resolution image synthesis.
In *Proceedings of the IEEE/CVF conference on computer vision and pattern recognition*, pp. 12873–12883, 2021.
- Gat et al. (2024)
Itai Gat, Tal Remez, Neta Shaul, Felix Kreuk, Ricky TQ Chen, Gabriel Synnaeve, Yossi Adi, and Yaron Lipman.
Discrete flow matching.
In *The Thirty-eighth Annual Conference on Neural Information Processing Systems*, 2024.
- Gushchin et al. (2023a)
Nikita Gushchin, Alexander Kolesov, Alexander Korotin, Dmitry Vetrov, and Evgeny Burnaev.
Entropic neural optimal transport via diffusion processes.
In *Advances in Neural Information Processing Systems*, 2023a.
- Gushchin et al. (2023b)
Nikita Gushchin, Alexander Kolesov, Petr Mokrov, Polina Karpikova, Andrey Spiridonov, Evgeny Burnaev, and Alexander Korotin.
Building the bridge of Schrödinger: A continuous entropic optimal transport benchmark.
In *Thirty-seventh Conference on Neural Information Processing Systems Datasets and Benchmarks Track*, 2023b.
- Gushchin et al. (2024a)
Nikita Gushchin, Sergei Kholkin, Evgeny Burnaev, and Alexander Korotin.
Light and optimal Schrödinger bridge matching.
In *Forty-first International Conference on Machine Learning*, 2024a.
- Gushchin et al. (2024b)
Nikita Gushchin, Daniil Selikhanovych, Sergei Kholkin, Evgeny Burnaev, and Alexander Korotin.
Adversarial Schrödinger bridge matching.
In *The Thirty-eighth Annual Conference on Neural Information Processing Systems*, 2024b.
URL [https://openreview.net/forum?id=L3Knnigicu](https://openreview.net/forum?id=L3Knnigicu).
- Heusel et al. (2017)
Martin Heusel, Hubert Ramsauer, Thomas Unterthiner, Bernhard Nessler, and Sepp Hochreiter.
GANs trained by a two time-scale update rule converge to a local Nash equilibrium.
In *Advances in neural information processing systems*, pp. 6626–6637, 2017.
- Hoogeboom et al. (2021)
Emiel Hoogeboom, Didrik Nielsen, Priyank Jaini, Patrick Forré, and Max Welling.
Argmax flows and multinomial diffusion: Learning categorical distributions.
*Advances in Neural Information Processing Systems*, 34:12454–12465, 2021.
- Hütter & Rigollet (2021)
Jan-Christian Hütter and Philippe Rigollet.
Minimax estimation of smooth optimal transport maps.
2021.
- Kholkin et al. (2024)
Sergei Kholkin, Grigoriy Ksenofontov, David Li, Nikita Kornilov, Nikita Gushchin, Evgeny Burnaev, and Alexander Korotin.
Diffusion & adversarial Schrödinger bridges via iterative proportional Markovian fitting.
*arXiv preprint arXiv:2410.02601*, 2024.
- Kim et al. (2024)
Jun Hyeong Kim, Seonghwan Kim, Seokhyun Moon, Hyeongwoo Kim, Jeheon Woo, and Woo Youn Kim.
Discrete diffusion Schrödinger bridge matching for graph transformation.
*arXiv preprint arXiv:2410.01500*, 2024.
- Korotin et al. (2024)
Alexander Korotin, Nikita Gushchin, and Evgeny Burnaev.
Light Schrödinger bridge.
In *The Twelfth International Conference on Learning Representations*, 2024.
- Ksenofontov & Korotin (2025)
Grigoriy Ksenofontov and Alexander Korotin.
Categorical Schrödinger bridge matching.
In *Forty-second International Conference on Machine Learning*, 2025.
URL [https://openreview.net/forum?id=RBly0nOr2h](https://openreview.net/forum?id=RBly0nOr2h).
- Léonard (2013)
Christian Léonard.
A survey of the Schrödinger problem and some of its connections with optimal transport.
*arXiv preprint arXiv:1308.0215*, 2013.
- Lou et al. (2023)
Aaron Lou, Chenlin Meng, and Stefano Ermon.
Discrete diffusion language modeling by estimating the ratios of the data distribution.
2023.
- Luo et al. (2024)
Xiaoshan Luo, Zhenyu Wang, Jian Lv, Lei Wang, Yanchao Wang, and Yanming Ma.
CrystalFlow: A flow-based generative model for crystalline materials.
*arXiv preprint arXiv:2412.11693*, 2024.
- Manole et al. (2024)
Tudor Manole, Sivaraman Balakrishnan, Jonathan Niles-Weed, and Larry Wasserman.
Plugin estimation of smooth optimal transport maps.
*The Annals of Statistics*, 52(3):966–998, 2024.
- Mokrov et al. (2024)
Petr Mokrov, Alexander Korotin, Alexander Kolesov, Nikita Gushchin, and Evgeny Burnaev.
Energy-guided entropic neural optimal transport.
In *The Twelfth International Conference on Learning Representations*, 2024.
URL [https://openreview.net/forum?id=d6tUsZeVs7](https://openreview.net/forum?id=d6tUsZeVs7).
- Peluchetti (2024)
Stefano Peluchetti.
Bm^2: Coupled Schrödinger bridge matching.
*arXiv preprint arXiv:2409.09376*, 2024.
- Pooladian & Niles-Weed (2021)
Aram-Alexandre Pooladian and Jonathan Niles-Weed.
Entropic estimation of optimal transport maps.
*arXiv preprint arXiv:2109.12004*, 2021.
- Qin et al. (2024)
Yiming Qin, Manuel Madeira, Dorina Thanou, and Pascal Frossard.
DeFoG: Discrete flow matching for graph generation.
*arXiv preprint arXiv:2410.04263*, 2024.
- Sahoo et al. (2024)
Subham Sahoo, Marianne Arriola, Yair Schiff, Aaron Gokaslan, Edgar Marroquin, Justin Chiu, Alexander Rush, and Volodymyr Kuleshov.
Simple and effective masked diffusion language models.
*Advances in Neural Information Processing Systems*, 37:130136–130184, 2024.
- Schrödinger (1931)
Erwin Schrödinger.
*Über die Umkehrung der Naturgesetze*.
Verlag der Akademie der Wissenschaften in Kommission bei Walter De Gruyter u. Company, 1931.
- Schrödinger (1932)
Erwin Schrödinger.
Sur la théorie relativiste de l’électron et l’interprétation de la mécanique quantique.
In *Annales de l’institut Henri Poincaré*, volume 2, pp. 269–310, 1932.
- Shi et al. (2025)
Juntong Shi, Minkai Xu, Harper Hua, Hengrui Zhang, Stefano Ermon, and Jure Leskovec.
Tabdiff: a mixed-type diffusion model for tabular data generation.
In *The Thirteenth International Conference on Learning Representations*, 2025.
URL [https://openreview.net/forum?id=swvURjrt8z](https://openreview.net/forum?id=swvURjrt8z).
- Shi et al. (2023)
Yuyang Shi, Valentin De Bortoli, Andrew Campbell, and Arnaud Doucet.
Diffusion Schrödinger bridge matching.
In *Thirty-seventh Conference on Neural Information Processing Systems*, 2023.
URL [https://openreview.net/forum?id=qy07OHsJT5](https://openreview.net/forum?id=qy07OHsJT5).
- Van Den Oord et al. (2017)
Aaron Van Den Oord, Oriol Vinyals, et al.
Neural discrete representation learning.
*Advances in neural information processing systems*, 30, 2017.
- Vargas et al. (2021)
Francisco Vargas, Pierre Thodoroff, Austen Lamacraft, and Neil Lawrence.
Solving Schrödinger bridges via maximum likelihood.
*Entropy*, 23(9):1134, 2021.
- Vignac et al. (2022)
Clement Vignac, Igor Krawczuk, Antoine Siraudin, Bohan Wang, Volkan Cevher, and Pascal Frossard.
Digress: Discrete denoising diffusion for graph generation.
*arXiv preprint arXiv:2209.14734*, 2022.
- Zhang et al. (2024)
Hengrui Zhang, Jiani Zhang, Balasubramaniam Srinivasan, Zhengyuan Shen, Xiao Qin, Christos Faloutsos, Huzefa Rangwala, and George Karypis.
Mixed-type tabular data synthesis with score-based diffusion in latent space, 2024.
URL [https://arxiv.org/abs/2310.09656](https://arxiv.org/abs/2310.09656).

## Appendix A Proofs

###### Proof of Theorem 3.1 .

We start from the expression of the static EOT minimization problem in ([8](https://arxiv.org/html/2509.23348v1#S3.E8))

$$ $\begin{split}\min_{q\in\Pi(p_{0},p_{1})}&\text{KL}\left(q(x_{0},x_{1})\|q^{\text{ref}}(x_{0},x_{1})\right)=\\ &=\min_{q\in\Pi(p_{0},p_{1})}-H(q)-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})-\text{const}\\ &=\min_{q\in\Pi(p_{0},p_{1})}\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q(x_{0},x_{1})-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})-\text{const}\end{split}$ (23) $$

Noting that the joint distribution factorizes as $q(x_{0},x_{1})=q(x_{0})q(x_{1}|x_{0})=p_{0}(x_{0})q(x_{1}|x_{0})$, and enforcing the marginal constraints $\sum_{x_{0}}p_{0}(x_{0})q(x_{1}|x_{0})=p_{1}(x_{1})$ and $\sum_{x_{1}}q(x_{1}|x_{0})=1$ (equivalently $q(x_{0})=p_{0}(x_{0})$), the corresponding Lagrangian can be formulated as

$$ $\begin{split}\mathcal{L}(q)&=\sum_{x_{0},x_{1}}p_{0}(x_{0})q(x_{1}|x_{0})\log\left(p_{0}(x_{0})q(x_{1}|x_{0})\right)-\sum_{x_{0},x_{1}}p_{0}(x_{0})q(x_{1}|x_{0})\log q^{\text{ref}}(x_{1}|x_{0})+\\ &+\sum_{x_{1}}\lambda(x_{1})\left(\sum_{x_{0}}q(x_{1}|x_{0})p_{0}(x_{0})-p_{1}(x_{1})\right)+\sum_{x_{0}}\tau(x_{0})\left(\sum_{x_{1}}q(x_{1}|x_{0})-p_{0}(x_{0})\right)\\ &=\underbrace{\sum_{x_{0},x_{1}}p_{0}(x_{0})q(x_{1}|x_{0})\log p_{0}(x_{0}))}_{=\sum_{x_{0}}p_{0}(x_{0})\log p_{0}(x_{0}))}+\sum_{x_{0},x_{1}}p_{0}(x_{0})q(x_{1}|x_{0})\log q(x_{1}|x_{0})-\\ -&\sum_{x_{0},x_{1}}p_{0}(x_{0})q(x_{1}|x_{0})\log q^{\text{ref}}(x_{1}|x_{0})+\sum_{x_{1}}\lambda(x_{1})\left(\sum_{x_{0}}q(x_{1}|x_{0})p_{0}(x_{0})-p_{1}(x_{1})\right)\\ &+\sum_{x_{0}}\tau(x_{1})\left(\sum_{x_{1}}q(x_{1}|x_{0})-1\right)\end{split}$ (24) $$

where $\lambda(x_{1})$ and $\tau(x_{0})$ denote the Lagrange multipliers associated with the marginal constraints on $x_{1}$ and $x_{0}$, respectively. Taking the pointwise partial derivative of $\mathcal{L}(q)$ with respect to $q(x_{1}|x_{0})$ then yields

$$ $\frac{\partial\mathcal{L}}{\partial q}=p_{0}(x_{0})\left(\log q(x_{1}|x_{0})+1\right)-p_{0}(x_{0})\log q^{\text{ref}}(x_{0},x_{1})+\lambda(x_{1})p_{0}(x_{0})+\tau(x_{1})=0$ (25) $$

Therefore, the optimal process $q^{*}$ can be written as

$$ $q^{*}(x_{1}|x_{0})=\exp(-\lambda(x_{1})-1)q^{\text{ref}}(x_{1}|x_{0})p_{0}\exp\left(-\frac{\tau(x_{0})}{p_{0}(x_{0})}\right)$ (26) $$

Setting $v^{*}(x_{1})=\exp(-\lambda(x_{1})-1)$ concludes the proof.
∎

###### Proof of Proposition 3.1 .

Assuming the CP parameterization introduced in ([10](https://arxiv.org/html/2509.23348v1#S3.E10)), and further assuming that the reference process factorizes across dimensions as $q^{\text{ref}}(x_{1}|x_{0})=\prod_{d=1}^{D}q^{\text{ref}}(x_{1}^{d}|x_{0})$, the normalized conditional distribution $q^{*}(x_{1}|x_{0})$ in ([9](https://arxiv.org/html/2509.23348v1#S3.E9)) can be rewritten as

$$ $\begin{split}q^{*}(x_{1}|x_{0})&=\frac{1}{c(x_{0})}\left(\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}r_{k}^{d}[x_{1}^{d}]\right)\prod_{d=1}^{D}q^{\text{ref}}(x_{1}^{d}|x_{0})\\ &=\frac{1}{c(x_{0})}\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}r_{k}^{d}[x_{1}^{d}]\,q^{\text{ref}}(x_{1}^{d}|x_{0}),\end{split}$ (27) $$

where the reference factors can be merged with the rank-1 components because they are independent of the mixture index $k$ and factorize over dimensions. From here, it is possible to obtain the normalizing constant $c(x_{0})$ by summing over all possible values of $x_{1}\in{\mathcal{X}}=\mathbb{S}^{D}$, where $x_{1}^{d}\in\{0,\dots,S-1\}$. The normalizing constant can then be rewritten as

$$ $\begin{split}c(x_{0})&=\sum_{x_{1}\in\mathbb{S}^{D}}\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}r_{k}^{d}[x_{1}^{d}]\,q^{\text{ref}}(x_{1}^{d}|x_{0})\\ &=\sum_{k=1}^{K}\beta_{k}\sum_{x_{1}\in\mathbb{S}^{D}}\prod_{d=1}^{D}r_{k}^{d}[x_{1}^{d}]\,q^{\text{ref}}(x_{1}^{d}|x_{0})\\ &=\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}\sum_{x_{1}^{d}=0}^{S-1}r_{k}^{d}[x_{1}^{d}]\,q^{\text{ref}}(x_{1}^{d}|x_{0}),\end{split}$ (28) $$

where $\sum_{x_{1}\in\mathbb{S}^{D}}=\sum_{x_{1}^{1}=0}^{S-1}\sum_{x_{1}^{2}=0}^{S-1}\dots\sum_{x_{1}^{D}=0}^{S-1}$. The exchange between the product and the sum is valid here because the summation is separable across dimensions, i.e., each factor depends only on its corresponding coordinate $x_{1}^{d}$.
∎

###### Proof of Proposition 4.1 .

We start from the standard KL minimization problem from the LightSB paper (Korotin et al., 2024) and define it in discrete space.

$$ $\displaystyle\begin{split}\text{KL}\left(q^{*}\|q\right)&=\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\left(\frac{q^{*}(x_{0},x_{1})}{q(x_{0},x_{1})}\right)=\sum_{x_{0},x_{1}}q^{*}\log q^{*}(x_{0},x_{1})-\sum_{x_{0},x_{1}}q^{*}\log q(x_{0},x_{1})=\\ &=-H(q^{*})-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q(x_{0},x_{1})=-H(q^{*})-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\left(q(x_{0})q(x_{1}|x_{0})\right)\\ &=-H(q^{*})-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\underbrace{q(x_{0})}_{=p_{0}(x_{0})}-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q(x_{1}|x_{0})=\\ &=-H(q^{*})-\sum_{x_{0}}\log p_{0}(x_{0})\underbrace{\sum_{x_{1}}q^{*}(x_{0},x_{1})}_{=q^{*}(x_{0})=p_{0}(x_{0})}-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q(x_{1}|x_{0})\end{split}$ $$

Now using ([9](https://arxiv.org/html/2509.23348v1#S3.E9)) on $q(x_{1}|x_{0})$ we can get

$$ $\begin{split}\text{KL}\left(q^{*}\|q\right)&=-H(q^{*})-\sum_{x_{0}}\log p_{0}(x_{0})p_{0}(x_{0})-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\left(\frac{v^{*}(x_{1})}{c^{*}(x_{0})}q^{\text{ref}}(x_{1}|x_{0})\right)=\\ &=\underbrace{-H(q^{*})-\sum_{x_{0}}\log p_{0}(x_{0})p_{0}(x_{0})-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})}_{=-\mathcal{L^{*}}}-\\ &-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\left(\frac{v^{*}(x_{1})}{c^{*}(x_{0})}\right)=\\ &=-\mathcal{L}^{*}+\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log c^{*}(x_{0})-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log v^{*}(x_{1})=\\ &=\sum_{x_{0}}p_{0}^{*}(x_{0})\log c^{*}(x_{0})-\sum_{x_{1}}q^{*}(x_{1})\log v^{*}(x_{1})-\mathcal{L}^{*},\end{split}$ $$

That concludes the proof.
∎

###### Proof of expression 7 .

Let $Q$ be the transition matrix in ([4](https://arxiv.org/html/2509.23348v1#S2.E4)), rewritten as

$$ $\displaystyle Q$ $\displaystyle=(1-\gamma)I+\frac{\gamma}{S-1}(\mathbf{1}\mathbf{1}^{\top}-I)$ $\displaystyle=\left(1-\gamma\frac{S}{S-1}\right)I+\frac{\gamma}{S-1}\mathbf{1}\mathbf{1}^{\top},$ $$

where $I$ is the identity matrix and $\mathbf{1}\mathbf{1}^{\top}$ is the all-ones matrix. Let

$$ $a=1-\gamma\frac{S}{S-1},\quad b=\frac{\gamma}{S-1},$ $$

so that $Q=aI+b\mathbf{1}\mathbf{1}^{\top}$ and note that $a+Sb=1$. We compute $Q^{N+1}$ using the binomial expansion. Since $I$ and $\mathbf{1}\mathbf{1}^{\top}$ commute:

$$ $\displaystyle Q^{n}$ $\displaystyle=(aI+b\mathbf{1}\mathbf{1}^{\top})^{n}$ $\displaystyle=\sum_{k=0}^{n}\binom{n}{k}a^{n-k}b^{k}(\mathbf{1}\mathbf{1}^{\top})^{k}.$ $$

Using $(\mathbf{1}\mathbf{1}^{\top})^{k}=S^{k-1}\mathbf{1}\mathbf{1}^{\top}$ for $k\geq 1$ and separating the $k=0$ term:

$$ $\displaystyle Q^{n}$ $\displaystyle=a^{n}I+\sum_{k=1}^{n}\binom{n}{k}a^{n-k}b^{k}S^{k-1}\mathbf{1}\mathbf{1}^{\top}$ $\displaystyle=a^{n}I+\frac{1}{S}\left(\sum_{k=1}^{n}\binom{n}{k}a^{n-k}(bS)^{k}\right)\mathbf{1}\mathbf{1}^{\top}.$ $$

The binomial expansion gives:

$$ $(a+bS)^{n}=\sum_{k=0}^{n}\binom{n}{k}a^{n-k}(bS)^{k}=a^{n}+\sum_{k=1}^{n}\binom{n}{k}a^{n-k}(bS)^{k}.$ $$

Since $a+bS=1$, we have $(a+bS)^{n}=1$, so
$\sum_{k=1}^{n}\binom{n}{k}a^{n-k}(bS)^{k}=1-a^{n}$. Thus,

$$ $Q^{n}=a^{n}I+\frac{1-a^{n}}{S}\mathbf{1}\mathbf{1}^{\top}.$ $$

Substituting $n=N+1$ and $a=1-\gamma\frac{S}{S-1}$ yields

$$ $q^{\text{ref}}(x_{1}|x_{0})=Q^{N+1}=\left(1-\gamma\frac{S}{S-1}\right)^{N+1}I+\frac{1-\left(1-\gamma\frac{S}{S-1}\right)^{N+1}}{S}\mathbf{1}\mathbf{1}^{\top}.$ $$

This completes the proof.
∎

###### Proof of Proposition 4.2 .

$$ $\displaystyle\text{KL}\left(q(x_{0},x_{\text{in}},x_{1})\|q^{\text{SB}}(x_{0},x_{\text{in}},x_{1})\right)=$ $\displaystyle=\text{KL}\left(q(x_{0},x_{1})\|q^{\text{SB}}(x_{0},x_{1})\right)+\underbrace{\text{KL}\left(q^{\text{ref}}(x_{\text{in}}|x_{0},x_{1})\|q^{\text{ref}}(x_{\text{in}}|x_{0},x_{1})\right)}_{=0}=$ (29) $\displaystyle=\underbrace{\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q(x_{0},x_{1})}_{=-H(q(x_{0},x_{1}))}-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{SB}}(x_{0},x_{1})=$ $\displaystyle=-H(q(x_{0},x_{1}))-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log\frac{v^{\text{SB}}(x_{1})q^{\text{ref}}(x_{1}|x_{0})}{c^{\text{SB}}(x_{0})}=$ (30) $\displaystyle=-H(q(x_{0},x_{1}))-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log v^{\text{SB}}(x_{1})-$ $\displaystyle-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})+\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log c^{\text{SB}}(x_{0})=$ $\displaystyle=-H(q(x_{0},x_{1}))-\sum_{x_{1}}\log v^{\text{SB}}(x_{1})\underbrace{q(x_{1})}_{=p(x_{1})=q^{*}(x_{1})}\underbrace{\sum_{x_{0}}q(x_{0}|x_{1})}_{=1=\sum_{x_{0}}q^{*}(x_{0}|x_{1})}-$ (31) $\displaystyle-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})+\sum_{x_{0}}\log c^{\text{SB}}(x_{0})\underbrace{q(x_{0})}_{=p(x_{0})=q^{*}(x_{0})}\underbrace{\sum_{x_{1}}q(x_{1}|x_{0})}_{=1=\sum_{x_{1}}q^{*}(x_{1}|x_{0})}=$ (32) $\displaystyle=\underbrace{-H(q(x_{0},x_{1}))-\sum_{x_{0},x_{1}}q(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})}_{=C_{1}}-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\frac{v^{\text{SB}}(x_{1})}{c^{\text{SB}}(x_{0})}=$ $\displaystyle=C_{1}-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\frac{v^{\text{SB}}(x_{1})}{c^{\text{SB}}(x_{0})}-$ $\displaystyle-\underbrace{\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})+\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})}_{=0}=$ (33) $\displaystyle=-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\frac{v^{\text{SB}}(x_{1})q^{\text{ref}}(x_{1}|x_{0})}{c^{\text{SB}}(x_{0})}+\underbrace{C_{1}+\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{\text{ref}}(x_{1}|x_{0})}_{=C_{2}}=$ $\displaystyle=C_{2}-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{\text{SB}}(x_{0},x_{1})+$ $\displaystyle+\underbrace{\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{*}(x_{0},x_{1})-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{*}(x_{0},x_{1})}_{=0}=$ (34) $\displaystyle=\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log\frac{q^{*}(x_{0},x_{1})}{q^{\text{SB}}(x_{0},x_{1})}+\underbrace{C_{2}-\sum_{x_{0},x_{1}}q^{*}(x_{0},x_{1})\log q^{*}(x_{0},x_{1})}_{C_{3}}=$ $\displaystyle=\text{KL}\left(q^{*}(x_{0},x_{1})\|q^{\text{SB}}(x_{0},x_{1})\right)+C_{3}$ $$

In ([29](https://arxiv.org/html/2509.23348v1#A1.E29)), we use the disintegration of the KL divergence to transition from the dynamic to the static formulation. In ([30](https://arxiv.org/html/2509.23348v1#A1.E30)), we apply our parameterization from ([9](https://arxiv.org/html/2509.23348v1#S3.E9)). Next, in ([31](https://arxiv.org/html/2509.23348v1#A1.E31)) and ([32](https://arxiv.org/html/2509.23348v1#A1.E32)), we use the properties of the reciprocal process $q$, which has the true marginals at $t=0$ and $t=1$. In ([33](https://arxiv.org/html/2509.23348v1#A1.E33)), we add a zero term to introduce $q^{\text{ref}}(x_{1}|x_{0})$ with the expectation taken over the optimal coupling $q^{*}(x_{0},x_{1})$. Finally, in ([34](https://arxiv.org/html/2509.23348v1#A1.E34)), we obtain the entropy term, completing the expression for the desired KL divergence.
∎

###### Proof of Proposition 4.3 .

We first derive the transitional distributions of the SB by recalling its well-known characterization (Léonard, 2013, Prop. 4.2):

$$ $q^{\text{SB}}\left(x_{t_{n}}|x_{t_{n-1}}\right)=q^{\text{ref}}\left(x_{t_{n}}|x_{t_{n-1}}\right)\frac{\phi^{\text{SB}}_{t_{n}}(x_{t_{n}})}{\phi^{\text{SB}}_{t_{n-1}}(x_{t_{n-1}})},\qquad\phi^{\text{SB}}_{t_{n}}(x_{t_{n}})=\mathbb{E}_{q^{\mathrm{ref}}(x_{1}|x_{t_{n}})}\big[v^{\text{SB}}(x_{1})\big].$ $$

Using the CP parametrization of $v^{\text{SB}}$ from ([10](https://arxiv.org/html/2509.23348v1#S3.E10)) and exploiting the conditional independence of dimensions under $q^{\text{ref}}$, the scalar-valued functions $\phi_{t_{n}}$ can be written as:

$$ $\phi^{\text{SB}}_{t_{n}}(x_{t_{n}})=\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}\mathbb{E}_{q^{\text{ref}}(x_{1}^{d}|x_{t_{n}}^{d})}\!\big[r_{k}^{d}(x_{1}^{d})\big]=\sum_{k=1}^{K}\beta_{k}\prod_{d=1}^{D}\underbrace{\sum_{x_{1}^{d}=0}^{S-1}\big[\overline{Q}^{\text{ref}}_{N+1-n}\big]_{x_{t_{n}}^{d},x_{1}^{d}}\,r_{k}^{d}[x_{1}^{d}]}_{u_{k,t_{n}}^{d}[x_{t_{n}}^{d}]},$ $$

where $u_{k,t_{n}}^{d}$ satisfy the following recursive relation:

$$ $u_{k,t_{n}}^{d}[x_{t_{n}}^{d}]=\sum_{x_{t_{n+1}}^{d}=0}^{S-1}\big[Q^{\mathrm{ref}}\big]_{x_{t_{n}}^{d},x_{t_{n+1}}^{d}}u_{k,t_{n+1}}^{d}[x_{t_{n+1}}^{d}],\qquad u_{k,t_{1}}^{d}=r_{k}^{d}.$ $$

Thus, we obtain the following transition distributions:

$$ $q^{\text{SB}}(x_{t_{n}}|x_{t_{n-1}})\propto q^{\text{ref}}(x_{t_{n}}|x_{t_{n-1}})\sum_{k=1}^{K}\beta_{k}\prod_{j=1}^{D}u_{k,t_{n}}^{j}[x_{t_{n}}^{j}].$ (35) $$

To obtain the $d$-th marginal transition distribution, we marginalize over $x_{t_{n}}^{-d}\stackrel{{\scriptstyle\text{def}}}{{=}}\{x_{t_{n}}^{j}\}_{j\neq d}$ as follows:

$$ $q^{\mathrm{SB}}\big(x_{t_{n}}^{d}|x_{t_{n-1}}\big)\propto\sum_{x_{t_{n}}^{-d}}\Bigg(\prod_{j=1}^{D}[Q^{\text{ref}}]_{x^{j}_{t_{n-1}},x_{t_{n}}^{j}}\Bigg)\Bigg(\sum_{k=1}^{K}\beta_{k}\prod_{j=1}^{D}u_{k,t_{n}}^{j}[x_{t_{n}}^{j}]\Bigg)=\\ =[Q^{\text{ref}}]_{x_{t_{n-1}}^{d},x_{t_{n}}^{d}}\sum_{k=1}^{K}\beta_{k}u_{k,t_{n}}^{d}[x_{t_{n}}^{d}]\prod^{D}_{\begin{subarray}{c}j=1\\ j\neq d\end{subarray}}\underbrace{\sum_{x_{t_{n}}^{j}}[Q^{\text{ref}}]_{x_{t_{n-1}}^{j},x_{t_{n}}^{j}}u_{k,t_{n}}^{j}[x_{t_{n}}^{j}]}_{u_{k,t_{n-1}}^{j}[x_{t_{n-1}}^{j}]\ \text{(by recursion)}}.$ $$

Finally, we obtain the desired expression up to normalization:

$$ $q^{\mathrm{SB}}\big(x_{t_{n}}^{d}|x_{t_{n-1}}\big)\propto[Q^{\text{ref}}]_{x_{t_{n-1}}^{d},x_{t_{n}}^{d}}\sum_{k=1}^{K}\beta_{k}u_{k,t_{n}}^{d}[x_{t_{n}}^{d}]\prod^{D}_{\begin{subarray}{c}j=1\\ j\neq d\end{subarray}}u_{k,t_{n-1}}^{j}[x_{t_{n-1}}^{j}].$ $$

Now we derive the sampling procedure. Sampling from the SB transitional distributions is based on the following factorization:

$$ $q^{\text{SB}}(x_{t_{n}}|x_{t_{n-1}})=\sum_{k=1}^{K}p(k|x_{t_{n-1}})\prod_{d=1}^{D}q^{\text{SB}}(x_{t_{n}}^{d}|x_{t_{n-1}},k).$ $$

Using the full joint SB transition distribution ([35](https://arxiv.org/html/2509.23348v1#A1.E35)), the probability of $k$ is

$$ $p(k|x_{t_{n-1}})\propto\beta_{k}\prod_{d=1}^{D}u_{k,t_{n-1}}^{d}[x_{t_{n-1}}^{d}].$ $$

Using the marginal distributions conditioned on $k$, the factors with $j\neq d$ are independent of $x_{t_{n}}^{d}$ and absorbed into normalization, yielding

$$ $q^{\text{SB}}(x_{t_{n}}^{d}|x_{t_{n-1}},k)\propto[Q^{\text{ref}}]_{x_{t_{n-1}}^{d},x_{t_{n}}^{d}}\,u_{k,t_{n}}^{d}[x_{t_{n}}^{d}].$ $$

Thus, sampling proceeds by first drawing

$$ $k^{*}\sim p(k|x_{t_{n}}),$ $$

and then sampling each coordinate independently as

$$ $x_{t_{n+1}}^{d}\sim q^{\text{SB}}(\cdot|x_{t_{n}},k^{*})\propto[Q^{\text{ref}}]_{x_{t_{n}}^{d},\cdot}\,u_{k^{*},t_{n+1}}^{d}[\cdot],\qquad d=1,\dots,D.$ $$

∎

###### Proof of Theorem 3.1 .

###### Proof of Proposition 3.1 .

###### Proof of Proposition 4.1 .

###### Proof of expression 7 .

###### Proof of Proposition 4.2 .

###### Proof of Proposition 4.3 .

## Appendix B Methods Details

This section provides additional theoretical and implementation details complementing \wasyparagraph[4](https://arxiv.org/html/2509.23348v1#S4), focusing on the methods used to evaluate our benchmark pairs.

#### CSBM.

In practice, the D-IMF procedure is usually implemented bidirectionally: the Markovian projection is applied using both forward and backward representations (see “Notation”), which is the approach we adopt in our experiments. This design has two advantages. First, it mitigates error accumulation caused by imperfect model fitting, as shown in (De Bortoli et al., 2024, Appendix F). Second, it enables the use of alternative starting couplings, as proposed in (Kholkin et al., 2024).

Limitations. To reduce the computational overhead of evaluating the full probability state space of size $S^{D}$, the authors propose factorizing transition probabilities across dimensions, reducing the space to $D\times S$. However, this parametrization constitutes a key limitation of CSBM, as it introduces approximation error.

#### DLightSB.

We optimize over the logarithm of the mixture weights $\beta\in\mathbb{K}$ and the logarithm of the CP cores $r_{k}^{d}$. This allows computation of the log terms in the surrogate loss ([18](https://arxiv.org/html/2509.23348v1#S4.E18)) by stable log-sum-exp operations.

Limitations. (1) In spite of being numerically stable, the log-sum-exp operations allocate extra memory, which can become a bottleneck when applied repeatedly in high-dimensional settings. (2) The CP parameterization requires an impractically large number of components to capture complex data, making it infeasible under memory constraints. Additionally, we approximate the summations using Monte Carlo samples from $p_{0}$ and $p_{1}$.

#### DLightSB-M.

Limitations. The practical implementation requires storing or recomputing $u_{k,t}^{d}$ at each iteration, which scales as $\mathcal{O}(B\times S^{2}\times K)$ in memory and computation. This quickly becomes prohibitive for high-dimensional data, limiting scalability to small state spaces.

## Appendix C Experiment Details

This section provides detailed descriptions of all methods and their configurations.

#### Shared Aspects.

Across all experiments, we use the AdamW optimizer with fixed beta values of $0.95$ and $0.99$. For the high-dimensional Gaussian benchmark (\wasyparagraph[5.2](https://arxiv.org/html/2509.23348v1#S5.SS2.SSS0.Px1)). Notably, for diffusion-based methods, we fully sample the Markov chain, in contrast to Austin et al. (2021), which applies an argmax operation at the final timestep. To evaluate the methods on the high-dimensional Gaussian mixture benchmark (\wasyparagraph[5.2](https://arxiv.org/html/2509.23348v1#S5.SS2.SSS0.Px1)), we use $20\ 000$ samples. Conditional metrics are computed using 156 instances of $x_{0}$, with $1\ 000$ samples of $x_{1}$ generated for each $x_{0}$.

#### CSBM and α \alpha -CSBM.

For CSBM and $\alpha$-CSBM, we use the official implementation from Ksenofontov & Korotin (2025):

[https://github.com/gregkseno/csbm](https://github.com/gregkseno/csbm).

To stabilize training and improve final performance, we apply Exponential Moving Average (EMA) parameter updates with a decay rate of $0.999$, tuned consistently across all experiments. Unlike Austin et al. (2021), we omit the $L_{\text{simple}}$ loss during training. We employ a simple MLP with three hidden layers of size $[128,128,128]$ and ReLU activations. Time conditioning is implemented via an embedding layer of the same size as dimensions, $D$. Both methods are trained for $5$ D-IMF iterations, using $120\ 000$ gradient updates in the first iteration and $40\ 000$ in each subsequent iteration. For $\alpha$-CSBM, we use a learning rate of $10^{-3}$ and halve the batch size for training a single model, following De Bortoli et al. (2024). For CSBM, we use a learning rate of $10^{-4}$.

#### DLightSB and DLightSB-M.

For all benchmark experiments, both methods use $K=1000$ components initialized from data samples and are trained for $100\ 000$ gradient updates. The learning rate is set to $10^{-2}$ for both, with DLightSB-M using independent coupling ($q^{0}(x_{0},x_{1})=p_{0}(x_{0})p_{1}(x_{1})$).

#### Computational Resources and Training Time.

All high-dimensional Gaussian mixture benchmark experiments were conducted on 1 A100 GPU unless otherwise specified, with training times reported inclusive of evaluation. For $D=2$, training is relatively short: CSBM and $\alpha$-CSBM each complete within about $5$ hours, DLightSB-M within $4$ hours, and DLightSB in roughly $20$ minutes. For $D=64$, CSBM completes in under $14$ hours, $\alpha$-CSBM in under $9$ hours, DLightSB-M in just under 2 days (on 2 A100 GPUs), and DLightSB in under $7$ hours.

## Appendix D Additional Experiments

**Table 3: Shape Score metric ($\uparrow$) on the high-dimensional Gaussian mixture benchmark. The best-performing method is highlighted in bold, and the second is underlined. Color code threshold: red for $<0.7$, yellow for $[0.7,0.9)$, and green for $\geq 0.9$.**
|  |  |  | $D\!=\!2$ | $D\!=\!16$ | $D\!=\!64$ |  |  |  |  |  |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
|  |  |  | gaussian | uniform | gaussian | uniform | gaussian | uniform |  |  |  |  |  |  |
| Method | Loss | $N\!+\!1$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ |
| DLightSB | – | – | ${\mathbf{\pagecolor{mygreen}{0.975}}}$ | ${\mathbf{\pagecolor{mygreen}{0.969}}}$ | ${\mathbf{\pagecolor{mygreen}{0.969}}}$ | ${\mathbf{\pagecolor{mygreen}{0.980}}}$ | ${\mathbf{\pagecolor{mygreen}{0.973}}}$ | ${\mathbf{\pagecolor{mygreen}{0.976}}}$ | ${\mathbf{\pagecolor{mygreen}{0.968}}}$ | ${\mathbf{\pagecolor{mygreen}{0.975}}}$ | ${\mathbf{\pagecolor{mygreen}{0.971}}}$ | ${\mathbf{\pagecolor{mygreen}{0.971}}}$ | ${\mathbf{\pagecolor{mygreen}{0.974}}}$ | ${\mathbf{\pagecolor{mygreen}{0.970}}}$ |
| CSBM | KL | 16 | $\pagecolor{myyellow}{0.855}$ | $\pagecolor{myyellow}{0.739}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{myyellow}{0.893}$ | $\pagecolor{myyellow}{0.890}$ | $\pagecolor{myyellow}{0.807}$ | $\pagecolor{myyellow}{0.855}$ | $\pagecolor{myyellow}{0.806}$ | $\pagecolor{mygreen}{0.953}$ | $\pagecolor{mygreen}{0.934}$ | $\pagecolor{mygreen}{0.951}$ | $\pagecolor{mygreen}{0.937}$ |
| 64 | $\pagecolor{mygreen}{0.936}$ | $\pagecolor{myyellow}{0.893}$ | $\pagecolor{mygreen}{0.955}$ | $\pagecolor{mygreen}{0.952}$ | $\pagecolor{mygreen}{0.959}$ | $\pagecolor{mygreen}{0.934}$ | $\pagecolor{mygreen}{0.962}$ | $\pagecolor{mygreen}{0.940}$ | $\pagecolor{mygreen}{0.966}$ | $\pagecolor{mygreen}{0.967}$ | $\pagecolor{mygreen}{0.963}$ | $\pagecolor{mygreen}{0.969}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.726}$ | $\pagecolor{myyellow}{0.704}$ | $\pagecolor{myyellow}{0.814}$ | $\pagecolor{myyellow}{0.850}$ | $\pagecolor{myyellow}{0.852}$ | $\pagecolor{myyellow}{0.782}$ | $\pagecolor{myyellow}{0.845}$ | $\pagecolor{myyellow}{0.754}$ | $\pagecolor{mygreen}{0.935}$ | $\pagecolor{mygreen}{0.936}$ | $\pagecolor{mygreen}{0.913}$ | $\pagecolor{mygreen}{0.903}$ |  |
| 64 | $\pagecolor{myred}{0.449}$ | $\pagecolor{myyellow}{0.843}$ | $\pagecolor{myyellow}{0.783}$ | $\pagecolor{myyellow}{0.774}$ | $\pagecolor{myyellow}{0.879}$ | $\pagecolor{mygreen}{0.903}$ | $\pagecolor{mygreen}{0.918}$ | $\pagecolor{mygreen}{0.915}$ | $\pagecolor{myyellow}{0.860}$ | $\pagecolor{mygreen}{0.944}$ | $\pagecolor{myyellow}{0.881}$ | $\pagecolor{mygreen}{0.949}$ |  |  |
| $\alpha$-CSBM | KL | 16 | $\pagecolor{myyellow}{0.829}$ | $\pagecolor{myyellow}{0.749}$ | $\pagecolor{mygreen}{0.925}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{myyellow}{0.887}$ | $\pagecolor{myyellow}{0.836}$ | $\pagecolor{myyellow}{0.888}$ | $\pagecolor{myyellow}{0.827}$ | $\pagecolor{mygreen}{0.965}$ | $\pagecolor{mygreen}{0.968}$ | $\pagecolor{mygreen}{0.959}$ | $\pagecolor{mygreen}{0.965}$ |
| 64 | $\pagecolor{mygreen}{0.902}$ | $\pagecolor{mygreen}{0.900}$ | $\pagecolor{mygreen}{0.965}$ | $\pagecolor{mygreen}{0.961}$ | $\pagecolor{mygreen}{0.963}$ | $\pagecolor{mygreen}{0.955}$ | $\pagecolor{mygreen}{0.954}$ | $\pagecolor{mygreen}{0.963}$ | $\pagecolor{mygreen}{0.964}$ | $\pagecolor{mygreen}{0.960}$ | $\pagecolor{mygreen}{0.953}$ | $\pagecolor{mygreen}{0.961}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.810}$ | $\pagecolor{myyellow}{0.712}$ | $\pagecolor{myyellow}{0.841}$ | $\pagecolor{myyellow}{0.887}$ | $\pagecolor{myyellow}{0.877}$ | $\pagecolor{myyellow}{0.821}$ | $\pagecolor{myyellow}{0.854}$ | $\pagecolor{myyellow}{0.819}$ | $\pagecolor{mygreen}{0.951}$ | $\pagecolor{mygreen}{0.947}$ | $\pagecolor{mygreen}{0.912}$ | $\pagecolor{mygreen}{0.930}$ |  |
| 64 | $\pagecolor{mygreen}{0.909}$ | $\pagecolor{mygreen}{0.903}$ | $\pagecolor{myyellow}{0.867}$ | $\pagecolor{myyellow}{0.883}$ | $\pagecolor{mygreen}{0.934}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{myyellow}{0.883}$ | $\pagecolor{mygreen}{0.929}$ | $\pagecolor{myyellow}{0.895}$ | $\pagecolor{mygreen}{0.925}$ | $\pagecolor{myyellow}{0.878}$ | $\pagecolor{mygreen}{0.930}$ |  |  |
| DLightSB-M | KL | 16 | $\pagecolor{mygreen}{0.924}$ | $\pagecolor{mygreen}{0.952}$ | $\pagecolor{mygreen}{0.961}$ | $\pagecolor{mygreen}{0.960}$ | $\pagecolor{mygreen}{0.919}$ | $\pagecolor{mygreen}{0.931}$ | $\pagecolor{mygreen}{0.957}$ | $\pagecolor{mygreen}{0.948}$ | $\pagecolor{mygreen}{0.935}$ | $\pagecolor{mygreen}{0.921}$ | $\pagecolor{mygreen}{0.947}$ | $\pagecolor{mygreen}{0.905}$ |
| 64 | $\pagecolor{mygreen}{0.909}$ | $\pagecolor{mygreen}{0.951}$ | $\pagecolor{mygreen}{0.964}$ | $\pagecolor{mygreen}{0.964}$ | $\pagecolor{mygreen}{0.905}$ | $\pagecolor{mygreen}{0.949}$ | $\pagecolor{mygreen}{0.960}$ | $\pagecolor{mygreen}{0.962}$ | $\pagecolor{mygreen}{0.922}$ | $\pagecolor{mygreen}{0.937}$ | $\pagecolor{mygreen}{0.962}$ | $\pagecolor{mygreen}{0.941}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.787}$ | $\pagecolor{mygreen}{0.944}$ | $\pagecolor{myyellow}{0.870}$ | $\pagecolor{mygreen}{0.920}$ | $\pagecolor{myyellow}{0.743}$ | $\pagecolor{mygreen}{0.921}$ | $\pagecolor{mygreen}{0.944}$ | $\pagecolor{mygreen}{0.950}$ | $\pagecolor{myyellow}{0.723}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{myyellow}{0.890}$ | $\pagecolor{myyellow}{0.850}$ |  |
| 64 | $\pagecolor{myyellow}{0.712}$ | $\pagecolor{mygreen}{0.942}$ | $\pagecolor{myyellow}{0.886}$ | $\pagecolor{mygreen}{0.908}$ | $\pagecolor{myred}{0.686}$ | $\pagecolor{mygreen}{0.915}$ | $\pagecolor{mygreen}{0.950}$ | $\pagecolor{mygreen}{0.937}$ | $\pagecolor{myred}{0.639}$ | $\pagecolor{mygreen}{0.903}$ | $\pagecolor{myyellow}{0.731}$ | $\pagecolor{myyellow}{0.879}$ |  |  |

**Table 4: Trend Score ($\uparrow$) on the high-dimensional Gaussian mixture benchmark. The best-performing method is highlighted in bold, and the second is underlined. Color code threshold: red for $<0.7$, yellow for $[0.7,0.9)$, and green for $\geq 0.9$.**
|  |  |  | $D\!=\!2$ | $D\!=\!16$ | $D\!=\!64$ |  |  |  |  |  |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
|  |  |  | gaussian | uniform | gaussian | uniform | gaussian | uniform |  |  |  |  |  |  |
| Method | Loss | $N\!+\!1$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ | $0.02$ | $0.05$ | $0.005$ | $0.01$ |
| DLightSB | – | – | ${\mathbf{\pagecolor{mygreen}{0.949}}}$ | ${\mathbf{\pagecolor{mygreen}{0.952}}}$ | ${\mathbf{\pagecolor{mygreen}{0.950}}}$ | ${\mathbf{\pagecolor{mygreen}{0.960}}}$ | ${\mathbf{\pagecolor{mygreen}{0.914}}}$ | ${\mathbf{\pagecolor{mygreen}{0.943}}}$ | ${\mathbf{\pagecolor{mygreen}{0.933}}}$ | ${\mathbf{\pagecolor{mygreen}{0.939}}}$ | ${\mathbf{\pagecolor{myyellow}{0.888}}}$ | ${\mathbf{\pagecolor{mygreen}{0.909}}}$ | ${\mathbf{\pagecolor{mygreen}{0.907}}}$ | ${\mathbf{\pagecolor{mygreen}{0.906}}}$ |
| CSBM | KL | 16 | $\pagecolor{myyellow}{0.797}$ | $\pagecolor{myred}{0.662}$ | $\pagecolor{myyellow}{0.876}$ | $\pagecolor{myyellow}{0.854}$ | $\pagecolor{myyellow}{0.809}$ | $\pagecolor{myred}{0.691}$ | $\pagecolor{myyellow}{0.747}$ | $\pagecolor{myred}{0.670}$ | $\pagecolor{myyellow}{0.877}$ | $\pagecolor{myyellow}{0.868}$ | $\pagecolor{myyellow}{0.881}$ | $\pagecolor{myyellow}{0.867}$ |
| 64 | $\pagecolor{mygreen}{0.907}$ | $\pagecolor{myyellow}{0.856}$ | $\pagecolor{mygreen}{0.926}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{myyellow}{0.899}$ | $\pagecolor{myyellow}{0.884}$ | $\pagecolor{mygreen}{0.909}$ | $\pagecolor{myyellow}{0.883}$ | ${\mathbf{\pagecolor{myyellow}{0.888}}}$ | $\pagecolor{mygreen}{0.905}$ | $\pagecolor{myyellow}{0.897}$ | $\pagecolor{mygreen}{0.904}$ |  |  |
| MSE | 16 | $\pagecolor{myred}{0.620}$ | $\pagecolor{myred}{0.633}$ | $\pagecolor{myyellow}{0.739}$ | $\pagecolor{myyellow}{0.781}$ | $\pagecolor{myyellow}{0.743}$ | $\pagecolor{myred}{0.651}$ | $\pagecolor{myyellow}{0.735}$ | $\pagecolor{myred}{0.618}$ | $\pagecolor{myyellow}{0.861}$ | $\pagecolor{myyellow}{0.864}$ | $\pagecolor{myyellow}{0.840}$ | $\pagecolor{myyellow}{0.824}$ |  |
| 64 | $\pagecolor{myred}{0.337}$ | $\pagecolor{myyellow}{0.774}$ | $\pagecolor{myyellow}{0.705}$ | $\pagecolor{myyellow}{0.724}$ | $\pagecolor{myyellow}{0.809}$ | $\pagecolor{myyellow}{0.826}$ | $\pagecolor{myyellow}{0.845}$ | $\pagecolor{myyellow}{0.838}$ | $\pagecolor{myyellow}{0.776}$ | $\pagecolor{myyellow}{0.865}$ | $\pagecolor{myyellow}{0.807}$ | $\pagecolor{myyellow}{0.867}$ |  |  |
| $\alpha$-CSBM | KL | 16 | $\pagecolor{myyellow}{0.772}$ | $\pagecolor{myred}{0.662}$ | $\pagecolor{myyellow}{0.887}$ | $\pagecolor{myyellow}{0.868}$ | $\pagecolor{myyellow}{0.821}$ | $\pagecolor{myyellow}{0.740}$ | $\pagecolor{myyellow}{0.798}$ | $\pagecolor{myyellow}{0.721}$ | ${\mathbf{\pagecolor{myyellow}{0.888}}}$ | $\pagecolor{mygreen}{0.905}$ | $\pagecolor{myyellow}{0.895}$ | $\pagecolor{myyellow}{0.899}$ |
| 64 | $\pagecolor{myyellow}{0.872}$ | $\pagecolor{myyellow}{0.853}$ | $\pagecolor{mygreen}{0.929}$ | $\pagecolor{mygreen}{0.914}$ | $\pagecolor{mygreen}{0.907}$ | $\pagecolor{mygreen}{0.916}$ | $\pagecolor{mygreen}{0.905}$ | $\pagecolor{mygreen}{0.918}$ | $\pagecolor{myyellow}{0.886}$ | $\pagecolor{myyellow}{0.899}$ | $\pagecolor{myyellow}{0.891}$ | $\pagecolor{mygreen}{0.900}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.733}$ | $\pagecolor{myred}{0.621}$ | $\pagecolor{myyellow}{0.759}$ | $\pagecolor{myyellow}{0.822}$ | $\pagecolor{myyellow}{0.790}$ | $\pagecolor{myyellow}{0.714}$ | $\pagecolor{myyellow}{0.771}$ | $\pagecolor{myyellow}{0.715}$ | $\pagecolor{myyellow}{0.864}$ | $\pagecolor{myyellow}{0.857}$ | $\pagecolor{myyellow}{0.824}$ | $\pagecolor{myyellow}{0.829}$ |  |
| 64 | $\pagecolor{myyellow}{0.860}$ | $\pagecolor{myyellow}{0.854}$ | $\pagecolor{myyellow}{0.803}$ | $\pagecolor{myyellow}{0.811}$ | $\pagecolor{myyellow}{0.855}$ | $\pagecolor{myyellow}{0.846}$ | $\pagecolor{myyellow}{0.802}$ | $\pagecolor{myyellow}{0.855}$ | $\pagecolor{myyellow}{0.816}$ | $\pagecolor{myyellow}{0.825}$ | $\pagecolor{myyellow}{0.777}$ | $\pagecolor{myyellow}{0.820}$ |  |  |
| DLightSB-M | KL | 16 | $\pagecolor{myyellow}{0.874}$ | $\pagecolor{mygreen}{0.933}$ | $\pagecolor{mygreen}{0.934}$ | $\pagecolor{mygreen}{0.935}$ | $\pagecolor{myyellow}{0.762}$ | $\pagecolor{mygreen}{0.902}$ | $\pagecolor{mygreen}{0.909}$ | $\pagecolor{mygreen}{0.908}$ | $\pagecolor{myyellow}{0.842}$ | $\pagecolor{myyellow}{0.864}$ | $\pagecolor{myyellow}{0.874}$ | $\pagecolor{myred}{0.665}$ |
| 64 | $\pagecolor{myyellow}{0.857}$ | $\pagecolor{mygreen}{0.928}$ | $\pagecolor{mygreen}{0.935}$ | $\pagecolor{mygreen}{0.937}$ | $\pagecolor{myyellow}{0.747}$ | $\pagecolor{mygreen}{0.906}$ | $\pagecolor{mygreen}{0.909}$ | $\pagecolor{mygreen}{0.911}$ | $\pagecolor{myyellow}{0.828}$ | $\pagecolor{myyellow}{0.865}$ | $\pagecolor{myred}{0.650}$ | $\pagecolor{myyellow}{0.792}$ |  |  |
| MSE | 16 | $\pagecolor{myyellow}{0.703}$ | $\pagecolor{mygreen}{0.920}$ | $\pagecolor{myyellow}{0.821}$ | $\pagecolor{myyellow}{0.892}$ | $\pagecolor{myred}{0.572}$ | $\pagecolor{myyellow}{0.865}$ | $\pagecolor{myyellow}{0.888}$ | $\pagecolor{mygreen}{0.902}$ | $\pagecolor{myred}{0.577}$ | $\pagecolor{myyellow}{0.822}$ | $\pagecolor{myyellow}{0.760}$ | $\pagecolor{myred}{0.548}$ |  |
| 64 | $\pagecolor{myred}{0.629}$ | $\pagecolor{mygreen}{0.915}$ | $\pagecolor{myyellow}{0.838}$ | $\pagecolor{myyellow}{0.879}$ | $\pagecolor{myred}{0.511}$ | $\pagecolor{myyellow}{0.845}$ | $\pagecolor{myyellow}{0.889}$ | $\pagecolor{myyellow}{0.890}$ | $\pagecolor{myred}{0.470}$ | $\pagecolor{myyellow}{0.791}$ | $\pagecolor{myred}{0.497}$ | $\pagecolor{myred}{0.685}$ |  |  |

Figure: (a) Input and Target
Refer to caption: images/benchmark/pairs_d2_g002.png

Figure: (a) Input and Target
Refer to caption: images/benchmark/pairs_d2_g005.png

Figure: (a) Input and Target
Refer to caption: images/benchmark/pairs_d2_u0005.png

Figure: (a) Input and Target
Refer to caption: images/benchmark/pairs_d2_u001.png

## Peer Reviews

confidence: 4
contribution: 1
initial_rating: 2
presentation: 3
questions: #### **On the Benchmark Design and Evaluation Bias**

1.  The paper's primary "winning" solver, DLightSB, is based on a Canonical Polyadic (CP) decomposition. The benchmark's ground truth is *also* constructed using a CP decomposition (Prop 3.1).  How can we be sure that the superior results of DLightSB are not an artifact of it being given the correct inductive bias, a bias that competing methods like CSBM (using an MLP) do not have?

2.  The paper's results show that DLightSB "wins" on a benchmark built on "the same principle used by the DLightSB solver." Does this not confirm a circular evaluation? How would DLightSB perform on a benchmark constructed with a *different* non-factorizable structure, such as one based on a deep neural potential?

3.  The ground-truth benchmark is generated with a low rank ($K=4$), while the DLightSB solver is given a capacity of $K=1000$. Does this 250-fold capacity advantage not make the task trivial for a solver that already shares the problem's functional form? What insights can be gained from a benchmark that is this simple and architecturally-matched to the solver?

#### **On the Choice of Methods and Comparisons**

4.  The remark in Section 2.4 dismisses "tabular EOT" solvers as non-generative. This seems to overlook the entire class of modern EOT solvers that parameterize the dual potentials ($\phi, \psi$) with neural networks, which are both generative and generalizable. Why were these potential-fitting neural EOT methods—which would provide a true "apples-to-apples" comparison against DLightSB—not included as a baseline?

5.  The evaluation compares DLightSB (a *static, dual-potential solver*) against CSBM (a *dynamic, primal-path solver*). Is this not a "apples-to-oranges" comparison? DLightSB solves a simpler, static problem, while CSBM solves a much harder, dynamic one. How can we conclude DLightSB is a better "Schrödinger Bridge solver" from this comparison?

#### **On the Choice of Evaluation Metrics**

6.  The entire purpose of the benchmark is that the ground-truth potential $v^*$ is known. Why, then, were the evaluation metrics (Shape Score, Trend Score) based on comparing *output samples*? Why did the authors not use the most direct and rigorous metric: the error between the learned potential $v_\theta$ and the true potential $v*$?

#### **On the Theoretical Contributions and Framing**

7.  Theorem 3.1 appears to be a restatement of the well-known dual solution to the EOT problem, where one potential is set to 1 by construction. Can the authors clarify the theoretical novelty of this theorem beyond this restatement?

8.  The paper's main innovation appears to be the *engineering* choice of using a CP decomposition to make the benchmark tractable. However, this choice seems to *require* the ground-truth problem to be low-rank ($K=4$). Does this not create a fundamental trade-off where the benchmark can be either (a) computable but trivial, or (b) complex but not computable?

10. The paper's "winning" solvers (DLightSB, DLightSB-M) are admitted to have "severe memory constraints" and become "prohibitive for high-dimensional data." This seems to create a contradictory message. What is the value of a solver that only "wins" on a simple, biased benchmark and is not scalable to the very high-dimensional problems it claims to address?
soundness: 2
strengths: 1. Given a source marginal, and a scalar valued function (what is known in the literature as potentials) it is easy to obtain a corresponding target marginal and the solution SB distribution. The construction is simple and practical and serves as a good synthetic test bed.

2. The parametrization of the scalar valued potential using rank-1 Canonical Polyadic decompositions seems to be versatile (in theory) given that they can act as universal approximators and possess nice theoretical properties.
summary: The paper identifies a valid and significant gap in the machine learning literature: the absence of a standardized, ground-truth benchmark for evaluating Schrödinger Bridge (SB) solvers on high-dimensional discrete state spaces. The authors' stated contributions are:
1.  A method for constructing benchmark pairs of probability distributions $(p_0, p_1)$ on discrete spaces ($\mathcal{X} = \mathbb{S}^D$) for which the analytical SB solution $q^*$ is known by design.
2.  The use of a Canonical Polyadic (CP) decomposition to parameterize a key scalar function ($v^*$), rendering the benchmark construction and sampling computationally tractable in high dimensions.
3.  The introduction of three new or adapted solvers (DLightSB, DLightSB-M, and $\alpha$-CSBM) as "byproducts" of this framework.
4.  An empirical evaluation of these new solvers against the existing CSBM method, using the proposed benchmark, which concludes that the DLightSB and DLightSB-M methods demonstrate superior performance.

Despite the commendable goal, the paper's execution suffers from several fundamental weaknesses in its methodology, evaluation, and framing that undermine its primary conclusions and limit its utility as a general benchmark for the community.
weaknesses: ### 1. Stated Theoretical Contribution: Benchmark Construction

The paper's central theoretical claim is the introduction of a method for constructing ground-truth benchmark pairs for the Schrödinger Bridge (SB) problem on high-dimensional discrete spaces. This method yields pairs of probability distributions $(p_0, p_1)$ for which the analytical SB solution $q^*$ is known by design, enabling rigorous evaluation.

This construction is formally presented in **Theorem 3.1**, which states:

1.  One begins with a known source (initial) distribution $p_0 \in \mathcal{P}(\mathcal{X})$.

2.  One also defines a scalar-valued function $v^*: \mathcal{X} \rightarrow \mathbb{R}$.

3.  The "ground-truth" joint distribution $q*$  is then constructed such that its $x_0$ marginal is $p_0$ and its conditional distribution $q ^* (x| y)$ is proportional to the product of this new function $v*$ and the reference transition kernel $q^{ref}(x_1|x_0)$:

    $$q ^* (x_1|x_0) \propto v^*(x_1) q^{ref}(x_1|x_0)$$,

or 

   $$q ^* (x_1|x_0) =v^*(x_1) q^{ref}(x_1|x_0) / c(x_0).$$

4.  The target distribution $p_1$ is then *defined* as the second marginal of this $q ^* (x_0, x_1) = p_0(x_0) q ^* (x_1|x_0)$.
5.  The theorem concludes that the $q^*(x_0, x_1)$ constructed this way is the true static SB solution between the resulting $p_0$ and $p_1$.


### 2. Marginal Theoretical Novelty

The theoretical novelty of this construction (Theorem 3.1) is only marginal. It is a "reverse-engineered" benchmark based on a restatement of the well-established dual solution to the Entropic Optimal Transport (EOT) problem, which the paper itself shows is equivalent to the static SB problem (Section 2.3).

* **Standard EOT Duality:** The solution $\pi ^*$ to the EOT problem $\min_{\pi \in \Pi(p_0, p_1)} \text{KL}(\pi || R)$ is uniquely characterized by two potentials (Lagrange multipliers) $\phi$ and $\psi$ such that the optimal coupling has the form:

    $$\pi^*(x_0, x_1) \propto \phi(x_0) \psi(x_1) R(x_0, x_1)$$

* **Deriving Theorem 3.1 from first principles:** The paper's construction can be easily derived from the known result above.

    * The paper's *reference measure* $R(x_0, x_1)$ is $q^{ref}(x_0, x_1)$, which is assumed to have the $p_0$ marginal: $R(x_0, x_1) = p_0(x_0) q^{ref}(x_1|x_0)$.
    * Suppose we construct the optimal coupling as $q ^* (x_0, x_1) = \phi(x_0) v ^* (x_1) q^{ref}(x_0, x_1)$, and we want its first marginal to be $p_0$. Integrating out $x_1$ gives $p_0(x_0) = \phi(x_0) p_0(x_0) \sum_{x_1} v ^* (x_1) q^{ref}(x_1 | x_0)$, so canceling out $p_0$ terms gives $\phi(x_0) = 1/c(x_0) = 1/ \sum_{x_1} v ^* (x_1) q^{ref}(x_1 | x_0)$. This is exactly the paper's construction. 

Since one can easily deduce the construction in Theorem 3.1 from a very well-known first principle in SB/EOT literature, it doesn't appear as a significant contribution in this work. 

### 3. The "Single-Potential Problem" Simplification

The paper's methodology *chooses* to construct a benchmark for a simplified theoretical problem.

* The general, two-potential SB problem requires solving for *both* potentials $\phi(x_0)$ and $\psi(x_1)$ to satisfy the two marginal constraints $\pi_0 = p_0$ and $\pi_1 = p_1$.
* The paper's benchmark, the only unknown potential is for the output domain; the input domain's potential is then determined. (Note that the potential $\phi$ on the input space is determined by the potential $v ^*$ on the output domain and the reference transition kernel $q^{ref}(x_1 | x_0)$.)
* This theoretical simplification is crucial because the paper's "winning" solver, **DLightSB**, is *explicitly* a static, dual-problem solver designed to solve this very particular problem. 



### 4. Theoretical Limitation of the Chosen Method

The paper's actual innovation is not Theorem 3.1, but the methodology to make it computationally tractable. The theoretical constructions above are trivial in low dimensions, but in the paper's high-dimensional setting ($\mathcal{X} = \mathbb{S}^D$), the potential $v ^*$ and the normalization constants are intractable, requiring sums over $S^D$ states. 

The paper's *actual* contribution is the application of a Canonical Polyadic (CP) decomposition to parameterize $v^*$. As shown in **Proposition 3.1**, this factorization cleverly turns the intractable normalization sum into a product of sums, making the computation feasible. This is a practical, **engineering contribution**, not a new theory for constructing benchmarks.

This reliance on a CP parameterization for tractability imposes a severe *theoretical limitation* on the benchmark itself: it guarantees the ground-truth problem is trivially simple.

* The computational tractability of the CP method is only achieved if the rank $K$ is small.
* The authors, in their experiments, use **$K=4$** to generate the ground-truth data.
* This means the "challenging" benchmark problem is, by construction, an *extremely simple, low-rank function*.
* This defeats the entire purpose of a benchmark, which is to test a method's ability to solve a complex problem. The paper's theoretical framework *requires* the problem to be simple for the ground truth to be computable.


### 5. Fundamental Methodological Bias and "Apples-to-Oranges" Evaluation
The paper's central weakness is a circular and biased evaluation methodology. The primary conclusion—that DLightSB and DLightSB-M are the superior solvers—is an artifact of the benchmark's design rather than a demonstration of general algorithmic strength.

* **Architectural Circularity:** The tractability of the benchmark (Prop. 3.1) relies entirely on parameterizing the scalar function $v^*$ using a **CP decomposition**. The authors' "new" solvers, DLightSB and DLightSB-M, are *also* explicitly parameterized using the **exact same CP decomposition**. This creates an evaluation that is not a fair comparison of general SB solvers. Instead, it tests how well a CP-based model can fit a function already known to have a CP structure. The authors explicitly acknowledge this circularity: "We attribute this to the benchmark pairs being built on the same principle used by the DLightSB solver".

* **"Apples-to-Oranges" Comparison:** The benchmark problem itself represents a *simplified* version of the full SB problem. The authors' "generative" construction (Theorem 3.1) fixes $p_0$ and the reference dynamics, which effectively reduces the problem to solving for a *single* static potential, $v^*$. The paper then compares two fundamentally different classes of algorithms:
    * **DLightSB (The "Apple"):** A **static, dual-problem solver** whose entire purpose is to fit this single potential $v_\theta$.
    * **CSBM (The "Orange"):** A **dynamic, primal-problem solver** based on Iterative Markovian Fitting (D-IMF), which is designed to learn the *full time-dependent transition probabilities* $m(x_{t_n} | x_{t_{n-1}})$.

This evaluation is inherently unfair. A specialized static solver (DLightSB) that is custom-built for this simplified problem (and given the correct architecture) is predictably going to outperform a general-purpose dynamic solver (CSBM) that is designed for a much harder task.

### 6. Triviality of the Benchmark Task
The reliance on the CP parameterization for tractability forces the ground-truth benchmark to be exceptionally simple, defeating the purpose of a benchmark, which should test a challenging problem.

* The ground-truth benchmark is generated using a CP decomposition with a rank of **$K=4$**.
* The authors' "winning" solvers, DLightSB and DLightSB-M, are then initialized with a capacity of **$K=1000$** components.

The task is thus to fit a simple $K=4$ mixture model using a massively over-parameterized ($K=1000$) model that shares the *identical* functional form. The 250-fold capacity advantage, combined with the correct inductive bias, makes the problem trivial for the DLightSB solvers. A benchmark that is this simple and architecturally-biased provides little insight into how these solvers would perform on complex, real-world problems.

### 7. Flawed Evaluation Metrics and Exclusion of Relevant Competitors
The paper's evaluation protocol is indirect and avoids the most rigorous tests of its own contributions by excluding the most relevant class of competing algorithms.

* **Exclusion of Neural EOT Solvers:** The paper's "Remark" dismisses all "tabular EOT" solvers as "non-generative." This is a strawman argument that ignores the modern, active field of research that *does* learn generalizable, generative couplings by parameterizing the dual potentials ($\phi, \psi$) with neural networks. A fair and informative benchmark would have *required* a comparison between **DLightSB (CP-potential)** and a **Neural EOT Solver (MLP-potential)**. This "apples-to-apples" comparison would have directly tested the paper's only true innovation: the claim that a CP decomposition is a good inductive bias for this problem.

* **Inadequate "Black-Box" Metrics:** By excluding other potential-fitting methods, the paper is forced into the "apples-to-oranges" (DLightSB vs. CSBM) comparison. Because these models learn different objects (potentials vs. dynamics), a direct comparison is impossible. The authors are thus forced to use weak, indirect metrics: **Shape Score** and **Trend Score**. These metrics are inadequate for two reasons:
    1.  **They Don't Measure the *Real* Ground Truth:** The *entire point* of the benchmark is that the true potential $v*$ is known. The most rigorous metric would be a direct comparison of the learned potential $v_\theta$ to $v*$. The paper avoids this.

    2.  **They Don't Guarantee Optimality:** These metrics only measure low-order statistics of the *output samples*. They confirm that the model can approximate the *final marginal* $p_1$, but they provide no guarantee that the model has learned the *true, optimal transport coupling* $q ^*$.

---

confidence: 4
contribution: 1
initial_rating: 4
presentation: 3
questions: * Is there a table that comprehensively summarizes the basic statistics of the GMM benchmarks? (D={2, 16, 64})
soundness: 3
strengths: Based on my reading, the paper has the following strengths. 
* Solid theoretical investigation and a well-grounded evaluation strategy for discrete SB problems. 
* An actual benchmark that will benefit the SB community. 
* Comprehensive empirical experiments and new discrete SB algorithms
summary: This paper presents a benchmark for SB methods on discrete spaces. The contributions of this work seem to be two-fold: (1) a theoretical investigation for constructing a benchmark for discrete SB problems and (2) an actual benchmark along with initial results, including three new implementations of SB algorithms. The authors constructed a practical benchmark based on Gaussian mixtures. The implemented DLightSB algorithms (discrete variants) achieved the best results.
weaknesses: I believe this work in its current format lacks novelty in many regards.

* The theoretical arguments presented in this paper can be seen as either well-known or a reformulation of continuous SB work. I generally view expanding discrete theoretical work to continuous spaces as more challenging and practical, not the other way around.
* The proposed benchmark setup relies heavily on Gaussian mixtures, which were originally intended for representing distributions in continuous spaces, and this tarnishes one of the central contributions of the paper (focusing on the new discrete problem), as the statistics will share many characteristics with continuous problems.
* It is widely known in the SB community that LightSB algorithms excel in Gaussian mixture EOT benchmarks, since LightSB actually has a parametric prior of Gaussian mixture models (see Korotin et al. and Gushchin et al.). Therefore, I am somewhat certain that the authors were highly knowledgeable of this information before the experimental design and engineered DLightSB and DLightSB-M to show similarly good performance on their platform. In other words, these results contain very similar implications to the two LightSB papers, which makes the novelty of this paper weak, since not much new information is given beyond the references.
* The sheer number of benchmarks is lacking, and four important algorithms (DLightSB, DLightSB-M, and α-CSBM) can be generally considered as reformulations from the work of a single research group. The authors are encouraged to put some more effort into implementing other discrete formulations of competitive SB algorithms.

---

confidence: 3
contribution: 3
initial_rating: 6
presentation: 3
questions: Why do CSBM algorithms have better performance on the benchmark with higher dimensionality? I would expect that performance would degrade with higher dimensions given a fixed model capacity but perhaps because of the low dimensional structure in the benchmark this is not the case?

Wouldn't it be possible to create a metric that doesn't just take into account single marginals or pairs of marginals by looking at distribution metrics? For example you could take a given starting point x0, sample p(x1 | x0) from the learned SB and the true SB then compute a distributional metric between those two sample distributions. Then you could average these values over different draws of x0.
soundness: 3
strengths: The use of the mixture of factorized functions for the scalar tilting function looks novel to me and is an interesting way to construct a non-trivial problem for which the schrodinger bridge is known. The same factorization technique is then also used to construct the Light Schrodinger algorithms which show very good performance compared to prior bridge matching approaches. It is important to have these Light Schrodinger bridge baseline algorithms for further comparison during algorithm development as they serve as useful reference points.
summary: In this paper, an analytic benchmark for discrete optimal transport is developed. In order to assess capabilities of current algorithms for solving the schrodinger bridge problem on discrete state spaces, there is a need for a benchmark for which the ground truth solution is known. The authors construct this analytic solution by writing the optimal transition distribution as a scalar function multiplying the reference transition distribution. Then a specific mixture of factorized forms is proposed for the scalar function such that the true schrodinger bridge transition distribution becomes analytically calculable. The authors propose two new light schrodinger bridge algorithms for discrete state spaces and compare them with existing methods on their benchmark.
weaknesses: The benchmark seems to already be too easy because the DLightSB method is already almost saturating the performance metrics. This means the benchmark is already losing its discriminative power between the best methods. Ideally, the benchmark would contain sub-tasks for which all methods struggle. Could the task be made more difficult by increasing the number of mixture components in the scalar function?

In the introduction in line 47 you state that one motivation for introducing this benchmark is to check whether methods actually solve the SB problem separating out 'true algorithmic performance from artifacts of specific parameterizations'. However, then on L431 you hypothesize that DLightSB performs the best on this benchmark because it is built on the same factorization principle as the benchmark itself. This would then seem to go against the stated goal of the benchmark to try and not be influenced by specific parameterizations of the method. Perhaps the benchmark needs to be made more complex as to not be 'defeated' by the DLightSB algorithm.

The differences between your two proposed methods DLightSB and DLightSB-M could have been investigated more. You hypothesize that DLightSB-M works worse due to error accummulation. This could have been verified or falsified with a further experiment, for example, changing the number of discretization time steps in the diffusion process to see if the discrepancy between the two methods gets smaller.

Equations (19) and (20) are quite hard to parse due to the circular definition of the set \mathcal{S}. The set itself is used within its own definition in the argmax and makes the whole equation hard to understand.

Overall, I think the contribution of the rank-1 Canonical Polyadic decomposition for constructing both ground truth SBs and for learning LightSBs is a worthy contribution to the conference notwithstanding the limitations in evaluations of the proposed methods.

---

confidence: 4
contribution: 2
initial_rating: 4
presentation: 2
questions: 1. Could you provide an additional section (like section 4.4 in CSBM) to verify the scalability?

2. i didn't get line 352: "be generalized to an arbitrary reference process qref, thereby enabling the application of the optimal projection in discrete space settings." If some algorithm accuracy is compromised (like low rank), does it still guarantee optimal projection? it is a bit confusing. 

3. I haven't fully checked CSBM yet. Is the algorithm simulation-free or not? To me, low-rank may be interesting but the simulation-free property is more important for the scalability.
soundness: 3
strengths: 1. Providing a mixture example to evaluate the accuracy of different approaches is helpful to the community.

2. The light (low rank) and alpha smoothing do have some potential to accelerate the algorithm.

3. IMF is more appealing and interesting than the IPF style of algorithms.
summary: The authors proposed a mixture-like benchmark for SB on discrete spaces to evaluate the effectiveness of different methods. The authors also proposed a light version of DSB and a smoothing version of the CSBM algorithm, α-CSBM,  to enhance the training stability and efficiency. Interesting comparisons are made on toy examples to show the strength and weakness.
weaknesses: 1. Even the "expensive" version of the CSBM algorithm is evaluated on real-world image experiments. The light and alpha-smoothed algorithm is only evaluated on the mixture benchmark. This is far from sufficient.

2. The writing is not good. The current draft is a mixture of lots of components, like a baby mixture benchmark, 2 light algorithms (without sufficient discussions on theoretical properties or empirical evaluations). The scope is ambitious, but for each contribution itself, it is a bit minor.

3. Although low0


Decision: Accept (Poster)