Learning turbulent flows with generative models for super resolution and sparse flow reconstruction

Domain: chemistry
Published: 2026-03-09
Editorial outcome: minor_revision

Abstract
Neural operators are promising surrogates for dynamical systems but when trained with standard L2 losses they tend to oversmooth fine-scale turbulent structures. Here, we show that combining operator learning with generative modeling overcomes this limitation. We consider three practical turbulent-flow challenges where conventional neural operators fail: spatio-temporal super-resolution, forecasting, and sparse flow reconstruction. For Schlieren jet super-resolution, an adversarially trained neural operator (adv-NO) reduces the energy-spectrum error by 15× while preserving sharp gradients at neural operator-like inference cost. For 3D homogeneous isotropic turbulence, adv-NO trained on only 160 timesteps from a single trajectory forecasts accurately for five eddy-turnover times and offers 114× wall-clock speed-up at inference than the baseline diffusion-based forecasters, enabling near-real-time rollouts. For reconstructing cylinder wake flows from highly sparse Particle Tracking Velocimetry-like inputs, a conditional generative model infers full 3D velocity and pressure fields with correct phase alignment and statistics. These advances enable accurate reconstruction and forecasting at low compute cost, bringing near-real-time analysis and control within reach in experimental and computational fluid mechanics.

Introduction
In this work, we explore how generative modeling and operator learning can be applied synergistically to address three central challenges in turbulence: spatio-temporal super-resolution, forecasting with limited data, and reconstruction from sparse observations. Our approach combines adversarial training with operator learning for super-resolution and forecasting.

We further leverage generative models to reconstruct turbulent flows from extremely sparse Particle Tracking Velocimetry (PTV)-like velocity measurements, achieving phase alignment and matching flow statistics even when trained with as few as 150 snapshots. To place these contributions in context, we first recall how turbulence modeling and simulation have evolved over the past decades.

In the early 1970s, Howard W. Emmons argued that direct numerical simulation (DNS) of turbulence would never be computable with standard Fourier methods, because the quadratic nonlinearity in Navier-Stokes would inevitably alias high-wavenumber energy back into the resolved band1. That pessimism was overturned when Patterson & Orszag introduced the simple 3/2-rule: zero-pad the Fourier modes to 1.5N in each direction, multiply in physical space, then truncate back to N, thereby removing all aliasing and enabling the first practical DNS of homogeneous isotropic turbulence on a (2π)3 grid2,3 in 1972. These ideas were later extended to complex-geometry flows using the spectral element method, see e.g.,4 in 1993.

However, simulating turbulent flows remains computationally expensive due to the need for fine temporal discretizations in DNS as well as relatively cheaper large-eddy simulations (LES). These small timesteps are typically required to ensure numerical stability and accurately resolve the nonlinear dynamics of turbulence. In contrast, neural operators (NOs)5,6 offer a promising alternative for modeling turbulent flows by effectively learning non-linear mappings across large time intervals, potentially bypassing restrictive time-stepping constraints.

Neural networks, and more recently neural operators, have gained traction in scientific computing due to their ability to approximate complex input-output functional mappings with fast inference once trained offline. They have been applied across a range of domains-from fluid dynamics and weather forecasting to materials science and quantum chemistry. Theoretical results like the universal operator approximation theorem by Chen and Chen7 have further accelerated interest by guaranteeing that neural operators can learn mappings between infinite-dimensional function spaces. Several architectures have been proposed, including the Deep Operator Network (DeepONet)5, which uses branch and trunk networks to learn functionals; the Fourier and Laplace Neural Operators (FNO & LNO)6,8, which use global convolutions; and the Convolutional Neural Operator (CNO)9, which replaces them with local regular convolutions. Despite their success, these models exhibit spectral bias10-that is, they struggle to capture high-frequency components, which limits their accuracy for systems with rich multiscale dynamics, such as turbulence.

The architecture of the neural operator plays a critical role in how well it handles high-frequency content. U-shaped convolutional architectures, such as CNO and Time-Conditioned UNet (TC-UNet)11,12, preserve the fine-scale details better than global architectures like FNO or DeepONet12,13. However, even these models struggle when applied to turbulent flows, where the energy spectrum follows a power-law decay across a wide range of wavenumbers. This difficulty becomes more pronounced at higher Reynolds numbers, where the inertial range extends deep into the high-frequency regime. In such cases, even expensive UNet-based operators often fail to model the dynamics associated with small-scale structures. In this work, we examine the root cause of this failure-not in the architecture, but in the way these models are trained.

Conventionally, neural operators are trained by minimizing a Euclidean loss, typically L1 or L2 error between the predicted and true states. While effective in reducing global error, this loss function implicitly favors low-frequency components, which dominate the energy spectrum in turbulent flows. As a result, the model can achieve low loss values while ignoring high-frequency content, leading to over-smooth predictions. Even architectures that are capable of representing high frequencies may fail to learn them under MSE-like training. Addressing this requires rethinking how we train neural operators.

Several recent works have targeted spectral bias in neural networks/operators by introducing clever techniques to penalize spectral content. They include phase-shift networks, multi-scale DNNs, and adaptive activations to better balance learning across low and high frequency modes14. A similar multi-scale branching strategy within the Fourier Neural Operator15 has also been investigated. Chakraborty et al.16 introduced the Binned Spectral Power loss, a frequency-domain loss that bins spectral energy and adaptively weights errors across bands to explicitly penalize mismatches in energy distribution. Khodakarami et al. introduced high-frequency scaling (HFS)17, where they decomposed the latent feature maps into non-overlapping patches, separating the low and high frequency components, which helped them achieve more accurate and localized representations of multiscale flows. Physics-informed neural networks18 also suffer from spectral bias19 because, again, the L2 norm of the PDE residuals is typically minimized. Recent remedies include Fourier feature embeddings20, architectural modifications21, and improved optimizers22,23. There have also been efforts to leverage recent advances in generative models for modeling multiscale systems like turbulence.

At its core, generative models are trained to synthesize new samples that fall within the training distribution, often conditioned on partial or low-fidelity information. By learning the full probabilistic structure of the data, they naturally mitigate the low-frequency bias that burdens deterministic operators. Generative adversarial networks (GANs)24, diffusion models25, and, more recently, flow matching26 techniques have all proven effective at learning turbulent flow dynamics. In the broader image-to-image translation literature, pix2pix-style conditional GANs demonstrate that combining pixelwise L1 or L2 reconstruction losses with an adversarial objective yields sharper and more realistic outputs than regression losses alone, an idea we want to adapt to turbulent flow fields27. Yousif et al.28 used a GAN to spatially super-resolve slices of the instantaneous velocity field in channel flow with longitudinal ribs from their downsampled counterparts. Li et al.29 employed a GAN to infer one instantaneous velocity component from another in rotating turbulent flow, and achieved closer statistical alignment with DNS than CNNs and EPOD. GAN-based workflows have also been used for spatially super-resolving isotropic30, wall31,32, and reactive turbulent flows33. Similarly, diffusion model-based frameworks34,35 have been applied to two-dimensional Kolmogorov flows. Although the diffusion models can estimate high-resolution states with good spectral fidelity, the iterative denoising makes them more expensive compared to GANs. So it would be useful to evaluate the trade-off between inference cost and spectral fidelity across different generative methods. Moreover, many of these studies focus on spatial super-resolution, whereas temporal super-resolution is also vital in turbulence, especially in settings like Schlieren systems, where achieving higher frame rates requires expensive high-speed cameras.

To simulate the temporal evolution of turbulent flows, Du et al.36 introduced the Conditional Neural Field Latent Diffusion (CoNFiLD) model, which integrates conditional neural-field encoding with latent diffusion to generatively synthesize spatiotemporal turbulence. While CoNFiLD produces time-coherent trajectories with strong spectral and statistical fidelity, its inference relies on iterative latent-diffusion sampling and is not formulated as a supervised, direct next-state forecaster mapping u_t → u_{t+∆t}. Oommen et al.13 trained a NO to directly forecast the future states, and used a diffusion model conditioned on the NO output to recover the fine-scale turbulent features that were smeared out due to the spectral bias of NO. Molinaro et al.37 introduced GenCFD by directly employing a diffusion model conditioned on the initial condition and physical lead time to generate time-evolving turbulent flow fields, achieving improved mean-field reconstruction and energy-spectrum alignment on several test cases compared to standard NOs. However, during inference, each timestep of the generated trajectory requires iterative denoising-based sampling using the diffusion model. While this is faster than DNS and can be further accelerated by batch-based parallelization in modern GPUs, avoiding iterative denoising is still desirable.

