Finite-Time Convergence Analysis of ODE-based Generative Models for Stochastic Interpolants
Authors: Yuhao Liu
Yu Chen
Rui Hu
Longbo Huang
Venue: ICLR 2026 Poster (2026)
Keywords: Stochastic Interpolants
Diffusion Models
Deterministic Sampler
Primary area: generative models

Abstract: Stochastic interpolants offer a robust framework for continuously transforming samples between arbitrary data distributions via ordinary or stochastic differential equations (ODEs/SDEs), holding significant promise for generative modeling. While previous studies have analyzed the finite-time convergence rate of discrete-time implementations for SDEs, the ODE counterpart remains largely unexplored. In this work, we bridge this gap by presenting a rigorous finite-time convergence analysis of numerical implementations for ODEs in the framework of stochastic interpolants. We establish novel discrete-time total variation error bounds for two widely used numerical solvers: the first-order forward Euler method and the second-order Heun's method. Our analysis also yields optimized iteration complexity results and step size schedules that enhance computational efficiency. Notably, when specialized to the diffusion model setting, our theoretical guarantees for the second-order method improve upon prior results in terms of both smoothness requirements and dimensional dependence. Our theoretical findings are corroborated by numerical and image generation experiments, which validate the derived error bounds and complexity analyses.

Published as a conference paper at ICLR 2026
FINITE-TIME CONVERGENCE ANALYSIS OF ODE- BASED GENERATIVE MODELS FOR STOCHASTIC IN- TERPOLANTS
Yuhao Liu Yu Chen Rui Hu Longbo Huang"
III, Tsinghua University
[liuyuhao21, chenyu23, hu-r2uuurraiiiii tsinghua. edu. cn longbohuangets iighua edu.
ABSTRACT
Stochastic interpolants offer a robust framework for continuously transforming samples between arbitrary data distributions via ordinary or stochastic differen- tial equations (ODEs/SDEss, holding significant promise for generative modeling. While previous studies have analyzed the finite-time convergence rate of discrete- time implementations for SDEs, the ODE counterpart remains largely unexplored. In this work, we bridge this gap by presenting a rigorous finite-time convergence analysis of numerical implementations for ODEs in the framework of stochastic interpolants. We establish novel discrete-time total variation error bounds for widely used numerical solvers: the first-order forward Fuler method and the second-order Heun's 's method. Our analysis also yields optimized iteration complex- ity results and step size schedules that enhance computational efficiency. Notably, when specialized to the diffusion model setting, our theoretical guarantees for the second-order method improve upon prior results in terms of both smoothness re- quirements and dimensional dependence. Our theoretical findings are corroborated by numerical and image generation experiments, which validate the derived error bounds and complexity analyses.
INTRODUCTION
Stochastic interpolants (Albergo & Vanden-Eijnden, 2023 Albergo et al. 2023) provide a powerful framework for constructing generative models by learning deterministic stochastic transformations that continuously map samples from an initial distribution to a target distribution Pl, governed by ordinary or stochastic differential equations (ODEs/SDES). To determine the differential equation, the approach constructs stochastic interpolations between Pu and #1 samples, then estimates a mean velocity field from these paths. With learned approximation of the velocity field, one can build generative model by solving the approximated ODE or SDE. This novel framework unifies flow matching (Lipman et al. 2023) and score-based diffusion (Song & Ermon 2020 Ho et al. 2020, Song et al. (021). offering significant design flexibility through its various choices of the initial distribution and generalized interpolation formulation, making it an important subject for generative
modeling (Chen et al. 2024b Albergo et al. 2024a Chen et al. 2024a).
Recent theoretical work has focused on understanding the convergence properties of stochastic interpolants. For transformations based stochastic differential equations (SDEs). theoretical guarantees on generation error bounds have been established for both continuous-tima (Albergo et al. 2023) and discrete-time (Liu al. 2025) scenarios. Continuous-time analysis focuses the exact solutions of the underlying equations, while discrete-time analysis further considers the numerical approximation used in practice since the equations usually cannot be analytically solved. However, for their deterministic counterparts, ODE-based transformations, the analysis has been limited to the continuous-time setting (Albergo & Vanden-Eijnden, 2023 Benton et al. 2024b), despite their widespread use in practical applications. Therefore, we address this critical gap by investigating
"Corresponding author.

---

Published as a conference paper at ICLR 2026
the following question: What is the iteration complexity for numerical ODE solvers to ensure accurate generation with stochastic interpolants?
This problem is important yet technically challenging. Although there have been insightful analyses for diffusion models and related approaches that can transform Gaussian to target distributions (c.g..
al. 2025b 2024b Huang al. 2025), they cannot be directly extended to stochastic interpolants due to the general data-to-data transformation structures which make the interpolant non-Markovian. While Liu al. 2025) recently established finite-time error bounds for the SDE case, their methods face inherent challenges in ODE cases due to the singular behavior in the degenerated process.
This work bridges the gaps by presenting the first systematic analysis of discretization schemes for stochastic interpolant ODEs. Specifically, we provide a complete theoretical characterization of both first-order (forward Euler) and second-order (Heun's 's) methods applied stochastic interpolants (Theorem 4.5 and 5.4). We explicitly quantified the discretization error introduced by numerical solvers under mild assumptions, and then derive dd² and O 0( iteration complexity bounds, respectively, for achieving distribution approximation error on d-dimensional data.
key technical innovation is the introduction of carefully designed continuous-ting interpolations. These interpolations reformulate discrete-time iterations into ODE forms comparable to the true process, which enables the use of continuous-ting analytical tools to quantify the approximation error by comparing the induced drifts and divergences with those of the true process.
Notably, our results establish a general analysis framework valid for every stochastic interpolant process satisfying mild regularity assumptions, including the famous surrince--resservie diffusion (Song et al., 2021). Furthermore, our dedicated analysis on the Heun's 's method achieves state-of- the-art iteration complexity results on diffusion models (see Section [.. and Appendix [B). This advancement arises from detailed error analysis, where carefully partition the drift and divergence errors into tractable expectation terms, yielding an improved error bound (check Appendix
B. for further comparison).
The contributions of this paper are summarized below:
We provide the first finite-time analysis for stochastic interpolant ODEs. We derive general TV error bounds for both first-order (Euler's) s) and second-orde (Heun's) s) methods, which are applicable to all stochastic interpolant process satisfying mild regularity assumptions.
We explicitly quantify the discretization error introduced by numerical ODE solvers, es- tablishing Oc=-ld2) and 0(e-1/2,3// iteration complexity bounds for both methods, respectively. Especially, when degenerating to the variance-preservie diffusion model, our result for Heun's 's method improves the existing analyses on diffusion models.
We employ carefully designed continuous-time interpolation to unify the comparison be- tween the discrete-time approximation and continuous- target processes, allowing us to quantify the approximation error by comparing the induced drifts and divergences. We also corroborate our findings with numerical experiments and validations on image generation.
RELATED WORKS
2.1
STOCHASTIC INTERPOLANTS ANALYSIS
The stochastic interpolant framework originates from continuous- normalizing flows, providing a principled approach for constructing data-to-dal generative models (Albergo & Vanden- Fijnden 2023). Subsequent work (Albergo et al., 2023) extended this framework through the incorporation of Gaussian perturbations, enabling learnable score functions and facilitating the application of stochastic differential equations for data transformation tasks.
Regarding ODE-based formulations, Albergo & Vanden-Eijnden 2023 established Wasserstein error
bounds for velocity field estimation under Lipschitz continuity assumptions. Benton et al. (2024b) advanced these results by considering time-dependent Lipschitz constants, thereby obtaining tighter error bounds. Their analysis further examined the control of Lipschitz constants specifically for linear interpolants. Both studies focused exclusively continuous-time settings.
2

---

Published as a conference paper at ICLR 2026
In the SDE context, Albergo al. (2023) derived continuous-time Kullback-Leible (KL) divergence bounds between target and estimated distributions, expressed in terms of the mean squared error of drift estimation. Liu et al. (2025) made progress by establishing finite-time error bounds for the Marrr--aruyama discretization scheme, representing the first discrete-ting analysis within the stochastic interpolant framework. They also investigated the impact of schedule selection.
2.2 ODE-BASED DIFFUSION MODELS ANALYSIS
Recent theoretical work has made significant progress in analyzing the convergence properties of probability-flow ODE diffusion models. Chen et al. (2023b) provided the first polynomial convergence bounds for DDIM-type samplers, Chen ct al. (2023a) developed a provably efficient sampling algorithm incorporating corrector steps within the ODE solver framework, though this approach introduces additional stochasticity to the process.
Recent theoretical analyses purely deterministic generation typically adopt two principal ap- proaches. The first directly examines discrete-time density evolution, where Li et al. (2024b) established a foundational framework for controlling the total variation (TV) distance between target and estimated distributions. Subsequent works (Li et al. 2024a 2025a Liang al.) 2025) extended this framework to derive tighter bounds and accelerated convergence guarantees. The second ap- proach considers equivalent continuous timee processes through partial differential equation (PDE) analysis. While both Huang et al. (2025) and (Li et al.) 2025b) employed PDE techniques to bound TV error, the latter achieved superior results through more refined error decomposition methods.
Current theoretical understanding indicates that first-order methods only require (d/=) iterations to achieve e-precision in TVV distance (Li et al. 2024c). However, existing analyses of higher-order methods remain either non-tight reliant overly restrictive assumptions, highlighting the need for more comprehensive theoretical frameworks for ODE-based approaches. See Appendix B and Table
D for a comparison between our results and existing works.
BACKGROUND ON STOCHASTIC INTERPOLANTS
Consider two probability distributions Po. P1 A (IRd). The stochastic interpolant framework
(Albergo Vanden-Eijnden, 2023; Albergo al. 2023) provides general approach for constructing
learnable mapping between Pg and 01: by identifying learnable vector field b(t,z: such that the solution of ODE (see Appendix A for SDE cases)
dX, = b(t, X,
satisfies X1 ~ P1 under initial condition X, ~ Po-
Specifically, z(t, is defined the average velocity b(t,z) Eiic/2, æ] of stochastic interpolant (#))ee[o.1] between Po and pr. given by:
I(t, 2t. x1) + (t)2,
where
(Io, F1) with marginals xo ~ Po. 01 ~ Pl. and NO, Ia) is an independent noise.
I is C2-smooth interpolation satisfying I(0,, Zo, I1) zo and [11, 20. 21) Z1.
7(1) > 0 for E (0, 1), so (()= serves to smooth the marginal distributions.
handle both paired and unpaired data, given (æo, 21) are dependent independent, respectively. Let p(l) denote the law of xe and p((,,) )) its density. Typically, we choose 70) >1) 0 so that p(0) Po and p(1) = Pr.
Notably, the coupling
Under this definition of bt,, the solution X, of ODE satisfies X, ~ p(t) for t E [o, 1]. In practice, since b(t, æ) is generally
analytically computable, usually use neural networks to obtain an estimate b(t, (see Appendix A for more details). By applying numerical methods to solve the ODE with b(t, x:), we can then transform samples from Po to P1.
Remark: The widely used diffusion model can be viewed special case of stochastic interpolants. For example, consider the linear interpolant 3; t)xr tan + 22 f)2 with vn Po
3

