FD-Bench: A Modular and Fair Benchmark for Data-driven Fluid Simulation

## Abstract

Abstract Data-driven modeling of fluid dynamics has advanced rapidly with neural PDE solvers, yet a fair and strong benchmark remains fragmented due to the absence of unified PDE datasets and standardized evaluation protocols. Although architectural innovations are abundant, fair assessment is further impeded by the lack of clear disentanglement between spatial, temporal and loss modules.
In this paper, we introduce FD-Bench , the first fair, modular, comprehensive and reproducible benchmark for data-driven fluid simulation. FD-Bench systematically evaluates nearly 85 baseline models across 10 representative flow scenarios under a unified experimental setup. It provides four key contributions: (1) a modular design enabling fair comparisons across spatial, temporal, and loss function modules; (2) the first systematic framework for direct comparison with traditional numerical solvers; (3) fine-grained generalization analysis across resolutions, initial conditions, and temporal windows; and (4) a user-friendly, extensible codebase to support future research. Through rigorous empirical studies, FD-Bench establishes the most comprehensive leaderboard to date, resolving long-standing issues in reproducibility and comparability, and laying a foundation for robust evaluation of future data-driven fluid models. The code is open-sourced at https://anonymous.4open.science/r/FD-Bench-15BC .

## 1 Introduction

The accurate simulation of complex fluid dynamical systems governed by partial differential equations (PDEs) is a cornerstone of scientific and engineering applications, ranging from aerodynamics [78, 21, 75], chemical engineering [16, 15, 136], biology [130, 106, 93], and environmental science [116, 79, 3, 42, 48, 84, 133]. Traditional numerical solvers have achieved remarkable success in offering high-fidelity solutions grounded in well-established mathematical formulations and rigorous convergence guarantees.
Recent advancements in machine learning have introduced numerous neural solvers for data-driven fluid simulation [107]. Data-driven approaches enable efficient prediction of complex dynamics with reduced computational cost. These methods have shown strong potential in capturing nonlinear behaviors beyond the reach of traditional solvers.

Despite rapid progress, the field faces significant challenges in establishing consistent evaluation standards.
Researchers in both the field of AI and traditional Computational Fluid Dynamics (CFD) often face three pressing questions. First, how effective are current neural solvers, and which methods demonstrate the best performance? Second, under specific conditions, how do neural solvers compare to conventional numerical solvers in balancing accuracy, efficiency, and generalization?
Third, and perhaps most critically, how well do these neural solvers generalize when applied to more practical and realistic scenarios? Therefore, as recently emphasized by [6], recent analysis reveals persistent challenges in making fair comparisons to prior methods. Advancing the field will require stronger and more consistent benchmark problems.

