We introduce a scalable Gaussian process (GP) framework with deep product kernels for data-driven learning of parametrized spatio-temporal fields over fixed or parameter-dependent domains. The proposed framework learns a continuous representation, enabling predictions at arbitrary spatio-temporal coordinates, independent of the training data resolution. We leverage Kronecker matrix algebra to formulate a computationally efficient training procedure with complexity that scales nearly linearly with the total number of spatio-temporal grid points. A key feature of our approach is the efficient computation of the posterior variance at essentially the same computational cost as the posterior mean (exactly for Cartesian grids and via rigorous bounds for unstructured grids), thereby enabling scalable uncertainty quantification. Numerical studies on a range of benchmark problems demonstrate that the proposed method achieves accuracy competitive with operator learning methods such as Fourier neural operators and deep operator networks. On the one-dimensional unsteady Burgers' equation, our method surpasses the accuracy of projection-based reduced-order models. These results establish the proposed framework as an effective tool for data-driven surrogate modeling, particularly when uncertainty estimates are required for downstream tasks.
Figures & tables
Figure 1: Illustration of the geometric transformation: mapping a parametrized physical domain Ω(μΩ) (left) to a fixed reference domain Ω~ (right) to enable tensor-product grid embedding.
Figure 2: 1D Burgers’ problem. Predicted solutions using the proposed GP framework with DPK and Matérn- 5/2 based product structured kernel at times t=3.5,7.0,10.5, and 14 are shown.
Method
Relative error
μtest(1)
μtest(2)
Deep-LSPG ( p=20 )
0.0009
0.001
Proposed GP (DPK-Matérn- 5/2 )
0.0033
0.0029
Deep-Galerkin ( p=20 )
0.015
0.013
Proposed GP (Matérn- 5/2 )
0.0211
0.0078
POD-Galerkin ( p=20 )
0.032
0.03
Table 1: 1D Burgers’ problem. Relative test errors for the test snapshots corresponding to the parameters μtest(1)=(4.3,0.021) and μtest(2)=(5.15,0.0285) are presented. Errors for the physics-based methods (POD-Galerkin, POD-LSPG, Deep-Galerkin, and Deep-LSPG) for reduced dimension p=20 are obtained from Figure 3 in [ 25 ] .
Figure 3: 1D Burgers’ problem. For test parameter μtest(1)=[4.3,0.021] , predicted solutions and confidence bounds ( ±2σ ) using the proposed GP framework with DPK (left) and Matérn- 5/2 based product kernel (right) at times t=3.5,7.0,10.5, and 14 are shown.
Figure 4: 1D Burgers’ problem. For test parameter μtest(2)=[5.15,0.0285] , predicted solutions and confidence bounds ( ±2σ ) using the proposed GP framework with DPK (left) and Matérn- 5/2 based product kernel (right) at times t=3.5,7.0,10.5, and 14 are shown.
Method
Relative error
training
testing
Proposed GP (DPK)
0.0181
0.0326
Geo-FNO (O-mesh)
0.0344
0.0363
FNO interpolation
0.0314
0.0508
DeepONet
0.0528
0.0965
Table 2: Elasticity problem. Relative errors for the proposed GP framework compared with ML-based operator methods, including Geo-FNO [ 26 ] and DeepONet [ 30 ] . Best testing results are highlighted in bold.
Figure 5: Elasticity problem. Comparison of the stress field predicted by the proposed GP framework (DPK) with the high-fidelity solution for a test geometry in the reference domain.
Method
Airfoil
Pipe
training
testing
training
testing
Geo-FNO
0.0134
0.0138
0.0047
0.0067
Proposed GP (DPK)
0.0051
0.0142
0.0062
0.0131
Proposed GP (Matérn- 5/2 )
0.0065
0.0316
0.0067
0.0171
FNO interpolation
0.0287
0.0421
0.0083
0.0151
UNet interpolation
0.0039
0.0519
0.0109
0.0182
Table 3: Airfoil and pipe problems. Relative errors comparing the proposed GP framework with ML-based operator methods, including Geo-FNO [ 26 ] and FNO [ 27 ] . Best testing results are highlighted in bold.
Appendix figures & tables4 assets
Supplementary material from the paper’s appendix.
Appendix
Problem
Learning rate
Maximum training iterations
1D unsteady Burgers’
10−2
1000
Hyper-elastic problem
1.5×10−2
450
2D Transonic flow
10−2†
368
Navier-Stokes Equation (Pipe)
10−2†
900
Appendix
Table 4: Learning rate and maximum number of training iterations for all the problems. Learning rate with (⋅)† denotes the initial learning rate. For the transonic flow and pipe flow test problems, we employed a step learning rate decay scheduler with a step size of 100 and a decay factor of 0.8.
Figure 6: The subdomains and the boundary conditions for the linear elastic block problem [ 15 ] are shown.
Figure 7: Elastic block problem. Convergence plot comparing the mean relative test errors corresponding to 1000 test snapshots for the proposed GP framework (using the deep product kernel (DPK)) with POD-GPR and physics-based POD-Galerkin with reduced dimension p=20 .
Figure 8: Elastic block problem. Box plot of relative test errors on 1000 test parameters for models trained on an increasing number of snapshots ( N=50,100,200,500,1000,2000,5000 ).
Neural operators provide fast surrogates for PDEs but their deterministic predictions limit their use in tasks requiring uncertainty quantification (UQ), especially under geometric variability. Existing approaches primarily model uncertainty in network parameters, largely overlooking the geometry-aware representations learned by the operator itself. We propose REEF-GP (Residual on Embedded Features Gaussian Process), a post-hoc UQ framework that fits a GP to the residuals of a frozen neural operator whose internal embeddings define the kernel feature space. Rather than learning a separate feature map, REEF-GP adapts the operator's intrinsic coordinate-feature representations to construct geometry-aware uncertainties. To ensure stability and scalability on unstructured domains, REEF-GP incorporates spectral-normalized projections, heteroscedastic geometry-aware noise, and efficient subset-based training that avoids restrictive low-rank approximations. Across five PDE benchmarks with varying geometries, REEF-GP preserves predictive accuracy while achieving calibrated uncertainty estimates competitive with deep ensembles but at a fraction of their cost. Our approach remains robust under geometric distribution shift, with uncertainty concentrating in physically meaningful regions (e.g., shock fronts). Our results demonstrate that accurate and scalable post-hoc UQ for neural operators can be achieved directly in their learned feature space, offering a practical alternative to parameter-centric approaches.
Oriol Vendrell-Gallart, Nima Negarandeh, Ramin Bostanabad
Department of Mechanical and Aerospace Engineering University of California, Irvine Irvine, CA 92617
We present a theoretically grounded Gaussian process framework that leverages neural feature maps to construct expressive kernels. We show that the learned feature map can be interpreted as an optimal low-rank approximation to a Gram matrix derived from an implied RKHS, from which we establish consistency of the GP posterior. We further analyse the spectral properties of the induced kernels and introduce product feature-map kernels to address oversmoothing. This simple yet powerful approach enables fast, scalable, and accurate exact GP inference with minimal upfront work. The flexibility of kernel design supports seamless application to both regression and classification tasks across diverse data modalities, including tabular inputs and structured domains such as images. On benchmark datasets, this approach surpasses pre-existing methods in terms of accuracy and training and prediction efficiency.
Gaussian Processes (GPs) are a powerful tool for modelling and quantifying uncertainty in functional relationships. However, they require practitioners to make a number of design decisions, such as the choice of the kernel and the observation model. Suboptimal choices can produce misspecified models that do not capture the underlying data generating process. We introduce Predictively Oriented Gaussian Processes (PrO-GPs), which treat predictive uncertainty as the primary inferential target and provide a robust alternative to standard GPs. Although direct computation of a PrO posterior for nonparametric models is intractable, we derive a reduced formulation and practical sampling scheme for efficient computation. Through synthetic and real data experiments, we show that PrO-GPs produce better calibrated predictive distributions under model misspecification compared to standard GP approaches.