The generative models discussed above learn the functional prior of turbulent flows by training on ensembles of simulated trajectories. During inference, they can sample from the posterior distribution of flow fields conditioned on available observations. However, in practice, previous works like CoNFiLD and GenCFD require on the order of 10^2 to 10^3 volumetric training trajectories to effectively learn these priors. In turbulence, generating 10^2 to 10^3 volumetric trajectories using high-fidelity DNS incurs a substantial computational burden, which displaces the practical advantage of faster downstream inference. Leveraging NOs with stronger inductive biases, enhanced with generative modeling principles to mitigate spectral bias, can make turbulence surrogates substantially more data efficient.

Beyond forecasting, generative models have also been adopted for reconstructing full turbulent flow fields from sparse measurements. Parikh et al.38 used flow matching to reconstruct full turbulent fields using only sparse wall measurements in near-wall turbulence. Similarly, Huang et al.39 and Buzzy et al.40 demonstrated that diffusion models conditioned on partial measurements can reconstruct full functions. Chakraborty et al.41 also used a conditional diffusion model to reconstruct ERA5 atmospheric data from spatially-unstructured IGRA radiosonde observations. In general, recent studies suggest diffusion models40–42 can serve as generative priors for assimilating sparse observations to reconstruct physically plausible fields. These advances suggest a natural extension to experimental fluid mechanics, where full turbulent flow fields could be reconstructed from sparse PTV measurements, providing a bridge between limited experimental data and high-fidelity flow representations.

Building on these developments, the major contributions of our work are as follows.

1. We propose an adversarially trained neural operator (adv-NO) that bridges operator learning and generative modeling for spatio-temporal super-resolution and forecasting, leading to 10× lower energy-spectrum error than vanilla NO and 114× speed-up compared to current generative model benchmarks.

2. We adapt a conditional diffusion model for “zero-shot” flow reconstruction (with no fine-tuning to unseen masking styles within the same geometry) from, for example, sparse PTV-style velocity measurements and achieve phase alignment in the observed regions while recovering correct global statistics, as setting directly relevant to experimental fluid mechanics.

3. We demonstrate both forecasting and flow reconstruction under severely limited data, training on only a single trajectory for each task, using 160 snapshots for forecasting and 150 for reconstruction.

Additionally, we perform detailed comparisons and analyses. (a) Compared to diffusion models, GANs, and VAEs, adv-NO offers a favorable balance between inference cost and spectral fidelity. (b) Physics-informed (divergence-penalized) variants fail to mitigate spectral bias. (c) Parseval-Plancherel analysis explain why conventional L2-trained NOs underrepresent high frequencies (spectral bias). (d) Fourier analysis of learned kernels reveal that adv-NO naturally develops a high-pass filtering behavior, mitigating spectral bias. A detailed outline of this study is illustrated in Figure 1. We present our results in the next section, followed by the Discussion and Methods sections.

Results
Overview of adversarially trained neural operators (adv-NO)

Conventional neural operators (NOs) are typically trained by minimizing a pointwise regression loss (e.g., an L1 or an L2 discrepancy between the target field u_tar and the predicted field u_f). In adv-NO, we retain this reconstruction term but augment it with a perceptual feature loss and an adversarial term, yielding the composite objective shown in Equation 8. During training, a discriminator D is optimized to distinguish real samples from generated ones (via L_D), while the surrogate operator is optimized to reduce the reconstruction/perceptual losses and to produce outputs that are scored as more realistic by D (via L_G). At inference time, only the trained surrogate operator is used; the discriminator is not needed. We first evaluate adv-NO on spatio-temporal super-resolution, where the goal is to recover high-resolution flow trajectories from coarse, sparsely sampled observations.

Spatio-temporal super-resolution of Schlieren visualizations

High-speed Schlieren imaging is widely used to visualize compressible flows43,44. Here, we consider impinging jet experiments to capture the interaction between a compressible air jet and a flat plate using Schlieren imaging, which visualizes gradients in fluid density. The flow is under-expanded with a nozzle exit Mach number of 1, and the jet is directed at ambient temperature onto the surface. A standard Z-type Schlieren setup is used to record the resulting shock structures and flow features. Further details can be found in 45,46. Here, we train our surrogate NO to perform spatio-temporal super-resolution, reconstructing high-resolution, high-frame rate (HRHF) states from low-resolution, low-frame rate (LRLF) Schlieren visualizations. Unlike typical computer vision super-resolution, which focuses solely on spatial upscaling, we simultaneously recover fine spatial details and intermediate temporal states.

The dataset consists of 1000 snapshots with a resolution of [128, 256] spanning the spatial domain [−0.8, 5.6] × [−1.45, 1.75]. Each snapshot is separated by a time interval of τ = 4.76 µs. The first 800 timesteps are used for training, the next 100 for validation and the last 100 for testing. The original dataset is the HRHF representation (u_HRHF(x,t)). The LRLF state (u_LRLF(x,t)) is created by subsampling the HRHF by a factor of 8 along space and 4 along time as shown in Figure 2 a). We train the surrogate NO, TC-UNet11, to learn the mapping from [u_LRLF(x,t+τ), u_LRLF(x,t+5τ)] → [u_HRHF(x,t+τ), u_HRHF(x,t+2τ), u_HRHF(x,t+3τ), u_HRHF(x,t+4τ), u_HRHF(x,t+5τ)].

As discussed earlier, the conventional NO trained to minimize the Euclidean error suffers from spectral bias and subsequently estimates over-smoothed states. We compare predictions of NO with other generative algorithms often used in the computer vision community for image super-resolution problems. Specifically, we train variational autoencoder (VAE), generative adversarial network (GAN), and diffusion model (DM) to super-resolve the over-smoothed NO predictions (NO+Gen) as shown in Figure 2 b). We also consider an NO that was adversarially trained against a discriminator (adv-NO). In Figure 2 c) and d), we compare the energy spectra and predicted density gradients, respectively, at t+3τ. Figure 2 e), compares the error in the log of energy spectrum against the per-sample inference cost in terms of giga-floating-point operations (GFLOPs). Additionally, we report the normalized-root-mean-squared-error (nRMSE) of the field and log-energy-spectrum, number of trainable parameters, and compute cost in terms of per-sample GFLOPs, inference time and peak VRAM usage in Table 1.

Adv-NO, NO+GAN, and NO+DM achieve the best energy spectrum, followed by NO+VAE. Adv-NO offers the best accuracy-compute trade-off, recovering spectra and sharp gradients at NO-like inference cost. In Figure 2 e), NO+GAN appears more computationally demanding because its generator super-resolves NO outputs using residual-in-residual dense blocks (RRDB)47, which incur higher compute than conventional residual blocks. Similarly, NO+DM performs iterative denoising-based super-resolution of NO predictions, adding 32 forward passes per sample, and thus further increases compute.

While adv-NO and NO+Gen models can recover high-frequency content, the fine scales of their predictions are not perfectly phase-aligned with the ground truth, yielding slightly higher field errors than NO, as shown in Table 1. We emphasize that this pointwise NRMSE is dominated by low-wavenumber content and is sensitive to small phase offsets; thus, a smoother L2-trained NO can achieve a low field error while underestimating high-frequency energy. To more directly compare instantaneous structure beyond global scalar metrics, Figure S2 compares the horizontal centreline profiles of the states predicted by different models, illustrating that models with comparable NRMSE can nonetheless differ significantly in how well they reconstruct the high-wavenumber structures. For this reason, we additionally report the energy-spectrum error, which more directly quantifies multiscale physical fidelity. We also compare super-resolved states estimated using NO and adv-NO with bicubic interpolation in Supplementary section S2. While these results demonstrate how generative models can improve spatial and temporal fidelity within known sequences, we now turn to the harder challenge of forecasting turbulent flow evolution beyond the observed training window.

Forecasting with limited training data

Training a data-driven surrogate to forecast turbulent flows is not notoriously difficult because generating high-fidelity volumetric, time-resolved ground truth is extremely expensive, and experimental measurements of full 3D fields are rarely available. In practice, one is often limited to a handful of snapshots from a direct numerical simulation. With that constraint, we train our neural operator to forecast three-dimensional forced homogeneous isotropic turbulence (HIT) using only 160 time-steps from a single simulated trajectory for training.

HIT is a classical benchmark because, unlike shear dominated flows, it isolates the non-linear energy cascade of small-scale eddies under statistically stationary conditions. We generated our HIT dataset with a spectral-element solver nekRS(48) at a Taylor Reynolds number Re_λ = 90, where Re_λ = u_rms λ/ν characterizes the ratio of inertial to viscous scales (λ is the Taylor microscale, ν the kinematic viscosity). The flow state was saved every t = τ (= 0.04) (normalized time units). The eddy-turnover time49 t_E was computed as u_rms^2 / ε, where ε = 2νS:S is the mean dissipation rate and S = 1/2(∇u + ∇u^T) is the rate-of-strain tensor. Using this single 160-step trajectory (which represents roughly 22 eddy-turnover timescales) we try to forecast the next 36 steps corresponding to five eddy-turnover timescales (see Training and Implementation Details in the Methods section for architecture sizes and training configurations). To effectively train the surrogate to learn the evolution dynamics under limited training data, we split the available single trajectory of 160 snapshots of u = [u, v, w, P] to multiple sub-trajectories. Specifically, the surrogate NO and its variants are trained to learn the mapping - [u(x, t−3τ), ···, u(x, t)] → [u(x, t+τ), ···, u(x, t+4τ)].

