From Numerical Simulators of PDEs to Neural Emulators and Back
Felix Koehler
cs.LG
2026-08-25
A TUM thesis shows neural PDE emulators can beat the solvers that produced their training data on rollout. APEBench covers 46 PDEs; PRDP cuts solver iterations by up to 86%.
Numerical PDE solvers get expensive once they sit inside a design loop that needs thousands of evaluations. Neural emulators train on solver-generated trajectories and promise cheap rollouts. The usual story treats the solver as a teacher and the network as a student that cannot beat the training signal.
This TUM dissertation by Felix Koehler argues the two sides are closer than that story allows. Convolutional layers look like finite-difference stencils, Fourier Neural Operators look like spectral methods, UNets look like multigrid. Solver error and network error can be read in the same Fourier multipliers. Once the solver's five roles in the learning pipeline are pulled apart, three published results follow: APEBench (NeurIPS 2024), PRDP (ICLR 2025), and Neural Emulator Superiority (NeurIPS 2025).
The five roles are: the data generator that produces training trajectories, the high-fidelity reference used at evaluation, the differentiable hybrid component that gradients must flow through, the classical method the architecture structurally resembles (this one need not run at all), and the baseline on the accuracy-speed trade-off. In the default setup the first, second, and fifth collapse onto one solver, so a benchmark cannot tell whether the network learned the physics or the numerical artifacts.
APEBench generates data on the fly with Exponax, a JAX pseudo-spectral solver, across 46 semi-linear PDE configurations on periodic uniform grids, split into linear, nonlinear, and reaction-diffusion families. Five architectures are compared: ConvNet, ResNet, UNet, Dilated ResNet, and FNO. Training follows a (T, B) taxonomy: T is the unrolled horizon, B is how often the prediction is aligned to a reference. Diverted chain sets B=1 and T>1, so the model rolls out autoregressively but is pulled back to the solver at every step. Default seed counts are 50 in 1D and 20 in 2D/3D.
PRDP treats differentiable physics as bilevel optimization: the outer loop learns network weights, the inner loop is an iterative linear solver. Networks typically stall around \(10^{-2}\) nRMSE and cannot use residual tolerances of \(10^{-5}\). Training starts near one inner iteration, bumps \(K\) by \(\Delta K\) when validation plateaus, and stops when a bump does nothing, calling that \(K{\max}\). Early cheap iterations are progressive-refinement savings; stopping below \(K{\epsilon}\) is incomplete-convergence savings.
The superiority analysis derives closed-form Fourier multipliers for linear advection, diffusion, and Poisson. A linear ansatz with the same structure as an explicit scheme is fit to implicit or unconverged solver data, then compared with the analytical multiplier. The superiority ratio \(\xi^{[n]}\) is the emulator's error against a high-fidelity reference, divided by the training solver's error against the same reference. Below 1, the student has beaten the teacher.
APEBench has no universally best architecture. ResNets are almost never the single winner and almost never collapse. FNOs win on problems with global spatial structure, notably low-Reynolds Navier-Stokes Kolmogorov flow, and lose on reaction-diffusion where energy sits in high-frequency patterns beyond their active modes; local convolutions win there. In higher dimensions UNets beat dilated ResNets, likely because the coarsening hierarchy covers the receptive field more evenly than axis-aligned dilations. On 1D advection, raising difficulty \(\gamma1\) (a CFL-like number) kills local networks once the receptive field no longer covers the domain of dependence; FNOs barely feel that cut. Local convolutions converge with as few as five training trajectories; FNOs need more. The full suite costs about 900 GPU-hours on 8× RTX 2080 Ti.
| Contribution | Number | Baseline |
| PRDP cumulative iterations | up to −86% | fully converged to \(\epsilon=10^{-5}\) |
| PRDP 2D NS training time | −62% | fully converged physics |
| APEBench coverage | 46 PDEs / 5 architectures | 50 seeds in 1D, 20 in 2D/3D |
| Emulator nRMSE ceiling | \(\approx 10^{-2}\) (best \(\approx 10^{-3}\)) | solvers can reach machine precision |
| Local-conv sample need | 5 trajectories to converge | FNOs need more |
PRDP matches fully converged physics on final network accuracy, within statistical noise. Cumulative solver iterations drop by up to 86%, in one breakdown 72% incomplete-convergence savings plus 14% progressive-refinement savings. The hardest case, a correction-based neural-hybrid solver for 2D Navier-Stokes, cuts wall-clock training time by 62%. The paper's own headline also quotes a 78% training-time reduction as the upper end. Tests cover Poisson inverse problems, 1D/2D heat, Burgers, and that NS hybrid. Explicit time steppers and spectral linear solves have no iterative linear solver in the graph, so PRDP does not apply.
Superiority splits in two. Autoregressive superiority: train and test on the same initial-condition distribution, \(\xi^{[1]} \ge 1\), then \(\xi^{[n]} < 1\) later in the rollout. ConvNets, ResNets, FNOs, UNets, and Transformers almost all get there, even when trained on an explicit upwind scheme inside its stability limit. State-space superiority: train on a single Fourier mode, test on others, and \(\xi^{[1]} < 1\) already at step one. That mainly appears for local receptive fields. FNOs and Transformers do not get it under the same distribution shift. Advection and diffusion show forward superiority (test modes above the training mode); Poisson with an unconverged Richardson iterator shows backward superiority, because low frequencies are where the iterator is slowest.
Emulator accuracy typically caps near \(10^{-2}\) nRMSE, with the best cases around \(10^{-3}\). Classical solvers can hit machine precision on their own discrete system. That is not the same as hitting the continuous PDE.
A PDE surrogate evaluated only on a frozen dataset such as PDEBench collapses generator, reference, and baseline into one opaque trajectory. A better emulator can then be penalized for not copying solver artifacts. Training data need not be generated to machine precision: current models cannot use that extra bit. Differentiable-physics training need not drive every inner solve to \(10^{-5}\). Architecture choice can follow what a good numerical method for that PDE looks like, matching receptive field to domain of dependence the way a CFL condition matches a stencil.
None of this is a reason to throw the solver away. Semi-linear equations on a periodic box are far from engineering meshes. Discretization error, inductive bias, and optimizer noise do not require periodicity; periodicity is what made them separable.
The thesis is explicit about the largest constraint: semi-linear PDEs, periodic boundaries, uniform Cartesian grids. Complex geometry is out. Dirichlet and Neumann conditions are mostly absent from the main experiments. High-Reynolds flow and Euler are excluded by the pseudo-spectral setup. Full-field nRMSE of \(10^{-2}\) does not imply usable stresses, lift, or drag. Each emulator is trained for one solver and one set of discretization choices, not across PDEs or coefficients. APEBench's main architecture sweep had no Transformer; the superiority paper added one, in a different experimental suite. PRDP is validated on iterative linear solvers for sparse discrete operators; Newton-type nonlinear solvers are described as straightforward in principle and were not run. Superiority is not a guarantee. It depends on the PDE, the error structure, the architecture, the objective, and the time horizon. The closed-form proofs stop at linear ansätze on linear PDEs; the nonlinear claims are empirical.