Green-Routed Neural Operators:\\Physics Determines Where the Network Reads

Chenhao Si; Ming Yan

cs.LGOct 4, 2026

Green-Routed Neural Operators:\Physics Determines Where the Network Reads

Authors: Chenhao Si, Ming Yan

Organizations: School of Data Science The Chinese University of Hong Kong, Shenzhen Shenzhen, China

Abstract

We identify a mismatch between the physical role of transport fields in many PDEs and their usual role in neural operators: PDEs use them to select read coordinates, whereas neural operators typically treat them only as input values. We address this mismatch with the Green-Routed Neural Operator (GRNO), which uses the governing equation to determine where latent features are sampled. A parameter-free equation adapter evaluates the diagnostic relation and constructs a departure map whose values are the read coordinates. A multiscale encoder-decoder combines centered and routed reads of latent features to learn the complete finite-time update. Across five two- and three-dimensional PDE systems, GRNO achieves the lowest mean final relative $L^2$ error on four under 40-step autoregressive evaluation and remains competitive on Keller-Segel. Fixed-weight route interventions reveal strong dependence on direction and spatial alignment in four systems, with weak dependence in Keller-Segel. In independently trained ablations, GRNO achieves lower mean errors than variants that supply the transport field only as an input feature, substitute a learned displacement for the equation-specified route, or apply the route with a spatial misalignment, across all five systems. It also substantially outperforms directly advecting the physical state and learning the remaining update, indicating that equation-specified read coordinates provide an effective structural prior for long-horizon PDE forecasting.

Figures & tables

Panel of Figure 1 — Figure 1: Equation-specified routing and GRNO. (A) The current state $q^{n}$ determines a diagnostic field $\psi[q^{n}]$ and a diagnostic carrier $\bm{v}[q^{n}]$ . Backtracing the carrier gives the departure map $\bm{\Xi}[q^{n}]$ , which selects the routed physical sample $q^{n}\!\circ\!\bm{\Xi}[q^{n}]$ ; $\mathcal{C}[q^{n}]$ is the remainder containing all other finite-time effects. GRNO explicitly realizes the carrier-to-coordinate-to-sample chain by using the carrier to construct read coordinates. (B) The parameter-free equation adapter evaluates the known diagnostic relation, supplies the carrier as an encoder feature, and uses a backward trace to construct the displacement $\bm{d}^{n}$ and departure map $\bm{\Xi}^{n}$ . (C) A local multiscale encoder produces features $s_{\ell}$ . At each scale, GRNO combines the centered read $s_{\ell}(\bm{x})$ with the routed read $r_{\ell}(\bm{x})=\mathcal{I}_{\ell}[s_{\ell}](\bm{x}+\bm{d}_{\ell}(\bm{x}))$ . The decoder fuses the two reads, predicts the complete state increment, and recomputes the diagnostic field, diagnostic carrier, and departure map from every predicted state during autoregressive rollout.

Model	Kolmogorov	SQG	Vlasov–Poisson	Keller–Segel	Porous convection
ReViT [ 20 ]	$27.84\pm 1.18$	$3.25\pm 0.03$	$4.86\pm 0.90$	$1.23\pm 0.02$	$19.35\pm 1.27$
Flowers [ 7 ]	$14.54\pm 5.39$	$13.12\pm 0.50$	$0.30\pm 0.01$	$\mathbf{0.32\pm 0.09}$	$16.51\pm 2.71$
NIPS-AR [ 17 ]	$134.25\pm 12.96$	$51.65\pm 0.06$	$19.33\pm 6.10$	$1.46\pm 0.39$	$69.33\pm 32.39$
RieszNO [ 21 ]	$39.01\pm 0.48$	$8.07\pm 0.63$	$0.49\pm 0.05$	$0.34\pm 0.01$	$34.84\pm 1.77$
GRNO (ours)	$\mathbf{7.80\pm 0.56}$	$\mathbf{2.34\pm 0.24}$	$\mathbf{0.16\pm 0.01}$	$0.37\pm 0.02$	$\mathbf{3.34\pm 0.16}$

Table 1: Long-horizon predictive accuracy. Final relative

L^{2}

error (

\times 10^{-2}

) at

H=40

. Entries are mean

\pm