In DNS, the state is advanced with a stable time-step ∆t_CFL, but saving full 3D fields at every ∆t_CFL is prohibitively memory intensive. Thus, trajectories are stored at a much coarser interval τ ≫ ∆t_CFL. Learning u(t) → u(t+τ) from such snapshots therefore represents a coarse-grained evolution operator, for which the reduced dynamics is generally non-Markovian and depends on history through an effective memory term (Mori-Zwanzig)50–52. We thus provide a short history window and predict a bundled horizon (4τ), balancing reduced error accumulation against the weak correlation at large lead times11.

During inference, the states beyond t+4τ are estimated autoregressively, by feeding the surrogate model’s prediction as its own input. In this experiment, we compare the standard NO, adv-NO, and a physics-informed variant of the NO that is constrained to satisfy the mass-conservation by penalizing the divergence of the velocity field vector computed using spectral derivatives. We call it NO+∇·u in this study. The DNS was performed using the spectral element method53, which satisfies the conservation laws at the quadrature points of each spectral element in a weak sense. For training NO and its variants, the DNS solution was interpolated using the solver’s polynomial basis from quadrature points onto an equi-spaced grid, where the divergence-free property is no longer satisfied. To address this issue, we enforce incompressibility by solving for φ from ∇^2φ = ∇·u, and using u − ∇φ as the updated solution. The results of this study are provided in Figure 3.

We compare the kinetic energy (KE) at ∆t = 2t_E (= 14τ) predicted by NO, adv-NO and NO+∇·u against the reference DNS simulated using the spectral element method in Figure 3 b). The velocity and pressure states are separately visualized in Figure S3. The KE and the [u, v, w, P] states predicted by adv-NO achieve better phase alignment with DNS compared to NO and NO+∇·u. The time-averaged energy spectrum compared in Figure 3 c) further highlight the over-smoothing/spectral bias suffered by NO. The energy spectra comparison between NO and NO+∇·u suggests that penalizing the divergence does not correct the spectral bias of a conventional L2-trained NO. This result is explainable because NO+∇·u also minimizes an L2 norm of the mass-conservation PDE residual, and therefore, struggles to minimize the fine-scale features that arise in the PDE residual of the NO prediction due to spectral bias. That being said, enforcing incompressibility does help NO+∇·u achieve a better phase-alignment than NO. In Figure 3 d), we investigate how the field error evolves across the entire forecast horizon of over five eddy-turnover timescales (which corresponds to 36τ). For the field error, we take the RMSE of the normalized velocity and pressure states at each timestep. The time-evolution of kinetic energy and energy spectrum is further visualized in Figure S4. The adv-NO is able to preserve the Kolmogorov scaling throughout the forecast horizon and maintains phase alignment well until 4t_E. We also observed that explicit time-conditioning (TC) is important for stabilizing long autoregressive rollouts in adv-NO. As shown in Figure S5, adv-NO with TC maintains coherent large-scale evolution over multiple eddy-turnover times, whereas removing TC leads to temporal stagnation of the low-wavenumber structures and the emergence of spurious high-wavenumber artifacts. We interpret this as TC regularizing the influence of adversarial training, which can otherwise over-emphasize sharp gradients/high-wavenumber structures.

We also benchmark adv-NO against GenCFD37, a conditional score-based diffusion model that conditions on the “lead time" ∆t and generates forecasts through iterative denoising. In our low-data setting, a one-shot GenCFD mapping from u(x,t) to [u(x,t+τ), ..., u(x,t+36τ)] did not converge because of the weak input-output correlation between states separated by a ∆t = 36τ (∼ 5t_E). We therefore adopt the same 4τ-ahead setup used for adv-NO and rollout autoregressively to 5t_E. Both models have ≈ 6M parameters and are trained for 48 GPU-hours. Except for the above-mentioned changes, the default training/inference settings were adopted unchanged. GenCFD is unable to reconstruct the fine-scale turbulent structures and outputs noisy predictions under the data-scarce setting, as shown in Figure 4. Forecasting 36 timesteps with GenCFD takes 228 s on an NVIDIA-H100 GPU, because of the iterative denoising of each forecasted timestep, making it over a hundred times slower than adv-NO. A detailed error-cost analysis is provided in Table 2.

We further assess small-scale structure via the invariants of the velocity-gradient tensor A = ∇u, defined as the coefficients of its characteristic polynomial λ^3 + Pλ^2 + Qλ + R = 0, with P = −tr A (P=0 for incompressible flow), Q = −1/2 tr(A^2), and R = −det A. Writing A = S + Ω with strain S = (A+A^T)/2 and rotation Ω = (A−A^T)/2, one has Q = 1/2 (||Ω||^2 - ||S||^2), so Q > 0 indicates rotation-dominated regions and Q < 0 strain-dominated regions. We normalize Q by ⟨A_ij A_ij⟩ and R by ⟨A_ij A_ij⟩^3/2. We compare the joint PDF p(Q,R) by overlaying isocontours (Fig. 4 b), which exhibits the typical tear-drop shape in the Q-R plane. Adv-NO aligns best with the DNS.

To further probe scalability beyond the moderate-Re regime, we additionally train NO and adv-NO on 3D HIT at Re_λ = 150 and evaluate a rollout at t = 2t_E (Figure 5). The same trends persist: adv-NO better retains small-scale content, with improved agreement in the high-wavenumber portion of E(k) and sharper midplane velocity/pressure structures, whereas NO exhibits stronger over-smoothing. Together, these results suggest that the accuracy gains of adv-NO extend to a more turbulent setting without modifying the overall workflow.

Flow reconstruction from partial observations

Unlike forecasting, where we predict future flow states from a history of complete snapshots, many experimental scenarios involve the inverse problem of reconstructing the full flow field from only sparse and partial measurements. In experimental fluid mechanics, the complete flow field is rarely accessible. Techniques such as particle tracking velocimetry (PTV)54 provide velocity measurements only at the locations of seeded particles, resulting in inherently sparse data. These measurements are often restricted to a smaller subdomain or a single plane of the flow55. Direct pressure measurements are even more challenging. Here, we aim to reconstruct the full three-dimensional velocity and pressure fields from sparse PTV-like velocity observations, ensuring phase alignment near the measurement region while also matching the global flow statistics, under very limited number of training samples.

Unlike Tasks 1 and 2, where we learn a forward operator surrogate that maps between fully observed flow fields for spatio-temporal super-resolution or time advancement, the sparse-measurements setting is an ill-conditioned inverse reconstruction problem. Formally, given partial measurements, the goal is to infer the full 3D state by learning and sampling from the posterior distribution p(full flow | partial observations), rather than approximating a forward evolution/super-resolution operator. To achieve this, we train a conditional generative model to sample directly from this posterior. In this study, the model is trained on high-fidelity DNS data of the cylinder wake at Re = U_∞ D / ν = 11,000, where U_∞ is the inflow velocity, D is the diameter of the cylinder and ν is the kinematic viscosity, generated by using the spectral element code nekRS(48) in the domain [−12D, 20D] × [−10D, 10D] × [0, 9.6D], where the cylinder axis is placed on x=0, y=0. Nonetheless, the training data is collected in the domain [−4D, 4D] × [−4D, 4D] × [0, 2D]. During training, the model learns to infer missing flow information by conditioning on partially observed inputs.

Specifically, we generate a binary mask for each of the four physical variables [u, v, w, P]. The entries are randomly set to zero (missing) or one (observed) according to a mask fraction sampled uniformly from U(0,1). Half of each mini-batch uses independent random masking per channel, while the other half uses a single random mask replicated across all channels. This mask is multiplied with the flow field to produce the masked input, which we then concatenate with the mask itself. This concatenated tensor (with 4+4 channels) serves as the conditioning input to the model. (See Algorithm 1 in SI). Sampling mask fractions that range from 0 to 1, enables the training to span the full spectrum between deterministic reconstruction (mask fraction=0, no missing data) and purely generative reconstruction (mask fraction=1, all data missing). This variability enables the model to get exposed to the entire continuum of partial observations.