To date, the three aforementioned questions remain unanswered, largely due to the inherent infinite-dimensional continuity
and variability of fluid dynamics, which makes definitive analysis challenging.
Significant limitations remain in three key areas: the formulation and accessibility of simulation data, the design of neural solver architectures, and the standardization of evaluation protocols, as illustrated in Figure [1](https://arxiv.org/html/2505.20349v1#S1.F1).
(1) Lack of unified PDE systems. From a data perspective, the diversity of benchmark settings, ranging in governing equations, domain geometries, spatial resolutions, temporal discretizations, and flow conditions, hinders fair comparison across methods. This variability often results in fragmented and non-reproducible findings. As shown in Table [1](https://arxiv.org/html/2505.20349v1#S1.T1), no two benchmarks employ the same dataset, making it difficult to draw reliable conclusions about model performance.

Figure: Figure 1: Limitations identified in three key areas.
Refer to caption: x3.png

(2) Entangled spatiotemporal modeling. Capturing fluid dynamics requires architectures that jointly model spatial structures and temporal evolution while maintaining physical consistency. However, many neural solvers entangle components from diverse domains without clear attribution, and lack systematic ablations, making it difficult to isolate the source of performance gains. For example, FNO [58] learns spatiotemporal features by stacking the time dimension in the frequency domain, while MP-PDE [8] introduces temporal bundling to enhance temporal evolution. These differing design choices make it nearly impossible to directly compare their individual contributions.
(3) Lack of standardized evaluation protocol. From an evaluation standpoint, the dynamical nature of fluid systems, characterized by sensitivity to initial conditions and parameter settings, makes it challenging to derive robust conclusions or establish reproducible baselines. Table [1](https://arxiv.org/html/2505.20349v1#S1.T1) presents a comparison of key aspects from our systematically constructed evaluation. This reveals three primary limitations:
(3.1) Inconsistent experimental configurations. It is common for different studies to adopt distinct experimental setups, including variations in discretization, flow field conditions, forecasting windows, and other hyper-parameters. Moreover, evaluation metrics are often not standardized, further limiting comparability across methods.
(3.2) Limited benchmarking against traditional solvers. Many neural solvers lack rigorous comparisons with classical numerical solvers, particularly under settings where accuracy, stability, and computational efficiency must be jointly evaluated.
(3.3) Absence of systematic generalization analysis. Few studies examine how models generalize across different spatial and temporal resolutions, boundary and initial conditions, or parameter regimes. This limits our understanding of whether models can extend to real-world or industrial-scale applications.

To address these limitations, this paper proposes a fair and comprehensive benchmark, FD-Bench, which systematically re-implements and evaluates 85 baseline models within a consistent experimental setup across 10 representative fluid flow scenarios.
FD-Bench brings four key innovations to the field.
First, it enables modular and fair comparisons by isolating spatial, temporal, and loss function modules, thereby eliminating confounding factors introduced by implementation-specific details. We construct the most rigorous and comprehensive leaderboard to date, providing a strong foundation for fair and reproducible evaluation in the field.
Second, FD-Bench establishes the first systematic framework enabling direct comparison with traditional numerical solvers, achieving a favorable balance between accuracy and computational efficiency.
Third, FD-Bench enables fine-grained generalization analysis, revealing how variations in initial conditions, spatial resolution, temporal windows, and model capacity affect performance in terms of accuracy, computational efficiency, and robustness across diverse physical regimes.
Finally, FD-Bench provides an easy-to-use and reproducible codebase, significantly lowering the barrier to entry for both researchers and practitioners. Simply specify the dataset and input a prompt with the required module combination to execute the codes. The modular design ensures that any baseline method can be reproduced or extended to your own design with minimal configuration.

A comprehensive comparison with existing benchmarks is presented in Table [1](https://arxiv.org/html/2505.20349v1#S1.T1), highlighting the unique scope, modularity, and evaluation depth of FD-Bench.
This work not only resolves long-standing comparability issues but also delivers a flexible toolkit for developing next-generation fluid dynamical system models.

**Table 1: A comprehensive comparison between our proposed FD-Bench and existing benchmarks. Each column represents a key evaluation dimension, where “comp." stands for “comparison". # Baselines indicates the number of existing methods reproduced and evaluated within each benchmark.**
| Benchmark | 1. Different<br>data topology | 2. Multiple<br>time steps | 3. Generalization<br>study on conditions | 4. Numerical<br>method comp. | 5. Temporal<br>modular comp. | 6. Spacial<br>modular comp. | 7. Loss<br>modular comp. | # Baselines |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
| PDEBench [100] | ✗ | ✓ | ✗ | ✓ | ✗ | ✗ | ✗ | 3 |
| CFDBench [70] | ✗ | ✓ | ✓ | ✗ | ✗ | ✗ | ✗ | 9 |
| PINNacle [32] | ✗ | ✓ | ✓ | ✗ | ✗ | ✗ | ✗ | 12 |
| DiffBench [46] | ✗ | ✓ | ✗ | ✗ | ✗ | ✗ | ✗ | 6 |
| Well [77] | ✗ | ✓ | ✗ | ✗ | ✗ | ✗ | ✗ | 4 |
| ML-PDE [72] | ✗ | ✗ | ✗ | ✓ | ✗ | ✗ | ✗ | 10 |
| PDENNEval [115] | ✗ | ✓ | ✓ | ✗ | ✗ | ✗ | ✓ | 12 |
| Posiden [35] | ✗ | ✓ | ✓ | ✗ | ✗ | ✓ | ✗ | 7 |
| FD-Bench | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | 85 |

## 2 FD-Bench

This benchmark aims to investigate four key questions in fluid simulation:

Q1. Does an optimal neural architecture exist? Modeling complex fluid dynamical systems often suffers from the “coupling curse,” [51, 123] where spatio-temporal flow conditions are inextricably intertwined, making it difficult to pinpoint the true source of performance gains.
Traditional benchmarks treat each baseline method as a monolithic whole, without probing the sources of its contributions and consequently lack a principled design.
Our FD-Bench systematically decouples fluid dynamical system modeling into three orthogonal dimensions: spatial representation, temporal representation, and loss function, shown in Figure [2](https://arxiv.org/html/2505.20349v1#S3.F2).
By isolating each axis, our design (1) breaks the coupling curse, enabling direct attribution of model behavior to specific design choices; (2) supports composable innovation, allowing, for example, a graph-based spatial module to pair seamlessly with an implicit neural representation temporal solver; and (3) yields taxonomy-driven insights by aligning 85 published baselines within a unified schema. More details are in Section [3](https://arxiv.org/html/2505.20349v1#S3).

Q2. Can neural solvers truly replace traditional numerical solvers?
Existing works [72, 115] often neglect comparison with traditional numerical solvers or offer unfair comparisons. To obtain a fair comparison of both computational efficiency and predictive accuracy, we propose benchmarking neural solvers against traditional numerical methods operating on coarser discretizations or using lower-order time-integration schemes (e.g., Euler’s method instead of fourth-order Runge–Kutta) that produce equivalent error. Thus, we extend the idea in [72] to compare the rollout error of the FNO [58] to that of these coarse-grid solvers. We perform experiments on three subsets with varying parameters, evaluate their long-horizon rollout performance, and ensure a fair comparison under identical system environments.

Q3. Which discretization scheme better supports their representative neural solvers?
Neural solvers employ one of two discretization schemes: Eulerian or Lagrangian. In Eulerian approaches (grid- or mesh-based) [98, 97], the spatial domain is discretized by a fixed arrangement of nodes. In Lagrangian (particle-based) methods [76], discretization is realized through material points that move with the local deformation of the continuum. Existing works only compare models with the same discretization schemes.
We investigate which discretization scheme most effectively supports its representative neural solvers in capturing core fluid dynamic behavior. To this end, we generate three subsets to train their representative neural solvers, and compare their rollout performance.

Q4. Whether knowledge learned from different systems boost generalization? As a data-driven approach, a neural solver is naturally influenced by the training data. Therefore, in practical applications, evaluating a neural solver’s performance in generalization scenarios is of greater significance. However, existing benchmarks offer very limited exploration in this regard. In Section [5.5](https://arxiv.org/html/2505.20349v1#S5.SS5), we provide a detailed investigation into how different modular designs perform under varying initial conditions, resolutions, model parameters, and longer rollout horizons. Specifically, we examine the zero-shot generalization ability of the models to unseen initial conditions, their convergence behavior when trained at different resolutions, their scalability with increased parameter size, and their ability to extrapolate in zero-shot rollout over extended time steps.

## 3 Decoupling in Modular Design

Figure: Figure 2: A schematic illustration of common approaches for each key module in data-driven neural PDE solvers. Note that self-attention is also applicable to temporal representation.
Refer to caption: x4.png

### 3.1 Spatial Representation

A key challenge is how to encode the continuous spatial field
$\displaystyle u\colon\Omega\times[0,T]\to\mathbb{R}^{d}$.
Let
$\Omega\subset\mathbb{R}^{d}$ denote the physical domain,
and let $u(\bm{x},t)\in\mathbb{R}^{c}$ be the physical variable (e.g., velocity, pressure), so that $u(\bm{x},t)\in\mathbb{R}^{c}$ represents the value of the field at spatial location $\bm{x}\in\Omega$ and time $t$. $c\in\mathbb{N}$ is the number of channels (for example, three velocity components and pressure). We collect $N$ spatial samples
$\{\bm{x}_{i}\}_{i=1}^{N}\subset\Omega$
(e.g., grid points, mesh nodes, or particles) and denote $u(\cdot,t)=\bigl{[}u(\bm{x}_{1},t),\,u(\bm{x}_{2},t),\,\dots,\,u(\bm{x}_{N},t)
\bigr{]}^{\top}\in\mathbb{R}^{N\times c}$.
Spatial representation is any map that encodes the continuous field $u(\bm{x},t)$ at time $t$
into a finite-dimensional feature $z$ as ${\mathcal{S}}_{\theta}:\bigl{\{}u(\cdot,t)\bigr{\}}\;\mapsto\;z\;\in\;\mathbb{
R}^{D}$.

Here, we consider the commonly used spatial representation paradigms in current neural operators.
Frequency-Mixing.
It transforms the spatial field to the spectral domain using frequency transform [57, 55, 134] (e.g., Discrete Fourier Transform), yielding spectral coefficients $\hat{u}(\mathbf{k},t)=\mathcal{F}(u)$, where $\mathbf{k}=(k_{1},\dots,k_{d})$ denotes the discrete frequency. A learnable spectral filter $\phi_{\theta}(\mathbf{k})$ is applied by element-wise multiplication to the Fourier coefficients $\hat{u}(\mathbf{k},t)=\mathcal{F}(u(\cdot,t))$, thereby efficiently capturing long-range correlations via global mixing in physical space while decoupling spectral modes to learn large-scale coherent structures.

$$ $\mathcal{S}_{\mathrm{Fourier}}(u(\cdot,t))=\mathcal{F}^{-1}(\phi_{\theta}\odot \mathcal{F}(u(\cdot,t)))\in\mathbb{R}^{N\times c}$ (1) $$

Self-attention. It allows each spatial location to attend to all others by computing weighted combinations of feature vectors based on learned similarity scores [120, 23], thereby capturing global interactions adaptively.
Given learnable projections $\mathbf{W}_{Q},\mathbf{W}_{K}\in\mathbb{R}^{c\times d_{k}}$ and $\mathbf{W}_{V}\in\mathbb{R}^{c\times d_{v}}$:

$$ $\mathcal{S}_{\mathrm{SA}}\left(u\left(\cdot,t\right)\right)=\mathrm{softmax} \bigl{(}\tfrac{1}{\sqrt{d_{k}}}\,(u(\cdot,t)\,\mathbf{W}_{Q})\,(u(\cdot,t)\, \mathbf{W}_{K})^{\top}\bigr{)}\,(u(\cdot,t)\,\mathbf{W}_{V})\;\in\;\mathbb{R}^ {N\times d_{v}}$ (2) $$

Convolution.
Convolution [86] aggregates information from a local neighborhood using a shared kernel, leveraging spatial locality and translation equivariance. Let $\mathcal{H}\subset\mathbb{Z}^{d}$ be the set of stencil offsets. Let $c_{\mathrm{in}}=c$ and $c_{\mathrm{out}}$ be the number of input and output channels. For each offset $h\in\mathcal{H}$, the kernel weight is $\mathbf{W}_{h}\in\mathbb{R}^{c_{\mathrm{out}}\times c_{\mathrm{in}}}$. Then, for each spatial index $i$:

$$ $\mathcal{S}_{\mathrm{Conv}}(u(\cdot,t))_{i}=\sum_{h\in\mathcal{H}}\mathbf{W}_{ h}\;u\bigl{(}\bm{x}_{i+h},t\bigr{)}\;\in\;\mathbb{R}^{c_{\mathrm{out}}}$ (3) $$

Graph.
Graph convolutions [8, 62, 108] generalize standard convolutional filtering to irregular domains by aggregating neighbor features weighted by graph connectivity, enabling modeling on meshes or point clouds. A typical form is:

$$ $\mathcal{S}_{\mathrm{Graph}}(u(\cdot,t))=\mathrm{softmax}(\tilde{D}^{-\tfrac{1 }{2}}\tilde{A}\tilde{D}^{-\tfrac{1}{2}}\,u(\cdot,t)\,W)\in\mathbb{R}^{N\times D}$ (4) $$

where $A$ is the binary adjacency, $\tilde{D}^{-1/2}\tilde{A}\tilde{D}^{-1/2}$ symmetrically normalizes the adjacency to control spectral properties, and $\sigma(\cdot)$ is a pointwise nonlinearity.

Reduced-Order Modeling.
Typical work like Proper Orthogonal Decomposition (POD) [87, 117] approximates $u(\cdot,t)$ by projecting onto the $K$ leading modes $\{\phi_{k}\}$ found via the covariance operator
$C(\bm{x},\bm{x}^{\prime})=\frac{1}{M}\sum_{i=1}^{M}u(\bm{x},t_{i})\,u(\bm{x}^{
\prime},t_{i})$,
and solving
$\int_{\Omega}C(\bm{x},\bm{x}^{\prime})\phi_{k}(\bm{x}^{\prime})\,\mathrm{d}\bm
{x}^{\prime}=\lambda_{k}\phi_{k}(\bm{x})$
under the inner product $\langle f,g\rangle=\int_{\Omega}f(\bm{x})g(\bm{x})\,\mathrm{d}\bm{x}$. The projection coefficients
$\alpha_{k}(t)=\langle u(\cdot,t),\phi_{k}\rangle$
give:

$$ $\mathcal{S}_{\mathrm{ROM}}(u(\cdot,t))=[\alpha_{1}(t),\dots,\alpha_{K}(t)]^{ \top}\in\mathbb{R}^{K}$ (5) $$

where $M$ is the number of snapshots, $N$ the spatial samples, and $K\ll N$ the retained modes.
Implicit Neural Representation.
It yields a continuous, differentiable approximation [17] of the field with sub-grid resolution and compact storage. The continuous field $u(\cdot,t)$ is modeled by implementing a coordinate-to-value mapping with sampled 2D coordinates $\{({x}_{i},{y}_{i})\}_{i=1}^{N}$ as:

$$ $\mathcal{S}_{\mathrm{INR}}\bigl{(}u(\cdot,t)\bigr{)}=\arg\min_{\theta\in \mathbb{R}^{D}}\frac{1}{N}\sum_{i=1}^{N}\bigl{\|}f_{\theta}(({x}_{i},{y}_{i})) -u(({x}_{i},{y}_{i})),t)\bigr{\|}^{2}$ (6) $$

### 3.2 Temporal Representation

Temporal representation encodes the evolution of the fluid field over time, providing essential features for predicting the future state of the flow field. Following the definition in Section [3.1](https://arxiv.org/html/2505.20349v1#S3.SS1), a discrete sequence of $M$ successive snapshots $\{u(t_{j})\}_{j=1}^{M}\;\subset\;\mathbb{R}^{N\times c}$ are learned with a mapping $\mathcal{T}_{\theta}\bigl{(}\{u(t_{j})\}_{j=1}^{M}\bigr{)}=h\;\in\;\mathbb{R}^
{D}$,
where $\theta$ parameterizes the temporal encoder and $h$ is a $D$-dimensional temporal feature.

Autoregression.
Autoregressive approaches [58, 56] decompose the sequence model into successive one-step predictions using a sliding history window of length $k$. If the input window at step $j$ is $\bigl{[}u(t_{j-k}),\,u(t_{j-k+1}),\,\dots,\,u(t_{j-1})\bigr{]}\;\in\;\mathbb{R
}^{k\times N\times c}$.
Then for $j$ from $k+1$ to $M$:

$$ $\hat{u}(t_{j+1})=\mathcal{T}_{\mathrm{AR}}\bigl{(}\bigl{[}u(t_{j-k+1}),\,\dots ,\,u(t_{j-1}),\,\hat{u}(t_{j})\bigr{]}\bigr{)}$ (7) $$

so that each new prediction is appended to the window while the oldest state is discarded.

Next-Step Prediction.
This simpler variant predicts a single step ahead from the most recent state and is often combined with a rollout mechanism [14, 22] to generate long-horizon forecasts as:

$$ $\hat{u}(t_{j})=\mathcal{T}_{\mathrm{Next}}\bigl{(}u(t_{j-1})\bigr{)},j=\{1, \dots,M\}$ (8) $$

Temporal Bundling.
It is inspired by MP-PDE [8], which stacks a fixed window of $k$ frames and applies a feed-forward or convolutional encoder to predict $k$ future steps. Thus a new vector $\mathbf{d}=\mathcal{T}_{\mathrm{TB}}(\mathbf{d}^{0})=(\mathbf{d}^{1},\cdots,
\mathbf{d}^{k})$ is used to update the solution:

$$ $u\bigl{(}t_{k+j}\bigr{)}=u\bigl{(}t_{k}\bigr{)}+\bigl{(}t_{k+j}-t_{k}\bigr{)} \,\mathbf{d}^{j},j=\{1,\dots,k\}$ (9) $$

Self-attention. It parallels the spatial representation mechanism; however, whereas attention scores are originally computed over spatial patches, they are here computed along the temporal dimension.
Neural ODE.
Neural ordinary differential equations [18, 99] provide a continuous-time framework for modeling temporal evolution in a latent space [28, 112]. Specifically, we define a hidden state $h(t)\in\mathbb{R}^{D}$ that evolves according to the parameterized dynamics:
$\frac{\mathrm{d}h}{\mathrm{d}t}=\mathcal{T}_{\mathrm{ODE}}\bigl{(}h(t),u(t_{0:
t-1}),t\bigr{)},h(t_{0})=\zeta(u({t_{0}})),u(t)=\xi(h(t))$,
where $\zeta(\cdot)$ and $\xi(\cdot)$ are parameterized to encode the dynamics of the hidden state. The sequence of latent states $\{h(t_{j})\}_{j=1}^{M}$ is obtained by solving ODE with a black-box integrator, and the final state:

$$ $h(t_{1}),\cdots,h(t_{j})=\mathrm{ODESolver}(\mathcal{T}_{\mathrm{ODE}},h(t_{0} ),(t_{0},\cdots,t_{j}))$ (10) $$

This continuous formulation naturally captures irregular time sampling and allows for arbitrary-horizon forecasting via the same learned dynamics.

### 3.3 Loss Function

In data-driven fluid simulation, the loss function choice determines not only numerical accuracy but also the ability to capture uncertainty and respect physical laws. We describe below four paradigms, each with its mathematical definition and a discussion of practical benefits and limitations.

Physical Variable Prediction Loss.
The mean squared error (MSE) [11, 30, 20, 22] between predicted and ground truth fields is commonly used in neural solvers given by:

$$ $\mathcal{L}_{\mathrm{MSE}}=\frac{1}{N}\sum_{i=1}^{N}\bigl{\|}\hat{u}(\bm{x}_{i },t)-u(\bm{x}_{i},t)\bigr{\|}_{2}^{2}$ (11) $$

This loss is straightforward to implement and optimize using standard gradient-based methods. It enforces pointwise fidelity and converges rapidly when the solution is unimodal.

Diffusion Denoising Loss.
Under the diffusion framework [36, 96], one perturbs the true field $u_{0}$ by a noise schedule
$u_{\tau}=\sqrt{\bar{\alpha}_{\tau}}\,u_{0}+\sqrt{1-\bar{\alpha}_{\tau}}\,
\varepsilon,\varepsilon\sim\mathcal{N}(0,I)$,
and trains a neural network $\varepsilon_{\theta}$ to estimate the added noise. The corresponding loss is:

$$ $\mathcal{L}_{\mathrm{Noise}}=\mathbb{E}_{u_{0},\varepsilon,{\tau}}\bigl{\|} \varepsilon-\varepsilon_{\theta}(u_{\tau},{\tau})\bigr{\|}_{2}^{2}$ (12) $$

This approach captures the full distribution of flow fields, enabling stochastic sampling and modeling of complex outcomes. Its drawbacks include significantly higher computational cost due to multiple diffusion time steps and the need for careful tuning of noise schedules and solver parameters.

Flow Matching Loss.
Flow matching [63] directly learns the instantaneous vector field by minimizing the discrepancy between a predicted velocity and a numerical approximation of the true flow with a simple (invertible) affine map $\psi$:

$$ $\mathcal{L}_{\mathrm{Flow}}=\mathbb{E}_{u,\tau}\Bigl{\|}v_{\theta}\bigl{(}\psi _{\tau}(u_{0}),\tau\bigr{)}-\frac{d}{d\tau}\psi_{\tau}(u_{0})\Bigr{\|}_{2}^{2}$ (13) $$

$\psi_{\tau}(u)$ can be $\mu_{\tau}(u_{1})+\sigma_{\tau}(u_{1})u,\tau\in[0,1]$. By training $v_{\theta}$ to match true dynamics, this loss facilitates stable multi-step and continuous forecasting, often producing sharper trajectories.

Physics-Informed Residual Loss.
Physics-informed neural networks [83, 113] enforce the governing PDE by penalizing the residual of its differential operator. Denoting the learned field $u_{\theta}$, the residual is $\mathcal{R}[u_{\theta}](x,t)=\partial_{t}u_{\theta}+\mathcal{N}(u_{\theta})-g(
x,t)$,
and the loss reads:

$$ $\mathcal{L}_{\mathrm{Residual}}=\mathbb{E}_{x,t}\bigl{\|}\mathcal{R}[u_{\theta }](x,t)\bigr{\|}_{2}^{2}$ (14) $$

Here $\mathcal{N}$ represents spatial operators (e.g., advection) and $g$ is any source terms. It enforces exact satisfaction of physical laws in the continuous limit and often generalizes well beyond training data.

## 4 Datasets

Figure: Table 2: Statistics of flows. The “shape” column represents, in order: the trajectory sample number, time steps, feature channels, and the flow field resolution.
Refer to caption: x5.png

### 4.1 Data Selection Criteria

To ensure a comprehensive evaluation, our benchmark incorporates datasets that rigorously satisfy the following requirements. (1) Diverse PDE types: Inclusion of PDE families with real-world physical significance, such as the N-S equations for fluid dynamics and Burgers’ equation for shock waves. This ensures the universality of tested methods across mathematical formulations.
(2) Heterogeneous data distributions: Coverage of varying initial conditions (e.g., smooth, discontinuous, stochastic) and boundary conditions to rigorously evaluate generalization capabilities under distribution shifts.
(3) Multi-scale field settings: Datasets spanning coarse-to-fine spatial resolutions and temporal discretizations, coupled with diverse numerical discretization schemes to assess resolution scalability and discretization invariance.
(4) Complex temporal dynamics: Inclusion of both short-term transient behaviors (e.g., rapid instabilities) and long-term equilibrium states, as well as systems with stiff dynamics or temporal discontinuities, to test robustness in temporal representation.
(5) Irregular geometries: Domains with non-Cartesian boundaries (e.g., curved and fractured).

### 4.2 Data Collection & Generation

For modular comparison, we collect three representative flow problems from existing benchmarks, namely Compressible N-S  [100], Diffusion-Reaction [100], and Kolmogorov Flow [100] for fair comparison. In addition, we generate a Stochastic N-S subset, driven by incompressible flow with stochastic forcing and initialized with 10 distinct initial conditions, to further enhance the diversity of our experiments. For more details, please refer to Appendix [F.1](https://arxiv.org/html/2505.20349v1#A6.SS1).

For comparison with traditional numerical solvers, we generate three subsets with varying parameters using a pseudo-spectral solver [10]. To simulate Incompressible N-S in vorticity form, we sample five Reynolds numbers: 50, 200, 500, 1k, and 2k. To simulate Burgers, we use five kinematic viscosities: $5\times 10^{-4}$, $1\times 10^{-3}$, $5\times 10^{-3}$, $1\times 10^{-3}$, and $5\times 10^{-2}$. Finally, to simulate Advection, we consider five spatially varying advection speeds. Details on the solver implementation and the PDE formulations are provided in Appendix [F.2](https://arxiv.org/html/2505.20349v1#A6.SS2).

For comparison on different discretizations, we use Smoothed Particle Hydrodynamics (SPH) [25] to generate three Lagrangian particle subsets, Taylor–Green vortex (Re = 100), Lid-driven cavity flow (Re = 100), and Reverse Poiseuille flow (Re = 10), each governed by the compressible N-S equations. Particle positions are recorded every 100 time steps, and neural models are trained to predict the flow evolution over these intervals. We convert Lagrangian particle datasets into Eulerian form. For grid representations, we generate uniformly spaced nodes and aggregate particle velocity information onto these nodes. For mesh representations, we randomly sample Eulerian nodes and perform the same velocity-aggregation procedure. In both cases, aggregation is performed using a quintic smoothing kernel. See Appendix [F.3](https://arxiv.org/html/2505.20349v1#A6.SS3) for additional details.

## 5 Experiments

### 5.1 Experimental Setup

Evaluation Metrics.
To ensure a comprehensive and objective evaluation, we employ a diverse suite of metrics. First, we assess the global performance using the root-mean-squared-error (RMSE) and its normalized version (nRMSE). However, these measures do not capture local performance nuances. We further incorporate additional metrics that focus on specific failure modes, the RMSE in Fourier space (fRMSE) evaluated separately in low, middle, and high-frequency regions. Moreover, we extend our assessment to include efficiency evaluations, such as computational time cost and memory consumption. More details can be seen in Appendix [G](https://arxiv.org/html/2505.20349v1#A7).

Experimental Setting.
We create a brand-new code base covering different data topologies, resolutions, time steps, and initial conditions. We use PyTorch to implement experiments on 8 $\times$ NVIDIA A6000 GPUs. For different experiments, we utilize different optimizers with different strategies, including weight decay, learning rates, and schedulers obtained
from grid search for all experiments.

**Table 3: Evaluation of different modular designs on prediction performance and computational efficiency across three representative fluid dynamics datasets. Note that the parameter counts, memory usage, and GFLOPs are averaged across the three tasks.**
|  |  | Compressible N-S | Diffusion-Reaction | Kolmogorov Flow | Efficiency |  |  |  |  |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
| Model | Param | RMSE$\downarrow$ | nRMSE$\downarrow$ | fRMSE$\downarrow$ | RMSE$\downarrow$ | nRMSE$\downarrow$ | fRMSE$\downarrow$ | RMSE$\downarrow$ | nRMSE$\downarrow$ | fRMSE$\downarrow$ | Mem | GFLOPs |
| $\mathcal{X}$ + Next-step + Physical variable |  |  |  |  |  |  |  |  |  |  |  |  |
| $\mathcal{X}$ = Graph rep | 3.13M | 0.2648 | 0.3013 | 0.0097 | 0.0128 | 0.0251 | 0.0011 | 0.0709 | 0.0240 | 0.0044 | 0.24G | 21.65 |
| $\mathcal{X}$ = Latent rep | 89.0M | 0.1043 | 0.1542 | 0.0088 | 0.0523 | 0.0208 | 0.0019 | 0.0151 | 0.0151 | 0.0012 | 1.52G | 81.96 |
| $\mathcal{X}$ = Fourier rep | 72.0M | 0.0787 | 0.1018 | 0.0083 | 0.0079 | 0.0152 | 0.0008 | 0.0044 | 0.0045 | 0.0004 | 2.43G | 38.99 |
| $\mathcal{X}$ = Self-atten rep | 81.9M | 0.0571 | 0.0695 | 0.0035 | 0.0067 | 0.0151 | 0.0007 | 0.0030 | 0.0030 | 0.0002 | 1.73G | 176.4 |
| $\mathcal{X}$ = Conv rep | 7.63M | 0.1325 | 0.1528 | 0.0074 | 0.0123 | 0.0274 | 0.0017 | 0.0675 | 0.0227 | 0.0041 | 0.42G | 21.49 |
| Self-attention representation + $\mathcal{Y}$ + Physical variable |  |  |  |  |  |  |  |  |  |  |  |  |
| $\mathcal{Y}$ = Self-attention | 36.4M | 0.2306 | 0.3085 | 0.0273 | 0.0427 | 0.1065 | 0.0042 | 0.0129 | 0.0133 | 0.0019 | 1.45G | 195.4 |
| $\mathcal{Y}$ = Neural ODE | 37.1M | 0.3357 | 0.4125 | 0.0391 | 0.0131 | 0.0295 | 0.0010 | 0.0289 | 0.0294 | 0.0031 | 1.46G | 214.1 |
| $\mathcal{Y}$ = Autoregression | 24.9M | 0.3943 | 0.4576 | 0.0478 | 0.0785 | 0.1398 | 0.0087 | 0.0182 | 0.0198 | 0.0026 | 1.24G | 87.58 |
| $\mathcal{Y}$ = Next + Rollout | 36.6M | 0.2563 | 0.3177 | 0.0309 | 0.0454 | 0.1382 | 0.0041 | 0.0085 | 0.0086 | 0.0008 | 0.82G | 79.15 |
| $\mathcal{Y}$ = Temporal Bund | 36.3M | 0.1357 | 0.1729 | 0.0122 | 0.0081 | 0.0205 | 0.0006 | 0.0049 | 0.0051 | 0.0005 | 1.71G | 97.33 |
| Fourier representation + $\mathcal{Y}$ + Physical variable |  |  |  |  |  |  |  |  |  |  |  |  |
| $\mathcal{Y}$ = Neural ODE | 10.1M | 0.2169 | 0.2882 | 0.0290 | 0.0120 | 0.0371 | 0.0008 | 0.0543 | 0.0587 | 0.0017 | 0.46G | 56.84 |
| $\mathcal{Y}$ = Autoregression | 18.0M | 0.2032 | 0.2676 | 0.0291 | 0.0734 | 0.2089 | 0.0113 | 0.0849 | 0.0871 | 0.0095 | 1.69G | 16.61 |
| $\mathcal{Y}$ = Next + Rollout | 11.1M | 0.1952 | 0.2586 | 0.0283 | 0.0180 | 0.0525 | 0.0016 | 0.0083 | 0.0085 | 0.0008 | 0.71G | 5.37 |
| $\mathcal{Y}$ = Temporal bund | 10.1M | 0.1879 | 0.2497 | 0.0190 | 0.0093 | 0.0192 | 0.0009 | 0.0076 | 0.0075 | 0.0008 | 0.59G | 3.97 |
| $\mathcal{X}$ + Next-step + Noise |  |  |  |  |  |  |  |  |  |  |  |  |
| $\mathcal{X}$ = Fourier rep | 81.2M | 0.3389 | 0.3652 | 0.0963 | 0.0819 | 0.0941 | 0.0247 | 0.1160 | 0.1129 | 0.0433 | 2.52G | 36.17 |
| $\mathcal{X}$ = Self-atten rep | 68.2M | 0.3209 | 0.3369 | 0.0765 | 0.0314 | 0.0607 | 0.0055 | 0.0900 | 0.0858 | 0.0201 | 1.41G | 51.44 |
| $\mathcal{X}$ = Conv rep | 74.2M | 0.1934 | 0.3329 | 0.0244 | 0.0425 | 0.0798 | 0.0077 | 0.1258 | 0.1125 | 0.0481 | 1.35G | 314.2 |
| $\mathcal{X}$ = Latent rep | 75.9M | 0.3518 | 0.3762 | 0.1006 | 0.0759 | 0.0962 | 0.0238 | 0.1345 | 0.1283 | 0.0559 | 1.07G | 80.72 |
| Self-attention representation + Next-step + $\mathcal{Z}$ |  |  |  |  |  |  |  |  |  |  |  |  |
| $\mathcal{Z}$ = Physical var | 36.5M | 0.0982 | 0.0871 | 0.0082 | 0.0108 | 0.0231 | 0.0011 | 0.0059 | 0.0058 | 0.0006 | 0.93G | 79.15 |
| $\mathcal{Z}$ = Noise | 68.2M | 0.3209 | 0.3369 | 0.0765 | 0.0314 | 0.0607 | 0.0055 | 0.0900 | 0.0858 | 0.0201 | 1.41G | 51.44 |
| $\mathcal{Z}$ = Flow | 42.6M | 0.1232 | 0.1494 | 0.0118 | 0.0629 | 0.1039 | 0.0098 | 0.0609 | 0.0613 | 0.0050 | 1.15G | 20.19 |
| $\mathcal{Z}$ = PDE Residual | 20.6M | 0.1373 | 0.1751 | 0.0120 | 0.1368 | 0.1652 | 0.0115 | 0.0313 | 0.0311 | 0.0032 | 0.45G | 14.22 |

### 5.2 Does an optimal neural architecture exist?

We begin by emphasizing that the comparative experiments are conducted on three representative 128$\times$128 2D flow problems to ensure fairness in evaluation. For a fair comparison, we present results under approximately consistent parameter sizes or computational costs (measured in GFlops).

To investigate the effect of spatial representations, we fix the prediction modes to y=next-step, and compare the performance of different methods. As shown in Table [3](https://arxiv.org/html/2505.20349v1#S5.T3), the self-attention-based representation consistently achieves the best results across multiple tasks, followed by the Fourier-based representation. It is worth noting that due to the high computational cost of training convolutional models, we were unable to evaluate their performance under larger parameter budgets. Additionally, for the graph-based model, we apply sampling strategies to accommodate the specific data topology before training. More analysis details (e.g., results on Stochastic N-S) are in Appendix [D](https://arxiv.org/html/2505.20349v1#A4).

To investigate temporal evolution, we fix the spatial encoder to either x=Self-attention or Fourier. As shown in Table [3](https://arxiv.org/html/2505.20349v1#S5.T3), temporal bundling consistently achieves the best performance across tasks, while maintaining relatively low computational overhead. This highlights the importance of developing modules that explicitly capture the temporal dynamics of fluid evolution. Notably, neural ODEs approach the performance of temporal bundling on several tasks, suggesting a promising direction for future research. More generation analysis is in Section [5.5](https://arxiv.org/html/2505.20349v1#S5.SS5)

Lastly, with respect to the loss function, physical variable prediction remains the mainstream approach. However, generative models such as diffusion-based and flow-based methods also demonstrate competitive performance, indicating promising directions for future investigation.

### 5.3 Can neural solvers truly replace traditional numerical solvers?

Figure: Figure 4: We compare FNO against traditional solver (i.e., pseudo-spectral solver) operating at lower resolutions on N-S equation across five Reynolds numbers, spanning laminar to turbulent regimes. Solver $x$ implies a pseudo-spectral solver operating at $x\times x$ resolution.
Refer to caption: extracted/6477461/fig/Solver/ns/RE50.png

Figure: Figure 5: Comparison of rollout performance between FNO (grid data), MeshGraphNets (mesh data), and GNS (particle data). We transform particle and mesh predictions to grid data to evaluate.
Refer to caption: extracted/6477461/fig/discretization/eulerian_mse/TGV.png

Figure: Figure 6: Runtime of FNO against pseudo-spectral solver operating at lower resolutions on the incompressible N-S. Solver $x$ implies pseudo-spectral solver operating at $x\times x$ resolution.
Refer to caption: x6.png

We present the results for incompressible N-S in Figure [4](https://arxiv.org/html/2505.20349v1#S5.F4).
Additional results are in Appendix [C](https://arxiv.org/html/2505.20349v1#A3). Across all benchmarks, FNO achieves substantially lower prediction error, especially at longer rollout steps. In terms of runtime, which is shown in Figure [6](https://arxiv.org/html/2505.20349v1#S5.F6), FNO outperforms both the high-fidelity simulator and coarse-grid solvers, delivering $10\times$ to $48\times$ speedups in the Burgers’ equation experiments. Data generated at $256\times 256$ is considered ground truth. We believe that neural solvers remain a powerful PDE simulator that could revolutionize scientific discoveries. All traditional methods, including those on coarser grids, employ adaptive time stepping governed by the Courant–Friedrichs–Lewy (CFL) condition to maintain numerical stability. Incorrect predictions may lead numerical solvers at lower resolution to take longer times. We refer the readers to Appendix [F.2](https://arxiv.org/html/2505.20349v1#A6.SS2) for solvers and Appendix [G](https://arxiv.org/html/2505.20349v1#A7) for FNO details.

### 5.4 Which discretization scheme better supports their representative neural solvers?

Figure [5](https://arxiv.org/html/2505.20349v1#S5.F5) compares the rollout MSE of three neural PDE solvers under different discretization schemes: the FNO [58] on uniform grids, MeshGraphNets (MGN) [80] on meshes, and the GNS [90] on particle representations. Across all experiments, the Eulerian discretizations achieve significantly lower rollout MSE than the Lagrangian scheme. We conclude that models employing Eulerian schemes more accurately capture fundamental fluid dynamics. In contrast, Lagrangian representations must model fine-grained details, which impedes their ability to learn robust representations and renders them susceptible to error accumulation over extended rollout time steps. We refer the readers to Appendix [G](https://arxiv.org/html/2505.20349v1#A7) for model and training hyperparameters.

Figure: Figure 7: Extrapolation evaluation on zero-shot generalization across diverse initial conditions, resolution generalization across spatial scales, and long-horizon generalization. All evaluations are performed on the Compressible N-S system using different methods for next-step prediction.
Refer to caption: x7.png

### 5.5 Whether knowledge learned from different systems can boost generalization?

Generalization on different initial conditions.
We begin by training all models on a dataset of the compressible N-S equations with Mach number $M=0.1$, and both shear and bulk viscosities set to $\eta=\zeta=10^{-8}$. After training, we directly evaluate zero-shot generalization on data generated from out-of-distribution initial conditions. As shown in the first panel of Figure [7](https://arxiv.org/html/2505.20349v1#S5.F7), all methods exhibit varying degrees of performance degradation when extrapolating to more complex turbulent regimes. Notably, spatial representations based on Fourier transforms and self-attention demonstrate superior robustness under such distribution shifts.

Generalization on different resolutions.
We train models using various methods on both the Compressible N-S and Diffusion-Reaction subsets across different spatial resolutions, and evaluate them on their corresponding test resolutions. As shown in the middle panel of Figure [7](https://arxiv.org/html/2505.20349v1#S5.F7), self-attention and latent-based models benefit from higher-resolution training, leveraging richer spatial information to improve performance. In contrast, Fourier-based spatial representations exhibit stronger performance at lower resolutions, suggesting superior inductive bias for coarse-scale generalization.

Generalization on long-horizon rollout.
We evaluate the long-horizon generalization ability of different methods under the setting y = next-step, by rolling out predictions over extended time steps. As shown in the right panel of Figure [7](https://arxiv.org/html/2505.20349v1#S5.F7), all methods experience some degree of error accumulation as the rollout length increases. However, the spatial representation based on self-attention remains the most stable, exhibiting the lowest performance degradation over time.

Generalization on scaling ability. Please refer to Appendix [E](https://arxiv.org/html/2505.20349v1#A5) for comparisons of scaling ability.

More visualization analysis can be seen in Appendix [H](https://arxiv.org/html/2505.20349v1#A8).

## 6 Promising Future & Conclusion

Our results show that (1) Self-attention is particularly effective for spatial encoding, and its performance can likely be improved further via techniques such as sparse attention, receptive field modulation, or transformer regularization. Meanwhile, frequency-based methods offer strong performance at minimal cost, highlighting the promise of embedding mathematical priors into neural solvers. (2) Simple yet effective temporal bundling strategies achieve strong performance under constrained compute budgets. This suggests a valuable research direction: accelerating temporal evolution modules like neural ODE through lightweight designs, potentially yielding a more favorable trade-off between predictive fidelity and efficiency in real-world applications. (3) While neural solvers have demonstrated promising capabilities compared to traditional numerical methods, our experiments suggest that the performance gains in terms of inference speed are less substantial than those often reported in the literature. This highlights the need for future research to go beyond improving predictive accuracy and to prioritize algorithmic and architectural strategies that effectively reduce inference latency without compromising stability. (4) We also observe that Lagrangian-based neural simulators tend to underperform in capturing essential fluid-dynamic structures. Their emphasis on simulating high-resolution local interactions leads to severe error accumulation and propagation across particles over time. To mitigate this, future directions may involve hierarchical or multi-scale architectures, such as U-Net-based designs, that can dynamically balance fine-grained resolution with global structural coherence, thus improving both stability and generalization.

Overall, this paper introduces a fair and insightful benchmark designed to systematically evaluate the performance of neural PDE solvers for fluid simulation, providing both practical comparisons and conceptual guidance for future research.

## References

- [1]
Benedikt Alkin, Andreas Fürst, Simon Lucas Schmid, Lukas Gruber, Markus Holzleitner, and Johannes Brandstetter.
Universal physics transformers: A framework for efficiently scaling neural operators.
In The Thirty-eighth Annual Conference on Neural Information Processing Systems, 2024.
- [2]
Jan-Hendrik Bastek, WaiChing Sun, and Dennis Kochmann.
Physics-informed diffusion models.
In The Thirteenth International Conference on Learning Representations, 2025.
- [3]
Kaifeng Bi, Lingxi Xie, Hengheng Zhang, Xin Chen, Xiaotao Gu, and Qi Tian.
Accurate medium-range global weather forecasting with 3d neural networks.
Nature, 619(7970):533–538, 2023.
- [4]
Suresh Bishnoi, Ravinder Bhattoo, Jayadeva Jayadeva, Sayan Ranu, and N M Anoop Krishnan.
Enhancing the inductive biases of graph neural ODE for modeling physical systems.
In The Eleventh International Conference on Learning Representations, 2023.
- [5]
Oussama Boussif, Yoshua Bengio, Loubna Benabbou, and Dan Assouline.
Magnet: Mesh agnostic neural pde solver.
Advances in Neural Information Processing Systems, 35:31972–31985, 2022.
- [6]
Johannes Brandstetter.
Envisioning better benchmarks for machine learning pde solvers.
Nature Machine Intelligence, 7(1):2–3, 2025.
- [7]
Johannes Brandstetter, Rianne van den Berg, et al.
Clifford neural layers for PDE modeling.
In The Eleventh International Conference on Learning Representations, 2023.
- [8]
Johannes Brandstetter, Daniel E. Worrall, and Max Welling.
Message passing neural PDE solvers.
In International Conference on Learning Representations, 2022.
- [9]
Andrey Bryutkin, Jiahao Huang, Zhongying Deng, Guang Yang, Carola-Bibiane Schönlieb, and Angelica I Aviles-Rivero.
HAMLET: Graph transformer neural operator for partial differential equations.
In International Conference on Machine Learning, 2024.
- [10]
Claudio Canuto, M. Yousuff Hussaini, Alfio Quarteroni, and Thomas A. Zang.
Spectral Methods in Fluid Dynamics.
Springer-Verlag, Berlin, Heidelberg, 1988.
- [11]
Qianying Cao, Somdatta Goswami, and George Em Karniadakis.
Lno: Laplace neural operator for solving differential equations.
Arxiv, 2023.
- [12]
Shuhao Cao.
Choose a transformer: Fourier or galerkin.
In A. Beygelzimer, Y. Dauphin, P. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, 2021.
- [13]
Shuhao Cao, Francesco Brarda, Ruipeng Li, and Yuanzhe Xi.
Spectral-refiner: Accurate fine-tuning of spatiotemporal fourier neural operator for turbulent flows.
In The Thirteenth International Conference on Learning Representations, 2025.
- [14]
Yadi Cao, Menglei Chai, Minchen Li, and Chenfanfu Jiang.
Efficient learning of mesh-based physical simulation with bi-stride multi-scale graph neural network.
In International Conference on Machine Learning, pages 3541–3558. PMLR, 2023.
- [15]
Yadi Cao, Xuan Gao, and Ri Li.
A liquid plug moving in an annular pipe–heat transfer analysis.
International Journal of Heat and Mass Transfer, 139:1065–1076, 2019.
- [16]
Yadi Cao and Ri Li.
A liquid plug moving in an annular pipe—flow analysis.
Physics of Fluids, 30(9), 2018.
- [17]
Honglin Chen, Rundi Wu, Eitan Grinspun, Changxi Zheng, and Peter Yichen Chen.
Implicit neural spatial representations for time-dependent pdes.
In International Conference on Machine Learning, pages 5162–5177. PMLR, 2023.
- [18]
Ricky TQ Chen, Yulia Rubanova, Jesse Bettencourt, and David K Duvenaud.
Neural ordinary differential equations.
Advances in neural information processing systems, 31, 2018.
- [19]
Zhao Chen, Yang Liu, and Hao Sun.
Physics-informed learning of governing equations from scarce data.
Nature communications, 12(1):6136, 2021.
- [20]
Abouzar Choubineh, Jie Chen, David A Wood, Frans Coenen, and Fei Ma.
Fourier neural operator for fluid flow in small-shape 2d simulated porous media dataset.
Algorithms, 16(1):24, 2023.
- [21]
Zhiwen Deng, Jing Wang, Hongsheng Liu, Hairun Xie, BoKai Li, Miao Zhang, Tingmeng Jia, Yi Zhang, Zidong Wang, and Bin Dong.
Prediction of transonic flow over supercritical airfoils using geometric-encoding and deep-learning strategies.
Physics of Fluids, 35(7), 2023.
- [22]
Meire Fortunato, Tobias Pfaff, Peter Wirnsberger, Alexander Pritzel, and Peter Battaglia.
Multiscale meshgraphnets.
Arxiv, 2022.
- [23]
Nicholas Geneva and Nicholas Zabaras.
Transformers for modeling physical systems.
Neural Networks, 146:272–289, 2022.
- [24]
Paul Ghanem, Ahmet Demirkaya, Tales Imbiriba, Alireza Ramezani, Zachary Danziger, and Deniz Erdogmus.
Learning physics informed neural odes with partial measurements.
In Proceedings of the AAAI Conference on Artificial Intelligence, volume 39, pages 16799–16807, 2025.
- [25]
R.A. Gingold and J.J. Monaghan.
Smoothed particle hydrodynamics: theory and application to non-spherical stars.
Monthly Notices of the Royal Astronomical Society, 181(3):375–389, 1977.
- [26]
John Guibas, Morteza Mardani, Zongyi Li, Andrew Tao, Anima Anandkumar, and Bryan Catanzaro.
Efficient token mixing for transformers via adaptive fourier neural operators.
In International conference on learning representations, 2021.
- [27]
Gaurav Gupta, Xiongye Xiao, and Paul Bogdan.
Multiwavelet-based operator learning for differential equations.
In A. Beygelzimer, Y. Dauphin, P. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, 2021.
- [28]
Kaiqiao Han, Yi Yang, Zijie Huang, Xuan Kan, Ying Guo, Yang Yang, Lifang He, Liang Zhan, Yizhou Sun, Wei Wang, et al.
Brainode: Dynamic brain signal analysis via graph-aided neural ordinary differential equations.
In 2024 IEEE EMBS International Conference on Biomedical and Health Informatics (BHI), pages 1–8. IEEE, 2024.
- [29]
XU HAN, Han Gao, Tobias Pfaff, Jian-Xun Wang, and Liping Liu.
Predicting physics in mesh-reduced space with temporal attention.
In International Conference on Learning Representations, 2022.
- [30]
Zhongkai Hao, Chang Su, Songming Liu, Julius Berner, Chengyang Ying, Hang Su, Anima Anandkumar, Jian Song, and Jun Zhu.
Dpot: Auto-regressive denoising operator transformer for large-scale pde pre-training.
ICML, 2024.
- [31]
Zhongkai Hao, Zhengyi Wang, Hang Su, Chengyang Ying, Yinpeng Dong, Songming Liu, Ze Cheng, Jian Song, and Jun Zhu.
Gnot: A general neural operator transformer for operator learning.
In International Conference on Machine Learning, pages 12556–12569. PMLR, 2023.
- [32]
Zhongkai Hao, Jiachen Yao, Chang Su, Hang Su, Ziao Wang, Fanzhi Lu, Zeyu Xia, Yichi Zhang, Songming Liu, Lu Lu, et al.
Pinnacle: A comprehensive benchmark of physics-informed neural networks for solving pdes.
Arxiv, 2023.
- [33]
Xiaodong He, Yinan Wang, and Juan Li.
Flow completion network: Inferring the fluid dynamics from incomplete flow information using graph neural networks.
Physics of Fluids, 34(8), 2022.
- [34]
Jacob Helwig, Xuan Zhang, Cong Fu, Jerry Kurtin, Stephan Wojtowytsch, and Shuiwang Ji.
Group equivariant fourier neural operators for partial differential equations.
ICML, 2023.
- [35]
Maximilian Herde, Bogdan Raonic, Tobias Rohner, Roger Käppeli, Roberto Molinaro, Emmanuel de Bézenac, and Siddhartha Mishra.
Poseidon: Efficient foundation models for pdes.
Advances in Neural Information Processing Systems, 37:72525–72624, 2024.
- [36]
Jonathan Ho, Ajay Jain, and Pieter Abbeel.
Denoising diffusion probabilistic models.
Advances in neural information processing systems, 33:6840–6851, 2020.
- [37]
Jiahe Huang, Guandao Yang, Zichen Wang, and Jeong Joon Park.
DiffusionPDE: Generative PDE-solving under partial observation.
In The Thirty-eighth Annual Conference on Neural Information Processing Systems, 2024.
- [38]
Xiang Huang, Zhanhong Ye, Hongsheng Liu, Shi Ji, Zidong Wang, Kang Yang, Yang Li, Min Wang, Haotian Chu, Fan Yu, et al.
Meta-auto-decoder for solving parametric partial differential equations.
Advances in Neural Information Processing Systems, 35:23426–23438, 2022.
- [39]
Zijie Huang, Yizhou Sun, and Wei Wang.
Learning continuous system dynamics from irregularly-sampled partial observations.
In Advances in Neural Information Processing Systems, 2020.
- [40]
Zijie Huang, Wanjia Zhao, Jingdong Gao, Ziniu Hu, Xiao Luo, Yadi Cao, Yuanzhou Chen, Yizhou Sun, and Wei Wang.
Physics-informed regularization for domain-agnostic dynamical system modeling, 2024.
- [41]
Jinsung Jeon, Hyundong Jin, Jonghyun Choi, Sanghyun Hong, Dongeun Lee, Kookjin Lee, and Noseong Park.
Pac-fno: Parallel-structured all-component fourier neural operators for recognizing low-quality images.
Arxiv, 2024.
- [42]
Zhongyi Jiang, Min Zhu, and Lu Lu.
Fourier-mionet: Fourier-enhanced multiple-input neural operators for multiphase modeling of geological carbon sequestration.
Reliability Engineering & System Safety, 251:110392, 2024.
- [43]
Pengzhan Jin, Shuai Meng, and Lu Lu.
Mionet: Learning multiple-input operators via tensor product.
SIAM Journal on Scientific Computing, 44(6):A3490–A3514, 2022.
- [44]
Xiaowei Jin, Shengze Cai, Hui Li, and George Em Karniadakis.
Nsfnets (navier-stokes flow nets): Physics-informed neural networks for the incompressible navier-stokes equations.
Journal of Computational Physics, 426:109951, 2021.
- [45]
G.A. Klaasen and W.C. Troy.
Stationary wave solutions of a system of reaction-diffusion equations derived from the fitzhugh–nagumo equations.
SIAM Journal on Applied Mathematics, 44(1):96–110, 1984.
- [46]
Georg Kohl, Li-Wei Chen, and Nils Thuerey.
Benchmarking autoregressive conditional diffusion models for turbulent flow simulation.
arXiv preprint arXiv:2309.01745, 2023.
- [47]
Katiana Kontolati, Somdatta Goswami, George Em Karniadakis, and Michael D Shields.
Learning in latent spaces improves the predictive accuracy of deep neural operators.
Arxiv, 2023.
- [48]
Remi Lam, Alvaro Sanchez-Gonzalez, et al.
Learning skillful medium-range global weather forecasting.
Science, 382(6677):1416–1421, 2023.
- [49]
Jae Yong Lee, SungWoong CHO, and Hyung Ju Hwang.
HyperdeepONet: learning operator with complex target function space using the limited resources via hypernetwork.
In The Eleventh International Conference on Learning Representations, 2023.
- [50]
Sangseung Lee and Donghyun You.
Data-driven prediction of unsteady flow over a circular cylinder using deep learning.
Journal of Fluid Mechanics, 879:217–254, 2019.
- [51]
Jiequan Li.
Two-stage fourth order: temporal-spatial coupling in computational fluid dynamics (cfd).
Advances in Aerodynamics, 1:1–36, 2019.
- [52]
Tianyi Li, Luca Biferale, Fabio Bonaccorso, Martino Andrea Scarpolini, and Michele Buzzicotti.
Synthetic lagrangian turbulence by generative diffusion models.
Nature Machine Intelligence, pages 1–11, 2024.
- [53]
Zijie Li and Amir Barati Farimani.
Graph neural network-accelerated lagrangian fluid simulation.
Computers & Graphics, 103:201–211, 2022.
- [54]
Zijie Li, Dule Shu, and Amir Barati Farimani.
Scalable transformer for pde surrogate modeling.
Advances in Neural Information Processing Systems, 36, 2024.
- [55]
Zongyi Li, Daniel Zhengyu Huang, Burigede Liu, and Anima Anandkumar.
Fourier neural operator with learned deformations for pdes on general geometries.
Journal of Machine Learning Research, 24(388):1–26, 2023.
- [56]
Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Andrew Stuart, Kaushik Bhattacharya, and Anima Anandkumar.
Multipole graph neural operator for parametric partial differential equations.
Advances in Neural Information Processing Systems, 33:6755–6766, 2020.
- [57]
Zongyi Li, Nikola Borislavov Kovachki, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew Stuart, Anima Anandkumar, et al.
Fourier neural operator for parametric partial differential equations.
In International Conference on Learning Representations.
- [58]
Zongyi Li, Nikola Borislavov Kovachki, et al.
Fourier neural operator for parametric partial differential equations.
In International Conference on Learning Representations, 2021.
- [59]
Zongyi Li, Nikola Borislavov Kovachki, et al.
Geometry-informed neural operator for large-scale 3d PDEs.
In Thirty-seventh Conference on Neural Information Processing Systems, 2023.
- [60]
Zongyi Li, Hongkai Zheng, Nikola Kovachki, David Jin, Haoxuan Chen, Burigede Liu, Kamyar Azizzadenesheli, and Anima Anandkumar.
Physics-informed neural operator for learning partial differential equations.
Arxiv, 2021.
- [61]
Guang Lin, Christian Moya, and Zecheng Zhang.
B-deeponet: An enhanced bayesian deeponet for solving noisy parametric pdes using accelerated replica exchange sgld.
Journal of Computational Physics, 473:111713, 2023.
- [62]
Mario Lino, Stathi Fotiadis, Anil A Bharath, and Chris D Cantwell.
Multi-scale rotation-equivariant graph neural networks for unsteady eulerian fluid dynamics.
Physics of Fluids, 34(8), 2022.
- [63]
Yaron Lipman, Ricky TQ Chen, Heli Ben-Hamu, Maximilian Nickel, and Matt Le.
Flow matching for generative modeling.
arXiv preprint arXiv:2210.02747, 2022.
- [64]
Jinxian Liu, Ye Chen, Bingbing Ni, Wei Ren, Zhenbo Yu, and Xiaoyang Huang.
Fast fluid simulation via dynamic multi-scale gridding.
In Proceedings of the AAAI Conference on Artificial Intelligence, volume 37, pages 1675–1682, 2023.
- [65]
Ning Liu, Siavash Jafarzadeh, and Yue Yu.
Domain agnostic fourier neural operators.
In Thirty-seventh Conference on Neural Information Processing Systems, 2023.
- [66]
Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and George Em Karniadakis.
Learning nonlinear operators via deeponet based on the universal approximation theorem of operators.
Nature machine intelligence, 3(3):218–229, 2021.
- [67]
Xiao Luo, Yiyang Gu, Huiyu Jiang, Hang Zhou, Jinsheng Huang, Wei Ju, Zhiping Xiao, Ming Zhang, and Yizhou Sun.
Pgode: Towards high-quality system dynamics modeling, 2024.
- [68]
Xiao Luo, Haixin Wang, Zijie Huang, Huiyu Jiang, Abhijeet Sadashiv Gangan, Song Jiang, and Yizhou Sun.
Care: Modeling interacting dynamics under temporal environmental variation.
In Thirty-seventh Conference on Neural Information Processing Systems, 2023.
- [69]
Xiao Luo, Jingyang Yuan, Zijie Huang, Huiyu Jiang, Yifang Qin, Wei Ju, Ming Zhang, and Yizhou Sun.
HOPE: High-order graph ODE for modeling interacting dynamics.
In Andreas Krause, Emma Brunskill, Kyunghyun Cho, Barbara Engelhardt, Sivan Sabato, and Jonathan Scarlett, editors, Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pages 23124–23139. PMLR, 23–29 Jul 2023.
- [70]
Yining Luo, Yingfa Chen, and Zhen Zhang.
Cfdbench: A comprehensive benchmark for machine learning methods in fluid dynamics.
Arxiv, 2023.
- [71]
Andrew J Majda, Andrea L Bertozzi, and A Ogawa.
Vorticity and incompressible flow. cambridge texts in applied mathematics.
Appl. Mech. Rev., 55(4):B77–B78, 2002.
- [72]
Nick McGreivy and Ammar Hakim.
Weak baselines and reporting biases lead to overoptimism in machine learning for fluid-related partial differential equations.
Nature Machine Intelligence, 6(10):1256–1269, 2024.
- [73]
R. Mikulevicius and B. L. Rozovskii.
Stochastic navier–stokes equations for turbulent flows.
SIAM Journal on Mathematical Analysis, 35(5):1250–1310, 2004.
- [74]
F. Moukalled, L. Mangani, and M. Darwish.
The Finite Volume Method in Computational Fluid Dynamics.
Springer, 1 edition, 2016.
- [75]
Bilal Mufti, Anindya Bhaduri, Sayan Ghosh, Liping Wang, and Dimitri N Mavris.
Shock wave prediction in transonic flow fields using domain-informed probabilistic deep learning.
Physics of Fluids, 36(1), 2024.
- [76]
Matthias Müller, David Charypar, and Markus Gross.
Particle-based fluid simulation for interactive applications.
In Proceedings of the 2003 ACM SIGGRAPH/Eurographics symposium on Computer animation, pages 154–159. Citeseer, 2003.
- [77]
Ruben Ohana, Michael McCabe, Lucas Meyer, Rudy Morel, Fruzsina Agocs, Miguel Beneitez, Marsha Berger, Blakesly Burkhart, Stuart Dalziel, Drummond Fielding, et al.
The well: a large-scale collection of diverse physics simulations for machine learning.
Advances in Neural Information Processing Systems, 37:44989–45037, 2024.
- [78]
Michael O’Connell, Guanya Shi, Xichen Shi, et al.
Neural-fly enables rapid learning for agile flight in strong winds.
Science Robotics, 7(66), 2022.
- [79]
Jaideep Pathak, Shashank Subramanian, et al.
Fourcastnet: A global data-driven high-resolution weather model using adaptive fourier neural operators.
Arxiv, 2022.
- [80]
Tobias Pfaff, Meire Fortunato, Alvaro Sanchez-Gonzalez, and Peter Battaglia.
Learning mesh-based simulation with graph networks.
In International Conference on Learning Representations, 2021.
- [81]
Yuan Qiu, Nolan Bridges, and Peng Chen.
Derivative-enhanced deep operator network.
Advances in Neural Information Processing Systems, 37:20945–20981, 2024.
- [82]
Md Ashiqur Rahman, Zachary E Ross, and Kamyar Azizzadenesheli.
U-NO: U-shaped neural operators.
Transactions on Machine Learning Research, 2023.
- [83]
Maziar Raissi, Paris Perdikaris, and George E Karniadakis.
Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations.
Journal of Computational physics, 378:686–707, 2019.
- [84]
Ellery Rajagopal, Anantha NS Babu, et al.
Evaluation of deep neural operator models toward ocean forecasting.
In OCEANS, pages 1–9. IEEE, 2023.
- [85]
Chengping Rao, Pu Ren, Qi Wang, Oral Buyukozturk, Hao Sun, and Yang Liu.
Encoding physics to learn reaction–diffusion processes.
Nature Machine Intelligence, 5(7):765–779, 2023.
- [86]
Bogdan Raonic, Roberto Molinaro, Tobias Rohner, Siddhartha Mishra, and Emmanuel de Bezenac.
Convolutional neural operators.
In ICLR 2023 Workshop on Physics for Machine Learning, 2023.
- [87]
Carlos JG Rojas, Andreas Dengel, and Mateus Dias Ribeiro.
Reduced-order model for fluid flows via neural ordinary differential equations.
arXiv preprint arXiv:2102.02248, 2021.
- [88]
Alexander Rudikov, Vladimir Fanaskov, et al.
Neural operators meet conjugate gradients: The FCG-NO method for efficient PDE solving.
In Forty-first International Conference on Machine Learning, 2024.
- [89]
Cristopher Salvi, Maud Lemercier, and Andris Gerasimovics.
Neural stochastic pdes: Resolution-invariant learning of continuous spatiotemporal dynamics, 2022.
- [90]
Alvaro Sanchez-Gonzalez, Jonathan Godwin, et al.
Learning to simulate complex physics with graph networks.
In International conference on machine learning, pages 8459–8468. PMLR, 2020.
- [91]
Jacob Seidman, Georgios Kissas, Paris Perdikaris, and George J Pappas.
Nomad: Nonlinear manifold decoders for operator learning.
Advances in Neural Information Processing Systems, 35:5601–5613, 2022.
- [92]
Yidi Shao, Chen Change Loy, and Bo Dai.
Transformer with implicit edges for particle-based physics simulation.
In European Conference on Computer Vision, pages 549–564. Springer, 2022.
- [93]
Hong Shen Wong, Wei Xuan Chan, Bing Huan Li, and Choon Hwai Yap.
Multiple case physics-informed neural network for biomedical tube flows.
Arxiv, 2023.
- [94]
Aleksei Sholokhov, Yuying Liu, Hassan Mansour, and Saleh Nabi.
Physics-informed neural ode (pinode): embedding physics into models using collocation points.
Scientific Reports, 13(1):10166, 2023.
- [95]
Aliaksandra Shysheya, Cristiana Diaconu, Federico Bergamin, Paris Perdikaris, José Miguel Hernández-Lobato, Richard E. Turner, and Emile Mathieu.
On conditional diffusion models for PDE simulations.
In The Thirty-eighth Annual Conference on Neural Information Processing Systems, 2024.
- [96]
Jiaming Song, Chenlin Meng, and Stefano Ermon.
Denoising diffusion implicit models.
arXiv preprint arXiv:2010.02502, 2020.
- [97]
Jos Stam.
Stable fluids.
In Proceedings of the 26th Annual Conference on Computer Graphics and Interactive Techniques, SIGGRAPH ’99, page 121–128, USA, 1999. ACM Press/Addison-Wesley Publishing Co.
- [98]
Shankar Subramaniam.
Lagrangian–eulerian methods for multiphase flows.
Progress in Energy and Combustion Science, 39(2-3):215–245, 2013.
- [99]
Fang Sun, Zijie Huang, Haixin Wang, Yadi Cao, Xiao Luo, Wei Wang, and Yizhou Sun.
Graph fourier neural odes: Bridging spatial and temporal multiscales in molecular dynamics.
arXiv preprint arXiv:2411.01600, 2024.
- [100]
Makoto Takamoto, Timothy Praditia, Raphael Leiteritz, Daniel MacKinlay, Francesco Alesiani, Dirk Pflüger, and Mathias Niepert.
Pdebench: An extensive benchmark for scientific machine learning.
Advances in Neural Information Processing Systems, 35:1596–1611, 2022.
- [101]
Maryam Toloubidokhti, Yubo Ye, Ryan Missel, Xiajun Jiang, Nilesh Kumar, Ruby Shrestha, and Linwei Wang.
Dats: Difficulty-aware task sampler for meta-learning physics-informed neural networks.
In The Twelfth International Conference on Learning Representations, 2023.
- [102]
F Arend Torres, Marcello Massimo Negri, Marco Inversi, Jonathan Aellen, and Volker Roth.
Lagrangian flow networks for conservation laws.
Arxiv, 2023.
- [103]
A. Tran, A. Mathews, et al.
Factorized fourier neural operators.
In ICLR, 2023.
- [104]
Mario Lino Valencia, Tobias Pfaff, and Nils Thuerey.
Learning distributions of complex fluid simulations with diffusion graph networks.
In The Thirteenth International Conference on Learning Representations, 2025.
- [105]
Simone Venturi and Tiernan Casey.
Svd perspectives for augmenting deeponet flexibility and interpretability.
Computer Methods in Applied Mechanics and Engineering, 403:115718, 2023.
- [106]
Paulien HM Voorter, Walter H Backes, et al.
Improving microstructural integrity, interstitial fluid, and blood microcirculation images from multi-b-value diffusion mri using physics-informed neural networks in cerebrovascular disease.
Magnetic Resonance in Medicine, 2023.
- [107]
Haixin Wang, Yadi Cao, Zijie Huang, Yuxuan Liu, Peiyan Hu, Xiao Luo, Zezheng Song, Wanjia Zhao, Jilin Liu, Jinan Sun, et al.
Recent advances on machine learning for computational fluid dynamics: A survey.
arXiv preprint arXiv:2408.12171, 2024.
- [108]
Haixin Wang, LI Jiaxin, Anubhav Dwivedi, Kentaro Hara, and Tailin Wu.
BENO: Boundary-embedded neural operators for elliptic PDEs.
In The Twelfth International Conference on Learning Representations, 2024.
- [109]
Rui Wang, Karthik Kashinath, Mustafa Mustafa, Adrian Albert, and Rose Yu.
Towards physics-informed deep learning for turbulent flow prediction.
In KDD, pages 1457–1466, 2020.
- [110]
Rui Wang, Robin Walters, and Rose Yu.
Incorporating symmetry into deep dynamics models for improved generalization.
Arxiv 2020.
- [111]
Sifan Wang, Hanwen Wang, and Paris Perdikaris.
Learning the solution operator of parametric partial differential equations with physics-informed deeponets.
Science advances, 7(40), 2021.
- [112]
Xiaoda Wang, Yuji Zhao, Kaiqiao Han, Xiao Luo, Sanne van Rooij, Jennifer Stevens, Lifang He, Liang Zhan, Yizhou Sun, Wei Wang, et al.
Conditional neural ode for longitudinal parkinson’s disease progression forecasting.
In Abstract in the Organization for Human Brain Mapping Annual Meeting, 2025.
- [113]
Yifan Wang and Linlin Zhong.
Nas-pinn: neural architecture search-guided physics-informed neural network for solving pdes.
Journal of Computational Physics, 496:112603, 2024.
- [114]
Wang, Rui and Walters, Robin and Yu, Rose.
Approximately equivariant networks for imperfectly symmetric dynamics.
In International Conference on Machine Learning, 2022.
- [115]
Ping Wei, Menghan Liu, Jianhuan Cen, Ziyang Zhou, Liao Chen, and Qingsong Zou.
Pdenneval: A comprehensive evaluation of neural network methods for solving pdes.
In Proceedings of the Thirty-Third International Joint Conference on Artificial Intelligence, pages 5181–5189, 2024.
- [116]
Gege Wen, Zongyi Li, Kamyar Azizzadenesheli, Anima Anandkumar, and Sally M Benson.
U-fno—an enhanced fourier neural operator-based deep-learning model for multiphase flow.
Advances in Water Resources, 163:104180, 2022.
- [117]
Christopher R Wentland, Karthik Duraisamy, and Cheng Huang.
Scalable projection-based reduced-order models for large multiscale fluid systems.
AIAA Journal, 61(10):4499–4523, 2023.
- [118]
Henning Wessels, Christian Weißenfels, and Peter Wriggers.
The neural particle method–an updated lagrangian physics informed neural network for computational fluid dynamics.
Computer Methods in Applied Mechanics and Engineering, 368:113127, 2020.
- [119]
Haixu Wu, Tengge Hu, Huakun Luo, Jianmin Wang, and Mingsheng Long.
Solving high-dimensional pdes with latent spectral models.
In International Conference on Machine Learning, 2023.
- [120]
Haixu Wu, Huakun Luo, et al.
Transolver: A fast transformer solver for PDEs on general geometries.
In Forty-first International Conference on Machine Learning, 2024.
- [121]
Hao Wu, Changhu Wang, Fan Xu, Jinbao Xue, Chong Chen, Xian-Sheng Hua, and Xiao Luo.
Pure: Prompt evolution with graph ode for out-of-distribution fluid dynamics modeling.
In A. Globerson, L. Mackey, D. Belgrave, A. Fan, U. Paquet, J. Tomczak, and C. Zhang, editors, Advances in Neural Information Processing Systems, volume 37, pages 104965–104994. Curran Associates, Inc., 2024.
- [122]
Hao Wu, Huiyuan Wang, Kun Wang, Weiyan Wang, ChanganYe, Yangyu Tao, Chong Chen, Xian-Sheng Hua, and Xiao Luo.
Prometheus: Out-of-distribution fluid dynamics modeling with disentangled graph ODE.
In Forty-first International Conference on Machine Learning, 2024.
- [123]
Hao Wu, Fan Xu, Yifan Duan, Ziwei Niu, Weiyan Wang, Gaofeng Lu, Kun Wang, Yuxuan Liang, and Yang Wang.
Spatio-temporal fluid dynamics modeling via physical-awareness and parameter diffusion guidance.
arXiv preprint arXiv:2403.13850, 2024.
- [124]
Tailin Wu, Takashi Maruyama, Qingqing Zhao, Gordon Wetzstein, and Jure Leskovec.
Learning controllable adaptive simulation for multi-resolution physics.
In The Eleventh International Conference on Learning Representations, 2023.
- [125]
Xiongye Xiao, Defu Cao, Ruochen Yang, Gaurav Gupta, Gengshuo Liu, Chenzhong Yin, Radu Balan, and Paul Bogdan.
Coupled multiwavelet operator learning for coupled differential equations.
In The Eleventh International Conference on Learning Representations, 2022.
- [126]
Zipeng Xiao, Zhongkai Hao, Bokai Lin, Zhijie Deng, and Hang Su.
Improved operator learning by orthogonal attention.
In Ruslan Salakhutdinov, Zico Kolter, Katherine Heller, Adrian Weller, Nuria Oliver, Jonathan Scarlett, and Felix Berkenkamp, editors, Proceedings of the 41st International Conference on Machine Learning, volume 235 of Proceedings of Machine Learning Research, pages 54288–54299. PMLR, 21–27 Jul 2024.
- [127]
Zipeng Xiao, Siqi Kou, Hao Zhongkai, Bokai Lin, and Zhijie Deng.
Amortized fourier neural operators.
Advances in Neural Information Processing Systems, 37:115001–115020, 2024.
- [128]
Wei Xiong, Xiaomeng Huang, Ziyang Zhang, Ruixuan Deng, Pei Sun, and Yang Tian.
Koopman neural operator as a mesh-free solver of non-linear partial differential equations.
Arxiv, 2023.
- [129]
Liu Yang, Siting Liu, Tingwei Meng, and Stanley J. Osher.
In-context operator learning with data prompts for differential equation problems.
Proceedings of the National Academy of Sciences, 120(39), 2023.
- [130]
Minglang Yin, Ehsan Ban, Bruno V Rego, Enrui Zhang, Cristina Cavinato, Jay D Humphrey, and George Em Karniadakis.
Simulating progressive intramural damage leading to aortic dissection using deeponet: an operator–regression neural network.
Journal of the Royal Society Interface, 19(187):20210670, 2022.
- [131]
Yuan Yin, Matthieu Kirchmeyer, Jean-Yves Franceschi, Alain Rakotomamonjy, and patrick gallinari.
Continuous PDE dynamics forecasting with implicit neural representations.
In The Eleventh International Conference on Learning Representations, 2023.
- [132]
Jingyang Yuan, Gongbo Sun, Zhiping Xiao, Hang Zhou, Xiao Luo, Junyu Luo, Yusheng Zhao, Wei Ju, and Ming Zhang.
EGODE: An event-attended graph ODE framework for modeling rigid dynamics.
In The Thirty-eighth Annual Conference on Neural Information Processing Systems, 2024.
- [133]
Qianru Zhang, Haixin Wang, Cheng Long, Liangcai Su, Xingwei He, Jianlong Chang, Tailin Wu, Hongzhi Yin, Siu-Ming Yiu, Qi Tian, et al.
A survey of generative techniques for spatial-temporal data mining.
arXiv preprint arXiv:2405.09592, 2024.
- [134]
Qianru Zhang, Peng Yang, Honggang Wen, Xinzhu Li, Haixin Wang, Fang Sun, Zezheng Song, Zhichen Lai, Rui Ma, Ruihua Han, et al.
Beyond the time domain: Recent advances on frequency transforms in time series analysis.
arXiv preprint arXiv:2504.07099, 2025.
- [135]
Rui Zhang, Qi Meng, Rongchan Zhu, Yue Wang, Wenlei Shi, Shihua Zhang, Zhi-Ming Ma, and Tie-Yan Liu.
Monte carlo neural operator for learning pdes via probabilistic representation.
Arxiv, 2023.
- [136]
Tianhan Zhang, Yuxiao Yi, Yifan Xu, Zhi X Chen, Yaoyu Zhang, E Weinan, and Zhi-Qin John Xu.
A multi-scale sampling method for accurate and robust deep neural network to predict combustion chemical kinetics.
Combustion and Flame, 245:112319, 2022.
- [137]
Qingqi Zhao, Xiaoxue Han, Ruichang Guo, and Cheng Chen.
A computationally efficient hybrid neural network architecture for porous media: Integrating cnns and gnns for improved permeability prediction.
Arxiv, 2023.
- [138]
Zhiyuan Zhao, Xueying Ding, and B. Aditya Prakash.
PINNsformer: A transformer-based framework for physics-informed neural networks.
In The Twelfth International Conference on Learning Representations, 2024.
- [139]
Su Zheng, Zhengqi Gao, Fan-Keng Sun, Duane Boning, Bei Yu, and Martin D Wong.
Improving neural ode training with temporal adaptive batch normalization.
Advances in Neural Information Processing Systems, 37:95875–95895, 2024.

## Appendix

## Table of Content

- •
Section [A](https://arxiv.org/html/2505.20349v1#A1)  Symbol Notation
- •
Section [B](https://arxiv.org/html/2505.20349v1#A2)  Decompositions Corresponding to Existing Works.
- •
Section [C](https://arxiv.org/html/2505.20349v1#A3)  Supplementary Results on Comparison with Traditional Solvers
- •
Section [D](https://arxiv.org/html/2505.20349v1#A4)  Supplementary Results on Different Modules
- •
Section [E](https://arxiv.org/html/2505.20349v1#A5)  Supplementary Results on Generalization Study
- •
Section [F](https://arxiv.org/html/2505.20349v1#A6)  More Dataset Details
- •
Section [G](https://arxiv.org/html/2505.20349v1#A7)  More Implementation Details
- •
Section [H](https://arxiv.org/html/2505.20349v1#A8)  Visualization Study

## Appendix A Notation

Table [4](https://arxiv.org/html/2505.20349v1#A1.T4) summarizes the key mathematical symbols and their meanings used throughout our benchmark framework and avoid symbol ambiguity across different methods, covering core components such as spatial representation, temporal modeling, and loss functions.

**Table 4: Unified notation across spatial, temporal, and loss components in our paper.**
| Symbol | Meaning |
| --- | --- |
| $\Omega$ | Spatial domain with dimension $d$, $\Omega\subset\mathbb{R}^{d}$ |
| $T$ | Temporal horizon, $t\in[0,T]$ |
| $u(\bm{x},t)$ | Field variable at location $\bm{x}$ and time $t$ with $c$ channels |
| $\{\bm{x}_{i}\}_{i=1}^{N}$ | Spatial locations of $N$ samples, $\bm{x}_{i}\in\Omega$ |
| $M$ | Number of temporal snapshots |
| $\mathcal{S}_{\theta},\mathcal{T}_{\theta}$ | Spatial and temporal encoder parameterized by $\theta$ |
| $z,h$ | Feature representation in spatial and temporal domain |
| $\mathcal{F},\mathcal{F}^{-1}$ | Fourier and inverse Fourier transforms |
| $\hat{u}(\mathbf{k},t)$ | Spectral coefficient at frequency $\mathbf{k}$ |
| $\phi_{\theta}(\mathbf{k})$ | Learnable spectral filter of frequency index $\mathbf{k}$ |
| $\odot$ | Elementwise (Hadamard) product |
| $\mathbf{W}_{Q},\mathbf{W}_{K},\mathbf{W}_{V}$ | Query, key and value projection matrices, $\in\mathbb{R}^{c\times d_{k}}$ |
| $d_{k},d_{k},d_{v}$ | Query/key and Value dimension in self-attention |
| $\mathcal{H}$ | Set of stencil offsets for convolution |
| $c_{\mathrm{in}},c_{\mathrm{out}}$ | Number of input and output channels to convolution |
| $\mathbf{W}_{h}$ | Convolution kernel at offset $h$, $\in\mathbb{R}^{c_{\mathrm{out}}\times c_{\mathrm{in}}}$ |
| $\tilde{A},\tilde{D}$ | Adjacency with self-loops, $\tilde{A}=A+I$, and degree matrix of $\tilde{A}$ |
| $K$ | Number of retained modes in ROM, $K\ll N$ |
| $C(\bm{x},\bm{x}^{\prime})$ | Covariance kernel between spatial locations |
| $\phi_{k}(\bm{x})$ | $k$-th POD basis function |
| $\lambda_{k}$ | Eigenvalue associated with $\phi_{k}$ |
| $\alpha_{k}(t)$ | POD projection coefficient at time $t$ |
| $\langle f,g\rangle$ | $L^{2}$ inner product over $\Omega$ |
| $\mathbf{d}^{j}$ | Predicted increment for bundling step $j$ |
| $h(t)$ | Latent state in Neural ODE |
| $\zeta(\cdot),\xi(\cdot)$ | Encoder/decoder for latent ODE state |
| $\varepsilon_{\theta}$ | Network estimating noise in diffusion models |
| $u_{\tau},\bar{\alpha}_{\tau}$ | Noisy version of field and noise schedule coefficient in diffusion |
| $v_{\theta},\psi_{\tau}(u)$ | Predicted velocity field and interpolated affine flow |
| $\mathcal{R}[u_{\theta}](x,t)$ | Residual of the PDE operator |
| $g(x,t)$ | Source term in PDE |

## Appendix B Modular Decomposition

**Table 5: We propose that the design of data-driven fluid simulation models can be analogized to a Taylor expansion. Design $=$ “Spatial representation” $+$ “Temporal representation” $+$ “Loss function” $+$ Additional Technique, where the additional item represents the key innovations introduced in each method.**
| # | Method | Publication | Spatial Repre. | Temporal Repre. | Loss | Additional Technique |
| --- | --- | --- | --- | --- | --- | --- |
| 1 | FNO [57] | ICLR 2021 | Fourier | Autoregression | Variable | Fourier integral operator |
| 2 | AFNO [26] | ICLR 2022 | Fourier | Autoregression | Variable | Mixing of tokens |
| 3 | Geo-FNO [55] | JMLR | Fourier | Autoregression | Variable | Data deformation |
| 4 | PINO [60] | ACM Data Sci. | Fourier | Autoregression | Variable & PDE loss | Efficient derivative compute |
| 5 | U-NO [82] | TMLR | Fourier | Autoregression | Variable | U-net architecture |
| 6 | F-FNO [103] | ICLR 2023 | Fourier | Autoregression | Variable | Factorize Fourier transform |
| 7 | CFNO [7] | ICLR 2023 | Fourier | Autoregression | Variable | Clifford algebras |
| 8 | CMWNO [125] | ICLR 2023 | Fourier | Autoregression | Variable | Multiwavelet decomposition |
| 9 | $G$-FNO [34] | ICML 2023 | Fourier | Autoregression | Variable | Group equivariant layers |
| 10 | GINO [59] | NeurIPS 2023 | Graph/Fourier | Autoregression | Variable | Discretization convergence |
| 11 | DAFNO [65] | NeurIPS 2023 | Fourier | Autoregression | Variable | Smoothed characteristic func |
| 12 | PAC-FNO [41] | ICLR 2024 | Fourier | Autoregression | Variable | Parallel structure |
| 13 | AM-FNO [127] | NeurIPS 2024 | Fourier | Autoregression | Variable | Amortized parameterization |
| 14 | DE-DON [81] | NeurIPS 2024 | MLP/Fourier | Autoregression | Variable & Derivative | Dimension reduction |
| 15 | MCNP [135] | TPAMI | Fourier | Autoregression | Neural monte carlo | Mathematical expectation |
| 16 | ST-FNO [13] | ICLR 2025 | Fourier | Time-depth Conv | Variable | Spectral fine-tuning |
| 17 | scOT [35] | NeurIPS 2024 | Self-atten | Temporal Bund | Variable | Large-scale pretraining |
| 18 | DiffPDE [37] | NeurIPS 2024 | Convolution | Autoregression | Noise | Guided Diffusion |
| 19 | Shysheya et al. [95] | NeurIPS 2024 | Convolution | Autoregression | Noise | Flexible pre- and post-training |
| 20 | Valencia et al. [104] | ICLR 2025 | Graph | Autoregression | Noise | Latent multi-scale GNN diffusion |
| 21 | Bastek et al. [2] | ICLR 2025 | Convolution | Autoregression | Noise | Physics informed diffusion |
| 22 | TA-BN [139] | NeurlPS 2024 | Convolution | Neural ODE | Variable | Temporal<br>adaptive BatchNorm |
| 23 | PINODE [94] | Sci. Rep. | MLP | Neural ODE | variable | Collocation Points |
| 24 | Ghanem et al. [24] | AAAI 2025 | MLP | Neural ODE | Variable | Physical-informed loss |
| 25 | PGODE[67] | ICML2024 | Self-atten | Neural ODE | Variable | Environment Parameter |
| 26 | Treat[40] | NeurlPS 2024 | Self-atten | Neural ODE | Vairable | Physical-informed regularization |
| 27 | Pure[121] | NeurlPS 2024 | Self-atten | Neural ODE | Vairable | Multi-view context information |
| 28 | EGODE[132] | NeurlPS 2024 | Self-atten | Neural ODE | Vairable | Event-attended information |
| 29 | Prometheus[122] | ICML 2024 | Self-atten | Neural ODE | Vairable | disentangled representations |
| 30 | HOPE[69] | ICML 2023 | Self-atten | Neural ODE | Vairable | High order information |
| 31 | LGODE[39] | NeurlPS 2020 | Self-atten | Neural ODE | Vairable | Temporal Self-attention |
| 32 | MPPDE[8] | ICLR 2022 | Graph | Temporal Bund | Variable | Model training stability |
| 33 | BENO[108] | ICLR 2024 | Graph | Autoregression | Variable | Complex geometry embedding |
| 34 | UPT[1] | NeurIPS 2024 | Self-atten | Autoregression | Variable | Scalability across discretization |
| 35 | Transolver[120] | ICML 2024 | Self-atten | Autoregression | Variable | Physics-Attention |
| 36 | OMO[126] | ICML 2024 | Self-atten | Autoregression | Variable | Orthogonal Attention |
| 37 | LSM[119] | ICML 2023 | Self-atten | Autoregression | Variable | Latent Propagation |
| 38 | GNOT[31] | ICML 2023 | Self-atten | Autoregression | Variable | Linear Attention |
| 39 | Galerkin-Trans[12] | NeurIPS 2021 | Self-atten | Autoregression | Variable | Petrov-Galerkin projection |
| 40 | MWT[27] | NeurIPS 2021 | Convolution | Autoregression | Variable | Multi-wavelet filters |
| 41 | CARE [68] | NeurIPS 2023 | Graph | ODE | Variable | Context acquirement |
| 42 | HAMLET [9] | ICML 2024 | Graph$+$SA | Next-step | Variable | Modular input encoders |
| 43 | LAMP [124] | ICLR 2023 | Graph | Next-step | Variable | GNN-based actor-critic for policy |
| 44 | DINo [131] | ICLR 2023 | Implicit | Neural ODE | Variable | Extrapolates at arbitrary loc |
| 45 | Zhao et al. [137] | Adv Water Resour | Graph Conv | Next-step | Variable | GNN Grad-CAM |
| 46 | FCN [33] | PHYS FLUIDS. | Graph | Next-step | Variable | Vortex force contribution |
| 47 | GNODE [4] | ICLR 2023 | Graph | Neural ODE | Variable | Encode the constraints explicitly |
| 48 | MAgNet [5] | NeurIPS 2022 | Graph | Next-step | Variable | Implicit neural representation |
| 49 | TIE [92] | ECCV 2022 | Self-atten | Next-step | Variable | Implicit Edges |
| 50 | Han et al. [29] | ICLR 2022 | Self-atten | Self-atten | Variable | Encoder-decoder structure |
| 51 | MGN [80] | ICLR 2021 | Graph | Next-step | Variable | Mesh graph representation |
| 52 | GNS [90] | ICML 2022 | Graph | Nest-step | Variable | Message passing layers |
| 53 | RSteer [114] | ICLR 2023 | Graph | Next-step | Variable | Approximately equivariant<br>networks |
| 54 | EquNet [110] | ICLR 2021 | Conv | Temporal Bund | Variable | Incorporate symmetry |
| 55 | TF-Net [109] | KDD 2020 | Conv | Next-step | Variable | Marry two simulation |
| 56 | DPUF [50] | J. Fluid Mech. | Conv | Next-step | Physical | Conservation of mass, momentum |
| 57 | Wessels et al. [118] | CoM Appl M | MLP | Next-step | Variable | Updated Lagrangian |
| 58 | FGN [53] | Comput Graph | Graph | Next-step | Variable | Node-focused and edge-focused |
| 59 | MCC [64] | AAAI 2023 | Conv | Next-step | Variable | Dynamic multi-scale gridding |
| 60 | LFlows [102] | ICLR 2024 | MLP | Temporal Bund | Flow | Satisfy the continuity equation |
| 61 | Li et al. [52] | Nat. Mach. Intell. | Conv | Next-step | Noise | Two different classes of DM |
| 62 | DeepONet [66] | Nat. Mach. Intell. | MLP | Next-step | Variable | Two sub-nets |
| 63 | PI-DeepONet [111] | Sci. Advection | MLP | Next-step | Physics | Soft penalty constraints for law |
| 64 | MIONet [43] | Siam J Sci Comput | MLP | Next-step | Physics | A low-rank approximation |
| 65 | Geneva et al. [23] | Neural Networks | Self-atten | Self-atten | Variable | Koopman dynamics |
| 66 | NOMAD [91] | NeurIPS 2022 | MLP | Next-step | Variable | Nonlinear decoder map |
| 67 | HyperDON [49] | ICLR 2023 | MLP | Next-step | Variable | Use a hypernetwork |
| 68 | B-DON [61] | J. Comput. Phy. | MLP | Next-step | Variable | Bayesian network |
| 69 | SVD-DON [105] | CoM Appl M | MLP | Next-step | Variable | Orthogonal decomposition |
| 70 | L-DON [47] | Arxiv 2023 | MLP | Next-step | Variable | Low-dimensional latent space |
| 71 | GNOT [31] | ICML 2023 | Self-atten | Next-step | Variable | Heterogeneous normalized atten |
| 72 | CNO [86] | NeurIPS 2023 | Conv | Autoregression | Variable | Adaptations for convolution |
| 73 | FactFormer [54] | NeurIPS 2023 | Self-atten | Self-atten | Variable | Axial factorized kernel |
| 74 | LNO [11] | Nat. Mach. Intell. | Laplace | Autoregression | Variable | Pole-residue relationship |
| 75 | KNO [128] | APL ML | Koopman | Autoregression | Variable | Approximate the Koopman operator |
| 76 | ICON [129] | PNAS | Self-atten | Self-atten | variable | In-context learning |
| 77 | Transolver [120] | ICML 2024 | Self-atten | Self-atten | Variable | Intrinsic physical states |
| 78 | FCG-NO [88] | ICML 2024 | Frequency | Autoregression | Variable | Flexible conjugate gradient |
| 79 | PINN [83] | J. Comput. Phy. | MLP | Next-step | Residual | - |
| 80 | PINN-SR [19] | Nat. Commun. | MLP | Next-step | Residual | Regularization loss |
| 81 | NSFnet [44] | J. Comput. Phys. | MLP | Next-step | Residual | For N-S Eq. |
| 82 | MAD [38] | NeurIPS 2022 | MLP | Next-step | Variable | Meta-learning |
| 83 | DATS [101] | ICLR 2024 | MLP | Next-step | Residual | Difficulty-aware task sampler |
| 84 | PINNsFormer [138] | ICLR 2024 | Self-atten | Next-step | Residual | Pseudo sequences |
| 85 | PeRCNN [85] | Nat. Mach. Intell. | Conv | Next-step | Variable | Forcibly encodes physics structure |

## Appendix C Supplementary Results on Comparison with Traditional Solvers

We present the rollout error of Burgers’ Equation in Figure [8](https://arxiv.org/html/2505.20349v1#A3.F8) and Advection Equation in Figure [8](https://arxiv.org/html/2505.20349v1#A3.F8). The runtime comparison is provided in Figure [10](https://arxiv.org/html/2505.20349v1#A3.F10). FNO not only achieves better predictive accuracy but also faster speed.

Figure: Figure 8: We compare FNO against pseudo-spectral solver operating at lower resolutions on Burgers’ equation across five kinematics viscosity. Solver$x$ implies pseudo-spectral solver operating at $x\times x$ resolution.
Refer to caption: extracted/6477461/fig/Solver/burgers/burgers1.png

Figure: Figure 9: We compare FNO against pseudo-spectral solver operating at lower resolutions on Advection equation across five avection speed types. Solver$x$ implies pseudo-spectral solver operating at $x\times x$ resolution.
Refer to caption: extracted/6477461/fig/Solver/advection/advection0.png

Figure: Figure 10: We compare the runtime of the FNO against pseudo-spectral solver operating at lower resolutions. Solver$x$ implies pseudo-spectral solver operating at $x\times x$ resolution. Data generated at $256\times 256$ are considered ground truth. Note that to ensure numerical stability, adaptive time-stepping governed by CFL condition is adopted even for solvers at low-resolution.
Refer to caption: x8.png

## Appendix D Supplementary Results on Different Modules

In this section, we provide additional experimental results on evaluation of different modular designs on prediction performance and computational efficiency across Stochastic N-S. The result is shown in Table [6](https://arxiv.org/html/2505.20349v1#A4.T6).

**Table 6: Evaluation of different modular designs on prediction performance and computational efficiency across Stochastic N-S.**
|  |  | Stochastic N-S |  |  |
| --- | --- | --- | --- | --- |
| Model | Param | RMSE$\downarrow$ | nRMSE$\downarrow$ | fRMSE$\downarrow$ |
| $\mathcal{X}$ + Next-step + Physical variable |  |  |  |  |
| $\mathcal{X}$ = Graph representaion | 3.13M | 0.0917 | 0.2256 | 0.0067 |
| $\mathcal{X}$ = Conv representation | 7.63M | 0.0438 | 0.0325 | 0.0042 |
| $\mathcal{X}$ = Fourier representation | 72.0M | 0.0398 | 0.0463 | 0.0050 |
| $\mathcal{X}$ = Latent representation | 89.0M | 0.0560 | 0.0392 | 0.0072 |
| $\mathcal{X}$ = Self-attention representation | 36.8M | 0.0500 | 0.0430 | 0.0047 |
| Self-attention representation + $\mathcal{Y}$ + Physical variable |  |  |  |  |
| $\mathcal{Y}$ = Self-attention | 36.4M | 0.1094 | 0.0740 | 0.0104 |
| $\mathcal{Y}$ = Temporal Bundling | 36.3M | 0.0998 | 0.0675 | 0.0089 |
| $\mathcal{Y}$ = Neural ODE | 37.1M | 0.0913 | 0.0617 | 0.0083 |

## Appendix E Supplementary Results on Generalization Study

#### Scaling Analysis.

We systematically assessed the impact of model capacity on predictive accuracy across five representative architectures by defining four scaling levels for each:

- •
Convolutional modules: number of layers $\in\{2,3,4,5\}$
- •
Fourier-based learning: hidden dimension $\in\{64,128,256,512\}$
- •
Latent learning: attention dimension $\in\{64,128,256,512\}$
- •
Self-attention learning: hidden size $\in\{96,192,384,768\}$
- •
Graph learning: hidden feature size $\in\{64,128,256,512\}$

As shown in Figure [11](https://arxiv.org/html/2505.20349v1#A5.F11), RMSE decreases monotonically from Scale 1 to Scale 4 for all five methods. The self-attention model achieves the largest absolute reduction (over $0.08$), demonstrating its superior ability to exploit increased representational capacity. The Fourier- and graph-based models follow, with improvements of approximately $0.06$ and $0.05$, respectively, indicating diminishing returns at higher dimensions. The convolutional baseline exhibits the smallest gain (about $0.03$), suggesting that depth alone yields limited expressivity for this task, while the latent flow model attains an intermediate reduction (about $0.04$). These findings underscore that transformer-style self-attention mechanisms benefit most from scaling, making them particularly well-suited for high-fidelity modeling in our PDE-driven benchmarks.

Figure: Figure 11: RMSE versus model capacity for five architectures at four scaling levels.
Refer to caption: x10.png

## Appendix F More Datasets Details

### F.1 Datasets for Modular Comparison

Compressible N-S Equations. The compressible N-S (CNS) equations govern the dynamics of fluid flows with variable density and internal energy. The compressible model accounts for variations in density $\rho$ and pressure $p$ due to thermodynamic effects, and captures phenomena such as shock waves, sound propagation, and high-speed turbulent mixing.

We utilize the 2D data from PDEBench [100]. The general form of the CNS equations in conservative form is given by:

$$ $\displaystyle\text{ Mass conservation (Continuity equation):}\quad\partial_{t} \rho+\nabla\cdot(\rho\textbf{v})=0,$ (15a) $\displaystyle\text{ Momentum conservation:}\quad\partial_{t}(\rho\textbf{v})+ \nabla\cdot(\rho\textbf{v}\otimes\textbf{v})=-\nabla p+\nabla\cdot\sigma,$ (15b) $\displaystyle\text{ Total energy conservation:}\quad\partial_{t}E+\nabla\cdot[ (E+p)\textbf{v}]=\nabla\cdot(\sigma\cdot\textbf{v})+\nabla\cdot(k\nabla T),$ (15c) $$

where $\rho$ is the mass density, $\textbf{v}\in\mathbb{R}^{d}$ is the velocity field, $p$ is the pressure, $E=\epsilon+\frac{1}{2}\rho\|\textbf{v}\|^{2}$ is the total energy per unit volume, $\epsilon$ is the internal energy, $T$ is the temperature, $k$ is the thermal conductivity, and $\sigma$ is the viscous stress tensor given by:

$$ $\sigma=\eta\left[\nabla\textbf{v}+(\nabla\textbf{v})^{T}\right]+\left(\zeta- \frac{2}{3}\eta\right)(\nabla\cdot\textbf{v})\mathbf{I},$ (16) $$

with $\eta$ and $\zeta$ denoting shear and bulk viscosities, respectively.
To close the system, an equation of state (EOS) is required. A common choice for an ideal gas is:

$$ $p=(\gamma-1)\epsilon,\quad\text{or}\quad p=\rho RT,$ (17) $$

where $\gamma$ is the adiabatic index and $R$ is the specific gas constant.

The CNS equations describe a wide range of fluid phenomena where compressibility effects are non-negligible. In particular, they can capture shock and rarefaction waves, supersonic and hypersonic flows, acoustic propagation, atmospheric and weather systems, combustion and detonation and astrophysical flows.

Diffusion-Reaction Flow.
We consider a 2D diffusion-reaction system involving two nonlinearly coupled variables: the activator $u=u(t,x,y)$ and the inhibitor $v=v(t,x,y)$. We utilize the data from PDEBench [100]. The governing equations are:

$$ $\displaystyle\partial_{t}u$ $\displaystyle=D_{u}(\partial_{xx}u+\partial_{yy}u)+R_{u},$ (18) $\displaystyle\partial_{t}v$ $\displaystyle=D_{v}(\partial_{xx}v+\partial_{yy}v)+R_{v}\,,$ (19) $$

where $D_{u}$ and $D_{v}$ denote the diffusion coefficients, and $R_{u}(u,v)$ and $R_{v}(u,v)$ represent the nonlinear reaction terms for the activator and inhibitor, respectively. The simulation is defined over the spatial domain $x,y\in(-1,1)$ and temporal range $t\in(0,5]$. This system is widely used to model biological pattern formation.
The reaction terms follow the FitzHugh–Nagumo model [45]:

$$ $\displaystyle R_{u}(u,v)$ $\displaystyle=u-u^{3}-k-v,$ (20) $\displaystyle R_{v}(u,v)$ $\displaystyle=u-v,$ (21) $$

with $k=5\times 10^{-3}$, and diffusion coefficients $D_{u}=1\times 10^{-3},D_{v}=5\times 10^{-3}$. The initial condition is sampled from a standard normal distribution: $u(0,x,y)\sim\mathcal{N}(0,1.0)$, for $x,y\in(-1,1)$. We provide simulation data discretized at high resolution $(N_{x}=N_{y}=512,N_{t}=501)$ and a downsampled version $(N_{x}=N_{y}=128,N_{t}=101)$ for training. Spatial discretization uses the finite volume method [74], and time integration is performed via a fourth-order Runge-Kutta scheme from the scipy library.

Kolmogorov Flow. It is a novel physical scenario and a 2D variant of the classical Kolmogorov flow [71]. We follow Poseidon [35] for data settings. This flow is governed by the incompressible N-S equations with an external forcing term. The governing equations take the form:

$$ $u_{t}+(u\cdot\nabla)u+\nabla p-\nu\Delta u=f,\quad{\rm div}~{}u=0,$ (22) $$

defined over the spatial domain $[0,1]^{2}$ with periodic boundary conditions imposed in both spatial directions. Here, $u=(u_{x},u_{y})$ denotes the velocity field, $p$ is the pressure, $\nu$ is the viscosity, and $f$ is the external forcing. The forcing term $f$ is smooth, spatially varying, and held constant over time. It is defined as:

$$ $f(x,y)=0.1\sin(2\pi(x+y)).$ (23) $$

This configuration introduces a diagonal sinusoidal excitation that drives the flow in a nontrivial manner, resulting in rich dynamics ideal for evaluating neural operator performance.
The initial conditions are based on the vorticity,
which is assumed to be constant along a uniform (square) partition of the underlying domain. Numerical simulations are conducted using the same discretization and time-stepping schemes adopted for other N-S-based flows.

This problem can be recast into the unified PDE framework by defining an augmented state variable $U=[u_{x},u_{y},f]$, and introducing the trivial evolution equation $\partial_{t}f=0$ to reflect the time-invariance of the forcing. The initial data is augmented accordingly using Equation ([23](https://arxiv.org/html/2505.20349v1#A6.E23)). The corresponding solution operator $\mathcal{S}(t,\cdot)$ maps the initial state $U^{0}_{x,y}$ to the solution at time $t$ as:

$$ $\mathcal{S}(t,U^{0}_{x,y})=[u_{x}(t),u_{y}(t),f]$ (24) $$

where $(u_{x},u_{y})$ evolve according to the forced N-S Equations ([22](https://arxiv.org/html/2505.20349v1#A6.E22)), and $f$ remains fixed over time.

Stochastic N-S Equation. We extend the existing 2D incompressible N-S equation by adding a stochastic forcing term, resulting in the 2D stochastic N-S equation, an emerging modeling class of fluid dynamics that aims to better handle uncertainty and study the mix-in behavior of fluid flows [73]. In particular, we present semilinear form over the space-time interval $(t,x)\in[0,T]\times[0,1]^{2}$, similar to the dataset in [89]:

$$ $\frac{\partial\mathbf{u}}{\partial t}=\nu\Delta\mathbf{u}+\mathbf{f}+\sigma\xi ,\quad\mathbf{u}(0,x)=\mathbf{u}_{0}(x),$ (25) $$

where $\xi=\dot{W}$ for $W$ a $Q$-Wiener process, colored in space and with scale $\sigma=0.05$ and viscosity $\nu=10^{-4}$. For our data generation, the initial condition is sampled according to $\mathbf{u}_{0}\sim\mathcal{N}(0,3^{3/2}(-\Delta+49I)^{-3})$ with periodic boundary conditions.

### F.2 Datasets for comparison with numerical solvers

2D Incompressible N-S Equation in Vorticity Form.
This pseudo-spectral N-S solver uses a semi-implicit, two-stage (Heun) time integrator and adaptive time stepping governed by CFL conditions. At each step, the vorticity field $\omega(x,y)$ is transformed to Fourier space, the stream-function $\psi$ is obtained via multiplication by the precomputed inverse Laplacian, and the velocity $(u,v)=(\partial_{y}\psi,,-\partial_{x}\psi)$ is recovered by spectral differentiation. Nonlinear advection is computed in physical space and re-transformed, while viscous diffusion is handled implicitly through Crank–Nicolson scaling. To suppress aliasing, high-wavenumber modes are zeroed via a $2/3$-rule mask.

We use a physical domain of $1\times 1$ ($L_{1}=L_{2}=1$) and a discretization grid size of $256\times 256$ ($s_{1}=s_{2}=256)$. The total time $T=20$, and we store vorticity data at each $\frac{1}{50}$ intervals. This means that the neural solvers make predictions with $\Delta t=\frac{1}{50}$. We generate 1200 sequences each with a sequence length of 1000. We use 1000 sequences for training, 100 sequences for validation, and 100 sequences for testing. The external forcing term is $f(x,y)=0.1(\sin(2\pi(x+y))+\cos(2\pi(x+y)))$. The initial vorticity is generated from a 2D gaussian random field.

Figure: Algorithm 1 Pseudo-Spectral Solver for 2D Incompressible N-S

2D Burgers’ Equation.
This solver advances the 2D Burgers equations as follows:

$$ $\displaystyle u_{t}+u\,u_{x}+v\,u_{y}$ $\displaystyle=\nu\!\left(u_{xx}+u_{yy}\right),$ (26) $\displaystyle v_{t}+u\,v_{x}+v\,v_{y}$ $\displaystyle=\nu\!\left(v_{xx}+v_{yy}\right).$ $$

using a pseudo-spectral RK4 integrator. Spatial derivatives are computed in Fourier space, while nonlinear advection is formed in physical space and then dealiased with a 2/3-rule mask. Viscous diffusion appears as a spectral multiplier $-\nu,K^{2}$ and is treated in the same explicit stages. Crucially, the time step $\Delta t$ is chosen adaptively via a CFL constraint.

We use a physical domain of $1\times 1$ ($L_{1}=L_{2}=1$) and a discretization grid size of $256\times 256$ ($s_{1}=s_{2}=256)$. The total time $T=20$, and we store velocity data at each $\frac{1}{50}$ intervals. We generate 1200 sequences each with a sequence length of 1000. We use 1000 sequences for training, 100 sequences for validation, and 100 sequences for testing. The initial velocity is generated from a 2D gaussian random field.

Figure: Algorithm 2 Pseudo-Spectral Solver for 2D Burgers Equation

ComputeRHS($\hat{u},\hat{v}$):
Compute nonlinear advection and diffusion in spectral space, i.e.

$$ $R_{u}=-\,\mathcal{F}\{u\,u_{x}+v\,u_{y}\}\,M\;-\;\nu\,K^{2}\,\hat{u},\quad R_{ v}=-\,\mathcal{F}\{u\,v_{x}+v\,v_{y}\}\,M\;-\;\nu\,K^{2}\,\hat{v}.$ $$

2D Advection Equation.
This solver advances the linear advection equation

$$ $\displaystyle u_{t}+v_{x}(x,y)\,u_{x}+v_{y}(x,y)\,u_{y}=0$ (27) $$

via a pseudo-spectral Runge–Kutta 4 scheme. At each step, the field is FFT-transformed, the spatial-dependent velocity is sampled, and an adaptive CFL-based time step is chosen to ensure stability. Nonlinear advection fluxes are computed in physical space, re-transformed, differentiated spectrally, and then combined in the RK4 stages; dealiasing enforces the 2/3-rule before the next step.

We use a physical domain of $1\times 1$ ($L_{1}=L_{2}=1$) and a discretization grid size of $256\times 256$ ($s_{1}=s_{2}=256)$. The total time $T=20$, and we store field data at each $\frac{1}{50}$ intervals. We generate 1200 sequences each with a sequence length of 1000. We use 1000 sequences for training, 100 sequences for validation, and 100 sequences for testing. The field velocity is generated from a 2D gaussian random field.

In the following, we specify the spatial dependent velocity.

Advection0:

$$ $\displaystyle v_{x}(x,y)$ $\displaystyle=\sin\!\bigl{(}2\pi y+t\bigr{)},$ (28) $\displaystyle v_{y}(x,y)$ $\displaystyle=\cos\!\bigl{(}2\pi x+t\bigr{)}.$ $$

Advection1:

$$ $\displaystyle v_{x}(x,y)$ $\displaystyle=U_{0}\frac{y}{L_{y}},$ (29) $\displaystyle v_{y}(x,y)$ $\displaystyle=0.$ $$

Advection2:

$$ $\displaystyle v_{x}(x,y)$ $\displaystyle=\sin\!\Bigl{(}\frac{\pi x}{L_{1}}\Bigr{)}\,\cos\!\Bigl{(}\frac{ \pi y}{L_{2}}\Bigr{)},$ (30) $\displaystyle v_{y}(x,y)$ $\displaystyle=-\cos\!\Bigl{(}\frac{\pi x}{L_{1}}\Bigr{)}\,\sin\!\Bigl{(}\frac{ \pi y}{L_{2}}\Bigr{)}.$ $$

Advection3:

$$ $\displaystyle v_{x}(x,y)$ $\displaystyle=-\bigl{(}y-\tfrac{L_{2}}{2}\bigr{)},$ (31) $\displaystyle v_{y}(x,y)$ $\displaystyle=x-\tfrac{L_{1}}{2}.$ $$

Advection4:

$$ $\displaystyle v_{x}(x,y,t)$ $\displaystyle=a\,x,$ (32) $\displaystyle v_{y}(x,y,t)$ $\displaystyle=-a\,y.$ $$

Figure: Algorithm 3 Pseudo-Spectral RK4 Solver for 2D Linear Advection

ComputeRHS($\hat{u},v_{x},v_{y}$):
Compute the spectral right-hand side of $u_{t}+\nabla\!\cdot\!(\mathbf{v}\,u)=0$ by
1. transforming $\hat{u}$ to $u$,
2. forming fluxes $u\,v_{x}$ and $u\,v_{y}$ and re-transforming,
3. differentiating in $x,y$ via $iK_{X},iK_{Y}$,
4. summing $-(\partial_{x}+\partial_{y})$ in spectral space.

### F.3 Datasets for comparison with discretization strategy

We employ Smoothed Particle Hydrodynamics (SPH) [25] to generate three Lagrangian particle subsets: the Taylor–Green vortex (TGV), the lid-driven cavity flow (LDC), and the reverse Poiseuille flow (RPF). Each subset is governed by the compressible N-S equations.
For Taylor–Green vortex, we set the Reynolds number $\mathrm{Re}=100$ and simulate on a $1\times 1$ periodic domain with spatial resolution $\Delta x=0.01$ and time step $\Delta t=4\times 10^{-4}$. We generate 100 sequences for training and 50 sequences for testing, each comprising 126 time steps, using $N=10,000$ particles. For Lid-driven cavity, we also use the same domain and discretization parameters ($\Delta x=0.01$, $\Delta t=4\times 10^{-4}$) at $\mathrm{Re}=100$, producing 10k frames for training and 5k frames for testing with $N=11,236$ particles. For Reverse Poiseuille flow, $\mathrm{Re}$ is set as 10 on a $1.12\times 1.12$ periodic domain with $\Delta x=1.25\times 10^{-2}$ and $\Delta t=5\times 10^{-4}$. We generate 20k training frames and 5k testing frames with $N=12,800$ particles. For all three datasets, particle positions are recorded every 100 time steps, so that the neural models are trained to predict the flow evolution over intervals of $100\Delta t$.

## Appendix G More Implementation Details

Modular Comparison. All of our experiments are implemented on PyTorch on 8 $\times$ NVIDIA A6000 GPUs. For different experiments, we utilize different optimizers with different strategies, including weight decay, learning rates, and scheduler obtained from grid search for all experiments. Please refer to the config folder in our anonymous GitHub codebase [https://anonymous.4open.science/r/FD-Bench-15BC](https://anonymous.4open.science/r/FD-Bench-15BC) for more details.

Traditional solvers. FNO uses width$=$12 and width$=$32. We train the model for 10 epochs using a cosine annealing learning rate with a starting learning rate of $1e-4$ and an end learning rate of $1e-6$. We use the Adam optimizer.

Discretization comparison. FNO uses width=32 and width=128. MeshGraphnets uses latent size=128 and message passing steps=10. GNS uses latent size=128 and message passing steps=10. All models are trained for 100 epochs using a cosine annealing learning rate with a start learning rate of $1e-4$ and an end learning rate of $1e-5$. We use Adam optimizer for all model training.

## Appendix H Visualization Study

In this section, we provide more visualizations for the simulated fluid and showcases as a supplement. We choose some specific methods and fluid types for demonstration. For example, we show the visualization results of Fourier-based methods on Compressible N-S in Figure [12](https://arxiv.org/html/2505.20349v1#A8.F12), [13](https://arxiv.org/html/2505.20349v1#A8.F13) and [14](https://arxiv.org/html/2505.20349v1#A8.F14) for prediction, ground truth and residual gap respectively. Besides, we also provide results of the Fourier-based method on Diffusion-Reaction in Figure [15](https://arxiv.org/html/2505.20349v1#A8.F15). And visualizations of three different methods on Stochastic N-S are shown in Figure [16](https://arxiv.org/html/2505.20349v1#A8.F16), [17](https://arxiv.org/html/2505.20349v1#A8.F17), and [18](https://arxiv.org/html/2505.20349v1#A8.F18).

Figure: Figure 12: Visualization of Fourier$+$Next-step$+$Variable on Compressible N-S.
Refer to caption: x11.png

Figure: Figure 13: Visualization of Fourier$+$Next-step$+$Variable on Compressible N-S.
Refer to caption: x12.png

Figure: Figure 14: Visualization of Fourier$+$Next-step$+$Variable on Compressible N-S.
Refer to caption: x13.png

Figure: Figure 15: Visualization of Fourier$+$Next-step$+$Variable on Diffusion-Reaction.
Refer to caption: x14.png

Figure: Figure 16: Visualization of Conv$+$Next-step$+$Variable on Stochastic N-S .
Refer to caption: x15.png

Figure: Figure 17: Visualization of Fourier$+$Next-step$+$Variable on Stochastic N-S .
Refer to caption: x16.png

Figure: Figure 18: Visualization of Latent$+$Next-step$+$Variable on Stochastic N-S .
Refer to caption: x17.png

## Appendix I Limitations and Broad Impact

Despite our efforts to cover a diverse set of architectures, our benchmark is limited by time constraints and thus does not yet include all possible representation paradigms (e.g., higher-order spectral methods, graph-wavelet embeddings, or learned Lagrangian bases). Expanding to these and other emerging modalities remains an important direction for future work.

Nonetheless, we have released our full evaluation pipeline as a modular, well-documented codebase that makes it straightforward to incorporate new model classes, data modalities, or loss functions. This design ensures that researchers can readily extend the benchmark, plugging in novel representations and comparators with minimal effort. In this way, our framework lays the groundwork for a truly fair, reproducible, and easy-to-use standard for comparing data-driven fluid simulation methods going forward.

## Peer Reviews

confidence: 2
contribution: 3
initial_rating: 6
presentation: 3
questions: See weaknesses.
soundness: 2
strengths: 1.This paper accurately identifies the core pain point of "fragmented evaluation" in the current field of neural PDE solvers, and the proposal of FD-Bench has clear practical significance.
2.This paper unifies and abstracts various existing methods into a quadruple of "spatial encoding + temporal encoding + loss function + tricks", and designs controlled variable experiments based on this to achieve "structure-performance" decoupling.
3.This paper provides an easy-to-use and reproducible codebase that addresses the long-standing lack of standardized evaluation protocols and helps to facilitate fair and reproducible evaluations in the future.
summary: The current field of data-driven fluid simulation lacks unified datasets and standardized evaluation protocols, and the highly coupled design of spatial, temporal, and loss modules makes fair comparisons between methods difficult. This paper proposes FD-Bench, which establishes a comprehensive leaderboard, addressing long-standing issues of repeatability and comparability, and laying the foundation for robust evaluation of future data-driven fluid models.
weaknesses: 1.This paper emphasizes comparisons across 89 baselines. However, Table 3 appears to extract components from existing methods, categorize them, and then compare the methods after reorganizing the components, which doesn't seem entirely consistent with the description.
2.FD-Bench’s primary metrics are RMSE, nRMSE, and frequency‑domain RMSE. These capture “fit” but do not directly assess physical consistency.

---

confidence: 5
contribution: 2
initial_rating: 4
presentation: 3
questions: I suggest these questions for consideration (with the understanding that many will be difficult to execute in the short response time, but with the hope that this will help the authors strengthen this good paper for the community)

*Geometry and BC realism*: Could the authors extend FD-Bench with non-rectilinear geometries (curved ducts, bluff-body wakes) and mesh-based Eulerian data (unstructured tri/tet, boundary-layer refinement) to stress spatial inductive biases outside regular grids? Can you introduce tasks with mixed boundary conditions (e.g., no-slip on walls, specified pressure at the outlet, and inflow profiles at the inlet) and report the per-boundary condition breakdown of errors? This would test if attention/Fourier advantages persist with wall physics.

*Generalization beyond periodic boxes*:  Several scenarios appear periodic or effectively box-confined. How do models fare under obstacle insertion or domain deformation (e.g., adding a cylinder/airfoil) without finetuning? 

*3D scaling/memory*: The efficiency section reports GFLOPs/memory on 2D. Do the memory/runtime advantages for Fourier vs. attention hold in 3D at comparable accuracy (e.g., (128^3), (256^3))? A 3D micro-suite could meaningfully extend the benchmark.

*Rollout step size & stability claims*: The paper argues neural solvers tolerate larger (\Delta t) with minor degradation compared to classical solvers constrained by CFL. Could you report absolute wall-clock vs. effective physical time advanced per step curves, plus failure cases where large (\Delta t) breaks stability/accuracy? This would help practitioners set safe schedules.

*Task difficulty calibration*: How sensitive are rankings to turbulence intensity (e.g., higher Reynolds for cavity/Poiseuille), shock strength in compressible cases, or stiff reaction rates? A difficulty-graded version of each task could probe robustness of the reported ordering.
soundness: 3
strengths: The modularization of neural PDE solvers into spatial/temporal/loss components, along with benchmarking controlled cross-comparisons of these choices, is a useful organizing principle that helps attribute where gains come from, rather than treating methods as monoliths. The inclusion of direct comparisons to classical solvers at coarser grids and lower-order schemes (with matched error targets) is a thoughtful step toward practicality.

The benchmark scope is broad (reported 89 baselines, 10 scenarios) with standardized tuning and clear metrics (RMSE/nRMSE, fRMSE bands, memory, GFLOPs). The paper transparently documents setups and seeds, and provides a public, modular codebase for reproduction.

The paper is well structured and clear to read. I enjoyed reading this paper!

A fair, reproducible reference point for data-driven CFD is valuable. The finding that Transformers dominate on accuracy while Fourier operators offer excellent cost-accuracy trade-offs, and that temporal bundling is a robust evolution strategy, can rationalize design choices in subsequent work. The numerical-solver comparisons and run-time/speed-up evidence increase the benchmark’s practical relevance.
summary: This paper introduces FD-Bench, a modular benchmark for data-driven fluid simulation that standardizes datasets, model components, and evaluation protocols. The core contribution is to decouple neural PDE solvers into three orthogonal modules (spatial representation, temporal representation, and loss) and to compare 89 baselines across 10 flow scenarios under a unified codebase and tuning budget. The benchmark also includes head-to-head comparisons with traditional numerical solvers, generalization tests (OOD initial conditions, resolution shifts, and rollout length), and discretization studies (Eulerian vs. mesh vs. particle). Empirically, the paper reports that self-attention/Transformers tend to yield the best accuracy, while Fourier-based operators are competitive and computationally efficient; temporal bundling is consistently strong among evolution strategies. The authors release an extensible codebase and claim improved fairness and reproducibility of comparisons for data-driven fluid simulation.
weaknesses: *Novelty of questions vs. synthesis*: The three headline questions (which neural architecture; can neural replace numerical; which discretization; how well do models generalize) are useful but not conceptually new; several prior surveys/benchmarks have articulated similar axes. FD-Bench’s novelty is primarily scope and modular rigor, not new task formulations. (This aligns with the authors’ own positioning against prior benchmarks in Table 1.)

*Insights largely confirm established intuitions*: The empirical takeaways (Transformers (self-attention) excel in accuracy, Fourier operators in efficiency/long-range mixing) are directionally consistent with recent literature and community experience (e.g., work by Mishra group at ETH, reviews from Brown, etc); FD-Bench consolidates rather than surprises here.

*Geometry and boundary-condition coverage feel limited for real-world deployment*: Despite listing heterogeneity and even “irregular geometries” among selection criteria, much of the benchmark centers on canonical, mostly rectangular/periodic domains (e.g., incompressible/Compressible N-S on uniform grids, Burgers, advection, Kolmogorov flow). The Lagrangian subsets (Taylor–Green, lid-driven cavity, reverse Poiseuille) add variety, but the dominant training/evaluation remains on boxes with simple/periodic BCs, and several Lagrangian cases are aggregated back to regular grids for evaluation. This limits conclusions about performance under complex geometries (curved walls, obstacles, multi-body systems), mixed/moving boundaries, and non-Cartesian meshes, which are typical in engineering applications.

*Boundary-condition diversity*: The protocol does not yet stress diverse BC types (Dirichlet/Neumann/Robin inflow–outflow mixes, no-slip vs. partial slip, pressure outlets), CFL-critical transient regimes near complex boundaries, or wall models and source terms relevant to real devices. This weakens claims about readiness for practical deployment, even if trends are robust on boxes.

*Limited 3D and multiphysics breadth*: The benchmark mix skews toward 2D and single-physics PDEs; many industrial use cases hinge on 3D, moving interfaces, or coupled physics (e.g., conjugate heat transfer, reacting flows). As a result, the reported rankings may not carry over.

---

confidence: 5
contribution: 2
initial_rating: 4
presentation: 2
questions: - How does your modular decomposition account for models that inherently couple spatial and temporal processing (e.g., neural ODEs or recurrent architectures)?

- Why were more physically meaningful metrics (e.g., enstrophy, divergence, spectral decay) not included in the evaluation?

- Can you provide evidence that the selected flow scenarios are sufficiently diverse to draw general conclusions about model performance?

- What steps were taken to ensure that the hyperparameter tuning for each model was fair and computationally feasible across 89 baselines?

Should you be able to satisfactorily address the points I've raised above, I will accordingly provide a positive rating.
soundness: 2
strengths: - The modular design is a step toward disentangling the contributions of different components in neural PDE solvers.

- The effort to unify and standardize evaluation across a large number of models is commendable.

- The inclusion of traditional numerical solvers as baselines is a valuable addition.

- The codebase is publicly released, which may facilitate future research.
summary: The paper introduces FD-Bench, a benchmark for data-driven fluid simulation aimed at addressing the lack of standardized evaluation in the field. The authors propose a modular decomposition of neural PDE solvers into spatial, temporal, and loss components, and evaluate 89 baseline models across 10 fluid flow scenarios. They also include comparisons with traditional numerical solvers and analyze generalization across resolutions, initial conditions, and rollout horizons.
weaknesses: - Lack of Theoretical or Algorithmic Novelty: The paper presents a benchmarking effort rather than a methodological advance. While useful, it does not introduce new models, theoretical insights, or algorithmic improvements. The decomposition into spatial/temporal/loss modules is intuitive but not novel, and the paper does not justify why this particular decomposition is the most meaningful or complete.

- Superficial Model Decomposition: The decomposition of 89 models into modular components is largely based on a post-hoc categorization rather than a unified theoretical framework. This risks oversimplifying architectural differences and may misrepresent the original contributions of the models being compared. For instance, grouping models under “self-attention” or “Fourier” ignores important nuances in their design and implementation.

- Limited Justification for Dataset Selection: While the authors curate 10 flow scenarios, the rationale for their selection is not deeply motivated from a fluid dynamics perspective. The paper does not adequately address whether these scenarios are representative of real-world challenges or if they cover the full spectrum of fluid behaviors (e.g., multi-phase flows, high Mach number regimes, complex geometries).

- Evaluation Metrics Lack Physical Interpretability: The reliance on RMSE and nRMSE as primary metrics is insufficient for fluid dynamics, where quantities like vorticity preservation, energy spectra, or divergence errors are often more informative. The paper does not justify why these standard fluid diagnostics are omitted, limiting the utility of the benchmark for the CFD community.

- Insufficient Discussion of Limitations in Neural vs. Traditional Solvers: The comparison with traditional solvers focuses on error and runtime but overlooks critical aspects such as stability, convergence guarantees, and physical consistency. Neural solvers are known to suffer from error accumulation and lack of long-term stability—issues that are not sufficiently addressed in the experiments or discussion.

- Generalization Claims Are Overstated: The generalization studies (e.g., zero-shot initial conditions, resolution shifts) are limited to narrow variations and do not convincingly demonstrate robustness to realistic distribution shifts. The paper does not test on truly out-of-distribution scenarios such as different geometries, boundary conditions, or physical parameters outside the training range.

- Reproducibility and Scalability Concerns: Although the code is released, the work does not provide sufficient details on computational budgets, hyperparameter tuning procedures, or the feasibility of reproducing all 89 models. The scalability of the benchmark to 3D or more complex PDE systems is also not demonstrated or discussed.


Decision: Reject