Paper deep dive
Cluster Attention Neural Operators for Solving Parametric Partial Differential Equations
Ming Zhong, Antonio Colanera, Gianluigi Rozza, Zhenya Yan
Intelligence
Status: not_run | Model: - | Prompt: - | Confidence: 0%
Entities (0)
Relation Signals (0)
No relation signals yet.
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Traditional simulations of parametric partial differential equations (PDEs) rely on repetitive computations for each parameter, which makes high-fidelity design impractical. Neural operators address this issue by learning solution operators, accelerating parameter-space mapping by orders of magnitude. Recent Transformer-based neural operators attempt to capture global dependencies, but often at the cost of quadratic attention complexity. Transolver resolves this problem by projecting physical states into a reduced slice space for attention computation. Although fast, this projection sacrifices fine spatial information. Moreover, by operating in this reduced space with shared weights across attention heads, it may constrain the model's flexibility, thereby limiting its capacity to capture complex phenomena. To address these issues, we propose the Cluster Attention Neural Operator (CANO), which reformulates attention via a novel cross-attention mechanism that dynamically clusters queries while preserving full-resolution keys and values. This avoids slice compression loss and removes weight-sharing limits. At the same time, the model remains fast without losing global interactions. Empirically, CANO achieves state-of-the-art performance across canonical PDE benchmarks, covering fluid and solid dynamics (e.g., Navier-Stokes, Airfoil, Plasticity), irregular unstructured geometries (e.g., Pipe Turbulence, Composites), and long-term temporal rollouts. Across solid deformation and turbulent flow benchmarks, CANO achieves lower errors than baselines and exhibits strong geometric adaptability and temporal consistency.
Tags
Links
- Source: https://arxiv.org/abs/2609.39914v1
- Canonical: https://arxiv.org/abs/2609.39914v1
Trouble viewing inline? Open PDF directly →
Full Text
72,709 characters extracted from source content.
Expand or collapse full text
Cluster Attention Neural Operators for Solving Parametric Partial Differential Equations Ming Zhong Affiliation: School of Advanced Interdisciplinary Sciences, University of Chinese Academy of Sciences, Beijing 100049, China Affiliation: State Key Laboratory of Mathematical Sciences, Academy of Mathematics and Systems Science, Chinese Academy of Sciences, Beijing 100190, China Affiliation: International School for Advanced Studies (SISSA), Trieste 34136, Italy Antonio Colanera Affiliation: International School for Advanced Studies (SISSA), Trieste 34136, Italy Gianluigi Rozza Affiliation: International School for Advanced Studies (SISSA), Trieste 34136, Italy Zhenya Yan Email: zyyan@mmrc.iss.ac.cn (Z. Yan) Affiliation: School of Advanced Interdisciplinary Sciences, University of Chinese Academy of Sciences, Beijing 100049, China Affiliation: State Key Laboratory of Mathematical Sciences, Academy of Mathematics and Systems Science, Chinese Academy of Sciences, Beijing 100190, China Affiliation: School of Mathematics and Information Science, Zhongyuan University of Technology, Zhengzhou 450007, China Abstract Traditional simulations of parametric partial differential equations (PDEs) rely on repetitive computations for each parameter, which makes high-fidelity design impractical. Neural operators address this issue by learning solution operators, accelerating parameter-space mapping by orders of magnitude. Recent Transformer-based neural operators attempt to capture global dependencies, but often at the cost of quadratic attention complexity. Transolver resolves this problem by projecting physical states into a reduced slice space for attention computation. Although fast, this projection sacrifices fine spatial information. Moreover, by operating in this reduced space with shared weights across attention heads, it may constrain the model’s flexibility, thereby limiting its capacity to capture complex phenomena. To address these issues, we propose the Cluster Attention Neural Operator (CANO), which reformulates attention via a novel cross-attention mechanism that dynamically clusters queries while preserving full-resolution keys and values. This avoids slice compression loss and removes weight-sharing limits. At the same time, the model remains fast without losing global interactions. Empirically, CANO achieves state-of-the-art performance across canonical PDE benchmarks, covering fluid and solid dynamics (e.g., Navier-Stokes, Airfoil, Plasticity), irregular unstructured geometries (e.g., Pipe Turbulence, Composites), and long-term temporal rollouts. Across solid deformation and turbulent flow benchmarks, CANO achieves lower errors than baselines and exhibits strong geometric adaptability and temporal consistency. keywordsParametric partial differential equation, Deep learning, Cluster attention neural operator, Transformers 1 Introduction Partial differential equations (PDEs) play a pivotal role in modeling physical phenomena and have became ubiquitous across science and engineering [2, 3, 4, 5, 6]. While these PDEs provide a theoretical framework for understanding natural phenomena, yet solving these equations numerically is often too slow for practical use [7, 8, 9, 10, 11]. Traditional numerical solvers, such as the finite element, finite difference, and spectral methods, typically rely on dense spatial discretizations [7, 8, 9, 10, 11]. However, the computational cost of standard methods grows steeply with grid resolution. As a result, high-fidelity simulations are often computationally infeasible. This leaves a significant gap between theoretical models and real-world engineering [12, 13]. Recently, scientific machine learning (SciML) has emerged as a promising alternative, replacing classical iterative solvers with trained neural networks [14, 15, 16, 17]. For example, the Deep Ritz method [18] solves PDEs via variational formulations, while Physics-Informed Neural Networks (PINNs) [19, 20, 21] penalize PDE residuals in the loss function. Other variants, such as PeRCNN [22], enforce physical laws directly through network architectures. However, these methods only solve one PDE instance at a time. Once the physical parameters or boundary conditions changed, the network must be re-trained from scratch, making it impractical for real-time applications[23, 24]. To address this issues, neural operators (NOs) [25, 26, 27, 28, 29] learn mappings between parameter or coefficient fields to solution fields. Instead of solving a single PDE instance, they approximate the underlying operator. Once trained, they provide zero-shot solution predictions for unseen configurations at a fraction of the classical computational cost. Existing neural operator architectures can be broadly categorized into MLP-based, spectral-based, and transformer-based methods. First, MLP-based neural operators, represented by Deep Operator Networks (DeepONet) [27], are fundamentally grounded in the universal approximation theorem for operators [30]. This framework has been extended through various architectures, including POD-DeepONet [31], PINN-DeepONet [32], MIONet [33], and Hybrid-preconditioned DeepONet [34]. In practice, however, these simple architectures struggle to resolve complex nonlinear behaviors. Second, spectral-based neural operators, such as the Fourier Neural Operator (FNO) [26, 25, 35, 36], exploit global convolutions in the wavenumbers domain to parameterize integral kernels, achieving high accuracy on structured grids. However, accommodating general unstructured domains significantly increases computational cost or requires lossy interpolations [35, 36]. Third, Transformer-based neural operators [37, 38, 40, 39, 41, 42] offer an alternative way to handle unstructured domains and irregular meshes. A standard self-attention block computes the output Y as [43, 44, 45] =Softmax(⊤d),Y=Softmax ( QK d )V, (1) where ,,∈ℝN×dQ,K,V ^N× d are the query, key, and value matrices, respectively, obtained from the projection of original features (such as function values or coordinates), where N is the number of tokens and d is the feature dimension. While this formulation captures long-range dependencies, it scales quadratically, (N2d)O(N^2d), making it impractical for fine physical discretizations. Linear attention variants reduce complexity to (Nd2)O(Nd^2) through matrix factorization [46, 47, 48, 49, 50]: =ϕ()(ϕ()⊤),Y=φ(Q) (φ(K) V ), (2) where ϕ(⋅)φ(·) denotes a feature mapping. Despite its efficiency, linear attention often suffers from overly smooth attention maps. While hierarchical or convolutional designs [51, 52, 53, 54, 55] mitigate this issue, they are strictly tailored for structured grids and do not readily apply to unstructured meshes. Another line of work tackles complex geometries by introducing latent tokens [56, 58, 57]. A notable example is Transolver [58], which learns a weight matrix to aggregate input features into a compact ”slice space”, performs attention within this reduced representation, and projects back. However, compressing physical points into slices loses local geometric details. In addition, Transolver computes attention only in this slice space and shares weights across attention heads, which limits its ability to capture complex dynamics. To overcome these challenges, we introduce the Cluster Attention Neural Operator (CANO)11 1 https://github.com/zhongming12342/CANO. Rather than projecting points into a reduced slice space, CANO computes cross-attention by clustering queries (Q) dynamically, while retaining keys (K) and values (V) in their original coordinate space. By clustering queries dynamically, CANO avoids compression loss in the slice space and removes the need for weight sharing across attention heads. As a result, the model maintains global receptive fields while scaling linearly with the number of points N, yielding strong generalization on complex physical systems. The main contributions of this work are twofold. First, we introduce Cluster Attention, which clusters queries while keeping keys and values at full resolution. This formulation prevents detail loss and scales linearly with the number of points N. Second, we validate CANO on 12 diverse PDE benchmarks, including complex geometries and long rollout horizons, where it achieves up to 41.67% lower error on solid mechanics and 33.73% on turbulent flows compared to previous models. The remainder of this paper is organized as follows. Section 2 briefly introduces the neural operator framework and related works. Section 3 presents the proposed CANO architecture. Section 4 reports extensive numerical experimental results on diverse benchmarks and model analysis, followed by concluding remarks in Section 5. 2 Related works In scientific machine learning, neural operators aim to learn the mappings between parameter spaces and corresponding PDE solution spaces [26, 27, 28, 29]. Rather than replacing classical solvers, these models act as fast surrogates that accelerate inference. From a mathematical perspective, we consider a generalized PDE formulated over a bounded spatial domain Ω⊂ℝdx ^d_x. Consider a physical system parameterized by an input field a(x)a(x), which can represent coefficients, initial/boundary conditions, or source terms. The governing equation can be written as (Lau)(x)=f(x),x∈Ω⊂ℝdx,u(x)=u0(x),x∈∂Ω, \ aligned (L_au)(x)=f(x),& x∈ ^d_x,\\ u(x)=u_0(x),& x∈∂ , aligned . (3) where f(x)f(x) denotes the source term. Let A and U be Banach spaces containing the input parameter field a and the target solution u:Ω→ℝu: , respectively. Given a parametric differential operator ℒa:→∗L_a:U ^*, the goal of operator learning is to approximate the underlying solution operator †:→G :A such that †(a)=uG (a)=u. Given a training dataset of input–solution pairs (a(i),u(i))i=1M\(a^(i),u^(i))\_i=1^M, we approximate †G using a parameterized neural operator θ:→G_θ:A with weights θ∈ℝpθ ^p. The optimal parameters θ∗θ^* are obtained by minimizing the empirical mean squared error in the Bochner norm [28]. In practice, the risk is approximated by the empirical risk over the training samples: minθ∈ℝp1M∑i=1M‖u(i)−θ(a(i))‖2. _θ ^p 1M _i=1^M \|u^(i)-G_θ(a^(i)) \|_U^2. (4) Within this framework, neural operator architectures such as the Graph Kernel Network (GKN) [62, 63] are closely motivated by Green’s functions. For a linear differential operator ℒaL_a with Green’s function G(x,y)G(x,y) satisfying ℒaG(x,⋅)=δxL_aG(x,·)= _x, the solution can be written as an integral convolution: u(x)=∫ΩG(x,y)f(y)y.u(x)= _ G(x,y)f(y)\,dy. (5) Analogously, iterative neural operators update a continuous hidden representation vt(x)v_t(x) through non-local kernel integration: vt+1(x)=σ(Wvt(x)+∫Ωκθ(x,y,a(x),a(y))vt(y)νx(y)).v_t+1(x)=σ (Wv_t(x)+ _ _θ (x,y,a(x),a(y) )v_t(y)\, _x(dy) ). (6) Here, an initial lifting layer P maps the local input (x,a(x))(x,a(x)) to the hidden state v0(x)=P(x,a(x))v_0(x)=P(x,a(x)). Each update layer combines a local linear transformation parameterized by weight matrix W with global kernel aggregation via the parameterized kernel κθ _θ, followed by an activation function σ and integration measure νx _x. After T iterative updates, a projection layer Q decodes the final hidden state to obtain the solution field u(x)=Q(vT(x))u(x)=Q(v_T(x)). This formulation allows the operator to capture both local features and long-range dependencies across the domain. To realize the operator learning framework described above, three of the most representative architectures are: the Fourier Neural Operator (FNO) [26, 25], the Deep Operator Network (DeepONet) [27], and Transformer-based neural operators [37, 38, 58]. The FNO [26, 25] accelerates non-local integration by assuming a translation-invariant kernel, κ(x,y)=κ(x−y)κ(x,y)=κ(x-y). Under this assumption, the integral operator reduces to a spatial convolution, which can be evaluated efficiently in the frequency domain via the convolution theorem: vt+1(x)=σ(ℱ−1(κ^θ⋅ℱ(vt))(x)+Wvt(x)),v_t+1(x)=σ (F^-1 ( κ_θ·F(v_t) )(x)+Wv_t(x) ), (7) where ℱF and ℱ−1F^-1 denote the forward and inverse Fourier transforms, and κ^θ κ_θ is a learnable Fourier multiplier. Truncating the spectrum to the lowest kmaxk_max modes allows FNO to run in (NlogN)O(N N) time using the Fast Fourier Transform (FFT). However, because standard FFTs require uniform Cartesian grids and periodic boundaries, standard FNO is difficult to apply directly to complex, unstructured meshes. An alternative approach is DeepONet [27, 33], which approximates the mapping G†G based on the universal approximation theorem for operators [30]. DeepONet decomposes the solution into a dot product of a branch network and a trunk network: Gθ(a)(x)=∑k=1mbk(a)tk(x),G_θ(a)(x)= _k=1^mb_k(a)t_k(x), (8) where bk(a)b_k(a) are outputs of the branch network encoding the input parameter function a, and tk(x)t_k(x) are outputs of the trunk network encoding the spatial coordinates x, m is a pre-specified parameter. While DeepONet offers significant mathematical generality, its performance is often constrained by the learning dynamics of its underlying multi-layer perceptrons (MLPs). Specifically, these architectures exhibit spectral bias [64], meaning they efficiently learn low-frequency macroscopic trends but struggle to resolve the sharp, high-frequency oscillations that are critical for accurate PDE modeling. To handle unstructured domains and capture complex long-range dependencies, attention mechanisms [43] have been increasingly adopted for PDE operators. Given a set of latent representations ∈ℝN×dX ^N× d, standard transformer blocks project these into query, key, and value matrices: =Q,=K,=V,Q=XW_Q, =XW_K, =XW_V, (9) where Q,K,V∈ℝd×dW_Q,W_K,W_V ^d× d. The standard Softmax attention computes the output as displayed in Eq. (1). While this formulation effectively models global pairwise interactions, it requires computing and storing the full N×N× N similarity matrix, with a quadratic computational cost of (N2d)O(N^2d). Linear attention methods [46, 47, 48, 49] reduce this complexity via a non-negative feature mapping function ϕ(⋅)φ(·) to decompose the kernel. By leveraging the associativity of matrix multiplication, the operation is rewritten as displayed in Eq. (2). This rearrangement evaluates the d×d× d context matrix first, avoiding the N×N× N attention matrix and reducing the complexity to (Nd2)O(Nd^2). However, this factorization leads to low-rank, smooth attention maps, which limits the model’s ability to capture high-frequency details. To avoid both the quadratic cost of softmax attention and the over-smoothing of linear attention, recent models compress physical tokens into a smaller set of latent representations. For example, Transolver [58] uses a physics-aware slicing mechanism to project the mesh into a reduced number of slice tokens. Specifically, given the input physical node features ∈ℝN×dX ^N× d, Transolver learns a normalized assignment matrix ∈ℝN×MA ^N× M to partition the domain into M slices, where M≪NM N. The model aggregates the high-dimensional node features into a compact ”slice space” via a projection: =⊤,∈ℝM×d.S=A X, ^M× d. (10) Standard Softmax attention is then executed exclusively on these M slice tokens to capture global interactions efficiently: ′=Softmax((Q)(K)⊤d)(V),S =Softmax ( (SW_Q)(SW_K) d )(SW_V), (11) where Q,K,V∈ℝd×dW_Q,W_K,W_V ^d× d are the learnable projection matrices. Note that in the implementation of Transolver22 2 https://github.com/thuml/Transolver/blob/main/PDE-Solving-StandardBenchmark/model/Physics_Attention.py, these matrices are defined as ∈ℝdhead×dhead ^d_head× d_head, forcing weight sharing across different attention heads. Here dheadd_head denotes the embedding dimension per head, which is typically calculated as d/Hd/H with H being the number of attention heads. Finally, a reverse projection (deslicing) is applied to broadcast the updated slice features back to the original N physical nodes: =′,∈ℝN×d.Y=AS , ^N× d. (12) Although this slice-and-broadcast design achieves (NMd)O(NMd) linear complexity, it has notable limitations. Specifically, constructing all queries, keys, and values from compressed slices discards local high-frequency details. In addition, Transolver shares attention weights across heads, which prevents individual heads from learning distinct physical patterns. To resolve these issues, our cluster attention mechanism preserves full-resolution keys and values, clustering only the queries dynamically. This allows unconstrained attention heads to capture physical interactions while maintaining linear (N)O(N) complexity. 3 Cluster Attention Neural Operator (CANO) As discussed in the previous sections, existing latent compression methods like Transolver cut fine-grained physical structures by projecting the input X into a highly compressed slice space. In such symmetric compression process, all three matrices including Q, K, and V, suffer from slice compression, which limits the model’s performance on complex problems. To resolve this issue, we propose the Cluster Attention Neural Operator (CANO). The key innovation of CANO is the asymmetric cross-attention mechanism. Instead of aggregating the raw input X, we selectively cluster only the query matrix Q, while preserving the keys K and values V at their full, uncompressed physical resolution. This achieves linear complexity (N)O(N) while retaining global context across the domain. Specifically, given the input representation ∈ℝN×dX ^N× d, we first project it into the standard query, key, and value spaces as done in Eq. (9). Unlike Transolver that learns a slice matrix from X, CANO derives a dynamic cluster assignment matrix ∈ℝN×KcC ^N× K_c from the query matrix Q, where KcK_c is the predefined number of clusters (Kc≪NK_c N). This matrix maps the N physical nodes to KcK_c cluster centroids, and it’s computed as =Softmax((Q)c),C=Softmax ((XW_Q)W_c ), (13) where the Softmax operator is applied across the cluster dimension, and c∈ℝd×KcW_c ^d× K_c is the learnable cluster projection matrix. Using C, we aggregate the query matrix into a compact set of clustered queries cQ_c: c=⊤,c∈ℝKc×d.Q_c=C Q, _c ^K_c× d. (14) With the clustered queries established, CANO executes a efficient cross-attention mechanism. The compact cQ_c attends to the full, uncompressed key matrix K to extract relevant features from the value matrix V: c=Softmax(c⊤d),c∈ℝKc×d.Z_c=Softmax ( Q_cK d )V, _c ^K_c× d. (15) This asymmetric operation is fundamentally different from Transolver. By avoiding the compression of K and V, the attention matrix inherently retains physical details. Finally, to reconstruct the updated representations for each individual physical node, we utilize the exact same assignment matrix C to dispatch (broadcast) the aggregated cluster features cZ_c back to the original spatial resolution: =c,∈ℝN×d.Y=CZ_c, ^N× d. (16) A comparison of different attention mechanisms (including vanilla self-attention, Transolver, and CANO) can be found in Figure 1. Figure 1: Comparison of the architecture of (a) Vanilla Self-Attention, (b) Transolver, and (c) CANO mechanisms. From a computational perspective, the calculation of the cross-attention matrix c⊤Q_cK and the subsequent value aggregation require (NKcd)O(NK_cd) operations. Since Kc≪NK_c N, CANO elegantly achieves linear complexity with respect to the sequence length N. Furthermore, because its mechanism, CANO inherently supports independent multi-head attention computations without the need to share weights across heads, avoiding the detail loss seen in Transolver. The theoretical connection between CANO and generalized linear attention [65] can be seen from Appendix A. The detailed algorithm flow is given in the algorithm 1, where FFN(⋅)FFN(·) and LayerNorm(⋅)LayerNorm(·) denote the feed-forward network and layer normalization, respectively. Algorithm 1 Cluster Attention Neural Operator (CANO) for solving PDEs Spatial mesh coordinates ∈ℝN×dx ^N× d_x , input function values a()∈ℝN×da(x) ^N× d_a , ground truth ugt()∈ℝN×du_gt(x) ^N× d_u Predicted solution field u^()∈ℝN×du u(x) ^N× d_u Initialize network parameters W , base learning rate α0 _0 , total epochs TepochsT_epochs ; epoch=1epoch=1 to TepochsT_epochsInput encoding: (0)←Encoder([,a()])∈ℝN×dh^(0) ([x,\,a(x)]) ^N× d ; Cluster Attention Blocks: =0l=0 to L−1L-1~(l)←(l)+Cluster-Attention(LayerNorm((l))) h^(l) ^(l)+Cluster-Attention(LayerNorm(H^(l))) ; (l+1)←~(l)+FFN(LayerNorm(~(l)))h^(l+1)← h^(l)+FFN(LayerNorm( h^(l))) ; Output decoding: ^()←Decoder((L))∈ℝN×du u(x) (h^(L)) ^N× d_u ; Loss Computation: ℒtrain←‖u^()−ugt()‖L2‖ugt()‖L2L_train← \| u(x)-u_gt(x)\|_L^2\|u_gt(x)\|_L^2 ; Backpropagate to compute parameter gradients ∇ℒtrain _WL_train ; Adjust current learning rate: αepoch←Schedule(epoch,α0,Tepochs) _epoch (epoch, _0,T_epochs) ; Update parameters via AdamW: ←AdamW(,∇ℒtrain,αepoch)W (W, _WL_train, _epoch) ; u^() u(x) 4 Numerical Experiment Results In this section, we provide a comprehensive evaluation of the proposed CANO from four aspects. First, we compare CANO with state-of-the-art baselines on standard PDE benchmarks [26, 35]. We then assess its ability to handle complex geometries through experiments on irregular-domain benchmarks [66]. Next, we investigate its temporal stability in challenging long-horizon autoregressive rollout settings [26, 73]. Finally, we analyze the impact of key hyperparameters on the performance of CANO. To validate the performances of our CANO approach, we compare it against different neural PDE solvers. This includes graph neural networks (e.g., GraphSAGE [67]), different neural operators, ranging from well known architectures like FNO [26] and DeepONet [27] to their variants (POD-DeepONet [31], Geo-FNO [35], F-FNO [68], U-FNO [69], and U-NO [70]), as well as recent innovations like LSM [71], NORM [66], and HPM [72]. Furthermore, we evaluate against state-of-the-art Transformer-based operators, including Galerkin [37], OFormer [38], FactFormer [39], GNOT [40], ONO [41], HT-Net [53], and Transolver [58]. As the strongest baseline, Transolver serves as our primary point of comparison throughout the experiments. For most baseline methods, we report the performance values provided in the corresponding original publications or in the comprehensive benchmark of the Transolver study. To enable a controlled comparison of different attention mechanisms, we additionally re-implemented the Galerkin Transformer and Transolver within the same training framework adopted for CANO. Table 1 outlines the experimental settings and hyperparameter choices, including cluster numbers (KcK_c) and transformer configurations, used for each benchmark. These settings have been selected empirically to provide a suitable trade-off between predictive accuracy and computational cost. For Transolver, the MLP expansion ratio (rmlpr_mlp in Table 1) was set to 2 for most benchmarks. 4.1 Standard PDE Benchmarks To ensure a fair comparison against all baseline methods, we have conducted experiments on six standard PDE benchmarks: Airfoil, Pipe, Plasticity, Navier-Stokes, Darcy flow, and Elasticity, as detailed in Table 2. Table 1: Unified training and model hyperparameters across benchmarks. All experiments utilize an initial learning rate of 10−310^-3. The model configuration is denoted as (L/H/d/rmlpr_mlp), representing layer, number of attention head, embedding dimension and MLP ratio in the transformer blocks. Benchmark Loss Function Epochs Batch Cluster Number L/H/d/rmlpr_mlp Standard benchmarks Airfoil Rel. L2L^2 500 4 32 8 / 8 / 128 / 2 Pipe Rel. L2L^2 500 4 16 8 / 4 / 128 / 2 Plasticity Rel. L2L^2 500 8 64 8 / 8 / 128 / 2 Navier-Stokes Rel. L2L^2 500 2 32 8 / 8 / 256 / 1 Darcy Rel. L2L^2 + 0.1L∇L_∇ 500 4 128 8 / 8 / 128 / 2 Elasticity Rel. L2L^2 500 1 128 8 / 8 / 128 / 2 Irregular domain benchmarks Irregular Darcy Rel. L2L^2 2000 16 128 4 / 4 / 64 /4 Pipe Turbulence Rel. L2L^2 2000 16 32 4 / 4 / 64 /4 Composite Rel. L2L^2 2000 16 32 4 / 4 / 64 / 4 Long-time rollout Navier-Stokes Rel. L2L^2 500 2 32 8 / 8 / 256 / 1 ICP Plasma Rel. L2L^2 1000 1 32 8 / 8 / 256 / 1 Heat Flow Rel. L2L^2 1000 1 16 8 / 8 / 256 / 1 Table 2: Overview of standard PDE benchmarks. Details on geometry type, spatial dimension, discretization size, dataset splits, and input features. Test Case Geometry Dim. Resolution Dataset Split Input Type Train Test Airfoil [35] Structured Mesh 2 221×51221× 51 1000 200 Boundary Shape Pipe [35] Structured Mesh 2 129×129129× 129 1000 200 Domain Shape Plasticity [35] Structured Mesh 2+1 101×31101× 31 900 80 External Force Navier-Stokes [26] Regular Grid 2+1 64×6464× 64 1000 200 Previous Vorticity Darcy [26] Regular Grid 2 85×8585× 85 1000 200 Coefficients Elasticity [35] Point Cloud 2 972 1000 200 Domain Shape This benchmark includes a diverse range of physical systems (e.g., fluid dynamics and solid mechanics), geometry types (regular grids, structured meshes, and point clouds), and input functions. The underlying physical equations and dataset details are introduced as in Appendix B.1. The quantitative results reported in Table 3 demonstrate that CANO consistently outperforms the considered baselines across all six standard PDE benchmarks. Table 3: Performance comparison of neural operators on standard benchmarks (CANO vs. baselines). All values represent the relative L2L^2 error. The best results are highlighted in bold, and the second-best are underlined. A slash (“/”) denotes benchmarks where the baseline is not applicable. Models marked with ∗ are reimplemented by us for fair comparison, while other results are taken from the original papers or the Transolver study. Model Structured Mesh Regular Grid Point Cloud Airfoil Pipe Plasticity Navier-Stokes Darcy Elasticity (×10−2)(× 10^-2) (×10−2)(× 10^-2) (×10−2)(× 10^-2) (×10−1)(× 10^-1) (×10−2)(× 10^-2) (×10−2)(× 10^-2) Classic models FNO [26] / / / 1.561.56 1.081.08 / DeepONet [27] 3.853.85 0.970.97 1.351.35 2.972.97 5.885.88 9.659.65 U-FNO [69] 2.692.69 0.560.56 0.390.39 2.232.23 1.831.83 2.392.39 Geo-FNO [35] 1.381.38 0.670.67 0.740.74 1.561.56 1.081.08 2.292.29 U-NO [70] 0.780.78 1.001.00 0.340.34 1.711.71 1.131.13 2.582.58 F-FNO [68] 0.780.78 0.700.70 0.470.47 2.322.32 0.770.77 2.632.63 LSM [71] 0.590.59 0.500.50 0.250.25 1.541.54 0.650.65 2.182.18 Transformer-based models Galerkin∗ [37] 0.640.64 0.440.44 0.170.17 1.041.04 1.971.97 0.840.84 HT-Net [53] 0.650.65 0.590.59 3.333.33 1.851.85 0.790.79 / OFormer [38] 1.831.83 1.681.68 0.170.17 1.711.71 1.241.24 1.831.83 GNOT [40] 0.760.76 0.470.47 3.363.36 1.381.38 1.051.05 0.860.86 FactFormer [39] 0.710.71 0.600.60 3.123.12 1.211.21 1.091.09 / ONO [41] 0.610.61 0.520.52 0.480.48 1.191.19 0.760.76 1.181.18 Transolver∗ [58] 0.51¯ 0.51 0.42¯ 0.42 0.12¯ 0.12 0.94¯ 0.94 0.52¯ 0.52 0.56¯ 0.56 CANO (ours) 0.430.43 0.330.33 0.070.07 0.660.66 0.410.41 0.460.46 Improvement 15.69% 21.43% 41.67% 29.79% 21.15% 17.86% Compared to Transolver, CANO consistently achieves lower relative L2L^2 errors across all cases. The largest improvements occur on the (2+1)-D Plasticity dataset, where the error drops by 41.67%, and on turbulent Navier–Stokes, with a 29.79% reduction. These gains confirm that preserving full-resolution keys and values helps resolve both large solid deformations and turbulent flow dynamics. Figure 2 provides a complementary view of the comparison between CANO and Transolver by showing the distribution of the relative L2L^2 error across test samples. Figure 2: Performance and distribution comparison: Relative L2L^2 error of CANO and Transolver on six standard PDE datasets. CANO exhibits lower median errors on all six datasets, consistently with the results reported in Table 3. The difference is particularly marked for Plasticity, where the two distributions are clearly separated. In addition, CANO generally shows narrower interquartile ranges and more compact error distributions, indicating improved consistency across different test cases. Figure 3 shows representative CANO predictions and absolute errors for three problems. On both the structured Airfoil mesh and the point-cloud Elasticity dataset, CANO captures the key spatial features with lower prediction errors. The Navier-Stokes example further shows that the model preserves the main flow structures over the autoregressive rollout up to timestep +10+10, with limited error accumulation. Figure 3: CANO predictions on representative PDE benchmarks. Steady-state physical fields and absolute errors for (a) Airfoil (structured mesh) and (b) Elasticity (point cloud). CANO accurately captures sharp gradients near complex boundaries. (c) Long-term temporal rollout predictions (from timestep +1 to +10) on the Navier-Stokes dataset. CANO effectively preserves intricate high-frequency eddies and structural integrity over time without severe error accumulation. To further examine the behavior of the proposed attention mechanism, Figure 4 compares the spatial organization learned by CANO and Transolver on the Airfoil benchmark. In particular, we visualize the latent assignment matrices extracted from the final layer of the first attention head: the cluster matrix ∈ℝN×KcC ^N× K_c for CANO and the slice matrix ∈ℝN×MA ^N× M for Transolver. In CANO, the assignment matrix C determines how query matrix Q are mapped to query clusters QcQ_c, while in Transolver, the slice matrix A determines the probability that the original token X belongs to the slice token S. As shown in Figure 4(b), these cluster assignments are highly localized, with several clusters concentrating around regions with strong flow gradients—specifically near the shockwave. This confirms that CANO adapts its latent clusters to physical features rather than simply partitioning the domain by geometric proximity. In contrast, the slice assignments in Transolver (Figure 4(c)) are smoother, showing little localized adaptation around flow discontinuities. While both methods capture the overall flow field (Figure 4(a)), CANO yields a significantly lower prediction errors. Figure 4: Visualization of the spatial assignments extracted from the final layer of the first head on the Airfoil dataset. (a) The corresponding ground truth, predictions and absolute errors. (b) The spatial activation maps of CANO’s cluster matrix C, showing allocations around discontinuities (e.g., shockwaves). (c) The slice allocation matrix of Transolver A, displaying relatively smooth regions. For neural operators, a low global error (e.g., relative L2L^2 norm) does not necessarily imply physical fidelity in the predicted fields. To evaluate spectral accuracy across spatial scales, we compute the turbulent kinetic energy (TKE) spectrum of the predicted voticity fields. To compute the one-dimensional (1D) kinetic energy spatial spectrum E(k)E(k) from a two-dimensional (2D) predicted vorticity field w()w(x) of the Navier-Stokes benchmark, we first obtain its Fourier coefficients in the wavenumber space: w^()=ℱw(), w(k)=F\w(x)\, (17) where =(kx,ky)k=(k_x,k_y) is the wavenumber vector and k=||=kx2+ky2k=|k|= k_x^2+k_y^2. By applying the incompressibility condition, the 2D spectral kinetic energy density is derived by scaling the vorticity power spectrum by the inverse square of the wavenumber (excluding the zero-frequency mean flow where ||=0|k|=0): E()=12|w^()|2||2.E(k)= 12 | w(k)|^2|k|^2. (18) Finally, the 1D radial energy spectrum is calculated by isotropically summing the discrete 2D energy density over concentric annuli of radius k and bin width Δk k: ()=∑k≤||<k+ΔkE(),E(k)= _k≤|k|<k+ kE(k), (19) which effectively quantifies the distribution of turbulent kinetic energy across varying spatial scales. Based on this formulation, Figure 5 visualizes the energy spectra E(k)E(k) of the ground truth and the CANO predicted fields, for the Navier-Stokes benchmarch, averaged over 200 test samples. Figure 5: Comparison of kinetic energy spectra (19) between CANO predictions and ground truth for the Navier-Stokes benchmark. Solid lines show the mean spectrum over 200 test samples; shaded regions indicate the mean ± 3 standard deviation. The spectral distribution predicted by CANO (dashed red line) matches well with the reference solution (solid blue line) across the entire wavenumber domain k, with the shaded regions (mean ± 3 standard deviation) largely overlapping. Even at small scales where the energy drops significantly, the spectrum predicted by CANO aligns well with the reference data. Table 4: Comparison of model parameters between Transolver and CANO across different standard cases. Benchmark Number of Parameters Transolver CANO (Ours) Airfoil 3.07M 2.44M Pipe 3.07M 2.49M Plasticity 3.11M 2.51M Navier-Stokes 12.28M 8.12M Darcy 3.09M 2.42M Elasticity 0.98M 1.33M To evaluate model efficiency, we compare the parameter counts of CANO and Transolver across all six examples (Table 4). Overall, CANO achieves higher accuracy while maintaining a more compact model size. In five out of the six benchmarks, CANO requires fewer parameters than Transolver. This difference is most noticeable on Navier–Stokes, where CANO reduces the parameter count from 12.28M to 8.12M (a 33.9% reduction). Although CANO uses slightly more parameters on the Elasticity dataset, the overall results show that it is more accurate while using fewer parameters. 4.2 Irregular domain benchmarks We next evaluate CANO on irregular geometries, testing its ability to handle complex physical domains. Table 5 summarizes the three irregular domain PDE examples [66] utilized in this section: Irregular Darcy, Pipe Turbulence, and Composite. Table 5: Overview of irregular domain PDE benchmarks. Details on geometry type, spatial dimension, node discretization, dataset splits, and input features. Test Case Geometry Dim. Nodes (N) Dataset Split Input Type Train Test Irregular Darcy [66] Point Cloud 2 2290 1000 200 Coefficients Pipe Turbulence [66] Point Cloud 2 2673 300 100 Previous State Composite [66] Point Cloud 3 8232 400 100 Temperature Unlike the previous benchmarks that mostly use regular grids or structured meshes, all geometries in this task are represented entirely as 2D or 3D point clouds, posing significant challenges. We investigate two 2D fluid problems: the Irregular Darcy flow (2,290 nodes) driven by varying the permeability coefficient field, and the Pipe Turbulence (2,673 nodes) requiring state predictions based on previous one. Furthermore, we test our evaluation to a challenging 3D solid mechanics task using the Composite material benchmark (8,232 nodes) driven by temperature profiles. These benchmarks are useful to validate the model’s capability to process irregular geometries. The underlying physical equations and dataset details are introduced as in Appendix B.2. A quantitative comparison with other baselines across these three cases is presented in Table 6. Table 6: Performance comparison of neural operators on Irregular Darcy, Pipe Turbulence, and Composite benchmarks (CANO vs. baselines). All values represent the relative L2L^2 error (×10−4)(× 10^-4). The best results are highlighted in bold, and the second-best are underlined. Model Irregular Darcy Pipe Turbulence Composite GraphSAGE [67] 6.73 23.60 20.90 DeepONet [27] 1.36 9.36 1.88 POD-DeepONet [31] 1.30 2.59 1.44 NORM [66] 1.05 1.01 1.00 HPM [72] 0.74 0.83 0.93 Transolver∗ [58] 0.89 1.01 1.61 CANO (ours) 0.71 0.55 0.89 Improvement 4.05% 33.73% 4.30% CANO consistently achieves the lowest relative L2L^2 error, outperforming baseline models ranging from foundational networks (GraphSAGE, DeepONet) to recent operators (NORM, HPM, and Transolver). Specifically, on the Pipe Turbulence dataset, CANO yields a 33.73% error reduction compared to the second-best baseline. Furthermore, while Transolver performs competitively on standard benchmarks, its error increases on these unstructured point clouds, particularly in the 3D Composite task. This observation suggests that the slicing mechanism struggles to effectively compress highly unstructured and non-uniform spatial information. In contrast, CANO’s cross-attention mechanism keeps keys and values at their full spatial resolution, maintaining stable predictive accuracy across different geometric discretizations. Figure 6 presents both the error distributions and the spatial distributions of the predictions. Figure 6: Performance and qualitative distribution on irregular domains. (a) Raincloud plots comparing the relative L2L^2 error distributions of CANO and Transolver. (b)-(d) Representative visual predictions of CANO for the (b) Irregular Darcy, (c) Pipe Turbulence, and (d) Composite benchmarks. Panel (a) displays raincloud plots comparing CANO and Transolver across all three datasets. CANO shifts the overall error distribution downward and suppresses extreme outliers in Irregular Darcy and Pipe Turbulence datasets. While some outlier errors occur in Composite dataset, the overall average error of CANO is still lower than that of Transolver. Panels (b–d) show qualitative predictions for the 2D irregular Darcy flow, Pipe Turbulence, and 3D Composite material, respectively. The predicted fields match the ground truth closely, and the corresponding absolute error maps confirm that CANO maintains low errors, which highlights its effectiveness on irregular physical geometries. 4.3 Long-time rollout benchmarks In practical physical simulations, neural operators are frequently required to extrapolate system dynamics far beyond the short temporal horizons observed during training. The challenge in such cases is to avoid severe error accumulation, numerical dissipation, and the eventual collapse of physical structures over time. The autoregressive stability and temporal generalization capabilities of our proposed model are evaluated through the three time-dependent datasets shown in Table 7: Navier-Stokes [26], Inductively Coupled Plasma (ICP) and Heat Transfer [73]. Specifically, the Navier-Stokes dataset models the evolution of chaotic fluid vorticity (ω) on a regular 64×6464× 64 grid. In contrast, the ICP Plasma and Heat Transfer datasets involve highly nonlinear multi-physics phenomena defined on point clouds, predicting complex state variables such as electron density, electron temperature, and velocity fields. Table 7: Overview of 2D long-time rollout benchmarks. Details on geometry type, node discretization, training/test window, dataset splits, and input features. Models are trained on short temporal windows (Tout=5T_out=5) but evaluated on significantly longer autoregressive rollouts (up to 15 or 50 steps) to test temporal stability. Test Case Geometry Nodes Train Window Test Window (ToutT_out) Data Split Input Type TinT_in ToutT_out Short Long Train Test Navier-Stokes [26] Regular Grid 64×6464× 64 5 5 5 15 1000 200 ω ICP Plasma [73] Point Cloud 3424 5 5 5 50 499 150 ne,Te,v,Tn_e,T_e,v,T Heat Transfer [73] Point Cloud 4887 (Avg.) 5 5 5 50 740 132 u,vu,v We use random-window training to enhance the model’s temporal stability across rollouts. During the training phase, all models can only see the historical of Tin=5T_in=5 steps, learning to predict the subsequent Tout=5T_out=5 steps. However, during the inference phase, the models are evaluated on two distinct scenarios: a short rollout matching the training horizon (Tout=5T_out=5), and a significantly long autoregressive rollout, extrapolating up to Tout=15T_out=15 steps for the Navier-Stokes dataset, and up to Tout=50T_out=50 steps for both the ICP Plasma and Heat Transfer datasets. Table 8 shows that CANO consistently outperforms Transolver throughout the autoregressive rollout. Although the relative L2L^2 error increases with the prediction horizon for both models, CANO maintains lower errors across all datasets and time horizons, indicating better long-term predictive stability. Table 8: Long-term predictive performance across different datasets. We evaluate the relative L2L^2 error (×10−1)(× 10^-1) at different future time steps (e.g., +5, +15, +50). The best results are highlighted in bold. Model Navier-Stokes ICP Plasma Heat Transfer +5 +15 +5 +50 +5 +50 Transolver 1.21 3.22 0.29 2.11 0.70 4.87 CANO 0.59 2.57 0.23 1.61 0.58 4.19 Figure 7 presents the corresponding error evolution and statistical distributions. The left column displays the relative L2L^2 error accumulation per step, where the errors for both models increases as the simulation advances. The middle and right columns show the error distributions at the short (Tout=5T_out=5) and long-time horizons. While the error distributions for both models are clustered at Tout=5T_out=5, the spread of the distribution increases for both models at the long-time horizons. In these extended rollouts, Transolver shows a more pronounced upper tail in the error distribution than CANO. Figure 7: Temporal error evolution and distribution analysis for long time rollouts. Panels (a), (b), and (c) correspond to the Navier-Stokes, ICP Plasma, and Heat Transfer datasets, respectively. Left: Per step relative L2L^2 error accumulation curves, showing CANO’s slower error growth rate. Middle and Right: Raincloud plots comparing the error distributions at the short (Tout=5T_out=5) and long (Tout=15T_out=15 or 5050) rollout stages. Figures 8 and 9 show the spatial distribution of the absolute error during the autoregressive rollout for the ICP Plasma dataset, considering the electron temperature TeT_e, and the Heat Transfer dataset, considering the velocity magnitude u2+v2 u^2+v^2. At the first several prediction steps, CANO and Transolver exhibit comparable error distributions. As the rollout progresses, however, the errors produced by Transolver become increasingly pronounced across the domain, while CANO maintains consistently lower error levels over time. Figure 8: Temporal error evolution on the ICP Plasma dataset. The visualization displays the absolute error of the predicted electron temperature TeT_e for Transolver and CANO at timesteps +1, +25, and +50 compared to the ground truth. Figure 9: Temporal error evolution on the Heat Transfer dataset. The visualization displays the absolute error of the predicted velocity magnitude (u2+v2 u^2+v^2) for Transolver and CANO at timesteps +1, +25, and +50 compared to the ground truth. Table 9: Relative L2L^2 error comparison across different training sample sizes. All values for Darcy Flow are scaled by 10−210^-2, and Navier-Stokes by 10−110^-1. Best results are highlighted in bold. Problem Training Number Transolver CANO (ours) Darcy Flow 200 2.172.17 1.801.80 400 0.980.98 0.720.72 600 0.710.71 0.540.54 800 0.670.67 0.460.46 1000 0.520.52 0.410.41 Navier-Stokes 200 3.763.76 3.093.09 400 3.143.14 1.961.96 600 2.872.87 1.011.01 800 2.492.49 0.900.90 1000 0.940.94 0.660.66 4.4 Model Analysis In this section, we analyze how CANO scales with training sample size and evaluate its sensitivity to the number of clusters (KcK_c). Training Sample Efficiency Table 9 evaluates the relative L2L^2 error on the Darcy Flow and Navier-Stokes datasets when the number of training samples is scaled from 200 to 1000. Both evaluated models exhibit a monotonic decrease in error as the data scale expands. Under identically restricted data regimes, CANO consistently maintains lower error bounds compared to the Transolver baseline. Notably, in the most data-scarce scenario (200 samples), CANO limits the error to 1.80×10−21.80× 10^-2 on Darcy Flow and 3.09×10−13.09× 10^-1 on Navier-Stokes, yielding a better performance than that of Transolver (2.17×10−22.17× 10^-2 and 3.76×10−13.76× 10^-1, respectively). Table 10: Relative L2L^2 errors of CANO with different cluster numbers (KcK_c) across various PDE datasets. The best results for each dataset are highlighted in bold. Cluster Numeber (KcK_c) Darcy Navier-Stokes Pipe Airfoil Elasticity Plasticity 8 5.385.38 9.059.05 3.833.83 4.964.96 5.085.08 0.880.88 16 4.834.83 7.677.67 3.323.32 4.774.77 4.784.78 0.780.78 32 4.734.73 6.586.58 3.553.55 4.294.29 4.724.72 0.810.81 64 4.364.36 6.366.36 3.483.48 4.884.88 4.814.81 0.710.71 128 4.084.08 5.605.60 3.823.82 4.384.38 4.624.62 0.770.77 Sensitivity to Cluster Number Table 10 shows the effect of the cluster number KcK_c across six datasets, with KcK_c varied from 8 to 128. The optimal choice of KcK_c depends on the specific physical system. For complex or turbulent flows (e.g., Navier–Stokes and Darcy Flow), the error systematically decreases with larger KcK_c, reaching the lowest value at Kc=128K_c=128. In contrast, for systems governed by localized dynamics (Pipe, Airfoil, and Plasticity), the lowest relative L2L^2 errors occur at moderate configurations (Kc=16K_c=16, 3232, and 6464, respectively). This indicates that aligning KcK_c with the spatial complexity of the PDE helps prevent representational redundancy. 5 Conclusions and Discussion In this work, we introduced the Cluster Attention Neural Operator (CANO) to address the information loss caused by latent space compression. CANO relies on an asymmetric cross-attention mechanism that clusters the query space while keeping keys and at full spatial resolution. This formulation allows the model to preserve fine-scale information while maintaining linear (N)O(N) computational complexity. Experiments across 12 PDE benchmarks, including regular and irregular geometries as well as long-horizon autoregressive prediction, show that CANO consistently improves predictive accuracy over the considered state-of-the-art baselines. Several directions remain open for future work. In the present work, the number of clusters KcK_c is fixed as a priori. An adaptive strategy that adjusts the cluster number according to solution complexity could provide greater performance. In addition, although CANO retains linear complexity, optimized implementations of the clustering and cross-attention operations will be important for extending the method to large-scale three-dimensional applications [59, 60, 61]. The results suggest that asymmetric cluster-based attention provides an effective alternative to symmetric latent compression for neural operator learning, particularly when fine-scale spatial information and long-term predictive stability are important. Overall, CANO offers a framework for scientific machine learning, enabling accurate and efficient modeling of complex physical systems. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Data availability The data that support the findings of this study are available within the article. Acknowledgments This work was partially supported by the National Key R&D Program of China (No. 2024YFA1013101), the National Natural Science Foundation of China (No. 12471242), and the China Scholarship Council. Appendix A Theoretical Connection to Linear Attention To motivate CANO mathematically, we analyze existing efficient neural operators through the unified framework of generalized linear attention [65]. To decouple the quadratic Softmax computation, the generalized linear attention framework employs independent feature mappings ϕ(⋅)φ(·) and ψ(⋅)ψ(·), and further incorporates an intermediate operation G to process the projected latent representations. The attention mechanism can thus be reformulated as [65]: =ϕ()∘G∘(ψ()⊤).Y=φ(Q) G (ψ(K) V ). (20) Revisiting Transolver from a Generalized Linear Attention Perspective. Recent analysis reveals that the Physics-Attention mechanism in Transolver [58] can be mathematically unified under this generalized linear attention paradigm. Specifically, Transolver compresses physical nodes into M latent slices using a slice-weight matrix, applies an intermediate self-attention operation G(⋅)G(·) within this slice space, and subsequently deslices the representations back to the physical domain. In the context of the generalized formulation Eq. (20), Transolver’s feature mappings are constructed as: ϕTransolver() _Transolver(Q) =Softmax()∈ℝN×M, =Softmax(XA) ^N× M, (21) ψTransolver()⊤ _Transolver(K) =Norm()⊤∈ℝM×N, =Norm (XA ) ^M× N, (22) where ∈ℝd×MA ^d× M is the learnable weight matrix for slicing. Notably, ϕTransolver _Transolver and ψTransolver _Transolver are derived from the exact same parameter A, differing only in their normalization strategies. This architectural choice introduces two structural limitations. First, the forced parameter-sharing between the slicing (ψ) and deslicing (ϕφ) mappings constrains the model’s performance (see Figure 4 ). Second, because the intermediate interaction G(⋅)G(·) operates strictly within the highly compressed slice space, the model’s ability to capture fine-grained dynamics is difficult by the information loss during the geometric compression. CANO as Dynamic Generalized Linear Attention. Our proposed CANO resolves these limitations. Instead of entirely decoupling the queries and keys through independent linear projections in Eq. (20), CANO dynamically clusters the queries while preserving the full-resolution keys, executing a cross-attention mechanism. Mathematically, by expanding our dispatching and aggregation equations, the full CANO computation is formulated as: =Softmax((⊤)⊤d).Y=C\,Softmax ( (C Q)K d )V. (23) Then we can establish an exact equivalence to the generalized linear attention given by Eq. (20) by defining the CANO feature mappings as: ϕCANO() _CANO(Q) ≜c=⊤(Q)∈ℝKc×d, _c=C (XW_Q) ^K_c× d, (24) ψCANO()⊤ _CANO(K) ≜Softmax(c(K)⊤d)∈ℝKc×N, ( Q_c(XW_K) d ) ^K_c× N, (25) Finally the CANO output can be written as: CANO=ϕCANO()(ψCANO()⊤).Y_CANO= _CANO(Q) ( _CANO(K) V ). (26) This highlights the main advantage of CANO’s architecture. First, CANO is asymmetric, effectively avoiding the structural issue of Transolver. The left mapping ϕCANO() _CANO(Q) acts as a dynamically computed cluster assignment derived strictly from the query space, while the right mapping ψCANO() _CANO(K) processes the interaction between the clustered queries cQ_c and the full resolution keys K. More importantly, unlike Transolver which compresses the entire input field, CANO preserves the keys and values in their original full resolution space. Appendix B Govering Equations In this part, we present the governing equations for the standard PDE benchmarks and the irregular domain benchmarks. B.1 Standard PDE Benchmarks Airfoil. This dataset models the transonic flow over an airfoil. Given the low air viscosity in this scenario, the viscous term is negligible, and the system is governed by the Euler equations (representing mass, momentum, and energy conservation, respectively): ∂ρf∂t+∇⋅(ρfU) ∂ _f∂ t+∇·( _fU) =0, =0, (27) ∂ρfU∂t+∇⋅(ρfUU+pI) ∂ _fU∂ t+∇·( _fU+pI) =0, =0, (28) ∂E∂t+∇⋅((E+p)U) ∂ E∂ t+∇·((E+p)U) =0, =0, (29) where ρf _f is the fluid density and E is the total energy. The data is generated on a structured mesh with a resolution of 221×51221× 51. The model learns to map the locations of these mesh points to the Mach number. Pipe. This benchmark focuses on incompressible flow through a pipe. The governing equations are deduced from the Navier-Stokes relations for Newtonian fluids: ∇⋅U ∇· U =0, =0, (30) ∂U∂t+U⋅∇U ∂ U∂ t+U·∇ U =f−1ρ∇p+ν∇2U. =f- 1ρ∇ p+ν∇^2U. (31) The dataset utilizes a structured mesh with a resolution of 129×129129× 129. We use the mesh structure as input to predict the horizontal fluid velocity within the pipe. Plasticity. Focusing on the plastic forging problem, this dataset simulates a plastic material impacted from above by an arbitrary-shaped die. The governing equation for solid material dynamics is given by the equilibrium equation: ρs∂2u∂t2+∇⋅σ=0, _s ∂^2u∂ t^2+∇·σ=0, (32) where ρs _s represents the solid density, u denotes the displacement vector of the material over time t, and σ is the stress tensor. The input is the die shape recorded on a 101×31101× 31 structured mesh, and the output is the deformation of each mesh point over the future 20 time steps. Navier-Stokes. This dataset simulates incompressible, viscous flow on a unit torus. Since the density of the fluid is constant, the energy conservation is independent, and the fluid dynamics are simplified to the vorticity formulation: ∇⋅U ∇· U =0, =0, (33) ∂ω∂t+U⋅∇ω ∂ω∂ t+U·∇ω =ν∇2ω+f, =ν∇^2ω+f, (34) where U=(u,v)U=(u,v) is the velocity vector, ω=∇×Uω=∇× U is the vorticity, and the viscosity is set to ν=10−5ν=10^-5. We predict the future 10 frames based on the past 10 frames on a 64×6464× 64 regular grid. Darcy Flow. This benchmark models fluid diffusion through porous media, governed by the 2D steady-state elliptic equation: −∇⋅(a(x)∇u(x))=f(x),-∇·(a(x)∇ u(x))=f(x), (35) where a(x)a(x) is the permeability coefficient field. This dataset tests the model’s ability to resolve high-frequency spatial discontinuities in the coefficient field a(x)a(x) on an 85×8585× 85 regular grid. Elasticity. These benchmarks estimate the inner stress of an incompressible material with an arbitrary central void under external tension. They share the same governing equilibrium equation Eq. (32) as the Plasticity dataset. The model maps the material structure to the inner stress. B.2 Irregular domain benchmarks Irregular Darcy. This 2D fluid problem is discretized with 2,290 unstructured nodes. It shares the same 2D steady-state elliptic equation as the standard Darcy flow Eq. (35). However, it challenges the model with highly irregular geometric boundaries. The task maps the varying permeability coefficient field a(x)a(x) to the pressure field u(x)u(x). Pipe Turbulence. Discretized with 2,673 nodes, this scenario requires state predictions of incompressible fluid dynamics within an unstructured pipe geometry. The governing equations are the as before Eq. (30). The model is tasked with predicting the future velocity states based on previous observations in the irregular domain. Composite. We scale our evaluation to a challenging 3D solid mechanics task using the Composite material benchmark (8,232 nodes). This dataset models the thermoelastic deformation of a 3D composite material driven by temperature profiles. The underlying physics is governed by the steady-state equilibrium equation coupled with thermal expansion: ∇⋅σ=0,whereσ=:ε−γΔT.∇·σ=0, σ=C: -γ TI. (36) Here, σ is the stress tensor, C is the stiffness tensor, ε is the strain tensor, and γΔTγ TI represents the thermal stress induced by the temperature change ΔT T. The model learns to map the input 3D temperature profiles directly to the structural stress fields. References (2) A. Friedman, Stochastic Differential Equations and Applications, Springer, 1975. (3) M. Braun, M. Golubitsky, Differential Equations and their Applications, Springer, 1983. (4) H. L. Smith, An Introduction to Delay Differential Equations with Applications to the Life Sciences, Springer, 2011. (5) J. D. Logan, Applied Partial Differential Equations, Springer, 2014. (6) G. F. Simmons, Differential Equations with Applications and Historical Notes, CRC Press, 2016. (7) A. Jaun, J. Hedin, T. Johnson, Numerical Methods for Partial Differential Equations, Swedish Netuniversity, 1999. (8) G. A. Evans, J. M. Blackledge, P. D. Yardley, Numerical Methods for Partial Differential Equations, Springer, 2000. (9) J. F. Wendt (Ed.), Computational Fluid Dynamics: An Introduction, Springer, 2008. (10) E. Tadmor, A review of numerical methods for nonlinear partial differential equations, Bull. Am. Math. Soc. 49 (4) (2012) 507–554. (11) W. F. Ames, Numerical Methods for Partial Differential Equations, Academic Press, 2014. (12) W. E, J. Han, A. Jentzen, A deep learning-based numerical methods for high-dimensional parabolic partial differential equations and backward stochastic differential equations, Commun. Math. Stat. 5 (2017) 349-380. (13) J. Jeon, Coupling of artificial intelligence and traditional numerical methods for accelerated and robust PDE solvers, JMST Adv. 7 (2025) 49-55. (14) J. Han, A. Jentzen, W. E, Solving high-dimensional partial differential equations using deep learning, PNAS 115 (34) (2018) 8505-8510. (15) J. Thiyagalingam, M. Shankar, G. Fox, T. Hey, Scientific machine learning benchmarks, Nat. Rev. Phys. 4 (2022) 413-420. (16) R. Pestourie, Y. Mroueh, C. Rackauckas, P. Das, S. G. Johnson, Physics-enhanced deep surrogates for partial differential equations, Nat. Mach. Intell. 5 (2023) 1458-1465. (17) S. L. Brunton, J. N. Kutz, Promising directions of machine learning for partial differential equations, Nat. Comput. Sci. 4 (2024) 483-494. (18) B. Yu, The deep Ritz method: a deep learning-based numerical algorithm for solving variational problems, Comm. Math. Stat. 6 (1) (2018) 1–12. (19) M. Raissi, P. Perdikaris, G. E. Karniadakis, Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations, J. Comput. Phys. 378 (2019) 686–707. (20) L. Lu, X. Meng, Z. Mao, G. E. Karniadakis, DeepXDE: A deep learning library for solving differential equations, SIAM Rev. 63 (1) (2021) 208–228. (21) G.E. Karniadakis, I.G. Kevrekidis, L. Lu, P. Perdikaris, S. Wang, L. Yang, Physics-informed machine learning, Nat. Rev. Phys. 3 (6) (2021) 422–440. (22) C. Rao, P. Ren, Q. Wang, O. Buyukozturk, H. Sun, Y. Liu, Encoding physics to learn reaction–diffusion processes, Nat. Mach. Intell. 5 (7) (2023) 765–779. (23) A. Krishnapriyan, A. Gholami, S. Zhe, R. Kirby, M. W. Mahoney, Characterizing possible failure modes in physics-informed neural networks, Adv. Neural Inform. Process. Syst. 34 (2021) 26548–26560. (24) S. Wang, X. Yu, P. Perdikaris, When and why PINNs fail to train: A neural tangent kernel perspective, J. Comput. Phys. 449 (2022) 110768. (25) N. Kovachki, S. Lanthaler, S. Mishra, On universal approximation and error bounds for Fourier neural operators, J. Mach. Learn. Res. 22 (2021) 1–76. (26) Z. Li, N. Kovachki, K. Azizzadenesheli, B. Liu, K. Bhattacharya, A. Stuart, A. Anandkumar, Fourier neural operator for parametric partial differential equations, arXiv preprint, arXiv:2010.08895 (2021). (27) L. Lu, P. Jin, G. Pang, Z. Zhang, G. E. Karniadakis, Learning nonlinear operators via DeepONet based on the universal approximation theorem of operators, Nat. Mach. Intell. 3 (3) (2021) 218–229. (28) N. Kovachki, Z. Li, B. Liu, K. Azizzadenesheli, K. Bhattacharya, A. Stuart, A. Anandkumar, Neural operator: Learning maps between function spaces with applications to PDEs, J. Mach. Learn. Res. 24 (2023) 1–97. (29) S. Lanthaler, Z. Li, A. M. Stuart, Nonlocality and nonlinearity implies universality in operator learning, Constr. Approx. 62 (2025) 261-303. (30) T. Chen, H. Chen, Universal approximation to nonlinear operators by neural networks with arbitrary activation functions and its application to dynamical systems, IEEE Trans. Neural Netw. 6 (4) (1995) 911-917. (31) L. Lu, X. Meng, S. Cai, Z. Mao, S. Goswami, Z. Zhang, G. E. Karniadakis, A comprehensive and fair comparison of two neural operators (with practical extensions) based on fair data, Comput. Methods Appl. Mech. Engrg. 393 (2022) 114778. (32) S. Wang, H. Wang, P. Perdikaris, Learning the solution operator of parametric partial differential equations with physics-informed DeepONets, Sci. Adv. 7 (40) (2021) eabi8605. (33) P. Jin, S. Meng, L. Lu, MIONet: Learning multiple-input operators via tensor product, SIAM J. Sci. Comput. 44 (6) (2022) A3490–A3514. (34) A. Kopanicˇ cáková, G. E. Karniadakis, DeepONet based preconditioning strategies for solving parametric linear systems of equations, SIAM J. Sci. Comput. 47 (1) (2025) C151–C181. (35) Z. Li, D. Z. Huang, B. Liu, A. Anandkumar, Fourier neural operator with learned deformations for PDEs on general geometries, J. Mach. Learn. Res. 24 (2023) 1–26. (36) Z. Li, N. Kovachki, C. Choy, B. Li, J. Kossaifi, S. Otta, A. Anandkumar, Geometry-informed neural operator for large-scale 3D PDEs, Adv. Neural Inf. Process. Syst. 36 (2023) 35836–35854. (37) S. Cao, Choose a transformer: Fourier or Galerkin, Adv. Neural Inf. Process. Syst. 34 (2021) 24924–24940. (38) Z. Li, K. Meidani, A. B. Farimani, Transformer for partial differential equations operator learning, arXiv preprint, arXiv:2205.13671 (2022). (39) Z. Li, D. Shu, A. Barati Farimani, Scalable transformer for PDE surrogate modeling, Adv. Neural Inf. Process. Syst. 36 (2023) 28010–28039. (40) Z. Hao, Z. Wang, H. Su, C. Ying, Y. Dong, S. Liu, J. Zhu, GNOT: A general neural operator transformer for operator learning, in: Proc. Int. Conf. Mach. Learn. 2023, p. 12556–12569. (41) Z. Xiao, Z. Hao, B. Lin, Z. Deng, H. Su, Improved operator learning by orthogonal attention, arXiv preprint, arXiv:2310.12487 (2023). (42) B. Shih, A. Peyvan, Z. Zhang, G. E. Karniadakis, Transformers as neural operators for solutions of differential equations with finite regularity, Comput. Methods Appl. Mech. Engrg. 434 (2025) 117560. (43) A. Vaswani, N. Shazeer, N. Parmar, J. Uszkoreit, L. Jones, A. N. Gomez, I. Polosukhin, Attention is all you need, Adv. Neural Inf. Process. Syst. 30 (2017). (44) J. Devlin, M. W. Chang, K. Lee, K. Toutanova, BERT: Pre-training of deep bidirectional transformers for language understanding, Proc. North Am. Chap. Assoc. Comput. Linguist, 2019, p. 4171–4186. (45) T. Brown, B. Mann, N. Ryder, M. Subbiah, J. D. Kaplan, P. Dhariwal, D. Amodei, Language models are few-shot learners, Adv. Neural Inf. Process. Syst. 33 (2020) 1877–1901. (46) A. Katharopoulos, A. Vyas, N. Pappas, F. Fleuret, Transformers are RNNs: Fast autoregressive transformers with linear attention, Proc. Int. Conf. Mach. Learn. 2020 5156–5165. (47) S. Wang, B. Z. Li, M. Khabsa, H. Fang, H. Ma, Linformer: Self attention with linear complexity, arXiv preprint, arXiv:2006.04768 (2020). (48) Z. Shen, M. Zhang, H. Zhao, S. Yi, H. Li, Efficient attention: Attention with linear complexities, Proc. IEEE/CVF Winter Conf. Appl. Comput. Vis. 2021, p. 3531–3539. (49) D. Han, T. Ye, Y. Han, Z. Xia, S. Pan, P. Wan, G. Huang, Agent attention: On the integration of softmax and linear attention, in: Proc. Eur. Conf. Comput. Vis. 2024, p. 124–140. (50) D. Han, Y. Pu, Z. Xia, Y. Han, X. Pan, X. Li, G. Huang, Bridging the divide: Reconsidering softmax and linear attention, Adv. Neural Inf. Process. Syst. 37 (2024) 79221–79245. (51) O. Ovadia, A. Kahana, P. Stinis, E. Turkel, D. Givoli, G. E. Karniadakis, VITO: Vision transformer-operator, Comput. Methods Appl. Mech. Engrg. 428 (2024) 117109. (52) S. Wang, J. H. Seidman, S. Sankaran, H. Wang, G. J. Pappas, P. Perdikaris, CVIT: Continuous vision transformer for operator learning, arXiv preprint, arXiv:2405.13998 (2024). (53) X. Liu, B. Xu, S. Cao, L. Zhang, Mitigating spectral bias for the multiscale operator learning, J. Comput. Phys. 506 (2024) 112944. (54) A. Dosovitskiy, An image is worth 16x16 words: Transformers for image recognition at scale, arXiv preprint, arXiv:2010.11929 (2020). (55) Z. Liu, Y. Lin, Y. Cao, H. Hu, Y. Wei, Z. Zhang, B. Guo, Swin Transformer: Hierarchical vision transformer using shifted windows, in: Proc. IEEE/CVF Int. Conf. Comput. Vis. 2021, p. 10012–10022. (56) T. Wang, C. Wang, Latent neural operator for solving forward and inverse PDE problems, Adv. Neural Inf. Process. Syst. 37 (2024) 33085–33107. (57) B. Alkin, A. Fürst, S. Schmid, L. Gruber, M. Holzleitner, J. Brandstetter, Universal physics transformers: A framework for efficiently scaling neural operators, Adv. Neural Inf. Process. Syst. 37 (2024) 25152–25194. (58) H. Wu, H. Luo, H. Wang, J. Wang, M. Long, Transolver: A fast transformer solver for PDEs on general geometries, arXiv preprint, arXiv:2402.02366 (2024). (59) H. Luo, H. Wu, H. Zhou, L. Xing, Y. Di, J. Wang, M. Long, Transolver++: An accurate neural solver for partial differential equations on million-scale geometries, arXiv preprint, arXiv:2502.02414 (2025). (60) M. A. Nabian, C. Liu, R. Ranade, S. Choudhry, X-meshgraphnet: Scalable multi-scale graph neural networks for physics simulation, arXiv preprint, arXiv:2411.17164 (2024). (61) M. Bleeker, M. Dorfer, T. Kronlachner, R. Sonnleitner, B. Alkin, J. Brandstetter, Neural computational fluid dynamics: deep learning on high-fidelity automotive aerodynamics simulations, arXiv preprint, arXiv:2502.09692 (2025). (62) A. Anandkumar, K. Azizzadenesheli, K. Bhattacharya, N. Kovachki, Z. Li, B. Liu, A. Stuart, Neural operator: Graph kernel network for partial differential equations, 2020, in: ICLR Workshop on Integration of Deep Neural Models and Differential Equations. (63) Z. Li, N. Kovachki, K. Azizzadenesheli, B. Liu, A. Stuart, K. Bhattacharya, A. Anandkumar, Multipole graph neural operator for parametric partial differential equations, Adv. Neural Inf. Process. Syst. 33 (2020) 6755–6766. (64) N. Rahaman, A. Baratin, D. Arpit, F. Draxler, M. Lin, F. Hamprecht, J.N. Bakas, A. Courville, On the spectral bias of neural networks, in: Int. Conf. Mach. Learn., 2019, p. 5301–5310. (65) W. Hu, S. Liu, P. Qiao, Z. Sun, Y. Dou, Transolver is a linear transformer: Revisiting physics attention through the lens of linear attention, Proc. AAAI Conf. Artif. Intell. 40 (1) (2026) 408–416. (66) G. Chen, X. Liu, Q. Meng, L. Chen, C. Liu, Y. Li, Learning neural operators on Riemannian manifolds, arXiv preprint, arXiv:2302.08166 (2023). (67) W. Hamilton, Z. Ying, J. Leskovec, Inductive representation learning on large graphs, Adv. Neural Inf. Process. Syst. 30 (2017) 1024–1034. (68) A. Tran, A. Mathews, L. Xie, C. S. Ong, Factorized Fourier neural operators, arXiv preprint, arXiv:2111.13802 (2021). (69) G. Wen, Z. Li, K. Azizzadenesheli, A. Anandkumar, S. M. Benson, U-FNO—An enhanced Fourier neural operator-based deep-learning model for multiphase flow, Adv. Water Resour. 163 (2022) 104180. (70) M. A. Rahman, Z. E. Ross, K. Azizzadenesheli, U-NO: U-shaped neural operators, arXiv preprint, arXiv:2204.11127 (2022). (71) H. Wu, T. Hu, H. Luo, J. Wang, M. Long, Solving high-dimensional PDEs with latent spectral models, arXiv preprint, arXiv:2301.12664 (2023). (72) X. Yue, Y. Yang, L. Zhu, Holistic physics solver: Learning PDEs in a unified spectral physical space, arXiv preprint, arXiv:2410.11382 (2024). (73) P. Liu, X. Ren, P. Wang, H. Yuan, Z. Hao, G. Chen, S. Cai, An efficient graph transformer operator for learning physical dynamics with manifolds embedding, arXiv preprint, arXiv:2512.10227 (2025).