Accordingly, we compare two conditional generative models: (i) a conditional GAN baseline (Real-ESRGAN adapted for 3D reconstruction) and (ii) a conditional score-based diffusion model (EDM). For clarity and consistency with the paper’s narrative, we reserve the term adv-NO for the adversarially-trained forward operator surrogates used in Tasks 1 and 2 (where adversarial losses are used to mitigate spectral bias and sharpen small scales in operator learning), and therefore refer to the Task 3 adversarial baseline as a conditional GAN. In this case, the generator is a UNet without time-conditioning as we are not learning in time. See Methods for details on how we adapt these generative models, originally designed for 2D image synthesis, to the volumetric reconstruction of 3D vector fields. Both models are conditioned identically on masked DNS data, allowing us to directly compare their ability to reconstruct full flow fields from sparse observations. Training, validation, and test datasets consist of 150, 18, and 17 snapshots, respectively. Again, our datasets are small to reflect the limited access to experiments and high-fidelity simulations. After training the models, we investigate the model’s ability to zero-shot reconstruct the velocity fields from sparse velocity observations in the domain (Test Setup 1), on a plane (Test Setup 2), in a sub-domain (Test Setup 3), and reconstruction of the full pressure field from sparse velocity observations (Test Setup 4), as shown in Figure 6. Here, “zero-shot” denotes applying a single pre-trained model without any fine-tuning to previously unseen masking styles/measurement layouts on new snapshots of the same cylinder-wake geometry (we do not consider generalization to unseen geometries in this study). The same random mask is used for each of the three velocity components (as shown in the region highlighted in green), and the pressure field is completely masked out to mimic the real PTV experiments. More details on the similarities and differences with PTV experiments56 are provided in Supplementary section S4.

In Test Setup 1, when the mask fraction is 0.5, both DM and GAN achieve good phase alignment with matching statistics. But GAN fails when the mask-fraction is increased to 0.99. The kinetic energy of the reconstructions, energy spectra, and intermittency for mask fractions = 0, 0.1, 0.3, 0.5, 0.7, 0.9, 0.99 are visualized in Figure S6. Both DM and GAN are able to successfully reconstruct the vortex street for all mask fractions (except GAN when mask fraction = 0.99). A detailed illustration of Test Setup 2 is provided in Figure S7. The DM achieves strong phase alignment on the measurement plane (z=0). The phase alignment gradually decreases on planes that are farther from the observation plane. The corresponding energy spectra and intermittency indicate that the DM recover turbulent statistics accurately. In contrast, GAN severely overestimates the higher velocity gradients, as seen in the intermittency plots (Figure 6) and fail to reconstruct the vortex street. We observe a very similar trend for DM and GAN in the context of phase alignment and turbulent statistics, in Test Setup 3 (see Figure S8). The pressure reconstructions corresponding to different mask fractions are illustrated in Figure S9. The diffusion model achieves strong phase alignment for all the mask fractions considered. The GAN reconstructions are suboptimal for mask fractions up to 0.9, and fail when the mask fraction is 0.99. Across all test setups, the conditional diffusion model demonstrates robust reconstruction capability from sparse observations, recovering both large and small-scale features even when the available data is severely limited.

Methods
Neural Operators (NO)

Neural operators G_θ are built to learn mappings between function spaces, with the goal of approximating the true underlying operator N. Formally, for an operator N : U → V defined between Banach spaces U and V, we train a surrogate operator G_θ to learn the mapping U → V. Specifically, in super-resolution we learn G_SR : U_LRLF → U_HRHF that maps a sparsely sampled/low-resolution field u_in = u_LRLF ∈ U_LRLF (defined only on a subset of (x,y,t)) to its dense/high-resolution counterpart u_tar = u_HRHF ∈ U_HRHF. In forecasting we learn G_FC that advances the state in time, taking as input u_in = u(·,t) or, more generally, a short history window u_in = [u(·,t−(m−1)∆t), ..., u(·,t)] and predicting the future field u_tar = u(·,t+∆t). Both tasks are trained by supervised loss minimization over pairs (u_in, u_tar):

θ* = argmin_θ E_{(uin,utar)} ||u_tar − G_θ(u_in)||_p, (5)

where for super-resolution (u_in, u_tar) = (u_LRLF, u_HRHF) and for forecasting (u_in, u_tar) = (u(·,t or history), u(·,t+∆t)). ||·||_p denotes an appropriate Lp norm (p is typically 1 or 2), which measures the discrepancy between the target u_tar and the neural operator prediction G_θ(u_in).

Here, we use a time-conditioned UNet11,12 as our neural operator, which takes u_in and ∆t as the input, and predicts u_tar. This makes the operator continuous in time, while the spatial side is still discretized on a grid. Our model is mesh-dependent, but it still qualifies as operator learning since the target is the underlying operator N. The lack of strict mesh-invariance could raise concerns, but in practice these architectures often outperform Fourier Neural Operator (FNO) and its variants11,12. One reason is the favorable inductive bias: convolutions naturally capture locality and translation symmetry in the data, and temporal conditioning makes the dynamics easier to model. This also makes the time-conditioned UNet especially effective in low-data regimes, which is often the case in turbulence, where DNS and experiments are expensive to generate/access.

Adversarially Trained Neural Operators (adv-NO)

The Lp objective in Equation 5 yields accurate low-wavenumber structure but under-penalizes errors at high wavenumbers, leading to over-smoothing (spectral bias) in both super-resolution and forecasting of multiscale systems like turbulence. We resolve this limitation while retaining the operator’s temporal conditioning and fast inference by training the same G_θ but adversarially against a discriminator.

Given an input u_in (LRLF snapshot for super-resolution or a history window for forecasting) and the target u_tar, let u_f = G_θ(u_in). Since u_f is estimated by the NO, we treat it as a “fake” state and u_tar as a “real” state. We adopt the relativistic average GAN (RaGAN) objective63, where a UNet-based discriminator estimates the relative realness between the real (u_tar) and fake (u_f) samples as shown below.

L_D = −E_{utar} [ log σ(D_ψ(u_tar) − E[D_ψ(u_f)]) ] − E_{uf} [ log (1 − σ(D_ψ(u_f) − E[D_ψ(u_tar)])) ], (6)

L_G = −E_{utar} [ log (1 − σ(D_ψ(u_tar) − E[D_ψ(u_f)])) ] − E_{uf} [ log σ(D_ψ(u_f) − E[D_ψ(u_tar)]) ], (7)

where σ(·) is the sigmoid and expectations are over minibatch samples.

During training, the discriminator parameters ψ are updated by minimizing L_D, which encourages it to assign higher relative scores to real samples and lower scores to generated ones. The surrogate operator was trained using a combination of L1 loss, perceptual loss, and adversarial loss (L_G) (Equation 7). L1 loss simply penalizes the difference between the u_tar and u_f, and guides the learning of the low-wavenumber structures. The perceptual loss is computed on intermediate feature maps from a pre-trained network. Let {φ}^L denote the first L convolutional layers of a pre-trained network. The perceptual loss is defined as ∑_{l=1}^L w_l ||φ^l(u_tar) − φ^l(u_f)||_1. For 2D data we use the first five layers of VGG-1964 with weights (w_1, ..., w_5) = (0.1, 0.1, 1, 1, 1). For 3D data we use the first three layers of 3D-ResNet10 from Med3D65 with weights (w_1, w_2, w_3) = (0.1, 0.1, 1). Adversarial loss L_G drives the time-conditioned neural operator to produce outputs that the discriminator considers more realistic than real samples on average. The surrogate operator’s parameters θ are updated by minimizing the composite loss term,

L_advNO = ||u_tar − u_f||_1 + ∑_{l=1}^L w_l ||φ^l(u_tar) − φ^l(u_f)||_1 + βL_G, (8)

where φ^l is the lth feature map from a pretrained network, weighed by w_l. β is a hyper-parameter coefficient that controls the contribution of L_G66 and was set to 0.1.

The first term in Equation 8 guides the large-scale accuracy. Adversarial loss L_G (third term), and the perceptual loss (second term) helps approximate the high-frequency turbulent structures without sacrificing the temporal coherence from the operator formulation. We alternate NO (θ ← argmin L_advNO) and discriminator (ψ ← argmin L_D) updates, progressively refining the NO by adversarially training it against an improving discriminator, during the course of training.

We also investigate a complementary strategy where we use generative models (VAE, DM, GAN) to super-resolve the smooth states estimated by the NO. The NO handles temporal evolution, while the generator restores fine-scale structure. Next, we overview the different generative models used in this work.

Variational Autoencoders (VAE)

We use an improved Vector-Quantized Variational Autoencoder (VQ-VAE)67, adapted for super-resolution tasks. The unique aspect of VQ-VAE is that the continuous Gaussian latent space of a standard VAE is replaced with a discrete codebook of learnable embeddings. The finite set of embeddings serves as a learnable dictionary of features, which helps the model to compress similar inputs into the same code, improving generalization and interpretability. Moreover, the mapping to integer code yields a compact representation which also helps with memory constraints. The training of VQ-VAE is governed by the loss L_VAE defined as,