---

Published as a conference paper at ICLR 2026
(N0, Ia) and 21 Pdata- In this case, the marginal distributions of as coincide with those of the pariance-preservinn (VP) diffusson model (see, e.g., Song et al. 2021).
foundational methods: the Euler's method (Section 4 and the Heun's ' method (SectionS). We derive finite-time TV error bounds for both methods (Theorem 4.5 and 5.4) respectively, and compare their convergence rates.
Below, conduct novel analysis for stochastic interpolants based
THE FIRST-ORDER RESULTS FROM THE FORWARD EULER METHOD
We start by analyzing the 3DE-based stochastic interpolants with the usage of the forward Euler method, an important first-order scheme for solving ordinary differential equations. We derive finite-time error bound for its numerical solution, which further allows us to analyze the iteration complexity.
4.1 ERROR BOUND FOR THE FORWARD EULER METHOD
To begin with, we formally define the iterations of the forward Euler method as follows. Consider time discretization given by the schedule {4 where to 6 < 02 IN. Given an
initial condition X the forward Euler method approximates the solution of ODE (*) itcratively by (8 Specifically, at iteration (A: + 1), we set
- = Xta + hy B(tt. Xtx).
(Euler)
where i(t, a)) is an estimator of the drift function b(t, vr), and we use h4 = 4++1 4 to denote the quantify and control the difference between the approximated terminal
step size. Our objective distribution p(tN) given by XIw and the true terminal distribution p(tx).
We make the following assumptions to rigorously analyze the forward Euler method.
Assumption 4.1 (Regularity). Let be the joint distribution of Po and p1. we assume that:
The fourth moment of the difference between the initial and goal distribution is bounded: E (xo.ar) [| =1|1*] 00.
The p-th order time derivatives (for p= 1,2) of T(t,, za, x1) and 7.2 are bounded w.r.t. t: There exist positive constants C₁, C. > O such that for all Fo. I1 E TRd and E {1,2).
d" dip
(1, S Cr| xill,
[2(0]]
SCy.
Assumption 4.1 ensures the regularity of both the stochastic interpolant process (z1)aε(0,1) and the function b(t, 2). For example, with the linear interpolant 3; = (1 t)zo + tæ + 2 - f)2, this assumption is met if Po and P1 both have finite fourth moments.
Assumption 4.2 (Drift Estimation Error). Let (x)) = [lō(tz, x)) b(ta, x) represent the error between the estimated drift coefficient b(tx, x)) and the true drift coefficient b(tx, x) at time ta and state z. We assume that
Σ R=0
hi
E mag <p((tt)
[21,)(
S EJnn < 00.
Assumption 4.2 is a standard condition for estimator quality, similar to those in (Benton et al. 2024al Liu et al. 2025). This assumption differs from continuous-time conditions by providing discrete- time formulation that focuses on the error at the time steps used the numerical method.
However, unlike SDE- based models, Assumption 4.2 alone is insufficient to fully control the distribu- tion error (Li et al. 2024b 2025b). Therefore, introduce further assumptions b(t,
4

---

Published as a conference paper at ICLR 2026
Assumption 4.3 (Jacobian Estimation Error), (2, = vb(te. Vb(tk, x)
we assume
hk R=O Sey
(22,(((e))) Ediv < co.
mp(tx)
Assumption 4.4 (Smoothness). b(t, x) is w.r.t. that for all k = =0,,, 00,,,..., ....,N and all n. E TR..
w.r.t. x. Furthermore, there is a constant L >> 0 such
(b(tes=)) S (28((z. L IV b(t,z)) V b(t,)))l S 1// 1 yll
Assumption 4.3 extends the estimation accuracy requirement of Assumption 4.2 to the Jacobian of bl,,x). Similar assumptions are used by Li et al. (2024a) and Li et al. 2025b Assumption 4.4 imposes uniform Lipschitz constant on both b(t, and its divergence. Its reasonableness is illustrated by cases where the data for Po and #1 are bounded in each dimension; in such cases, Lemma F.7 in the Appendix shows that L is of order (d))
With the preceding assumptions established, we are now ready to present the main theoretical result concerming the application of the forward Euler method.
Theorem 4.5. Under Assumptions 4.1
and 4.4 suppose the forward Euler method is initialized with Xto ~ ((to),, and the step sizes satisfy h2 S 21. Then,
TV(p(IN), TV(p(()), ((oo)) Ediv + Eani (a)/35(7, l0 12) + 12
d 7/ M
Discretization Error
Here, we have defined:
7/k
((1)
f.o 7 -t)dtt and M
max d, L, / E. 1 - xn [1
Theorem 4.5 provides a comprehensive upper bound on the total variation (TV) distance between the true target distribution p(t,) and its approximation p(t,v) obtained via the forward Euler method. This result quantifies the convergence behavior of the Euler scheme, and elucidates how several factors influence the approximation error: the initialization error (p((o)) (((x)), dimension the
distance between source and target distribution captured by E., x1 [['). the Lipschitz constant N-1
L, the estimation errors Edrift and Ediv, and the step size schedule { (hx))
The step sizes (hx N= specifically influence the last term discretization error in Theorem 4.5 which quantifies the influence of using the discrete-time numerical solver. The other terms are independent with the step size schedule, reflecting the inherent quality achievable with an infinite number of steps. Since the discretization error term is of order O(h2) for each step, summing them up will yield error bound proportional to the step size (see Section 4.2 for detailed analysis when the schedule is specified), so it approaches zero as the step sizes tend to zero. In addition, the O discretization error per step coincides with prior works for diffusion models with first-order samplers (Li et al. 2025b Huang al. 2025).
To avoid an unbounded right-hand side in Theorem 4.5 when (0) = or ((1) (1) = (which ensures p(0) Po and p(1) P1 in the stochastic interpolant definition), we follow Liu et al.) (2025) and simulate the process within a sub-interval [fo, [N] C (o, 1). This means our sampling starts from an estimation of p(to) rather than Po. and similarly aims for estimation of p(tx) instead of pj. When t close to 0 and fny close to 1, p(to) and p(tx) are close to Po and P1, respectively. In addition,
the initialization error can be made very small (e.g., when forcing I [(tg 2g. 21) .g)) This technique
also known as early stopping in the context of diffusion models (Song et al. 2021 Hongrui et al.)
2023 Benton et al. 2024a).
5

---

Published as a conference paper at ICLR 2026
Main Ideas of the Proof We sketch the main ideas of the proof of Theorem 4.5 with the detailed proof in Appendix D Since we are comparing a continuous-tin time process with a discrete- time estimator, the first key step is to bridge this gap. We do this by reinterpreting the discrete-time process {X J
-eefiid continuous process. Specifically, we define
X
Fb, -6 -o (X.,) == Xt. + (t tx tx) Xt.)
for & E [1xs2+t)) Under Assumption 4.4 and for sufficiently small step sizes hi, the mapping Fa-r:: : x + x (1 Lx 1x)6((k x)) a diffeomorphism on R&d, which allows us to write dX, - b(t,, Xe, )d b(4,, (X, By defining the effective drift b(t, x) == b(4x, F we obtain dx, b(L, X2), which is of the same form the ground truth ODE.
To quantify the similarity between the two processes, the main technical challenge is to show that the approximate drift b(t,z) closely matches the true drift b(t, x), as well as their respective divergences. Let (Xe) denote the ground truth process. We apply the following sequential approximation scheme to decompose the errors for both the drift and the divergence (see Appendix D.2 and D.3):
b(4, Xx) blle, X Re UI, Xt),
(Drift Error)
(() = tr Vo(te: (X.))
VF
= tr Vb(te. X Id] v b(t,
(Divergence Error)
Combined with the analysis on partial differential equations (see Lemma D.2 we can then control the overall TV distance by the controlled drift and divergence errors.
4.2 COMPLEXITY OF FORWARD EULER METHOD
Theorem 4.5 built on a general stochastic interpolant framework, which is available to control the convergence rates for every valid stochastic interpolants (i.e., valid / and 7) under necessary assumptions. Below, we provide a detailed analysis on the convergence complexities for some widely used stochastic interpolants.
Stochastic Interpolants with (t) = / /at(1- -) For this case, (t = at - -t) is the variance of Brownian bridge process, which is popular choice (Albergo et al. 2023 De Bortoli et al. 2024).
To optimize the schedule, we focus on minimizing the discretization error in Theorem 4.5
77k
dd
Fx 2M
The error is proportional to h/ which suggests setting hg &h to balance the error. Inspired by the tine-scheduling method built for SDE in Liu et al.) (2025). we construct the following schedule: define midpoint m at tm. 0.5 and select a step size scale parameter 0. The schedule (tx N then constructed by exponentially decaying on both sides:
{ 4 (1 1.)m-k
KS m:
to 1
(1-1)*
k
>> m.
Schedule i ensures hy T Oh)) and the total number of steps N N= log [o(l-()) By
))
setting 8 min(to, 1 1.]. Theorem 4.5 directly gives the error bound for this scheduling:
TV(p(t)),((ww)) & (V(pttu), p(to)) + Eaiv + Eann lo./"log//2(1// + L122
+ h dd (og((10)) + 2
Discretization Ervor
indicating that achieving e-T'V emor requires N = 0 [d2 (o + M² log(1/5)) steps.

---

Published as a conference paper at ICLR 2026
Variance--Preerrvnng Diffusion Models Consider stochastic interpolant defined follows:
sample 21 Pdata and (0, La), then set tz + v1 11--22 f2z. Then 2, (tr1,, (1 - La) when #1 is given, so the marginal distributions p(t) of the stochastic interpolant match those of the "ariance-preserving diffusion model (under change of parameter,
e.g., Albergo al. (023))
Similar to the previous case, we aim for h4 x 7² =
= 6 - t). This leads to the schedule:
(1 h). 0 s k S N;
where the scale parameter controls the step sizes.
Let ty be the early-stopping time, the number steps satisfies N = log 4). If we further assume data support bounded in each dimension, Theorem 4.5 yields that the forward Euler method with this scheduling requires O dd log )) iterations achieve e-error, matching the prior theoretical findings (Li al. 2024b 2025b 2024ch. Noting that Li et al. (2024c) reports
faster (d/e) complexity, their result relies on specific analysis that is restricted to Gaussian-to-data
scenarios. Theorem 4.5 enables the extension to the more general data-to-data cases.
THE SECOND-ORDER RESULTS FROM THE HEUN'S METHOD
This section focuses on deriving finite-time error bound and the corresponding iteration complexity for Heun's method, a widely used second-order approximation for solving ordinary differential equations (ODEs) (see, e.g.. Stoer et al. 1980). We will then compare the error and convergence rates of Heun's s method to those of the forward Euler method, highlighting the Heun's 's method as an acceleration of the latter.
5.1 ERROR BOUND FOR THE HEUN'S METHOD
Similarly, we first provide formal definition of the Heun's 's method. Given a sequence of discrete time steps (tx - let Xx, denote the estimated solution at time fg. Beginning with an initial condition X1. the method proceeds via the following iterative scheme:
Xie+s = 822 + h4 Xet),
(Heun)
Xixm = 8ez + hx
+ b(tx+1, X
Similar to the first-order case, aim to quantify and control the error between the estimated terminal distribution and the true terminal distribution. However, to accommodate this higher-order method, must modify some of our initial assumptions in Section 4
Assumption 5.1 (Regularity). Let v be the joint distribution of Po and pr. we assume that:
The sixth moment of the difference between the initial and goal distribution is bounded: E ( [| - 10] co.
The time derivatives [tt, x1) and (1) are bounded up to the third order: There exist positive constants Cr, C. > O such that for all xo, 21 E Rd and p E {1, 21,,33,
[0,0 x1)1 S Cyllær xill,
d"
[2(00]
s Cyy-
dtv
Assumption 5.2 (Drift Estimation Error). Let (0()),u(x) [1b(te,)) b(tR, x))l-
Then we assume
N-1
h&E [1, (1 S Egang x0.
Assumption 5.3 (Jacobian Estimation Error). Let (x) = vo(te.z) (6(t).2)
We assume
h (2t. 62,441 221//2 < Edio 00.
0=

---

Published as a conference paper at ICLR 2026
Assumptions 5.1 5.2, and 5.3 are similar to Assumptions 4.1 1
4.2 and 4.3 respectively, while they
are adapted for the higher-order Heum's method. Assumption 5.1 requires a stricter regularity since are analyzing higher-order integrator. Assumptions 5.2 and 5.3 incorporate additional terms, reflecting the fact that Heun's 's method performs two evaluations b(t, per step. is worth noting that Assumption 5.3 utilizes 2, instead of E2, which is slightly different with Assumption 5.3
With the preceding assumptions established, we now present our main theorem for Heun's 's method.
Theorem 5.4. Under Assumptions 5.7 5.2 5.3 and 4.4 if the equation is solved using Heun's method with an initial condition X ~ ((to), and provided that the step sizes satisfy h& S
min 12 Bl - zıll -1111)-1/ 1/3 d 222) then the TV error is bounded as:
(V(p(tN). & TV(p(to),P(ta) Eai + sarin (d to, 12) + 1 1122)
N d
+ + Md]
Discretization
Here, we have defined 7 == infre r(t), S(y,to,ty) = for Y (t)dt, and M ==
max { (d, L, E, [||00
[6]] 1/3
21
Similar Theorem 4.5 Theorem 5.4 bounds the TV error the target distribution. The primary distinction lies in the discretization error term, where the order the step sizes is enhanced from h to hg. reflecting the higher-order accuracy. The O error for each step aligns with prior analyses on diffusion models (see, e.g., Huang et al. 2025).
Main Ideas of the Proof The proof of Theorem 5.4 (see Appendix E shares similar ideas as the first-order case. central novelty lies the following construction of the continuous-tims interpolation for the discrete Heun updates:
X. = Geeene(
(t - tx) = 11x + (t- (t-tx) 2(tx+1 tx)
(t - tx)² b(tb+1, Fi Fic-ue (81x))
b(tk, Xex) +
2(tx+1 tx)
f-11-t b(tk) Xxx) Te-r 1- Fir Fiee (Xew)). This interpolation yields velocity di 1041-14 which linearly interpolates between the drift estimates and achieves higher-order accuracy. Unlike some alternatives (e.g., Huang al, (022)). our method only requires b(t, .)) at discrete timesteps.
d dx,
For small hy, the mapping Gt se is similarly a diffeomorphism, enabling the representation d Xx X.) := o,G -ve (G - (X.)). As in the first-order case, we establish that ō(t, ar) and its divergence closely approximate their ground-truth counterparts, ultimately enabling controlling the TV error. The linear velocity interpolation yields an O(h)) approximation error, improving the previous Ohh)) error for Euler's method. The detailed decomposition is omitted here for brevity.
at
5.2 COMPLEXITY OF HEUN'S s METHOD
Similar to the forward Euler case, Theorem 5.4 can be applied to all stochastic interpolants satisfying necessary assumptions. Below, we analyze the complexity for the same stochastic interpolants discussed in Section 4.2 with the second-order Heun's 's method replacing the forward Euler method.
Stochastic Interpolants with ((1) = Val(1 () The discretization error in Theorem 5.4 is pro-
portional to h Hence, balancing the error suggests h x 5/ again, Icading to the exponentially min(to, tN). Applying Theorem 5.4 shows that achieving
decaying schedule pp. Let
3/2 (1/5) M3/2 log(1 /8))) steps. This represents an s-VV error requires N O TE [ 3/2
reduction in complexity from (11/s) (for forward Fuler) 0(1/Vo).
8

---

Published as a conference paper at ICLR 2026
A
B
Hort D.kA
Po
lock
2.2
P1
(a) 2D density transforms, labeled task A/B/C, respectively.
(b) TV step size h.
(c) TV squared step size h².
Figure 1: Empirical verification of convergence rates for numerical methods.
Euler
Recl Scint
Pud J067
Heun
2-10. Hour
p-15,
(a) Empirical TV error v.s. h. (b) Error growth w.r.t. d (Fuler). (c) Error growth w.r.t d (Heun).
Figure 2: Empirical TV error for d-dimensional Gaussian mixtures.
SriianccePreserving Diffuuson Models Recall that the interpolant o, tæ + -iv 111--1- P2 shares the
same marginal distributions as puriance-preseevvie diffusion models. Using our proposed schedule 2
we derive a complexity of N = O 55-1/2
1/2 d/2 log 3/2 (1/6) + 3 11 log(1/8) where 5 = 1 tw
the carly-stopping time and M max L, E 1 /91/3). When the data support is bounded in each dimension, which is common scenario in practice, M is typically (d), which simplifies the complexity to N 0 3//2 log 3/2 //2(1118/)) (1/5) This (s-11/2 convergence rate consistent with Huang et al. 2025 Notably, our dependence the dimension is
prior works Li et al.
2024a
lower than the O 2 complexity reported by
Li al. 2025a) for second-order method.
While Huang al. 2025) also achieve comparable 0(2-1/2 complexity, their analysis relies on stricter assumptions. Specifically, they require uniform boundedness the third-order derivatives of the score estimator s(t, a:) with respect both and xr, whereas we only assume bounded second- order derivatives of b(t, z)) with respect to æ. This significant distinction, especially since the true drift function's 's I-derivative can be unbounded. Additionally, their analysis considers errors over the entire continuous interval, while ours focuses on the discrete time steps used by Heun's 's method. See Appendix for further comparison.
NUMERICAL EXPERIMENTS
This section presents experimental validation of our theories. Using the interpolant x, = (1- f) tæ + 12a((1 f) and schedule 11). we evaluate both forward Fuler and Heun's 's methods on datascts and d-dimensional Gaussian mixtures. Our experiments characterize the TV error growth w.r.t.: (i) step size scale and (ii) dimension
We evaluate framework on three 2D dataset pairs from Grathwohl et al. 2019). We have designed three different generation tasks, as shown in Figure la Using neural network to estimate b(t. x). compare the forward Euler and Heun's 's methods mitiialized p(to). Figures and and Ic show the empirical TV distances between p(tN) and p(tx). ). confirming our theoretical complexity analysis: the discretization error bounds (hh) (Fuler) and ( (Heun).

---

Published as a conference paper at ICLR 2026
For d-dimensional Gaussian mixtures where b(t, admits analytical solutions (see Albergo et al. 2023 we evaluate the empirical TV eπror without model training. Figure 2a demonstrates the error growth rate versus step size h, confirming
theoretical convergence rates. Figure ab [c further examines the dimensional dependence at fixed h. While our theory establishes ( (Fuler) and O(h2d²) (Heun) error bounds, empirical observations suggest linear growth in both cases. This theory-experiment gap indicate that the current bound could be further improved, warranting future investigation.
In addition, to further demonstrate the applicability of our theory, we evaluate both Euler's and Heun's methods on real-world image generation tasks, as detailed in Appendix C
7 CONCLUSIONS
This paper presented finite-time analysis of discrete-time numerical implementations for ODEs derived within the stochastic interpolants framework. We established total variation distance error bounds for both the first-order forward Buler and second-order Heun's methods, quantifying the discrepancy between true and approximated target distributions. Furthermore, we analyzed the iteration complexity for both methods, elucidating their convergence rates. Numerical experiments corroborated our theoretical convergence findings.
GRANT ACKNOWLEDGEMENT
This work is supported by the National Natural Science Foundation of China Grant 52494974.
REFERENCES
Michael S. Albergo, Nicholas M. Boffi, and Eric Vanden-Eijndee Stochastic interpolants: A unifying framework for flows and diffusions. 2023.
Michael Samuel Albergo and Eric Vanden- Building normalizing flows with stochastic interpolants. In The Eleventh International Conference on Learning Representations, 2023.
Michael Samuel Albergo, Nicholas Matthew Boffi, Michael Lindsey, and Eric Vanden-Fijnden Multimarginal generative modeling with stochastic interpolants. In The Twelfth International Conference on Learning Representations, 2024a.
Michael Samuel Albergo, Mark Goldstein, Nicholas Matthew Boffi, Rajesh Ranganath, and Eric Vanden-Eijnden. Stochastic interpolants with data-dependent couplings. In Forty-first International Conference on Machine Learning, 2024b.
Joe Benton, Valentin De Bortoli, Arnaud Doucet, and George Deligiannidis. Nearly d-linear con- vergence bounds for diffusion models via stochastic localization. In The Twelfth International Conference on Learning Representations, 2024a.
Joe Benton, George Deligiannidis, and Arnaud Doucet. Error bounds for flow matching methods. Transactions on Machine Learning Research, 2024b.
Sitan Chen, Sinho Chewi, Holden Lee, Yuanzhi Li, Jianfeng Lu, and Adil Salim. The probability flow ode is provably fast. In A. Oh, T. Naumann, A. Globerson, K. Saenko, M. Hardt, and S. Levine (eds.), Advances in Neural Information Processing Systems, 2023a.
Sitan Chen, Giannis Daras, and Alexandros G. Dimakis. Restoration-degradatii beyond linear diffusions: A non-asymptotic analysis for ddim-type samplers. In International Conference on Machine Learning, 2023b.
Yifan Chen, Mark Goldstein, Mengjian Hua, Michael Samuel Albergo, Nicholas Matthew Boffi, and Eric Vanden- Eijnden.. Probabilistic forecasting with stochastic interpolants and föllmer processes. In Forty-first International Conference Machine Learning, 2024a.
10

---

Published as a conference paper at ICLR 2026
Yu Chen, Marin Biloš, Sarthak Mittal, Wei Deng, Kashif Rasul, and Anderson Schneider. Recurrent interpolants for probabilistic time series prediction. In NeurlPS 2024 Third Table Representation Learning Workshop. 2024b.
Valentin De Bortoli, Iryna Korshunova, Andriy Mnih, and Arnaud Doucet. Schrodinger bridge flow for unpaired data translation. A. Globerson, L. Mackey, D. Belgrave, A. Fan, U. Paquet, J. Tomezak, and C. Zhang (eds.), Advances in Neural Information Processing Systems, 2024.
Will Grathwohl, Ricky T. Q. Chen, Jesse Bettencourt, and David Duvenaud. Scalable reversible generative models with free-form continuous dynamics. International Conference on Learning Representations, 2019.
Jonathan Ho, Ajay Jain, and Pieter Abbeel. Denoising diffusion probabilistic models. In Advances in Neural Information Processing Systems, 2020.
Chen Hongrui, Lee Holden, and Lu Jianfeng. Improved analysis of score-based generative model- ing: frienffre bounds under minimal smoothness assumptions. In Proceedings of the 40th International Conference Machine Learning, 2023.
Daniel Zhengyu Huang, Jiaoyang Huang, and Zhengjiang Lin. Convergence analysis of probability flow ode for score-based generative models. IEEE Transactions on Information Theory, 71(6), 2025.
Gen Li, Yu Huang, Timofey Efimov, Yuting Wei, Yuejie Chi. and Yuxin Chen. Accelerating convergence of score-based diffusion models, provably. In Forty-first International Conference on Machine Learning, 2024a.
Gen Li, Yuting Wei, Yuxin Chen, and Yuejie Chi. Towards faster non-asymptotic convergence for diffusion-based generative models, 2024b.
Gen Li, Yuting Wei, Yuejie Chi, and Yuxin Chen. A sharp convergence theory for the probability flow odes of diffusion models, 2024c.
Gen Li. Yuchen Zhou, Yuting Wei, and Yuxin Chen. Faster diffusion models via higher-order approximation, 2025a.
Runjia Li, Qiwei Di, and Quanquan Gu. Unified convergence analysis for score-based diffusion models with deterministic samplers. In The Thirteenth International Conference on Learning Representations, 2025b.
Jiadong Liang, Zhihan Huang, and Yuxin Chen. Low-dimensional adaptation of diffusion models: Convergence in total variation, 2025.
Yaron Lipman, Ricky T. Q. Chen, Heli Ben-Hamu, Maximilian Nickel, and Matthew Le. Flow matching for generative modeling. The Eleventh International Conference on Learning Repre- sentations, 2023.
Yuhao Liu, Yu Chen, Rui Hu, and Longbo Huang. Finite-time analysis of discrete-ting time stochastic interpolants. In Forty-second International Conference on Machine Learning, 2025.
Yang Song and Stefano Ermon.
Improved techniques for training score-based generative models. In H. Larochelle. M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin (eds.), Advances in Neural Information Processing Systems, 2020.
Yang Song, Jascha Sohl-Dickstein, Diederik P Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based based generative modeling through stochastic differential equations. In International Conference on Learning Representations, 2021.
Josef Stoer, Roland Bulirsch, R Bartels, Walter Gautschi, and Christoph Witzgall. Introduction to numerical analysis, volume 1993. Springer, 1980.
Patrick von Platen, Suraj Patil, Anton Lozhkov. Pedro Cuenca, Nathan Lambert, Kashif Rasul, Mishig Davaadorj. Dhruv Nair, Sayak Paul, William Berman, Yiyi Xu, Steven Liu, and Thomas Wolf. Diffusers: State-of-the-art diffusion models. https:: /tttps://itaae thub www/huggiingaaceeee...oo.ooooeee di f fusers) 2022.
11

---

Published as a conference paper at ICLR 2026
USE OF LARGE LANGUAGE MODELS
Parts of this manuscript, including sections the main text and appendix, were polished and refined using large language model (LLM). The LLM was exclusively used for editing and proofreading to improve clarity, grammar, and flow. All research ideas, technical content, results, and analyses were generated solely by the authors. This specific statement was also generated with the assistance of an LLM.
NOTATIONS
Unless otherwise specified, we use the following notations throughout the paper:
For a vector x E R2,, ||=|| denotes its & norm.
For a matrix A, ||_|]] denotes the operator norm induced by the la norm, and ||-A]] r represents the Frobenius norm.
The notations d aB Øi and θ all refers to the (partial) derivative with respect to I.
For two vectors x, y E Rd, denotes their inner product, and æ @ y denotes their tensor product.
The notation V denotes the the gradient operator. For a scalar-valued function f. V ff is the gradient vector. For vector-valued function Vf is the Jacobian matrix. Unless specified, the gradient is always taken with respect to the spatial variable z, not the temporal variable t.
The divergence of function written V f 2kml Jet
8f
The Laplacian of function f is denoted by A = Σ1
8²1
km 0.7
We use the notation f S g to indicate that / = Og))
A MORE DETAILS ON STOCHASTIC INTERPOLANTS
Recall that the stochastic interpolant is a stochastic process defined by
Ft Itt, zo. 21) (t)
Although xy establishes connection between the two distributions, its computation depends both Zo and ø1, making it unsuitable as a generative model yet. To solve this problem, [Albergo al. (023) demonstrated that the marginal density p(t,.) satisfies the following transport equation (where the notation represents the gradient operator):
Oep(t, x) V (b(t,))((t,z)) O.
where the velocity field b(t, æ) is given by
b(t,z) E[itk = x] E[0,/ (t, zo, x1) + + = x].
This transport equation reveals that any process (X,) satisfying the initial condition Xo ~ po and solving the ODE
dX,; = Xt)
will share identical marginal distributions with the stochastic interpolant (zz))eto..)) i... X, p(t) for all & [0, 1]. Crucially, the temporal derivative of X, depends only on the current time and position X,. Consequently, when initialized with a sample xo ~ Po. solving (4) yields X,
Furthermore, for any non-negative function <()) the SDE
Albergo et al. 2023 shows that the solution to
dX2, [b(t. X,) + (tt)s(t, X,)]dt + /2e(t)dW,
also transforms samples from Pa into P1. where W, denotes the standard Wiener process and s(t, = V log p(t, x) = Efy- | or æ represents the well-known score function. Notably, the ODE in emerges special case of this SDE when e()) 0.
12

---

Published as a conference paper at ICLR 2026
In practice, estimator b(t, x) of the expected velocity field b(t, 2²) can be obtained by minimizing the quadratic loss:
b(t, di,
where Ha ~ ptt) denotes the stochastic interpolant process. This loss differs from the mean squared enror
MSE /
1 { b(t, 2x) dt
by a constant that is independent b, as seen by expanding the quadratic term. The resulting estimator b(1, x)) when used in + yields generative model that transports to P1-
Table 1: Comparison with Previous Results
Reference
Order
Smoothness
Support
Complexity
ei al. (2024b)
No Requirements
Bounded
Ô
Theorem
Li et al.
(2024c)
No Requirements
Bounded
Ô (
Theorem
Li ct al.
(2025b)
I--Lipschitz
Bounded
Ö ( LR'd
Theorem 5.5
Li al. (2024a)
2
No Requirements
Bounded
Theorem
Ve
Huang et al.
2025)
Bounded Derivatives
(LDd)t+1/p
Thcorem 3.10
p
up (p 1)-th Order"
Bounded
el7v
Li al. (2025a)
ut-2/K
K
No Requirements
Bounded
Theorem
el/K
This paper
Bounded Derivatives
Theorem 4.5
up to 2nd Order (Assumption 4.4)
(Data-to-Data
This paper
Bounded Derivatives up 2nd Order (Assumption 4.4)
Theorem 4.5
M3/2 o VE
(Data-to-Data)
"They assume that 10 0% ( (x) S L for 1 S t |al Sp 1, where denotes the estimated score function.
B COMPARISON WITH PREVIOUS WORKS
In this section, we present comparison between our results and prior convergence analyses of ODE-based diffusion models, summarized in Table The table reports the iteration complexities of existing methods alongside their respective assumptions on the data distribution.
mdx L, E, - = is typically of order (d) in practical settings, such as when the data support bounded. Under this assumption, the iteration complexity simplifies to L and O
For Theorem 4.5 and 5.4 the value
respectively. Compared to the complexity established by Li et al. (2024a; 2025aj, this bound exhibits improved dependence on the dimension when we apply second-order methods. Furthermore, relative to the result of Huang al. (2025), when restricted to second-order solvers, our analysis relaxes the smoothness assumptions: specifically, we reduce the requirement of third-order smoothness in both z and alone (Assumption 4.4).
J3/2
smoothness in
B. MORE TECHNICAL COMPARISON WITH LI ET AL. 2025BB AND HUANG ET AL. 2025
Our technical approach shares similarities with the works by Li et al. 2025b) and Huang et al. (2025), which also employ interpolations reformulate discrete- iterations into continuous-time
13

---

Published as a conference paper at ICLR 2026
processes. In this section, we highlight the key technical distinctions between analysis and prior works.
Comparison with Li al. (025b) Although our analysis relies on Lemma 3.2 from Li al. (2025b) (see Lemma D.2) to bound the total variation distance between the target and estimated distributions, there exist significant technical distinctions our proofs.
First, rather than controlling the density ratio as in Li et al. (2025b) to extract the estimation error, we directly decompose the drift and divergence errors via the triangle inequality in a more intuitive manner, where estimation eтror and discretization error are considered in a more unified approach. Furthermore, when bounding the discretization error via time derivatives, we derive bounds in expectation rather than uniform bounds. Through careful analysis, this leads to improved dependence on the dimension d and Lipschitz constant L.
Second, while the general framework Li et al. (2025b) (Theorem 4.3) can estimate distribution errors for a variety of deterministic samplers, it is primarily tailored to first-order algorithms, as their discretization eπrоr analysis only considers first-order approximations, In contrast, our error decomposition framework readily extends to the second-order Heun's method, by replacing the first-order approximation with second-order one.
Comparison with Huang et al. (2025) Huang et al. (2025) analyze general Runge-Kurte Kutta methods for probability flow ODEs in diffusion models, which include Heun's method as special case. Compared to their results for Heun's method, our analysis offers two key improvements: (1) we reduce the smoothness requirement in from third-order second-order, and (2) we eliminate the smoothness requirement in t.
For the first advantage, this advantage comes from Lemma 3.2 by Li et al. (0025b), see Lemma D.2 Huang et al. applied Gagliardo- Niiennbbrr inequality, which requires that higher-order derivative is bounded, while a direct partition into drift error and divergence error has less requirements on the smoothness.
Regarding the first improvement, our approach leverages Lemma 3.2 from Li al. (2025b) (Lemma D.2p. While Huang et al. (2025) rely on the Gagliardo-Nirenberg inequality, which demands bounded higher-order derivatives, our direct decomposition into drift and divergence errors relaxes the smoothness assumptions.
As for the second improvement, avoid explicit time-smoothness assumptions for estimator ö(t, x) by implicitly utilizing the smoothness of the true drift function ztt, This is achieved through different interpolation scheme for the discrete-time iterates. Let (x:)) == (t - Lx 16)4,()
denote the interpolated forward Fuler update. In Huang et al. (2025), Heun's 's method is interpolated
+ t -- tu Nett, + b(t, Fi, Fe
X X,,
2
whereas our analysis uses the alternative choice:
- (+ t. XI, + (x)) 2(1++1 [x)
(t - tx)2 b(tenl Fig- (8te)) B(lx, X.,) ) 2(1++1 (x)
Both interpolations are valid representations of Heun's method, but they differ in motivation. The interpolation in HHang et al. (2025) ensures that the Taylor expansion around = tx closely approximates the true dynamics, requiring control over time derivatives of b(t, x). In contrast, our method approximates the intermediate drift dt X, = b(t, X.) using a linear interpolation:
d X, t-+1 b(tk: X..) + t- tx b(tx+1+ X
dt
t++1
tu
tx+1 tu
This form does not involve time derivatives of b(t, x). As long the estimator is accurate at fx and = t 1. our bounds only rely on the smoothness of the underlying drift b(t, x)), which is already established without further assumptions.
14

---

Published as a conference paper at ICLR 2026
C ADDITIONAL VERIFICATION ON GAUSSIAN-TO-IMAGE GENERATION
This section presents additional experiments on image generation from Gaussian noise, employing the same definition of the stochastic interpolant introduced in the main text.
Recall that VP diffusion models, the noisy data at step s is given by
Jts dig + v 0z
where 2 ~ (0, La) is the Gaussian noise. To align this with the stochastic interpolant framework, define the time variable t v a₂ and construct the stochastic interpolant as
= f + t
with 21 Pdata and NO, Ia). The corresponding drift term can then be expressed as
6 - t2
=(t,z) = Ef
B B[z|are] [[] -
s(t,
In our experiments. we leverage pre-trained score functions from Huggingface (von Platen et al. 2022) and compare the performance of the forward Euler method and Heun's method using the drift b(t, x) defined above. The resulting generated images are displayed in Figure 3 and Figure 4 As illustrated, although both methods exhibit converging errors as the number of steps increases, Heun's second-order method achieves faster convergence, consistent with our theoretical predictions.
10 steps
20 steps
30 steps
10 steps
50 steps
70 steps
100 steps
Euler
Heun
Figure 3: Comparison of Euler's 's and Heun's 'su's method on LSUN-Churches dataset.
steps
10D steps
steps
steps
25 steps
100 stops
25 stepss
100 steps
Euler
Euler
Heun
Heur
(a) I.SUN-Churches.
(b) -UUN-edrroo.
Figure 4: More comparison between Euler's and Heun's 's method on image datasets.
D PROOF OF THEOREM 4.5
D.1 INTERPOLATION OF THE DISCRETE-TIME PROCESS
First, for the forward Euler solver, we apply the following interpolation to reformulate the process as continuous-time process:
X, (x..) X1, (1 (6-14)6662 11)6ble X..), Vt (4x,)
15

---

Published as a conference paper at ICLR 2026
That is,
dx,- - =btx, ) = Ôe Ft (xi-
To write b(1x, Xx.)) as a function of (1, X.), we need to show that Fig-st is a diffeomorphism from F itself.
Lemma D.I. Under Assumption 4.4 suppose that the step size ha S zt. then Ft, Pey-st is a diffeomor- phism, and
Vz E R4,
VF 22 (x))l 22
Proof of Lemma D.1
VF, t2 (z)) = +
+ (t tx. vb(
So,
VF = (t- x)11 S h4L 1 2 <.
The above inequality shows that the Jacobi matrix V (x) is invertible, so F is a local diffeomprhism. the same time, note that global inverse function theorem, it is a global diffeomorphism on Rd. Moreover, the matrix norm of the inverse of its Jacobi matrix satisfies [[
F (z) 2 - - ylll so by Hadamard's
< 2.
By Lemma D.1 we know that X, satisfies the ODE
dxx b(t, X.)dt = b(te,
Now, we introduce an existing result by Li et (2025b control the overall TV distance. Lemma D.2. (Lemma 3.2 by et al. 2025b) For two processes X, and X, solving dX, b(t, X,)dt and dX, = b(t, X, dt, respectively. Let p(t) and p(t) be their corresponding laws. Then:
d di
/ (V ((1,)))))) V (ll,,
TVV(((),
f Ja. (b(1,z)
where Ω4 = fx E Rd p(t, x)) > p(t, x)}
So according to Lemma D.2 TV(p(t)), ((tx)) (V ((to))
AN
Exixp(e)
b(t, b(t, x.)] dt
Exixp(t)
[1o(t,xe) b(t. V In p(t, Xx)I| dt.
In the following sections, we denote by (X1)eE(to)w)) the true solution of the original ODE. When t E [4x, fx+1) are given, we denote 2 = Fi (X,, for simplicity. The rest of the problem is to control the drift error b(t,X,) b(t,X;) and the divergence error V btt V b(t, Xi).
D.2 CONTROLLING THE DRIFT ERROR
By a simple triangle inequality,
[|b(t, X²) b(t, XIIl < 16(te. 2) - 8(ta. X(e)
|b((t, Xxx) - b((ttt Xox)ll
||b(1 X(a) 6(1,, (.)x
16

---

Published as a conference paper at ICLR 2026
Below we discuss the above terms respectively. For simplicity, the notation 61,6(z) [lb(tz, x) b(tz,=)) denote the error of b(t, 2:) (1,,
For the term A, by Assumption 4.4
l16(t;,2) 2) b(tx. Xxx )1 - L|| -Xeell
= L| + - (Fu
S Xt - Pe
where the last inequality uses Lemma introduce the following lemma:
1 D.1 and the assumption that hg S 2T
For the term C, we
Lemma D.3.
(2) [
(b(s, x.)1| ds (t (x
||, (X S (t
bx
Proof of Lemma D.3
Consider taking derivatives and then integrate, we obtain that
/' X.
| Fe (X,,)ll
ds ds Fi, - (Xx,)(X))))
ds
1/
/ (b(s - b(tz: Xen) ds
[BbssX,)
b(tx. Xx )1 ds tu (x,,e,)(
["" du
(b(u, X,)) du + (t- (0-,t))e,,xX.
d
(1 1x)
ds + (1- (1-1)))1,((()) (b(s, X.)
ds
Then, by Lemma D.3 term can be bounded below:
[ d
[lb(le, 2) b((2, xx.) S 2T(1 (x) 1 ds (b(s. x.)) ds + e1,x
1/') ds (b(s. x.)) do 81,x
where the last inequality uses hy 2T.
For the term B, it is just (X,,).
For the term C, apply the similar method term A:
[btt,, Xx.) - x.t, <
b(s, Xx) ds.
1/ ds
17

---

Published as a conference paper at ICLR 2026
Hence, by summing up the above results, we obtain
E Ex))
[lo(t, - b(t, X+)I
X+)l|
11 (3)
[1ō(t, X,) - b(t, xx) x1) 3/4
E Eximp())
E Ex (t)((1)))]1))
4/37 3/4
7() "d [V/']]
d
bis, b(s, xx)
ds + (X
ds
-(t) /2/2 xt()
(t (2)t)))/3 [ do b(s,, xx)
d
4/3
ds t1, )4/3
3/4
5 -n(t)-1d112 10/2
(t 1//
[
E
d ds
b(s, Xx) Xx)
1/33
+E (x. 4/3
3/4
(8)
S ~1(11-1 /// [e1, )x 8
+7(1)_1d1/1( - -xa) E [1:000 - 2211" 1/2 7/ 11/2 1=1/ + t/x 83/2)
5 2001) [sia (Xe, )& 1/2 + 1 2dE [[1=0 x1 11]] 1/2
119 + h ld.
Here, (a) applies Hölder's inequality with power and 4/3; (b) plugs in the results for terms A, B and C; (c) applies the Jensen's 'ss inequality to the integral; and (d) uses Lemma F.1 to obtain the bound.
Therefore,
lx Ex Exiep()) [1o(1, X.) - blt, X.)I - p(,, de
to
I 11/2
7(t) "1 (X ) dt
to
N
2
2d [[|-0 - x21111] 11/] + 3xp 4d
k=
lx
1/10 C'N
11/2
-(t)-2d1
EX 2 de
N
dI
[11000
x11 11/
+ + 408
7k
(b)
S d Eann l0 1N2 1/2
N-1
1
7k
dE |||50 - x 1/11:11/2 + + 182)
km
where
- pxx
S(7 to (N) =
7(t)
(3)
For the above derivation, the inequality (a) applies (b) uses Assumption 4.2
inequality, while the inequality
D.3 CONTROLLING THE DIVERGENCE ERROR
First, according to the chain rule of derivatives,
Vū(1, X.) = vi(le,2). VFr.-n(2)
18

---

Published as a conference paper at ICLR 2026
Therefore, by applying triangle inequality,
V b(t, - V -(4 x.) - (illtn1)) vo(tn,
(bb(inn)) X1c) ((==-
+ tr [o((ts Vb(te, X(x) (()
+ tr Xoet (1=) Ia)]
+ V b(tR, ((tR, V b(t,
Now, we deal with the above terms respectively. Similarly to the velocity error, we use the notation
(2) = Vb(tk. x) Vb(tk. 2)1F.-
First, since | (V) - Tall= S L.hg S 1/4 (er --al=
F I1(Ia Fu))(2)) -aa) 11 Ialls
1 + Σ (la (z))+ -all
Lha)* s 2Lh
So for the term A,
Vi((e; 2) Vi(tk, x )] VF, VF (2)--1(2)-11 |1]
< V b(tz, 2) V -b(te,, X ) + 21/1 Vb(lR, - -
< 2L 1/
X1al 1l S 3/2 |x, - (X )1
[_]]
d
< 22
(b(s, X.)) + 21 1/2 (2,,
ds
where the last inequality applies Lemma D.3 For the term B,
ir
Vō(te) Xxx) - Vo(tk, Xie)
(() I]
$ 2 Vo(tk, Xx.) - Vo(tk, Xax /
< 25 (Xte)
For the term C,
|tr [b(tz, ) (VFr [a)] s 2Lh2
Vb(tt. /
m..
For the term D, consider taking derivatives w.r.t. t:
d (V b(s, Xo)) ds.
b(th. X.x)) V Xtt <
ds
Therefore, combining Lemmas
F.14 F.13 and
F.
we can get
E V b(t, X.)
b(t, Xx)
$ he E [1120 - x 1 (Th d + 7p- 1,22
+ 3/ 268
L 1d2 1El - 211 + L 2d]
+ 1/2 [s
(X
19

---

Published as a conference paper at ICLR 2026
Then, taking the integral and applying Assumption 4.2 and 4.3
v b(t, b(t, xi) dt
Ex(p(t)
80
N
E |||-00 I1 1110///5 P 2d + "x- td df
A=0
T A d1/2
2 d E1 22112 1i
E El52a((Xx))
dl +
N
S Ediv + 1
Sarin
+
E
11:000
2z 2111% 11138
2/3
7x 2d 7k" '88
k==
L5k
d=/2
E
21
1
d
Combining Lemma D.2 the bound velocity error and the bound on divergence error, obtain the error bound given in Theorem
E PROOF OF THEOREM 5.4
E. INTERPOLATION THE DISCRETE-TIME PROCESS
The proof of Theorem shares the same basic ideas as Theorem In order to compare the discrete-time iterations with the true process, we first apply an interpolation on the Heun's 's method: -el (X.,)) The interpolation is defined as follows:
4.5
X, Gea-te(X
(a)) x + f"d
(((t),2) 1++1 - 14
(bltu+1) Fu-ter b(tk, w)) ds
tx
(1 (1-14) 1x)2
(1- (1-1x)²
=x + (1 - 11)
b(12, + 2(tx+1 - tx)
8(16+1, Fru (2)),
2(fx+1 [x)]
where recall that (x)) = x + (t la (u)((ts, x)) is the interpolation for the Euler's method.
Similarly, we want to show that X, admits the same ODE form as the true process, which relies on the following Lemma.
Lemma E.1. Suppose Assumption 4.4 holds, then if hy S aT. Gea-st is diffeomprhism, and for all a 12.
[V 22
(v) S2.
Proof. Similarly to the proof of Lemma D.1 we only need to show that [VGta+e(z)) Iall S 1 2-
(t 12
VG Ia +
(t tk)
vb(te,2)
2(t++1 - tx
(t )2
+ 2(1441 - -4a) Vb(tR+11 Feu-tess (x)) Fin-
(t - (t))
Id + - -tb) (t++
V)(te,
2(t - tx)
(t - tx
Vb(tu+1. Fec-tee (x)) (Is he vb(te,z)))
+
2(1x+1 - -ta)
20

---

Published as a conference paper at ICLR 2026
So,
1
[VG - Lall < h&L + 2 h h₂L) <
which completes the proof.
By Lemma E.1 if we write (Xt), then
d
X1 = b(t, = tx+1 b(tk: 2) + t- tk
b(tb+1, Fi Fic-teri (2)).
dt
th+1 - to
tx+1 - tx
Then by Lemma D.2
(((((t)),(tx)) TV(p(to),
E Ex(e((t)
V V Xx) V V-b Xe)l dt
AN
Ex(e((t)
[io(t,x)) - Xi)l V p(t, Xo)l| dt.
E.2 CONTROLLING THE DRIFT ERROR
Still, use the notation ε (x) = [bb(tz, x) - ((tk,x)); then we can control b(t, - b(t, in the following way:
X(, - X.)
B(tee 2)
(((tet), Fer
2ttx,2))
b(t, Xx)
tk+1 - tk
tk+1 -
bl(te, 2) - b(tx, ((x) X + t- tk
[[ttet1, F'x Pon- (2)) - b(tx41, Xtema
16+1 11
lu+1 Lu
A-bias
1++1 -
[ilk, Xe.) bllk, Ntn)
[ille1,, Xemra) Ule+11 Xemra)
+ t&+1 t&
fk+1 - tb
Brestimation error
+ b(tk:Xix) + t tk
(b(tu+1, X b(tx, Xxx)) b(t, X.).
tx+1 ty
Crdiscretization ertor
Below bound the terms respectively. First, for the bias error A,
lR+1 L
T.|l Neall + 1- l1
A 5
L Fea (2) X Xeota
tk+1 t₁z
tx+1
tk+1 -t t
Ll[z - -xaul
l++1 - 1x
11
+ t tk
I (1l -ta+s (z) - Fie-teus (X,,)ll + | - (Xe,) - Xe II)
d
< 2L
X + L(t tx)
b(s. X., )
ds + L(t tx)E1, (Xte)-
ds
The last step uses Lemma D.3
For ||2
Xx, |l. we observe that
|| Xx,II < |G (GTLL G (Gt, $ 2
Now, introduce the following lemma:
21

---

Published as a conference paper at ICLR 2026
Lemma E.2. When h2 t. suppose Assumption 4. holds,
|X Gu,)e(X(,)) 111 2 (t (Xt,) + (t ts (x))e1
+ (t 2 / 1e41
- b(s, X,) ds
ds
d2
+ tx ++x(1++1 - -tx)
b(s, X,)
ds.
ds²
Proof of Lemma E.2 First partition the error into several parts.
(Xt())) - Xell
(Xtx)) Xtx+ - (X, Xxx)ll
(6-1x) 1x)2
|b(tes Xte) - bllx, X(x)
(x)) 2(1x+1 (x)
(1)
+ (t- tx)²
X - b(tut1, (lxet))Xemmm) Xemue
2(1++1 - (4)
(t-ta 12x2)2
||olkAn Pen-sta (Xtx)) - Blla+1, X(x11)
+
2(1x+1 - -14)
(ii)
+ 1
Xex) + 8 tx
(b(tx+1, Xtex1) - b(tk, - b(s, X,) ds
t4+1 to
t/
(N)
We then need to provide upper bounds for the parts. First.
(i) + (ii) s 2(t (x)51u((()) + (t (N)ea))
For the third term, apply Lemma D.3
(I - (x)2 )2
(iii) <
T Pta (Xta)) - Xtatill
(t&+1 (x)
ptat
d
< (1 - L
b(8 Xa) ds + (1 (6-1x) Le
ds
For the last term, note that
tx
((the+1,(ee4)) -bbt, b(tk X.e)] [b(t,X,) b(tx, Xex )
tk+1 to
to
plats d
ds bbs,, '"d d bss,,,,,0
d
1++1 4x
44
- 44
pl1l [
blt, Xi) ds 1"' d b(ss, - di
ds blss,, Xs) dt
b(t Xi) ds
1x+1 tk
J1z ds
bx
plati
ds f d²
duz bu, X,)du + 14+1 1 1: ds S d²
duz bu, x, dd.
1++1
t
tx+1 - to
So.
f
|x+1
d²
(iv) S (1 (4) - (x)
b(s, Xx) ds.
ds
The proof is finished by directly combining the four parts.
22

---

Published as a conference paper at ICLR 2026
By Lemma E.2 the bias error
A S L(t - (
+ L(t - +1 (Xtet2)
d
L(t (x)
b(s, Xs) ) ds
ds
AN
pllll
d²
+ (t tu the
b(s, Xx) ds. do
For the estimation error B, clearly,
f E (Xww fk+1 - tk
t- tx
B
E (Xe.:).
tk+1 - ts
For the discretization error C, according to the proof of Lemma
E.2
we obtain that
ix+1
d²
||||| S (fx+1 - tx)
(ss
X,
ds.
ds
Combine the above three terms, and apply Lemmas F.13 and F.15 Then, we get:
Ex Ex,xp()
((((t, b(t, x.)l In x,)
[lb(t, (x)1 5/6
Ex(x((t) [t] X x1) 1/6
$ Ex(x
~0(1) 12/2 [B [[lsTa (X.,, (1x1)) 1/2 + E [[le1,e (X
m)13]] 1/2
d b(s, X3)
6/5
5/6
+ (d))-
(t t2=s))2/
ds
ds
(1 - (1) 1/5 1/1
d
6/5
5/6
b(s, Xz)
ds
ds
$ -t)) 12/2
[ [E ) )
(X )
)211/2]
+ 9/6 8/2///
11:00
x1 21110 1/2 + d
+ 19 2 dE
|1100~1111 16]
-
1/3
+ Idd
Integrate and write M- = max E
|||000 x1 191 L,, d). then
Ex.xo())
[ X.) b(l, x.)ll | In p(l, X.>l| di
N-1
lik
7kk
d8
+ M3337
k=0
E.3 CONTROLLING THE DIVERGENCE ERROR
Use the notation (x) Vb(te, 76(te, First,
tk41 Vilta, 2) (G(z-e(2))
V. X() V b(t,X,)
t++1 t
- tx
tx+1 to Fo(tt-41,FF--1 (2)) Fi.-Fuuu (2) (z)
Vb(1, (2)
23

---

Published as a conference paper at ICLR 2026
Then, by Lemma E.1
tk+1 vb(lk:2)
V b(1, V. b(1 xi)
-
tr
kR+1 4
t tk
vo(tk+I, Fex-teen (2)) -e2 (2)
tk+1 tk
Vb(l, V.)
|tr DDiis + Dest + Daisll
where the terms are defined as:
Doias
tk+1
vb(tu,2) Vo(tu, Xte)
1++1 - b4
t to
Vb(tk+1,Fuu-i Fu- (2)) V.,, +1 (z) V6(te+1; V +124-1241 (Xee)
tk+1 tk
(Xtx))
Dest
tk+1
Vi((n, X.,) Vl(le; X
1++1 L4
Vi(letI, ) V (Xt,) Vb(lx+I, X VG,, -1 (Xtn))
+ t++1 t
1++1
Vb(le, X V (Xte) 1- 44 Vb(le+1, Xtar VG staps (Xta) tx-1 ts
Dais
tx+1 tx
VG
Bias Error
For Daias, with Assumption 4.4
tu+1 1 13/2 13/2|1=-X11 Xxal
Urissi
1x+1 - 4x
+ t. t
Vo(te+I. Fux Fe (z)) - Vb(tx+1, (itt+t,,xxxn))) X ||V Fi. -ce+s (2)lle
tk+1 t+
t t
-1Vbi(l(( Xt F | P, m (2) VGrumann (Xtn)ll
4+1 - 4x
Vb((, X.)IIF Gr Vo
1 1 ///2 Xt. +
t t-th t 1//
Fix (2) Ximail
fx+ - to
1x+1 - 1x
t to
I I1(1,, 1/ nu vi(te,))
1
Ne Vō(ta: Xee)
tx+ tx
2
1
11 vo(lat1, Pen- (X F (Xtx)
2
+ vx)(t
-L3/2 Xtx - =1
24

---

Published as a conference paper at ICLR 2026
(6)
2 1|2 - Xte 1| (1 Vbll + (Xce) - X12+i
Lh VFLL-0441 (xx) Iall [1V((t2+11 x lsse
Lhk
Vb(tk, X Vb(ta+1,X
(1 + h4 Vb(t, x))lF)
(LEIA((O)) + h4 (
d
d²
h
b(s, Xx)
ds +
b(s X₂) ds
ds
ds2
3/2 hk
b(s, Xx) ds
ds
+ Lhf (52.
+ 22 (Xest1))
Lhtg
(bb(s Xz))
ds.
The inequality (a) applies Assumption 4.4 the inequality (b) rearranges the terms and applies triangle inequalities; the inequality (c) expands the differences into the form of integrals by
Lemma D.3 F.2 and we control the difference [vb(t,, Xt,) - vb(la+1, X by control-
ling Vb(lk. X
Notably, in the first term, when di S and h h2 'B 21 1/11/// S (p 1 is a small constant), we have lg [/V,t,, X.) S 1, so based on the assumptions on hg, Lemmas F.13
F.15 F.14 and Hölder's is inequality,
E |r [oius]
|1
x h
hk
E I0 xn 11 1/3
%x
1--1d d1/2 5 9-3)))
+ L 3/2 hi
E ||:00- - x1 1/110///2 1/2
d%
-2 + db/2gt0
+ LhE E [[ 21 110 1/1 di + dz
+ ///2// (EE5), (Xea)l ]1 + B|SI, (Xte41)l)
+ Lh (X, ) + BlEg (Xee41))) )l)
Estimation Error For Dest.
[r[[Desimate] ]1 S tk+1 22 (Ner) + t tv
-2 (Xtrur)
| (Xu)) )l
1411 tx
to+1 tu
< (Xtn((())) + 22 (Xtr.)).
Discretization Error For Dais note that || (x - Iall is of order O(h), we cannot use this approximation because are expecting the error of order O(h?). Thus we different approach to deal with this term. By taking the second-order derivative,
[rs]Dis]] < (1241 (x)
Vb(s, X.)G.
d
We now explicitly write the derivatives of VG below for further discussion:
VGr,- = la + (11 ((--a)) (1-1x) (x)
vi(le,
2hv
(t (1-ba) 12)2
Vo(te+2, Fu-mun (x)) (Ia he v)(t,,))),
+
2hp
25

---

Published as a conference paper at ICLR 2026
d t v)(
VG v
dt
hg
1 Vb(lktI, Frs tutts (2)) (La + he v)(tx.z)),
d²
1 [ Fuataai))) (Ia + he Vb(tz, 2)
vo
=
dt2
h*
So, by our assumption that h4 L 21 1.
d v( SL.
(GG=-( 21,
dt
F
VG Vau-i(Xue)
2 L (E2,(()))) (()
dt2
T,3/2
| (Xic) - Xeatille
h&
pt++1
d Vb(s, Xx)
+ hk
ds
ds
tkk
p
62
(Xa,)
+
62.+1
h&
h
The constants omitted by the notation
"S" above is uniform for all x (or xx,) ) and t. So,
pla+1
d
E [tr[Daia]] S ho
Vb(8,
Xs)
ds?
F
d
+ LE
Vb(8, X.)
ds
F
d² VG
1
+ E
Vb(s, Xr
(X
ds
ds?
P
hk
E
[e2, t, 121 1/2
+ |2 X )
+ 1,3/2 1.3/2/ (X,,, 991/22)
(d)(1)=2 + 22))
E [1100 x 11772/2)
+ 2 E [1|=r0- - x1 | 1/2 d حج 18/
+ d dZ 6
+ (L + d EE [1koo - x 11/ 1/3 di + 3/ d.)
Therefore, by adding Driss, Dess and Dais together and integrating with respect to f, we can obtain that
AN Eximp())
v b(t, V b(t,,, x(t] dt
N
h& E (xxx 1912/ 1/2 + )2 29))/2)
ha 2//2 (X 991//2 - -E ((Iaa(Xe)
N
h (daq-o 6 + M29-1).
E [1|zo - x1 19 1/3 T., d
where we denote M = max
Combining the above discussions, we can get the bounds given in Theorem 5.4
26

---

Published as a conference paper at ICLR 2026
F TECHNICAL LEMMAS
In this section, we use the notation 2; to denote the stochastic process
2t = [tt, Zo. 21) + +
defined in the stochastic interpolant, and use the notation X, to denote the process that solves the equation
d X, = X(t
dt
and satisfies X p(().
Recall that in stochastic interpolants, we have defined b(t, x) = (4, and s(t,z) v log 2) (t)]]l2 z]. For convenience, define v(t,z)
xn, + 5z/ar, = x]
Et zo, x1)|=,
F.1 THE CONCRETE FORM OF DERIVATIVES OF s(t, x) AND v(t, x)
Before introducing the lemmas, we first define the notation ft m1 21)
so that the conditional density p(zz|zo, æ1) exp(ft). The first lemma shows that how the jacobian matrices of v and can be expressed.
2x(()²
Lemma F.1.
Vz0(t,x) Cov(S,T, Vafilar =x); Vxs(t,z) = Cov(Vaft V Veftl = x)) - -7-2 la.
Proof.
Vou(t,z)= V. I (exp(f()..ta
eexp(
The conditions for Lebesgue dominated convergence theorem be easily checked, the order of differential and integral can be alternated. Then
V (t, x) (exp(()(( ( [exp(S()O.T @ VaJidv [xx((I))d + f (exp(si) @ Vexp(fi)V
J eexp(s))d
efxxp
Vov z). Here the notation @ represents the tensor product. Similarly,
(t,(t [ I ('ex(fi) & Vafide /cxp(fi)dv Vex & exx((f))Va..dd.
I /exp(fi)dr
[S exxp(())d
| 1
= Cov(VV,J VaJt V2
2
= Cov(V V.fi x) 7 2
2
In the above calculations, the key is to change the order of taking differential and integral, and rearrange the final formula into the form of conditional expectations. The higher order derivatives can be obtained similarly, but since the calculation is too long, will give the results without showing the detailed proof. Lemmas and E. are similar as Lemma F. while we consider the recond-order and third-order derivatives.
27

---

Published as a conference paper at ICLR 2026
Lemma F.2.
s(t,,))) = x], -et.0 [t, = -
- x(t, & (VaJ. - s(,, x)) 8 (V ah - stt, x))|,, = x| Note that v(t, x) and stt, are the conditional expectations of 0 and V vz respectively.
Lemma F.3.
=) = E[ u(t,z) & (Vzfi - s(t, 83 D4 = xj -T(Cov(0//,V.fikk 2) Cov( Ve Vefelz, - x)),
=]
- -T(Cov(VV.... Vzfil = x) )7
Here, the notation A represents the tensor product of n tensors A. The operator T is defined as
T: Rdxdddddd Raxdxdxd,
T(X) T1231 (X) + (X) + (X),
where
[TPipeppPa (X) Xtp fog tog
Similarly to the derivatives with respect x, the following lemmas (Lemmas F.4 calculates the derivatives of v, s and their Jacobian matrices with respect to the time t.
F.5 and F.6
Lemma F.4
0,o(1,2) x)) = =a
=
0,5(1, - x].
V x)
=x)
+ ((.0000-~ot)))000V__33t( & (0rf: -ma, )ar z). =
V.V,s(t,
f.fe/ = x)
+ + @ - st, x)) @ (0rft - ma,se) = x].
In the above formula, we use the notation my = E[V[2; x) for simplicity.
Lemma F.5.
0, = = + + RE((O,T - (0.
=i]
ma ) = =x], (t = 4 Cov Vo. x]
SiVVzfis(tta))) (8.f ma, (0.f. ma |x) =].
Lemma F.6.
V.Vav(t,)) x) = -
=
+ +Cov(01,88V.f/eee
+ +2E[(0 m(p @ (V.fi - st,,x)) 8 (0.f. - ma. ma. = x]
+ 2B[((a.1 v(l, -
ma, v.Ji) 8 (0.f: ma,j)) = x]
=[[[aI vtt. x)) mv 8 (0 mm² ( 12) = 2]
E[[0,T- - () (0., )8 ,3,
[[ (:t0::/0.00.000000.0..00..... 0.f = =)))
0. [t) = -00 18
0, Vafilzz = X))
+ +CovVV//..VV///lle
+ (. (V - st, x)) @ (0.fe m (x) x]
+
- -tt.f//) (0. - ma. 1.) @ (0.ft - mm. = x]
+ IE[( =s(,,4) @ (V - mv. @ (0 f: map is ) =x]
[V((V/f.- (0.): ma,,.) 11 x]
(CCoovVV......ft @ Cov(O.fe, 0. l = x))
28

---

Published as a conference paper at ICLR 2026
F.2 UPPER BOUNDS ON THE DERIVATIVES
We first provide uniform upper bound to explain why Assumption 4.4 is reasonable.
Lemma E.7. Suppose that for all t E [0, (t, 1 R) 1, i.e. the data is bounded. Then under Assumption I,
[V((t, *R²
[[v x) 5 91² IV(V b(t,)))] 2 RR'
Proof. According to Lemma F.1 since γ = (1), for any u E R by noticing that || @ ylll= ||>||| llyl.
x))( =
S - + IV - s(t, VV.. s(t, = x²]
Since ||0,7|| & R, V s(t,x) -- -502 E8- -0 xy x] = (1) B[T म. z)) 2
((t) -R, we have
[Vb(1, S "
Similarly,
b ( 5 (() (tt |I(V b(t,)))) &
Lemma F.7 justifies the choice of Lipschitz constants in Assumption 4.4 where we set T, and T,3/2 as the Lipschitz constants for b(t, x) and its spatial derivative, respectively. When R O(vd), as occurs when data are bounded each dimension, it follows that L Odd).
The rest of this section analyzes expectation- upper bounds, rather than uniform bounds with respect and
Lemma E8. For Gaussian random variable z ~N La), for any constant p 2 2.
E 17 C
where the constant C(p)= > 0 only depends on
Proof. First, by Jensen's inequality,
D/2
||=||-- (1//13)p//
dp/2
2;
p
so
[[|2|P]] <
/
Here, C(p) E/21/P < oo is a constant that only depends on p since #1 ~ N 1) is standard 1 1-dimensional Gaussian variable.
The following lemmas (Lemma F9 to Lemma E.12 provide upper bounds the moments of the time and spatial derivatives of o and s. Notably, for these lemmas, the expectation taken over I, p(t), where recall that x, is the stochastic interpolant.
Lemma F.9. Under Assumption
4.
for p 2 2,
E [1) [P) S B|| x x1111,,
& /P/
29

---

Published as a conference paper at ICLR 2026
Proof. By law of total probability and Jensen's 's inequality,
E 21))(1 - E ||0, & E || 1 12
E (t. - (t) " [|*] 2 (tt)
Lemma E.10. Under Assumption 4.1 for any p N 1, we have
E V x)) - (() -p [[ -2 [2r]
x1
[ S 7(1) 2rd
Proof.
[[Cov(a,1, V2filz; = x)|lr -
=
Then, by Jensen's inequalty and Cauchy-Schwarz inequality,
E V x)) E 1
(tt r x1 [2
Similarly,
+(, 2 7() 2p)0)22pp0 2010 E -I 22
Lemma F.11. Under Assumption 4.1 for any p 2
(1, x - B x0 - x 20]
p p/2 E [[|00
x1 1|001//8 7 dp
1/2
E [1]0,s(t, x;) s 3p 3p E I0 2 21 112001/2 7-0
Lemma F12. Under Assumption 4.1 for any p 2 1,
E | V o(1,2,)) 2 3p 3p/ E [| x1 ||00] 1/3 + - [1 - =1 2///3
E B SeVV.s(t,zz) 11 $ - Ap 2o + 30 3 [|:20 xz 111100114 1/3
The proofs for the above two lemmas are omitted as they almost repeat the proof of Lemma F.10
Now, with the previous lemmas, we are ready to provide upper bounds on the derevvvvisss of b(t, X,) (Lemmas F.13 and F.15) and its divergence V b(t, X,) (Lemmas F.14 and F.16). where X, the true solution the ODE dX, b(t, Xt) that satisfies X, ~ p(t). These results are later used to control the discretization error of numerical methods.
Lemma F.13. Under Assumption
for 2 1,
E
b(,, X.) di
d
S Io - x [|00] 2/3 7(1) 1/20 + 17(1)-3000)))
-30
Proof.
dl 1(L, Xx) = 0,b(1, X,) + V.)( (()
+ =,,t( X,) + 2 (i, + (t,((,
+ V X.)] b.t.
30

---

Published as a conference paper at ICLR 2026
By Jensen's inequality,
k= [x=1,:) nn11 Ek= la,l,, so
|
E dt
b(t, xx)
$ E [l0,eut X
+ +B X;)
X
V(t, + + x.)I (1((,,)))1
E [10,vv(1, x=)
+ -E| +Ell x + E [||s(t, x
[|V X,) + (t, X,)l 30/2 2/3
E [|b X1)l
E [1 - x [| 2/3 P1) + 711)
The inquality (a) uses Hölder's inequality, while the inequality (b) uses the results of previous lemmas.
Lemma F.14. Under Assumption 4.1
for p 1.
S E -10 21 21115p]2/8 2/3 y 22 + + 7- dpp
E
V
b(t, X.))
dt
Vb(l, x.)
2 E [1 [11=0=x1||00
|| 2/3 --20 p0 7 4p 6p
dt
Proof. Similarly to the proof of Lemma F.13
(V b(t, xi))
I E ||r (8t
+ E tr (v²
(o)(,,
For the first term on the right hand side, we can use Lemma any four vectors v, y, 2, w, have
F.12 For the second term, note that for
ttr [I @ 2 = = |tr [(> @ v)
(z w-]l S || ||yll || ||/]].
Hence, by the conditional expectation form in Lemma F22
E tr (t) x.)
E [
2/3 γ -20 d + + - 2p8
So.
E (V b(t,
S EE |
2/1330/2/3
2/3
y- 2p d 7 A2
The case for the Jacobian matrix is the same except where we use || @ vl the inequality for the trace.
|||||| ||y||
instead of
Lemma F.15. Under Assumption
5.1 for p > 1,
d b(t Xt)
E
E [|100 xz 13008 d' 2p d0012
dl
Proof.
d²
6(1,, X.) = 8,6(1, X,) + 20, VI(1 X.) 6(,, + Vb(1, a,b(1, X.) d1z
v X,) b6( 0² 02
bb x b(t, X,)
Consider Lemma F.I to F5 where the conditional expectation forms of the above terms are given, then use Hölder's inequality and Jensen's inequality obtain the upper bound in the Lemma. n
31

---

Published as a conference paper at ICLR 2026
Lemma F.16.
Under Assumption 5.7
d²
E
(V.
b(t, xi))
2 E || - 2111 2/1111122 (t) + ( (t)
dt2
d² xi)
E
Vb(l,
23 13 [1:00- - -211191//2 1/30)/28//7001 711= -3 + d'((1)-
di
Proof.
d²
(,), = + X.] + V2b(e, (,) X.)]
b(t, x,) X.) 00 Vib(t, X.) Vb(t, X,) b(t,
The rest is similar to the previous lemma.
32


Reviews:
[{"review_id": "Nvuk8VIKBD", "reviewer": "Reviewer_TgjG", "summary": "This paper investigates finite-time convergence of ODE samplers built from stochastic interpolants (SI). While previous studies have established finite-time bounds for SDEs, the ODE scenario has notably lacked discrete-time guarantees. The authors fill this gap by providing total variation (TV) distance error bounds and detailing the iteration complexity for both the forward Euler and Heun's methods under regularity and estimator-error assumptions.", "questions": "See the weaknesses.", "limitations": "", "rating": 4, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "- The paper addresses a significant theoretical gap for SI-based ODEs by providing finite-time TV bounds and iteration complexities for both Euler and Heun methods within the SI framework.\n- Mathematical proofs are solid and clear. Discrete-to-continuous interpolation yields a piecewise ODE, enabling drift/divergence comparisons to bound TV.\n- 2D tasks and Gaussian-mixture tests verify $O(h)$ (Euler) and $O(h^2)$ (Heun) discretization orders.", "weaknesses": "- The requirement of uniform Lipschitz on $\\hat{b}$ and its divergence (Assumption 4.4) seems kind of idealized. Is it possible to provide an empirical verification?\n- Experiment demos are limited to three 2D transformations and d-dim Gaussian mixtures; Is it possible to provide results on real-data benchmarks?\n- The theory in this paper establishes error bounds that scale with $d^2$ and $d^3$. However, the empirical findings indicate a roughly linear correlation, suggesting room for further refining on the bounds.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "A0ckGGrmyy", "reviewer": "Reviewer_Pr2N", "summary": "This paper gives convergence analysis of ODE methods for stochastic interpolants under general conditions, in particular, the first-order (forward) Euler method and second-order Heun's method.\n\nFor Euler's method, the requirements are L^2 error bound on the drift $\\epsilon_{drift}$, an expected Jacobian estimation error in Frobenius norm $\\epsilon_{div}$, and smoothness bounds (Frobenius norm bound in the Jacobian & 2nd derivative of the score and Lipschitzness of divergence). The theorem gives $O(d^2/\\epsilon)$ iteration complexity.\n\nFor Heun's method (with Jacobian estimation error in Frobenius norm squared), the theorem gives $O(d^{3/2}/\\epsilon^{1/2})$ iteration complexity.\n\nExperiments were run to obtain the empirical scaling of TV error with step size and dimension.", "questions": "It would be helpful to explain where the extra factor of $d$ come from compared to diffusion models. \n\n* p. 6: \"Divergencce\" -> \"Divergence\"\n* p. 8: \"folloiwing\" -> \"following\"\n* p. 8: Missing period after \"unbounded\"", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The paper is the first to give an analysis of the ODE discretization for stochastic interpolants under general conditions. This works for any interpolations with path $I$ and noise magnitude $\\gamma(t)$ satisfying some smoothness conditions, and estimates satisfying the assumptions. Although the analysis techniques do not seem particularly novel (in being similar to previous ODE flow analyses), it is nevertheless valuable to work out the precise bounds in new setting under carefully detailed general assumptions, as in this paper. In particular, the theory shows the benefits of a higher-order ODE integrator and under what smoothness/error assumptions we can expect the improvement.", "weaknesses": "There is an additional factor of $d$ in Euler's method compared to diffusion models, and it is unknown whether this is a limitation of the analysis, or necessary in this general setting. (The experiments suggest the true scaling is linear.)", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "p1Cnlme9Kj", "reviewer": "Reviewer_9RHo", "summary": "This paper presents a theoretical study of iteration complexity for discretization schemes of stochastic-interpolant ODEs. As representative schemes, we analyze the first-order forward Euler method and the second-order Heun’s method. Under certain assumptions, complexity bounds for achieving $epsilon$-TV error for each scheme are derived. The analysis proceeds via carefully constructed continuous-time interpolants, which allow us to import continuous-time techniques to quantify discretization error and evaluate approximation accuracy.", "questions": "(1) Definition of $d$ seems not appear in Introduction.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "**(S1) Important research theme.**\nInvestigating the iteration complexity of ODE solvers to ensure accurate generation with stochastic interpolants is an important research direction. Although diffusion models (in discretization form of neural ODE) are a major focus in machine learning, most evidence remains empirical. This work may provide the missing theoretical analyses underpinning diffusion models.\n\n**(S2) Theoretical contributions.**\nFor general stochastic-interpolant ODEs, we derive the iteration complexity required to achieve $\\epsilon$-TV error when using the first-order forward Euler method and the second-order Heun’s method.", "weaknesses": "**(W1) Positioning relative to prior work.**\nIn the Introduction (lines 048–049), you state that for ODE-based transformations, “the analysis has been limited to the continuous-time setting.” Meanwhile, Sec. 2.2 (lines 123–127) and Appendix B (Table 1) summarize prior results on the time needed for an ODE to reach $\\epsilon$-TV error. Am I correct that all of these results are derived in continuous-time setting and do not account for discretization error? If YES, please revise Sec. 2.2 to make this explicit.\n\n**(W2) Challenge level of the proposed analysis technique.**\nMy understanding is that the techniques in Sec. 4.1 (“Main ideas of the proof”)—which bridge the gap between a continuous-time process and a discrete-time estimator—are the key methodological contribution, enabling tools from continuous-time ODE analysis to be imported into the discretized setting. That said, it is somewhat difficult for me to gauge how technically challenging this bridge is. Could you additionally explain the main difficulties encountered in introducing this technique?", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "pzu2DnbQTV", "reviewer": "Reviewer_9Pnh", "summary": "This work provides the first finite-time convergence analysis for ODE-based stochastic interpolants, establishing total variation error bounds for the Euler and Heun methods. The results yield optimized iteration complexity and step-size schedules, offering improved smoothness and dimensional dependence compared to prior diffusion model analyses.", "questions": "In my understanding, many techniques are closely related and can be applied to previous analysis of diffusion ODEs. Do you think there are any properties specific to the stochastic interpolant problem (I.e., starting from another distribution instead of Gaussian noise)? For instance, can you provide some concrete examples to show the difference between the stochastic interpolant and the original diffusion process?", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "(1) The paper is the first to analyze the convergence of ODE-based stochastic interpolants. The analysis is rigorous, and the presentation is very clear. The convergence results show improvement compared with prior works.\n\n(2) It conducts synthetic experiments to exhibit the TV distance.", "weaknesses": "(1) There are some typos and ambiguities in the proof.\nLine 730: In the triangle inequality, the last term (term C) seems to be $b(t,X_t)$, instead of $b(t,X_{t_k})$. The remaining proof should be correct.\n\nLine 792: It takes me some time to verify the correctness of these inequalities. However, the reasoning skips too much explanation of why they work. I suggest adding a detailed description of why we need these inequalities and how they are proved, including the calculation of Holder coefficients.\n\nLine 927: The $G_{t_k \\rightarrow t_{k+1}}$ should be $F$. Moreover, the definition of $F$ should be restated to make it clearer.\n\n(2) At first glance, I feel that the proof technique is very similar to that in Li et al.2025(b) with the application of Lemma C.2. After a careful examination, I find that it contains some novel techniques that should be highlighted. Could you add an additional section (Maybe in the appendix) to discuss the novel techniques (not just bounds) different from those used in Li et al.2025(b), and why these calculations can lead to a better result? In my opinion, this might be very helpful to readers. Same as the analysis of Heun's method, for example, Huang et al. 2025.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)