ParaFlow: Parallel Sampling for Flow Matching Models

Conference: ECCV 2026
Paper: ECCV Official
Code: https://github.com/Jianrong-Lu/ParaFlow
Area: Image Generation
Keywords: Flow Matching, Parallel Sampling, Triangular Nonlinear Equations, Velocity Extrapolation, Sampling Acceleration

TL;DR

ParaFlow reformulates the autoregressive ODE sampling of Flow Matching models into a system of Triangular Nonlinear Equations (TNEs) and introduces a lightweight velocity extrapolation scheme via Newton's first-order divided difference, achieving up to a 4.3x wall-clock speedup on a single GPU with fewer total function evaluations (NFEs) than sequential Euler integration.

Background & Motivation

Flow Matching (FM) has established itself as a leading paradigm in continuous-time generative modeling, serving as the foundational engine for cutting-edge text-to-image and video synthesis architectures such as Stable Diffusion 3, Flux, and HunyuanVideo. By parameterizing a time-dependent vector field, Flow Matching maps simple Gaussian priors to complex empirical data distributions along ordinary differential equation (ODE) trajectories. However, generating samples requires numerically integrating this ODE step by step in an autoregressive fashion. This dependency forces systems to execute dozens or hundreds of strictly sequential forward evaluations through massive multi-billion parameter neural networks, posing a severe computational bottleneck for interactive and real-time generation.

Prior acceleration attempts have predominantly followed three distinct directions, each encumbered by notable trade-offs. Tailored high-order numerical solvers (such as DPM-Solver and Bespoke Solvers) reduce the total number of discretization steps, yet their execution remains strictly sequential. Model distillation and consistency methods (such as Latent Consistency Models and Sana-Sprint) enable few-shot generation, but require computationally intensive re-training and frequently degrade fine-grained visual details and text prompt alignment. Recent parallel sampling frameworks developed for diffusion models (such as ParaDiGMS and ParaSolver) introduce Picard iterations to predict and refine steps concurrently; however, they impose an exorbitant computational "tax", often doubling or tripling the total number of function evaluations (NFEs). This inflated workload necessitates multi-GPU clusters with high-bandwidth interconnects and leads to severe speed slowdowns when deployed on a standard single GPU.

The key insight of this paper lies in exploiting the geometric dynamics unique to Flow Matching models. Unlike the erratic and highly curved probability flows observed in score-based diffusion setups, Flow Matching models naturally optimize for straight-line probability paths, resulting in exceptionally smooth and locally near-linear trajectories. Core idea: recast autoregressive Flow Matching sampling into a system of Triangular Nonlinear Equations (TNEs), and leverage the high temporal smoothness of flow trajectories via Newton's first-order divided difference velocity extrapolation to achieve lossless parallel sampling on a single GPU with lower total NFEs than sequential baselines.

Method

Overall Architecture

ParaFlow transforms the stepwise autoregressive ODE integration process into a root-finding problem over a system of nonlinear algebraic equations defined over sampling trajectories. This formulation decouples step-level dependencies and enables parallel execution across GPU hardware threads. The complete inference framework proceeds through three integrated stages: first, the discretized ODE trajectory within a temporal window is converted into Triangular Nonlinear Equations (TNEs), allowing intermediate latent states to be solved simultaneously via fixed-point iteration (FPI); second, to eliminate the computational tax of evaluating full neural network forward passes across all parallel steps, the algorithm computes Newton's first-order divided difference from reference states to extrapolate future velocity fields; third, an adaptive sliding window paired with an error threshold criterion dynamically tracks convergence and strides through the full trajectory.

%%{init: {'flowchart': {'rankSpacing': 24, 'nodeSpacing': 28, 'padding': 6, 'wrappingWidth': 400}}}%%
flowchart TD
    A["Initial prior noise x0 ~ p0<br/>and initial velocity evaluation"] --> B["Triangular Nonlinear Equation Reformulation<br/>cast sequential ODE into parallel TNEs"]
    B --> C["Newton Divided Difference Velocity Extrapolation<br/>evaluate r steps and extrapolate remainder"]
    C --> D["Adaptive Sliding Window Mechanism<br/>residual threshold check and dynamic striding"]
    D -->|not reached final step tN| C
    D -->|reached final step tN| E["Output final generated image xtN"]

Key Designs

1. Triangular Nonlinear Equation Reformulation: converting autoregressive integration into parallel root finding