L_VAE = ||u_tar − VAE(G_θ(u_in))||_2^2 + ||γ[z_e(u_tar)] − e_k||_2^2 + β||z_e(u_tar) − γ[e_k]||_2^2, (9)

where VAE(G_θ(u_in)) is the state super-resolved by the VQ-VAE that takes NO predictions as its input. z_e(u_tar) is the continuous latent representation mapped by the encoder, e_k is the nearest codebook vector quantized from its continuous representation z_e(u_tar), γ[·] is stop gradient function and β is a hyperparameter balancing the encoder commitment. The stop-gradient function γ[·], prevents gradients from flowing through certain parts of the model, effectively treating them as constants during training. This is crucial in VQ-VAE because the quantized output is non-differentiable, as it’s constructed using nearest-neighbor lookup. In the loss term ||γ[z_e(u_tar)] − e_k||_2^2, the stop gradient function ensure that the encoder output z_e(x) is kept constant, so that only the codebook vectors e_k are updated. In contrast, in the commitment loss term β||z_e(u_tar) − γ[e_k]||_2^2, the codebook vectors are treated as constant, ensuring only encoder parameters are updated. We observed that VQ-VAEs in their original form struggle to capture the fine-scaled turbulent structures. This occurs due to the compression and quantization during the codebook construction. Small oscillations become indistinguishable during nearest-neighbor lookup and hence, the decoder struggles to reconstruct high-wavenumber features solely from these quantized coarse latents. To address these issues, we use UNet-style skip connections to propagate the intermediate features from the encoder through concatenation at the decoder. This facilitates reconstructions of fine edges and small oscillations, relaxing the codebook quantized space bottleneck. Importantly, this modification aligns well with our super-resolution objective, since the goal is not dimensionality reduction but accurate recovery of fine-scale turbulent features.

Diffusion Models (DM)

We also consider a conditional score-based diffusion model to refine fields via denoising score matching across multiple noise levels. Let u denote a field and u_cond the conditioning signal. The score network s_φ(u; σ, u_cond) estimates the gradient of the log-density of the Gaussian-perturbed conditional distribution,

s_φ(u; σ, u_cond) ≈ ∇_u log p_σ(u|u_cond),

where p_σ is the data distribution corrupted with N(0, σ^2I). Following EDM68, we sample σ from a log-normal schedule and minimize a preconditioned denoising score-matching objective over pairs (u_cond, u_tar):

min_φ E_{σ,ε,(ucond,utar)} [ ω(σ) ||s_φ(u_tar + σε; σ, u_cond) + ε/σ||_2^2 ],

where ε ∼ N(0, I) and ω(σ) is the EDM preconditioning weight. In other words, the score networks sees the target corrupted by Gaussian noise at level σ (i.e., u_tar + σε) and is trained to predict the corresponding noise-scaled direction −ε/σ that removes this corruption. For super-resolution, the conditioning signal u_cond = G_θ(u_in), the over-smoothed NO output and for flow reconstruction, u_cond is the sparse observation. We start sampling from a Gaussian noise at a large σ_0 and progressively lower σ → 0, updating the sample with the learned score at each step. We use EDM’s noise schedule and a Heun-style sampler for efficiency and stability. Training across a log-spaced range of σ yields multiscale fidelity, where large σ captures coarse structure and small σ restores high-wavenumber detail, but at the cost of iterative sampling.

Generative Adversarial Network (GAN)

The GAN implemented in our work is based on Real-ESRGAN66, which is an extension to the ESR-GAN47, designed to achieve enhanced super-resolution performance compared to previously developed GAN-based models. A deep convolutional upsampling network with residual-in-residual dense blocks (RRDB)47 is used as the generator, and a UNet architecture is used as the discriminator. The RRDB stacks three Residual Dense Blocks, each with five 3×3 convolutions. These nested convolutions substantially increase compute, explaining the higher FLOPs we report for NO+GAN in Table 1. For training, we follow the relativistic average GAN (RaGAN) objective described earlier, and include a VGG-19-based perceptual loss. RaGAN compares the relative realness of real vs. generated samples, which sharpens textures and stabilizes adversarial updates. A UNet discriminator with spectral normalization is used to increase discriminator capacity and stabilize training. Perceptual loss is computed on VGG-19 features before activation, together with an L1 content loss to learn large-scale structures, following the standard practice of defining feature-space losses69 on pretrained classification networks for super-resolution and image transformation70 tasks.

Extending from 2D to 3D: Memory-efficient network adaptations

For the spatio-temporal super-resolution problem, the perceptual loss in adv-NO and GAN was computed on the 2D feature representations projected by the pretrained layers of VGG-19, which is a 2D convolutional network as discussed earlier. For forecasting three-dimensional isotropic turbulence with adv-NO and reconstructing the three-dimensional cylinder wake with GAN, we utilize a pre-trained 3D ResNet from the Med3D65, which was originally trained on biomedical images of human organs. Specifically, we used the first three pre-trained layers of 3D-ResNet10 from65 with weights of {0.1, 0.1, 1}. However, we encountered out-of-memory errors. To make our implementation memory-efficient, we modified how we utilized the discriminator. Instead of sending the full 128^3 or 256×256×32 dimensional snapshots for HIT and cylinder cases, respectively, as the input to the UNet-based discriminator and pretrained 3D-ResNet10, we utilized randomly sampled smaller sub-domains of size 32^3 to compute the adversarial and perceptual losses. This strategy helped us avoid out-of-memory errors. Our 3D implementation of the diffusion model68 for the cylinder wake reconstruction also encountered out-of-memory errors. In this case, the culprit was the linear and full attention blocks in the neural network architecture, which we turned off to make the network more memory efficient.

Training and implementation details

All models used in the spatio-temporal super-resolution of 2D Schlieren images contained approximately two to four million trainable parameters, while those used for 3D forecasting and sparse flow reconstruction contained approximately six million parameters. We intentionally selected relatively small models to ensure efficient training on a single GPU and to enable effective learning in the low-data regime.