standard deviation over three training seeds. Lower is better; bold marks the lowest mean.

Panel of Figure 2 — Figure 2: Long-horizon GRNO accuracy and system-dependent route use. a , Relative $L^{2}$ error over 40-step free rollouts. GRNO gives the lowest error on four systems and remains competitive on Keller–Segel. b , Reference fields, GRNO predictions and absolute errors at $h=40$ for Kolmogorov flow and porous convection; pointwise errors remain localized around the recovered structures. c , Ratios of final error under the zero, reverse, and shift interventions relative to the correct route, with network weights fixed. In four systems, intervention errors are 5.0–124 times those obtained with the correct route; in Keller–Segel, the ratios remain near one. d , The error differences appear early and persist throughout the rollout. e , Space–time concentration along the selected porous-convection transect for the reference, Carrier feature only , and GRNO. In this selected window, Carrier feature only develops temporal banding, whereas GRNO follows the advancing front.

System	Correct	Zero	Reverse	Shift
Kolmogorov	$\mathbf{7.80\pm 0.56}$	$210.18\pm 0.17$	$196.34\pm 0.21$	$130.67\pm 0.63$
SQG	$\mathbf{2.34\pm 0.24}$	$11.77\pm 0.36$	$24.72\pm 0.27$	$27.58\pm 0.31$
Vlasov–Poisson	$\mathbf{0.16\pm 0.01}$	$6.71\pm 0.09$	$19.75\pm 0.21$	$20.31\pm 0.06$
Keller–Segel	$\mathbf{0.37\pm 0.02}$	$0.38\pm 0.02$	$0.39\pm 0.01$	$0.38\pm 0.01$
Porous convection	$\mathbf{3.34\pm 0.16}$	$41.08\pm 0.91$	$47.76\pm 1.17$	$44.42\pm 1.21$

Table 2: Fixed-weight route interventions. Final relative

L^{2}

error (

\times 10^{-2}

) at

H=40

. Entries are mean

\pm

sample standard deviation across the three seed means.

Variant	Kolmogorov	SQG	Vlasov–Poisson	Keller–Segel	Porous convection
Centered only	8.1659	4.2917	0.3158	0.3563	6.5352
Carrier feature only	9.3106	3.7392	0.2849	0.3788	6.2905
Route only	9.0460	3.3361	0.2256	0.3710	4.8331
Full (GRNO)	7.8023	2.3354	0.1635	0.3718	3.3426
Learned route	8.7553	2.8387	0.1782	0.3892	3.4775
Shifted equation-specified route	9.5599	3.3667	0.1913	0.3791	4.0076

Table 3: Independently trained ablations and route controls. Mean final relative

L^{2}

error (

\times 10^{-2}

) at

H=40

, for seed 0 and 12 prescribed test windows per system. The first four rows form the

2\times 2

factorial. Bold marks the lowest mean among evaluated variants.

Appendix figures & tables19 assets

Supplementary material from the paper’s appendix.

Appendix

Question	Control	Evidence
Does the complete model predict well?	Independently trained models share data, loss, update budget, checkpoint rule, and evaluator.	Main Table 1
Does the trained model use the route?	Keep weights and carrier features fixed; change only the departure displacement.	Figure 2 and Appendix I
Does equation-specified routing improve prediction after training?	Independently train models with different carrier features, coordinate constructions, and routing targets.	Main Table 3 and Appendix I

Appendix

Table 4: Claim–evidence map. This table separates evidence for predictive utility from evidence for the proposed routing mechanism.

System	Evolved state $q$	Diagnostic relation	Diagnostic carrier	Boundary rule
Kolmogorov	$\omega$	$\Delta\psi=\omega-\bar{\omega}$	$\bm{u}=\nabla^{\perp}\psi$	periodic $x,y$
SQG	$\theta$	$\psi=-(-\Delta)^{-1/2}\theta$	$\bm{u}=\nabla^{\perp}\psi$	periodic $x,y$
Vlasov–Poisson	$f(x,v)$	$\rho=\int f\,\mathrm{d}v$ , $\partial_{x}E=1-\rho$	$(a_{x},a_{v})=(v,-E)$	periodic $x$ ; zero outflow in $v$
Keller–Segel	$\rho$	$(-\Delta+\lambda^{2})c=\rho-\bar{\rho}$	$\bm{v}=\chi\nabla c$	periodic $x,y,z$
Porous convection	$C$	$\bm{u}=-\nabla p+C\bm{e}_{z}$ , $\nabla\!\cdot\bm{u}=0$	Darcy velocity $\bm{u}$	periodic $x,y$ ; bounded $z$