Standard Euler integration progresses autoregressively: \(x_{t_{i+1}} = x_{t_i} + h_i v(x_{t_i}, t_i, \theta)\). To eliminate this sequential barrier, ParaFlow treats the entire sequence of states \(\{\hat{x}_{t_0}, \dots, \hat{x}_{t_N}\}\) as unknown variables and formulates a global system of \(N+1\) equations \(\mathcal{F}(\hat{x}_{t_0}, \dots, \hat{x}_{t_N}) = 0\), where each state satisfies: $$ \hat{x}{t_i} = \hat{x}, t_j, \theta) $$ Because of its triangular dependence structure, the system possesses a unique root that strictly coincides with the ground-truth autoregressive ODE trajectory. ParaFlow applies fixed-point iteration to update all candidate states concurrently at iteration round } + \sum_{j=0}^{i-1} h_j \cdot v(\hat{x}_{t_j\(k\). Theoretically, the underlying update operator exhibits strict nilpotency, guaranteeing exact convergence to the autoregressive trajectory in at most \(N\) iterations. In practice, the system reaches high numerical precision in far fewer iterations.

2. Velocity Extrapolation via Newton's Divided Difference: estimating future velocities with negligible compute

Evaluating the complete neural network across all candidate steps during fixed-point iterations introduces a severe computational tax, multiplying total NFEs and eroding wall-clock gains. To overcome this limitation, ParaFlow capitalizes on the smooth, near-linear characteristics of Flow Matching vector fields. In each parallel execution window, the model explicitly evaluates the neural network for only \(r\) timesteps (with \(r=1\) for single-GPU setups), while the remaining \(p-r\) future velocities are approximated via Newton's first-order divided difference polynomial extrapolation: $$ D_1 = \frac{v_{t_{i+r-1}} - v_{t_{i+r-2}}}{t_{i+r-1} - t_{i+r-2}}, \quad v(\hat{x}{t) $$ Unlike naive zero-order replication, the first-order divided difference explicitly captures the local temporal acceleration trends of the vector field. The authors mathematically prove that under spatial Lipschitz continuity and bounded time-derivatives, the cumulative trajectory approximation error decays cubically at }}, t_{i+j}, \theta) \approx v_{t_{i+r-1}} + D_1 (t_{i+j} - t_{i+r-1\(\mathcal{O}(1/N^3)\) with respect to total steps \(N\). For practical schedules (\(N \ge 25\)), the denominator dominates, rendering approximation errors statistically imperceptible while reducing total NFEs below the sequential baseline.

3. Adaptive Sliding Window Mechanism: block-parallel propagation via residual convergence checks

Attempting to solve all \(N\) steps simultaneously across the entire interval \([0, 1]\) would incur excessive GPU memory overhead and accumulate early estimation noise. ParaFlow introduces a sliding window of size \(p \ll N\) (defaulting to \(p=8\)) to restrict parallel iterations to local subproblems. Newly entered window states are initialized using the last converged state from the preceding window. Following each iteration round, the framework computes the normalized mean squared residual between consecutive iterates: \(e_j = \frac{1}{D} \|\hat{x}^{(k)}_{t_{i+j+1}} - \hat{x}^{(k-1)}_{t_{i+j+1}}\|^2\). By comparing against a tolerance threshold \(\delta\), the algorithm locates the first unconverged index \(s\). When the first \(s\) contiguous states satisfy \(e_j \le \delta^2\), the window shifts forward by \(s\) steps, freezing the converged prefix. This strategy ensures adaptive numerical striding: advancing rapidly across linear trajectories while iterating carefully through high-curvature segments.

A Worked Example

Consider generating an image using FLUX with total steps \(N=50\), window size \(p=8\), evaluation count \(r=1\), and tolerance \(\delta=0.01\): 1. Window Initialization: For the initial window \(t \in [t_0, t_7]\), all 8 latent state guesses are initialized to the starting Gaussian noise \(x_{t_0}\). 2. Parallel Evaluation & Extrapolation: In the first iteration, only a single neural network forward pass is executed on the GPU for \(t_0\) to obtain \(v_{t_0}\). Using historical divided difference \(D_1\), the velocities for the remaining 7 steps (\(t_1 \dots t_7\)) are computed instantly via closed-form polynomial extrapolation. 3. Parallel State Update & Residual Check: All 8 intermediate states are updated in parallel using the cumulative step sums. The change from the prior iterate is evaluated against the threshold \(\delta^2 = 1 \times 10^{-4}\). If the first 3 steps (\(t_1, t_2, t_3\)) meet the convergence criterion while \(t_4\) does not, the stride is set to \(s=3\). 4. Window Shift: States \(t_1, t_2, t_3\) are locked as converged. The sliding window shifts forward by 3 steps to span up to \(t_{10}\), with new incoming slots inheriting the value of \(t_3\). Rather than requiring 50 sequential network evaluations, the entire sampling run completes in approximately 16 lightweight iterations.

Key Experimental Results

Main Results

Experiments were conducted on Ascend 910B hardware across two text-to-image Flow Matching models: Stable Diffusion 3 Medium and FLUX.1-dev. Evaluations were performed on 5,000 random prompts sampled from the FluxPrompting benchmark at \(1024 \times 1024\) resolution, measuring Fréchet Inception Distance (FID), CLIP score, wall-clock generation time per image, and overall speedup relative to the standard sequential Euler solver.

The quantitative comparison on Stable Diffusion 3 across different sampling budgets is summarized below:

Sampling Steps \(N\)	Method & Hardware Config	Iterations (Iters)↓	Evaluations (NFE)↓	CLIP Score↑	FID↓	Latency (s)↓	Speedup↑
100 steps	Euler baseline (1 GPU)	100.00	100.00	32.38	61.54	18.22	1.0×
	Euler + ParaFlow (\(r=100\), w/o approx, 8 GPUs)	26.47	211.78	32.37	61.74	5.55	3.3×
	Euler + ParaFlow (\(r=1\), w/ approx, 1 GPU)	28.34	28.34	32.17	61.96	5.20	3.5×
75 steps	Euler baseline (1 GPU)	75.00	75.00	32.37	62.14	13.57	1.0×
	Euler + ParaFlow (\(r=75\), w/o approx, 8 GPUs)	22.96	183.66	32.37	62.78	4.95	2.7×
	Euler + ParaFlow (\(r=1\), w/ approx, 1 GPU)	26.61	26.61	32.17	62.77	4.81	2.8×
50 steps	Euler baseline (1 GPU)	50.00	50.00	32.37	62.24	9.37	1.0×
	Euler + ParaFlow (\(r=50\), w/o approx, 8 GPUs)	19.94	159.49	32.36	62.79	4.43	2.1×
	Euler + ParaFlow (\(r=1\), w/ approx, 1 GPU)	23.29	23.29	32.20	62.65	4.42	2.1×
25 steps	Euler baseline (1 GPU)	25.00	25.00	32.31	62.83	5.17	1.0×
	Euler + ParaFlow (\(r=25\), w/o approx, 8 GPUs)	16.02	128.12	32.32	62.61	3.73	1.4×
	Euler + ParaFlow (\(r=1\), w/ approx, 1 GPU)	15.97	15.97	31.86	62.64	3.23	1.6×

On the significantly larger FLUX.1-dev model under a 100-step schedule, standard sequential Euler takes 87.94s, whereas ParaFlow (\(r=1\), 1 GPU) completes generation in 20.45s, achieving a 4.3× wall-clock speedup while slashing the required NFE from 100 to 22.93 and maintaining virtually identical FID (62.28 vs. 61.47).

Ablation Study

To rigorously evaluate the effect of the number of parallel neural network evaluations (\(r\)) per window, an ablation was conducted on FLUX with \(N=50\) steps, \(512 \times 512\) resolution, tolerance \(\delta=0.01\), and window size \(p=50\):

Method	Evaluation Count \(r\) (Hardware Config)	Iterations (Iters)↓	Evaluations (NFE)↓	CLIP Score↑	FID↓	Latency (s)↓	Speedup↑
Standard Euler	50 (1 GPU, sequential baseline)	50.00	50.00	32.13	59.96	23.80	1.0×
ParaFlow	50 (8 GPUs, no extrapolation)	14.76	472.44	32.14	59.94	18.31	1.3×
ParaFlow	25 (8 GPUs)	14.32	368.02	32.15	59.97	12.53	1.9×
ParaFlow	8 (8 GPUs)	13.58	116.34	32.17	59.94	8.50	2.8×
ParaFlow	4 (4 GPUs)	15.34	64.98	32.16	59.95	7.93	3.0×
ParaFlow	1 (1 GPU, full extrapolation)	16.00	16.00	32.14	59.96	7.44	3.2×

In addition, analyzing window size \(p\) demonstrates that \(p=8\) delivers the optimal performance trade-off, reaching a 2.8× speedup on FLUX 50-step generation with zero fidelity drop. Relaxing tolerance \(\delta\) from 0.005 to 0.01 allows practitioners to dial in additional speedups while maintaining high visual quality.

Key Findings

• Eliminating the computational tax of parallel solvers: As shown in the ablation study, reducing \(r\) from 50 to 1 slashes total NFE from 472.44 down to 16.00. Consequently, wall-clock latency drops from 18.31s on 8 GPUs to 7.44s on a single GPU, boosting speedup from 1.3× to 3.2×. This invalidates the common assumption that parallel sampling inherently requires massive multi-GPU concurrency to offset evaluation overhead.
• Cubic error decay preserves generation quality: Across all benchmarks on both SD3 and FLUX, FID and CLIP scores under ParaFlow remain indistinguishable from the sequential Euler baseline (FID deviations remain strictly \(<0.5\)). This confirms the theoretical \(\mathcal{O}(1/N^3)\) global error bound and demonstrates that first-order divided difference extrapolation introduces no perceptual artifacts.
• Amplified speedup on large-scale models and schedules: The efficiency advantage of ParaFlow grows substantially as model parameter size expands and step counts increase. On FLUX with 100 steps, single-GPU speedup reaches 4.3×, effectively bridging the throughput gap for large-scale Diffusion Transformers (DiTs).

Highlights & Insights

• Rooted in flow geometry rather than brute-force iteration: Rather than applying heavy diffusion Picard iterations, ParaFlow exploits the near-linear probability paths of Flow Matching models, substituting costly multi-step neural network passes with closed-form divided difference extrapolation.
• Single-GPU friendly and completely training-free: By eliminating the NFE penalty, ParaFlow is the first parallel generative sampler to achieve superior acceleration on a single consumer-grade or standard GPU without requiring model fine-tuning, architectural changes, or distributed communication infrastructure.
• Dynamic adaptive striding: The convergence-driven sliding window mechanism acts as a numerical wave front, advancing aggressively over smooth trajectory segments while automatically dedicating iterations to high-curvature regions.

Limitations & Future Work

• Diminishing returns in extreme few-step regimes: In ultra-low-step regimes (\(N \le 10\)), large step intervals reduce the local smoothness of the vector field, narrowing the wall-clock speedup margin to 1.4–1.6×.
• Exploration of higher-order polynomial extrapolation: The current framework relies on first-order divided difference; incorporating adaptive higher-order Hermite interpolation could further mitigate error accumulation under complex classifier-free guidance scales.
• Extension to spatio-temporal video flow models: While demonstrated on text-to-image synthesis, extending ParaFlow to video Flow Matching models (e.g., HunyuanVideo, Sana-Video) offers even greater acceleration potential due to the heavy autoregressive cost of 3D spatio-temporal representations.

Related Work & Insights

• vs ParaDiGMS / ParaSolver: ParaDiGMS and ParaSolver target diffusion models and suffer from a 2–3× NFE computational tax that requires multi-GPU clusters; ParaFlow is custom-tailored for Flow Matching dynamics and achieves lower total NFEs than sequential solvers on a single GPU.
• vs DPM-Solver / Bespoke Solvers: High-order ODE solvers optimize the step schedule but remain bound to sequential execution; ParaFlow provides an orthogonal acceleration dimension through step-level parallelism and can be integrated into advanced solver pipelines.
• vs Model Distillation (LCM / Consistency Models / Sana-Sprint): Distillation compresses generation to 1–4 steps but incurs heavy re-training costs and often sacrifices sample fidelity; ParaFlow is fully training-free, preserving 100% of the original model's visual quality and generalization capability.

Rating

• Novelty: ⭐⭐⭐⭐⭐ Formulates the first parallel sampling framework for Flow Matching models via Triangular Nonlinear Equations and introduces a zero-overhead divided difference velocity extrapolation scheme.
• Experimental Thoroughness: ⭐⭐⭐⭐⭐ Validated thoroughly across two modern generative models (SD3 and FLUX) across varying steps, resolutions, evaluation counts, and window sizes.
• Writing Quality: ⭐⭐⭐⭐⭐ Clear mathematical derivations, coherent motivation, and well-structured experimental analyses.
• Value: ⭐⭐⭐⭐⭐ Solves the longstanding computational tax problem of parallel samplers, providing an immediately deployable acceleration tool for resource-constrained environments.