All models were trained using the Adam optimizer71 with a learning rate of 10^-4, and training was performed until convergence, with a maximum wall-clock time of 48 hours per experiment. For the conventional NO trained solely to minimize the L2 loss, we additionally employed a linear warm-up followed by a cosine annealing learning rate scheduler12 to accelerate convergence. While this setup enabled faster training, the conventional NO consistently suffered from spectral bias, as discussed earlier. To mitigate this limitation, we further trained the NO adversarially (adv-NO) and also evaluated VAE-, GAN-, and diffusion-based frameworks for super-resolving the NO-predicted states. In general, adv-NO, VAE, GAN, and diffusion models exhibited slower convergence compared to conventional NOs but achieved superior reconstruction of fine-scale structures and turbulent dynamics as reported in Figure 2. All training and inference were conducted on a single NVIDIA H100 GPU using the PyTorch framework72.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "With this manuscript, Vivek Oommen et al. benchmark a unified framework combining neural operators (NOs) and generative modeling (GANs, diffusion models, VAEs) to address three fundamental challenges in turbulence modeling: (1) spatio-temporal super-resolution, (2) forecasting with limited data, and (3) sparse flow reconstruction. The authors introduce an adversarially trained neural operator (adv-NO) that mitigates spectral bias and achieves high-frequency fidelity with low inference cost. The study demonstrates impressive performance on Schlieren jet data, homogeneous isotropic turbulence (HIT), and cylinder wake reconstruction, showing large gains in spectral accuracy, phase alignment, and computational efficiency.\nIn my opinion, this manuscript reports on excellent and timely work addressing a pressing challenge in data-driven turbulence modeling — the over-smoothing and spectral bias in neural operator learning. The integration of generative modeling techniques with neural operators is elegant and well-justified theoretically but not completely novel. The authors provide clear evidence across diverse benchmarks that adversarial training restores fine-scale features while maintaining the operator’s low inference cost. The experiments appear to be carefully designed, and the theoretical analysis (via Parseval–Plancherel identity and kernel Fourier analysis) strengthens the contribution.\nI believe that the manuscript is well-written and technically rigorous, with significant breadth — spanning theoretical insight, numerical validation, and practical applications in turbulence. Therefore, I recommend acceptance with minor revisions and would like to encourage the authors to consider the following comments:\n(1) The benefits of combining L1/L2 and adversarial losses are well known from pix2pix-style GANs (Isola et al., 2016). I would like to encourage the authors to acknowledge this line of research.\n(2) Although I appreciate the idea to test the sparse-data reconstruction performance with a “PTV-like” dataset, the description of the methodology left me confused here.\nPlease clarify if the random mask was set for the three velocity components independently or had equal values (zero/one) – the latter would be the case for PTV measurements.\nThe fraction of available data was varied, however this metric’s relation to experimental conditions is difficult to assess.\nPlease consider estimating the tracer image density (in particles per pixel) based on some reasonable assumptions. This would make your findings much more accessible for experimentalists. You may also want to consider the recent article by Sciacchitano et al. (Exp. Fluids, 2025) on the topic of flow reconstruction based on PTV data.\n(3) The study focuses on moderate Re (e.g., Reλ = 90). A discussion or short test at higher Re (even if 2D) would strengthen the claim of scalability.\n(4) The impact of perceptual loss vs adversarial loss on spectral recovery could be quantified more explicitly. You could elaborate on the sensitivity to β (Eq. 8) to show robustness of adversarial tuning.\n(5) While the diffusion model performs best, the reasons for GAN instability at extreme sparsity could be expanded — perhaps with visualization of discriminator loss behavior.\n(6) Please clarify whether “zero-shot reconstruction” implies no fine-tuning or full generalization to unseen geometries.", "reviewer_2": "This manuscript provides a thorough characterization of a novel approach to operator learning in fluid mechanics. The authors introduce and demonstrate an “adv-NO”, which overcomes known issues with neural operators trained with L1/L2 error with respect to error in high-frequency content. The proposed method is applied to three settings: spatio-temporal super-resolution of Schlieren visualizations, time forecasting of homogeneous isotropic turbulence, and flow reconstruction from partial observations. The authors demonstrate distinct advantages for adv-NO in the first two of these settings, with some difficulties with zero-shot flow reconstruction relative to the diffusion model used for comparison. Overall, this manuscript presents a methodical and well-reasoned study with high-impact results.\nThere are a few minor revisions we might suggest based on our understanding of this work:\n* While the overall contributions of the work and importance of the method are made clear in the introduction, it may be beneficial to more explicitly introduce readers to how adv-NO training diverges from the Lp loss typically used for NOs. Obviously, GANs are generally familiar, but it may add a bit of clarity prior to the results section.\n* The field error difference between adv-NO and NO in Table 1 is surprising–while the difference in numerical error makes sense given that the NO is learning to minimize L2 loss, qualitatively from Figure 2d NO appears far less physical. This is attributed to imperfections in the phase alignment, but it would be interesting to hear more about this.\n* For the forecasting setting, it would be helpful to state how many eddy-turnover timescales are contained within the 160-step trajectory.\n* While not necessary for this paper, I would be curious to see how the adv-NO performed relative to the other models in the forecasting task when only given a single observation as input (e.g. u(x, t) -> u(x, t + tau)). While beyond the scope, I am curious how would you expect this time embedding impacts adv-NO relative to vanilla NO.\n* The sparse measurements setting could benefit from some additional clarity in the text regarding “adv-NO” vs. “GAN with a UNet generator that is not time conditioned”. It is stated that this is due to the fact that the reconstruction is done without learning the dynamics. However, given this difference, it may be useful to add text explicitly stating how this relates to the overall narrative of the paper with respect to adv-NO.\nThis is a well-presented and considerate manuscript, with results that have the potential to shift how the fluid mechanics community views and trains surrogate models. Other than the minor suggestions for clarification noted above, I would recommend accepting this work to Nature Communications.", "reviewer_3": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts."}, "author_response": "Point-by-point responses to the reviewers’ comments on NCOMMS-25-71696\nDear Reviewers,\nWewouldliketothankyouforthecriticalreadingofmanuscriptNo.NCOMMS-25-71696and\nthe constructive suggestions for improvement of our work.\nWe believe that we have addressed all comments from the reviewers, within the intended\nscope of this work, in the revised version of this manuscript. We noted how the reviewers’\nsuggestions have been implemented into the text by highlighting it in the Authors’ Reply in\nyellow. Changes in the revised manuscript have been highlighted in red font in the marked\nrevision. Please see the responses below for itemized details.\n–2– December 31, 2025\nReviewer 1: With this manuscript, Vivek Oommen et al. benchmark a unified\nframework combining neural operators (NOs) and generative modeling (GANs,\ndiffusion models, VAEs) to address three fundamental challenges in turbulence\nmodeling: (1) spatio-temporal super-resolution, (2) forecasting with limited\ndata, and (3) sparse flow reconstruction. The authors introduce an adversar-\nially trained neural operator (adv-NO) that mitigates spectral bias and achieves\nhigh-frequency fidelity with low inference cost. The study demonstrates im-\npressive performance on Schlieren jet data, homogeneous isotropic turbulence\n(HIT), and cylinder wake reconstruction, showing large gains in spectral accu-\nracy, phase alignment, and computational efficiency.\nIn my opinion, this manuscript reports on excellent and timely work addressing a\npressing challenge in data-driven turbulence modeling - the over-smoothing and\nspectral bias in neural operator learning. The integration of generative modeling\ntechniques with neural operators is elegant and well-justified theoretically but\nnot completely novel. The authors provide clear evidence across diverse bench-\nmarks that adversarial training restores fine-scale features while maintaining the\noperator’s low inference cost. The experiments appear to be carefully designed,\nand the theoretical analysis (via Parseval-Plancherel identity and kernel Fourier\nanalysis) strengthens the contribution.\nI believe that the manuscript is well-written and technically rigorous, with sig-\nnificant breadth - spanning theoretical insight, numerical validation, and practi-\ncal applications in turbulence. Therefore, I recommend acceptance with minor\nrevisions and would like to encourage the authors to consider the following com-\nments:\nWe thank the reviewer for the constructive and thorough review of our work. We have\naddressedallyourcommentsbelowandinthesubstantiallyrevisedversionofourmanuscript.\nR1.1. The benefits of combining L1/L2 and adversarial losses are well known\nfrom pix2pix-style GANs (Isola et al., 2016). I would like to encourage the au-\nthors to acknowledge this line of research.\nAuthors’ Reply 1.1: We thank the reviewer for reminding us to acknowledge the relevant\nprior works. We have now added citations to pix2pix conditional GANs [1], SRGAN [2], and\nthe work of Johnson et al. [3], which introduced perceptual losses for image transformation\nand super-resolution tasks.\nR1.2. Although I appreciate the idea to test the sparse-data reconstruction per-\nformance with a “PTV-like” dataset, the description of the methodology left me\nconfused here. Please clarify if the random mask was set for the three velocity\ncomponents independently or had equal values (zero/one) – the latter would\nbe the case for PTV measurements. The fraction of available data was varied,\nhowever this metric’s relation to experimental conditions is difficult to assess.\n–3– December 31, 2025\nPlease consider estimating the tracer image density (in particles per pixel) based\non some reasonable assumptions. This would make your findings much more ac-\ncessible for experimentalists. You may also want to consider the recent article\nby Sciacchitano et al. (Exp. Fluids, 2025) on the topic of flow reconstruction\nbased on PTV data.\nAuthors’ Reply 1.2: We thank the reviewer for this comment. For all the test setups\n(1 through 4), the same random mask was used for each of the three velocity components to\nimitate the real PTV setup. We have clarified this point in the flow reconstruction part of\nthe Results section.\nWe agree that providing the relation between mask fraction and experimental conditions can\nbenefit turbulence experimentalists. Here, our sparse reconstruction input data is different\nthan actual PTV experimental data in a sense that we directly use the sparse velocity\ninformation randomly taken from a uniform grid of pixels (same random mask repeatedly\nused for each velocity component), while the PTV measurements are generally at irregular\ngrids and require further interpolation to map them into a uniform grid with possible further\nresolution degradation.\nHowever, it is still possible to draw an analogy to the experimental setup. We report the\nsparseactivepixeldensitysimilartotheparticledensitygenerallyreportedbyparticles/mm3,\nand mean inter pixel spacing as the average inter particle spacing consistent with experi-\nmental particle resolution.\nGiven our domain size of [-4D, 4D], [-4D, 4D], and [0, 2D], where D is the cylinder diameter,\nand the grid resolution of 256 × 256 × 32, we can calculate active pixel density depending\non the mask fraction. For example, for mask fractions of 0.99, 0.9, and 0.5 (the fraction of\npixels masked out), the active pixel density is 163.84 pixels/D3, 1638.40 pixels/D3, and 8192\npixels/D3, respectively. For active inter-pixel spacing, we implemented the nearest neighbor\nalgorithm on the randomly masked data to find the closest active pixel to another active\npixel and calculate the average over all the active pixels. Here, active pixels refer to those not\nmasked in the random process. In particular, we used cKDTree algorithm of Scipy library\nin Python to find the minimum distance between two pairs of active pixels as shown below:\nd = min ∥r −r ∥ ,\ni i j 2\nj=1,...,N\nj̸=i\nwhere N is the number of active pixels, and r is the location of the active pixel in 3D. The\naverage inter pixel spacing is then calculated by:\nN\n1 (cid:88)\n¯\nd = d\ni\nN\ni=1\nFor example, for a cylinder of D = 1 cm, and a mask fraction of 0.99, the average inter pixel\nspacing is 1.03 mm. This value is comparable to the average inter-particle spacing ranging\nfrom 0.5 to 1.58 mm used in actual PTV studies [4]. The results for different mask fractions\n–4– December 31, 2025\nare reported in the following table. We have included this discussion in the Supplemen-\ntary Information, under a new subsection - Mimicking Particle Tracking Velocimetry (PTV)\nmeasurements.\nTable 1. Average inter pixel distance calculate based on nearest neighbor algorithm for\neach mask fraction. D is the cylinder diameter.\nMask fraction Average inter pixel spacing [D]\n0.99 0.103\n0.9 0.050\n0.7 0.036\n0.5 0.032\n0.3 0.031\n0.1 0.031\nR1.3. The study focuses on moderate Re (e.g., Re = 90). A discussion or short\nλ\ntest at higher Re (even if 2D) would strengthen the claim of scalability.\nAuthors’ Reply 1.3: We thank the reviewer for this suggestion. We have now added\nan additional higher-Re forecasting case at Re = 150, at the end of the HIT forecasting\nλ λ\nsubsection.\nR1.4. The impact of perceptual loss vs adversarial loss on spectral recovery\ncould be quantified more explicitly. You could elaborate on the sensitivity to β\n(Eq. 8) to show robustness of adversarial tuning.\nAuthors’ Reply 1.4: We thank the reviewer for their suggestion. We performed an ex-\nperiment systematically varying β ∈ {1,10−1,10−2,10−3,10−4} to investigate the robustness\nof adversarial training. We have now included an additional subsection in the Discussions\nsection - Sensitivity to adversarial weighting - to discuss the same.\nR1.5. While the diffusion model performs best, the reasons for GAN instability\nat extreme sparsity could be expanded - perhaps with visualization of discrimi-\nnator loss behavior.\nAuthors’ Reply 1.5: Wethankthereviewerfortheircomment. Wehavenowexpandedthe\ndiscussion on adversarial training under extreme sparsity in the first paragraph of Discussion\nsection. We also visualize the evolution of discriminator loss (L ), generator loss (L ), and\nD G\nthe discriminator score gap (|E[D(u )]−E[D(u )]|) with epochs under moderate and\nreal fake\nhigh observation sparsity in the Supplementary Figure S10.\nR1.6. Please clarify whether “zero-shot reconstruction” implies no fine-tuning\nor full generalization to unseen geometries.\n–5– December 31, 2025\nAuthors’ Reply 1.6: We thank the reviewer for their comment. We clarified that “zero-\nshot” refers to applying the same pre-trained model without any fine-tuning to previously\nunseen masking styles/measurement layouts on new snapshots within the same cylinder-\nwake geometry. We do not claim generalization to unseen geometries. We have modified\npoint two in the contributions, the flow reconstruction case in the Results section and the\ncaption of Figure 6, to make this clarification.\n–6– December 31, 2025\nReviewer 2: This manuscript provides a thorough characterization of a novel\napproach to operator learning in fluid mechanics. The authors introduce and\ndemonstrate an “adv-NO”, which overcomes known issues with neural operators\ntrained with L1/L2 error with respect to error in high-frequency content. The\nproposed method is applied to three settings: spatio-temporal super-resolution of\nSchlieren visualizations, time forecasting of homogeneous isotropic turbulence,\nand flow reconstruction from partial observations. The authors demonstrate\ndistinct advantages for adv-NO in the first two of these settings, with some\ndifficulties with zero-shot flow reconstruction relative to the diffusion model\nused for comparison. Overall, this manuscript presents a methodical and well-\nreasoned study with high-impact results.\nThere are a few minor revisions we might suggest based on our understanding\nof this work:\nWe appreciate and thank the reviewer for positive feedback and review of our work. Below\nare responses to your comments.\nR2.1. While the overall contributions of the work and importance of the method\nare made clear in the introduction, it may be beneficial to more explicitly intro-\nduce readers to how adv-NO training diverges from the Lp loss typically used\nfor NOs. Obviously, GANs are generally familiar, but it may add a bit of clarity\nprior to the results section.\nAuthors’ Reply 2.1: We agree with the reviewer that the sudden transition to spatio-\ntemporal super-resolution in the Results section reads abrupt. We added a short subsection,\n“Overview of adversarially trained neural operators (adv-NO)”, at the beginning of the\nResults section (immediately before the spatio-temporal super-resolution results) to briefly\ncontrast standard Lp-trained NOs with adv-NO training (reconstruction + perceptual + ad-\nversarial losses and alternating updates). This provides a smoother entry to spatio-temporal\nsuper-resolution as the first case for assessing small-scale recovery.\nR2.2. The field error difference between adv-NO and NO in Table 1 is surpris-\ning - while the difference in numerical error makes sense given that the NO is\nlearning to minimize L2 loss, qualitatively from Figure 2d NO appears far less\nphysical. This is attributed to imperfections in the phase alignment, but it would\nbe interesting to hear more about this.\nAuthors’ Reply 2.2: We thank the reviewer for their comment. We investigate this-\ndiscrepancyy further by comparing the centreline profiles of the impinging jet states in Sup-\nplementary Figure S2. We have also provided a more detailed explanation in the Results\nsection about spatio-temporal super-resolution.\nR2.3. For the forecasting setting, it would be helpful to state how many eddy-turnover timescales are contained within the 160-step trajectory.\n–7– December 31, 2025\nAuthors’ Reply 2.3: Thank you for the suggestion. We now explicitly mention that the\n160-step training trajectory contains roughly 22 eddy-turnover time scales in the forecasting\npart of the Results section.\nR2.4. While not necessary for this paper, I would be curious to see how the adv-\nNO performed relative to the other models in the forecasting task when only\ngiven a single observation as input (e.g. u(x,t) → u(x,t+τ)). While beyond the\nscope, I am curious how would you expect this time embedding impacts adv-NO\nrelative to vanilla NO.\nAuthors’ Reply 2.4: We thank the reviewer for their insightful thought. In DNS, the\nstate is advanced with a numerically stable time-step ∆t set by the CFL criterion, but\nCFL\nsaving full 3D fields at every ∆t is prohibitively expensive; thus, trajectories are typically\nCFL\nstored at a much coarser interval τ ≫ ∆t . As a result, learning a “single-step” map\nCFL\nu(t) (cid:55)→ u(t + τ) from saved snapshots corresponds to a coarse-grained evolution operator\nthat compresses many fine solver steps into one update. Under such coarse-graining, the re-\nduced dynamics is generally non-Markovian and depends on history through a memory term\n(Mori-Zwanzig) [5, 6, 7], which motivates supplying a short history window to approximate\nthese effective memory effects. Accordingly, we use a 4-snapshot history [u(t−3τ),...,u(t)]\n(spanning ≈ 0.5t ) to predict [u(t+τ),...,u(t+4τ)] and then roll out autoregressively.\nE\nWe also note that the “single-observation” setting interacts strongly with forecast horizon:\none-step targets (t+τ) require many more autoregressive steps (increasing error accumula-\ntion), while very large one-shot horizons suffer from weak input-output correlation [8]. We\nobserved that a one-shot mapping across ∆t = 36τ (≈ 5t ) fails to converge due to the low\nE\ncorrelation at such large lead times. We therefore adopt a bundled horizon (4τ) as a com-\npromise, and rollout autoregressively for estimating states further into the future. We have\nadded a brief clarification in the manuscript distinguishing ∆t from τ and motivating\nCFL\nthe use of short-history inputs and bundled horizons for forecasting from sparsely sampled\nturbulence data.\nTo address the second part of the reviewer comment R2.4, we conducted an experiment\nto compare adv-NO with and without time-conditioning (TC). In particular, we find that\nTC improves long-horizon rollout stability. Adv-NO with TC maintains coherent large-scale\nevolution across multiple eddy-turnover times, whereas removing TC leads to temporal stag-\nnation of the low-wavenumber structures and the emergence of spurious high-wavenumber\nartifacts during autoregressive rollout. We have now updated the forecasting part of the\nResults section, and visualized the impact of TC on adversarial training in Supplementary\nFigure S5.\nR2.5. The sparse measurements setting could benefit from some additional clar-\nity in the text regarding “adv-NO” vs. “GAN with a UNet generator that is not\ntime conditioned”. It is stated that this is due to the fact that the reconstruction\nis done without learning the dynamics. However, given this difference, it may\nbe useful to add text explicitly stating how this relates to the overall narrative\n–8– December 31, 2025\nof the paper with respect to adv-NO.\nAuthors’ Reply 2.5: We thank the reviewer for this important point and agree that the\nsparse-measurementssettingcouldbemadeclearer. Wewanttoclarifythat, adv-NOdenotes\nan adversarially trained neural operator, whose objective is to approximate the true operator\nN , which is the mapping from function space U to V (N : U (cid:55)→ V ). For task 1 (spatio-\ntemporal super-resolution), U represents space of low-resolution low-framerate states, while\nV represents the space of high-resolution high-framerate states. For task 2 (forecasting), U\nrepresents the space of functions in the short history window and V represents that in the\nprediction window.\nIn contrast, Task 3 is an ill-posed inverse reconstruction problem from sparse measure-\nments, where the objective is to sample from the conditional probability distribution of full\nfields given partial measurements, rather than to approximate a forward evolution/super-\nresolution operator. Accordingly, we directly use conditional generative models - GAN and\ndiffusion models for Task 3. We have revised the Flow Reconstruction part of the Results\nsection (and added a few clarifying sentences in the Discussion) to connect the terminol-\nogy to the paper’s central narrative: Tasks 1 and 2 use generative-model-inspired objectives\nto refine operator learning and mitigate spectral bias, while Task 3 addresses sparse\nreconstruction as an ill-posed inverse problem by directly learning the conditional posterior\np(full flow | partial observations) with conditional generative models (GANs and diffusion\nmodels).\nRemarks on code availability: We were unable to access the code through the\nDOI provided\nAuthors’ Reply: Thank you for reviewing our manuscript. For confidential peer review,\nthe editors and referees can access the code by following the directions given in the private\nZenodo link provided below.\n• Private Zenodo link (clickable): Click here\n• Full URL (copy/paste fallback): https://zenodo.org/records/17088765?preview=\n1&token=eyJhbGciOiJIUzUxMiJ9.eyJpZCI6IjZiZDg4MGY5LThlNDEtNDE0YS05NGE0L\nWYwYzQ4MGE0YTI2MCIsImRhdGEiOnt9LCJyYW5kb20iOiJkZmE0YTNiMTk3NWI5YTkwMzJ\nkMDE4MDU0NTAxMzY5MiJ9.S1G98LA98njOa2t5vT99OOQ9IqugXIbAQK9FdnmINxcXzQ8p\nnI2LWXSy4O90OSC8HYyjaxZQG3avGCCBsnkJ9A\nTheZenodorecordhasbeenpreparedforthesoftware(reservedDOI:10.5281/zenodo.17088765)\nand will be made public after the manuscript is published.\nThis is a well-presented and considerate manuscript, with results that have the\npotential to shift how the fluid mechanics community views and trains surrogate\nmodels. Other than the minor suggestions for clarification noted above, I would\nrecommend accepting this work to Nature Communications.\n–9– December 31, 2025\nAuthors’ Reply: We sincerely thank the reviewer for the positive assessment of our\nmanuscript and for the constructive suggestions. We have carefully addressed all the com-\nments in the revised manuscript, and we indicate the corresponding changes in our point-\nby-point responses above.\n–10– December 31, 2025\nReviewer 3: I co-reviewed this manuscript with one of the reviewers who pro-\nvided the listed reports. This is part of the Nature Communications initiative\nto facilitate training in peer review and to provide appropriate recognition for\nEarly Career Researchers who co-review manuscripts.\nRemarks on code availability: We were unable to access the code through the\nDOI provided\nAuthors’ Reply: Thank you for reviewing our manuscript. For confidential peer review,\nthe editors and referees can access the code by following the directions given in the private\nZenodo link provided below.\n• Private Zenodo link (clickable): Click here\n• Full URL (copy/paste fallback): https://zenodo.org/records/17088765?preview=\n1&token=eyJhbGciOiJIUzUxMiJ9.eyJpZCI6IjZiZDg4MGY5LThlNDEtNDE0YS05NGE0L\nWYwYzQ4MGE0YTI2MCIsImRhdGEiOnt9LCJyYW5kb20iOiJkZmE0YTNiMTk3NWI5YTkwMzJ\nkMDE4MDU0NTAxMzY5MiJ9.S1G98LA98njOa2t5vT99OOQ9IqugXIbAQK9FdnmINxcXzQ8p\nnI2LWXSy4O90OSC8HYyjaxZQG3avGCCBsnkJ9A\nTheZenodorecordhasbeenpreparedforthesoftware(reservedDOI:10.5281/zenodo.17088765)\nand will be made public after the manuscript is published.\nReferences\n[1] Phillip Isola, Jun-Yan Zhu, Tinghui Zhou, and Alexei A Efros. Image-to-image trans-\nlation with conditional adversarial networks. In Proceedings of the IEEE conference on\ncomputer vision and pattern recognition, pages 1125–1134, 2017.\n[2] Christian Ledig, Lucas Theis, Ferenc Husz´ar, Jose Caballero, Andrew Cunningham, Ale-\njandro Acosta, Andrew Aitken, Alykhan Tejani, Johannes Totz, Zehan Wang, et al.\nPhoto-realistic single image super-resolution using a generative adversarial network. In\nProceedings of the IEEE conference on computer vision and pattern recognition, pages\n4681–4690, 2017.\n[3] Justin Johnson, Alexandre Alahi, and Li Fei-Fei. Perceptual losses for real-time style\ntransfer and super-resolution. In European conference on computer vision, pages 694–\n711. Springer, 2016.\n[4] Andrea Sciacchitano, Benjamin Leclaire, and Andreas Schr¨oder. On the accuracy of\ndata assimilation algorithms for dense flow field reconstructions. Experiments in Fluids,\n66(2):42, 2025.\n[5] Hazime Mori. Transport, collective motion, and brownian motion. Progress of theoretical\nphysics, 33(3):423–455, 1965.\n–11– December 31, 2025\n[6] Robert Zwanzig. Memory effects in irreversible thermodynamics. Physical Review,\n124(4):983, 1961.\n[7] Saad Qadeer, Panos Stinis, and Hui Wan. Stabilizing pde–ml coupled system. arXiv\npreprint arXiv:2506.19274, 2025.\n[8] OdedOvadia, VivekOommen, AdarKahana, AhmadPeyvan, EliTurkel, andGeorgeEm\nKarniadakis. Real-time inference and extrapolation with time-conditioned unet: Appli-\ncations in hypersonic flows, incompressible flows, and global temperature forecasting.\nComputer Methods in Applied Mechanics and Engineering, 441:117982, 2025."}, {"round": 2, "reviewer_comments": {"reviewer_1": "I have reviewed the revised manuscript and find that the authors have adequately addressed all of my previous comments. The revisions have improved the clarity and quality of the paper, and I have no further concerns. I am satisfied with the current version of the manuscript and recommend it for publication.", "reviewer_2": "I appreciate the thoughtful and well-constructed responses to reviewer feedback and congratulate the authors on an excellent manuscript.\nI believe this is very strong work, and am happy to recommend it for publication in Nature Communications.", "reviewer_3": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts."}, "author_response": "Point-by-point responses to the reviewers’ comments on NCOMMS-25-71696A\nDear Reviewers,\nWewouldliketothankyouforthecriticalreadingofmanuscriptNo.NCOMMS-25-71696Aand\nthe constructive suggestions for improvement of our work.\nWe believe that we have addressed all comments from the reviewers, within the intended\nscope of this work, in the revised version of this manuscript. We noted how the reviewers’\nsuggestions have been implemented into the text by highlighting it in the Authors’ Reply in\nyellow. Changes in the revised manuscript have been highlighted in red font in the marked\nrevision. Please see the responses below for itemized details.\n–2– February 3, 2026\nReviewer 1: I have reviewed the revised manuscript and find that the authors\nhave adequately addressed all of my previous comments. The revisions have\nimproved the clarity and quality of the paper, and I have no further concerns.\nI am satisfied with the current version of the manuscript and recommend it for\npublication.\nAuthors’ Reply: Thank you for recommending our work for publication and taking the\ntime to review our revised manuscript. We are grateful for your constructive feedback, which\nhelped improve the clarity and overall quality of the paper.\nReviewer 2: I appreciate the thoughtful and well-constructed responses to re-\nviewer feedback and congratulate the authors on an excellent manuscript. I\nbelieve this is very strong work, and am happy to recommend it for publication\nin Nature Communications. Remarks on code availability: The code appears\ncomplete and well documented for reproducing the results shown.\nAuthors’Reply: Thankyouforyourthoughtfulreviewandforrecommendingourmanuscript\nfor publication. We appreciate your time in evaluating both the manuscript and the code.\nReviewer 3: I co-reviewed this manuscript with one of the reviewers who pro-\nvided the listed reports. This is part of the Nature Communications initiative to\nfacilitate training in peer review and to provide appropriate recognition for Early\nCareer Researchers who co-review manuscripts. (Remarks on code availability):\nThe code is well documented and organized.\nAuthors’ Reply: Thank you for your note regarding co-review as part of the Nature\nCommunications initiative. We appreciate your time in reviewing the manuscript and the\ncode."}]

Decision letter