Appendix

Table 5: Diagnostic carriers used in the experiments. Each row lists the evolved state, diagnostic relation, diagnostic carrier that determines the read coordinates, and boundary rule.

System	Grid	Train	Validation	Test	Stored $\Delta t$
Kolmogorov	$128^{2}$	$64\times 51$	$8\times 101$	$12\times 101$	$0.10$
SQG	$128^{2}$	$64\times 51$	$8\times 101$	$12\times 101$	$0.10$
Vlasov–Poisson	$128_{x}\times 128_{v}$	$64\times 101$	$8\times 121$	$12\times 161$	$0.10$
Keller–Segel	$64^{3}$	$64\times 51$	$8\times 101$	$12\times 101$	$0.25$
Porous convection	$64^{3}$	$64\times 101$	$8\times 101$	$12\times 101$	$0.10$

Appendix

Table 6: Dataset contract. Shapes list trajectories

\times

stored frames; every frame also contains the grid shown in the second column. Splits are fixed before training.

Coordinate type	Coordinate handling	Sample outside domain
Periodic	wrap modulo domain length	wrapped value
Bounded wall	clamp to the nearest grid index	value at the nearest boundary grid point
Phase-space outflow	mask each out-of-range index	zero contribution

Appendix

Table 7: Boundary behavior of routed interpolation. The rule is applied identically at every latent scale.

Quantity	2D / phase space	3D
Input channels	5	6
Base width	48	24
Stage widths	$(48,96,144,192)$	$(24,48,72,96)$
Encoder depths	$(2,2,3,4)$	$(2,2,3,4)$
Decoder blocks per scale	1	1
Local depthwise kernel	$5\times 5$	axial $5$ -tap

Appendix

Table 8: Dimension-specific GRNO configuration. “Axial” in 3D denotes the average of three depthwise kernels, one per spatial axis.

Setting	2D / phase space	3D
Epochs	100	100
Windows per epoch	2,048	512
Micro-batch / accumulation	$4/1$	$1/4$
Effective batch	4	4
Optimizer updates per epoch	512	128
Training rollout $K$	4	4

Appendix

Table 9: Matched training and evaluation protocol. Four 3D micro-batches are accumulated to obtain an effective batch size of four.

Model	Retained interaction mechanism	Protocol adaptation
ReViT	canonical local frames and rotationally equivariant transformer blocks	Physical-feature adapter constructs the vector/scalar layout from each task; one-frame residual decoder.
Flowers	learned multihead coordinate warps and latent resampling	The same state, context, boundary sampler, and autoregressive evaluator are used; read coordinates are learned from forecast data.
NIPS-AR	channel-independent Fourier kernel and linear attention	Released spatial block placed in a recurrent one-frame wrapper; reported as an autoregressive adaptation.
RieszNO / RNO-3D	Riesz directional representation plus spectral mixing	Released 2D construction; direct axial 3D extension is explicitly marked as a dimensional adaptation.

Appendix

Table 10: Baseline roles and adaptations. For dimensional or autoregressive adaptations, we retain each published model’s interaction mechanism and modify only the input and output interfaces required by the common one-frame protocol. The table summarizes the evaluated implementations.

Task	Additional records
Kolmogorov / SQG	correlation, spectral error, variance error, and kinetic-energy or scalar-variance error
Vlasov–Poisson	total-mass drift, negative-density fraction, charge-density error, and electric-energy error
Keller–Segel	mass drift, negative-mass fraction, peak and variance error, and screened free-energy error
Porous convection	anomaly error, plume-mask IoU, top dissolution-flux error, plume-depth error, scalar-dissipation error, variance error, and fractions below zero or above one

Appendix

Table 11: Task-specific physical metrics recorded by the common evaluator. These metrics are reported in addition to relative

L^{2}

error.

