Neural surrogates for computational fluid dynamics (CFD) are trained and ranked on datasets that each vary one design family in one flow domain. Such datasets cannot test the least that an aerodynamic foundation model must do: predict CFD-quality incompressible flow in seconds for new shapes in new enclosures. On DrivAerML, larger models and longer training barely lower Transolver's error; doubling the training cases cuts it by about 12% up to 400 cases, exhausting the dataset. We introduce AeroZoo, a benchmark and training corpus of 3,815 time-mean solutions to the unsteady Reynolds-averaged Navier--Stokes (RANS) equations with a transition model, across seven geometry families, eight enclosure configurations and Reynolds numbers from $10^3$ to $10^8$. Each sample poses a well-posed PDE problem in the form a solver receives it: the boundary meshes grouped by condition type, the prescribed values and the fluid properties. The 2,816 samples with a stable time mean form the canonical benchmark: two tasks, engineering metrics and 22 frozen splits, disjoint in geometry, that hold out geometry families, enclosures, boundary-condition layouts and Reynolds ranges. On the random split, boundary-condition tokens cut GeoTransolver's median velocity error from 76.9% to 8.4% of the inlet speed in single-seed runs (uniform inlet flow: 11.0%). ISLA, a transformer for boundary-driven PDEs, reaches 5.3% over four seeds, about 4 times the change from a finer mesh. On the same test samples, a held-out geometry family barely changes ISLA's error relative to uniform flow, and a held-out enclosure or boundary-condition layout raises it 2.4 to 9.2 times. The data and the scoring code will be released.
A foundational principle of physics is that all coordinate systems are equally valid: the behavior of a PDE-governed system does not depend on the origin or orientation of the frame, nor on the units or mesh used to write it down, and neither should a surrogate model. Yet the neural surrogates with the lowest errors on today's aerodynamics benchmarks lack these invariances and fail under simple transformations of the same problem: for models trained on PDEs of cars in one pose (i.e., translation and rotation) and from one mesher, rotating the full PDE raises the surface-pressure error from 4% to more than 100%, and re-sampling its boundary raises it to about 30%. We introduce ISLA (Invariant Slice Attention), a slice-attention transformer for steady, boundary-driven PDEs in which no learned operation receives a coordinate and every sum over the boundary is a quadrature. Positions enter only as distances and relative angles to learned slice anchors, lengths are measured in a unit computed from the boundary, and interior points query the boundary through a tree-accelerated learned kernel without influencing one another. With every boundary cell a token and an exact boundary read, ISLA is equivariant to rigid motions and reflections, invariant to the units of its declared data, and consistent under mesh refinement, while at every setting each query is independent of all other queries. On AeroZoo, 2,816 flows that vary in shape, size, speed, and boundary conditions, its mean reference-scaled velocity and pressure errors are at least 3.1 and 3.3 times lower than any baseline's. On AirfRANS its errors are second only to those of a boundary-integral model (GLOBE) and at least 1.6 times lower than those of Transolver and GeoTransolver.
We introduce GLOBE, a nonlinear neural surrogate for boundary-driven, homogeneous, weakly nonlinear PDEs that draws inductive bias from boundary-element methods and equivariant ML, representing solutions as superpositions of learnable Green's-function-like kernels evaluated from boundary faces to targets. The architecture is translation-, rotation-, and parity-equivariant; discretization-invariant; and units-invariant via rigorous nondimensionalization. On AirFRANS (steady incompressible RANS over NACA airfoils), a compact ~117k-parameter GLOBE model has lower mean-squared error than published values for the strongest reported ML surrogates on both the interpolation and low-data splits, in single-run comparisons under the released AirFRANS normalization. The model evaluates fields at arbitrary points, handles non-watertight meshes, and generalizes under Reynolds number and angle-of-attack distribution shifts. We further introduce a hierarchical Barnes-Hut acceleration that exploits the kernel's far-field decay to reduce evaluation cost from quadratic to near-linear in the number of boundary elements, with a tunable opening angle that recovers exact dense evaluation in the zero-angle limit. Quantitative validation focuses on 2D AirFRANS; on 3D DrivAerML cars, GLOBE predicts surface and volume fields on the full meshes, and we compare it with a matched Transolver baseline. Together, these results show that boundary-element-inspired inductive biases can substantially improve the accuracy of ML surrogates for boundary-driven PDEs while scaling to large 3D meshes.