Mode	Operation on displacement $\bm{d}$	Question tested
Correct	$d$ from the equation adapter	reference computation
Zero	$d\leftarrow 0$	Is the carrier used as a coordinate rather than only as an input feature?
Reverse	$d\leftarrow-d$	Does route orientation matter?
Shift	deterministic spatial roll of $d$	Does alignment between route and current state matter?

Appendix

Table 12: Fixed-weight route interventions. Network weights, construction of carrier features, centered branches, and prescribed test windows are held fixed.

Figure 3: Long-horizon qualitative comparison on Kolmogorov flow. The first row shows the reference vorticity field at ten horizons spanning the 40-step free rollout. The second and third rows show the pointwise absolute errors of GRNO and Flowers, respectively. The error scale is shared by both models and all displayed horizons.

Figure 4: Long-horizon qualitative comparison on dissipative SQG. The first row shows the reference scalar field at ten horizons spanning the 40-step free rollout. The second and third rows show the pointwise absolute errors of GRNO and ReViT, respectively. The error scale is shared by both models and all displayed horizons.

Figure 5: Long-horizon qualitative comparison on Vlasov–Poisson. The first row shows the reference phase-space distribution at ten horizons spanning the 40-step free rollout. The second and third rows show the pointwise absolute errors of GRNO and Flowers, respectively. The error scale is shared by both models and all displayed horizons.

Figure 6: Long-horizon qualitative comparison on Keller–Segel. The first row shows the reference density field on the fixed visualization slice at ten horizons spanning the 40-step free rollout. The second and third rows show the pointwise absolute errors of GRNO and Flowers, respectively. The error scale is shared by both models and all displayed horizons.

Figure 7: Long-horizon qualitative comparison on porous convection. The first row shows the reference concentration field on the fixed visualization slice at ten horizons spanning the 40-step free rollout. The second and third rows show the pointwise absolute errors of GRNO and Flowers, respectively. The error scale is shared by both models and all displayed horizons.

Figure 8: Route interventions on Kolmogorov flow at $h=40$ . The top row shows the reference field and predictions using the correct, zero, reversed, and shifted routes. The bottom row shows the corresponding pointwise absolute errors. All predictions use the same fitted GRNO checkpoint; only the departure route is changed.

Figure 9: Route interventions on dissipative SQG at $h=40$ . The top row shows the reference field and predictions using the correct, zero, reversed, and shifted routes. The bottom row shows the corresponding pointwise absolute errors. All predictions use the same fitted GRNO checkpoint; only the departure route is changed.

Figure 10: Route interventions on Vlasov–Poisson at $h=40$ . The top row shows the reference phase-space distribution and predictions under the Correct , Zero , Reverse , and Shift interventions. The bottom row shows the corresponding pointwise absolute errors. All predictions use the same fitted GRNO checkpoint; only the departure route is changed.

Figure 11: Route interventions on Keller–Segel at $h=40$ . The top row shows the reference field and predictions on the fixed visualization slice using the correct, zero, reversed, and shifted routes. The bottom row shows the corresponding pointwise absolute errors. All predictions use the same fitted GRNO checkpoint; only the departure route is changed.

Figure 12: Route interventions on porous convection at $h=40$ . The top row shows the reference concentration field and predictions on the fixed visualization slice using the correct, zero, reversed, and shifted routes. The bottom row shows the corresponding pointwise absolute errors. All predictions use the same fitted GRNO checkpoint; only the departure route is changed.

Explore similar work

May 25, 2026cs.LG

Autoregression-Free Neural Operators for Time-Dependent PDEs

Neural operators learn mappings from function-dependent inputs to solutions, providing an effective framework for solving partial differential equations (PDEs). For time-dependent PDEs, existing methods typically perform long-horizon prediction through autoregressive rollout directly in high-dimensional physical field spaces, where each predicted state is recursively fed back as the input for the next step. Although effective for short-term prediction, this autoregressive rollout and the lack of continuous-time modeling lead to progressive error accumulation over long-horizon rollouts. In this work, we propose Autoregression-Free Neural Operators (AFNO), which map the time evolution of PDEs into a latent space and model continuous-time vector fields within it. AFNO uses flow matching to learn the latent vector field, thereby enabling continuous evolution over extended horizons, avoiding autoregressive rollout and capturing dynamics under varying parameter configurations through explicit conditioning on physical parameters. Theoretical analysis and extensive experiments on six PDEs demonstrate that AFNO improves long-horizon prediction stability and consistently reduces rollout errors compared with the baselines.

Jiaquan Zhang, Caiyan Qin, Haoyu Bian +7

School of Computer Science and Engineering, University of Electronic Science and Technology of China, Chengdu, Sichuan 611731, China · School of Robotics and Advanced Manufacture, Harbin Institute of Technology, Shenzhen, Guangdong 518055, China · School of Mathematical Sciences, Capital Normal University, Beijing 100048, China +3

Neural Operators Neural Partial Differential Equation Solvers Partial Differential Equations Long-Horizon Forecasting Autoregressive Rollout

Ask PDF

May 10, 2026cs.LG

HS-FNO: History-Space Fourier Neural Operator for Non-Markovian Partial Differential Equations

Neural operators provide fast surrogate models for time-dependent partial differential equations, but their standard autoregressive use usually assumes that the instantaneous field

u(t,\cdot)

is a complete state. This assumption fails for delay equations, distributed-memory systems, and other non-Markovian dynamics: two trajectories may agree at time

t

and nevertheless have different futures because their histories differ. We introduce the History-Space Fourier Neural Operator (HS-FNO), a neural operator for delay and memory-driven PDEs formulated on the lifted state

u_t(θ,x)=u(t+θ,x)

θ\in[-τ,0]

. The key computational step is to decompose one history-state update into a learned predictor for the newly exposed future slice and an exact shift-append transport for the portion of the history window already known from the previous state. This avoids learning deterministic history coordinates, reduces the learned output dimension, and enforces the natural discrete history update. We test HS-FNO on five benchmark families covering delayed reaction--diffusion, spatial epidemiology, nonlocal neural-field dynamics, delayed waves, and distributed-memory closures. Across ten random seeds, HS-FNO attains the lowest aggregate one-step, history-space, and rollout errors among the principal baselines. The largest gain occurs in autoregressive prediction, where aggregate rollout error decreases from

0.241

0.188

, and

0.185

for current-state, lag-stack, and unconstrained history-to-history operators, respectively, to

0.094

. The same model uses fewer parameters than unconstrained history prediction. These results indicate that enforcing the discrete shift structure of history-state evolution is an effective inductive bias for non-Markovian PDE surrogate modeling.

Lennon J. Shikhman

College of Computing, Georgia Institute of Technology Atlanta, GA, United States

Fourier Neural Operators Neural Operators Partial Differential Equations Surrogate Models Neural Fields History

Ask PDF

Jul 31, 2026cs.AI

Geometry-aware Incremental Neural Operator for Long-Horizon PDE prediction

Neural operators have shown strong potential for learning solution operators of partial differential equations (PDEs). However, long-horizon autoregressive prediction remains challenging: local errors accumulate as spectral inconsistency, phase misalignment, or mean drift. Existing methods mainly improve state representations and operator backbones, while leaving the repeatedly applied latent transition increment weakly structured, allowing spectral errors and unstable channel couplings to accumulate during rollout. To address these issues, we propose a geometry-aware incremental neural operator (GeoIncNO) for stable long-horizon PDE prediction. GeoIncNO predicts latent increments for residual advancement and uses lightweight low-rank projectors to regulate channel coupling within active frequency bands derived from the increment spectral energy distribution. To reduce physical-space reconstruction errors, GeoIncNO further introduces a mean--fluctuation decoupled reconstruction mechanism, where stable mean structures and dynamic fluctuations are fused separately, and phase correction is applied only to the zero-mean fluctuation component. Extensive experiments on six PDE benchmarks, covering 1D, 2D, and 3D dynamical systems, show that GeoIncNO achieves consistently strong prediction accuracy, improved rollout stability, and better spectral fidelity compared with competitive neural-operator baselines.

Jiaquan Zhang, Shuxu Chen, Haifan Meng +6

School of Information and Software Engineering, UESTC · Electronics and Information Convergence Engineering, KHU · School of Computer Science and Engineering, UESTC +4

Neural Partial Differential Equation Solvers Long-Horizon Forecasting Partial Differential Equations Fluctuation-Dissipation Theorem Geometry-Aware Long-Horizon Task Planning

Ask PDF