Neural Operator induced Gaussian Process framework for probabilistic solution of parametric partial differential equations
Sawan Kumar Rajdip Nayek Souvik Chakraborty

Identifying non-Hermitian critical points with quantum metric
    Dan-Wei Zhang
    April 29, 2024
=============================================================




The study of neural operators has paved the way for the development of efficient approaches for solving partial differential equations (PDEs) compared with traditional methods. However, most of the existing neural operators lack the capability to provide uncertainty measures for their predictions, a crucial aspect, especially in data-driven scenarios with limited available data. In this work, we propose a novel Neural Operator-induced Gaussian Process (NOGaP), which exploits the probabilistic characteristics of Gaussian Processes (GPs) while leveraging the learning prowess of operator learning. The proposed framework leads to improved prediction accuracy and offers a quantifiable measure of uncertainty. The proposed framework is extensively evaluated through experiments on various PDE examples, including Burger's equation, Darcy flow, non-homogeneous Poisson, and wave-advection equations. Furthermore, a comparative study with state-of-the-art operator learning algorithms is presented to highlight the advantages of NOGaP. The results demonstrate superior accuracy and expected uncertainty characteristics, suggesting the promising potential of the proposed framework. 







§ INTRODUCTION

Neural Operators (NO)<cit.> have recently emerged as a viable alternative to traditional numerical methods for solving partial differential equations. Some recent applications of NO in computational mechanics include weather modeling <cit.>, medical diagnosis using elastography <cit.>, fracture propagation <cit.>, and fault diagnosis <cit.>, among others. However, a major bottleneck in the widespread adoption of NO in critical engineering systems (e.g., nuclear reactors) resides in the fact that they fail to quantify the uncertainty in the trained model. This has been recognized by the research community, and some efforts in this direction can be found in <cit.>.
An alternative to the neural operators is to employ probabilistic machine learning algorithms such as the Gaussian Process (GP) for solving partial differential equations. Unlike NO, GP is inherently probabilistic and hence, can seamlessly quantify the uncertainty in the model. Recent developments on the use of GP for solving PDEs include <cit.>. Unfortunately, the developments on GP and NO are carried out somewhat in parallel. To that end, the objective of this paper is to investigate how NO and GP can benefit from each other, and if a fusion of the two can result in a superior model. 


Neural operators learn and map two infinite-dimensional function spaces. Some neural operators that have been studied in the past and gained popularity are Wavelet Neural Operator (WNO) <cit.>, Physics Informed WNO <cit.>, DeepONet <cit.>, Fourier Neural Operator (FNO) <cit.>, PCA-Net <cit.>, Graph Neural Operators (GNO) <cit.> and many more. Compared to classical PDE solvers, neural operators have significant computational advantages; once trained, they can make predictions for different initial and boundary conditions efficiently. The performance of both WNO and FNO is promising when compared with the state-of-the-art operator learning methods in terms of prediction accuracy. However, FNO utilizes the Fast Fourier Transform (FFT), resulting in a loss of spatial information. Consequently, frequency-localized basis functions in FFT pose challenges in learning complex domain problems. To address this variants such as geometric-FNO <cit.> and geometry-informed neural operators <cit.> have been proposed.
Alternatively, WNO utilizes wavelet transform which retains space-frequency information and can handle signals with spikes and discontinuity, thus enabling superior learning of the patterns in images. Some recent improvements of the classical WNO include <cit.>. However, as previously stated, both WNO and FNO are deterministic and cannot yield the predictive uncertainty that is essential for decision-making in safety-critical systems. 


Data in real-world scenarios are accompanied by inherent noise and sources of uncertainty. These sources include errors stemming from the limitation of precision of the measuring instrument, incorrect measurement of observations, missing data, etc. Additionally, models may be built upon simplifying assumptions that may not fully capture the complexity of the underlying phenomena, which also contributes to the uncertainty. This is often addressed using probabilistic frameworks like Gaussian Process Regression (GPR) <cit.> and Bayesian neural operator frameworks <cit.>. GPR is a non-parametric machine learning model, and by design, learns distribution over the functions. Similar to neural networks and neural operators, GPR is also a popular choice among researchers for solving partial differential equations. This includes the development of physics-informed GPR for linear and nonlinear systems <cit.>. However, GPR suffers from the curse of dimensionality and does not scale with an increase in training samples; hence, GPR is less popular as compared to neural network-based approaches.


The developments on NO and GP are mostly carried out in parallel. To that end,
the objective of this paper is to propose a novel architecture referred to as the Neural Operator induced Gaussian Process (NOGaP) that attempts to bridge the gap between the NO and GP-based approaches. NOGaP ingeniously combines the strengths of Gaussian process regression with NO. Specifically, the high-level idea is to use NO as a mean function within the GP framework. We hypothesize that with NO as the mean function, NOGaP will be able to exploit the accuracy of NO and the probabilistic nature of the GP. Without loss of generality, we use the recently proposed WNO as the mean function with the proposed approach. Additionally, to address the challenge associated with the scalability of GP, we propose to use the Kronecker -product-based approach proposed in <cit.>.

The proposed framework has the following key features: 
 
     
  * Accuracy:
     NOGaP, a fusion of NO and GP, is anticipated to offer superior accuracy compared to both NO and GP individually. This has been empirically illustrated in the numerical section.  
     
  * In-built uncertainty quantification capability: NOGaP possesses built-in probabilistic capabilities, allowing it to assess the epistemic uncertainty due to limited and noisy datasets. This feature is crucial, especially for ensuring the reliability of the model in critical engineering systems. 
     
  * Scalability: Leveraging the Kronecker product property, NOGaP exhibits scalability, enabling its efficient handling of large datasets. This marks a departure from conventional GPs, which often encounter scalability issues due to the curse of dimensionality.

    
    
    
    
 

 
The remainder of the paper is arranged as follows: Section <ref> elaborates on the mathematical formulation of the proposed framework. Section <ref> explores the different numerical examples on which the NOGaP is tested along with the results. Finally, the paper is concluded in the last section <ref> with results, observations and future work. 



§ NEURAL OPERATOR INDUCED GAUSSIAN PROCESS (NOGAP)

Let us define the hypothesis space, where 𝒵⊂ℝ^d represents the input space and 𝒴⊂ℝ^o denotes the output space. Suppose we have N distinct pairs of observations with inputs, 𝐙 = (x_1, …, x_N,a_1, …, a_N), and a corresponding set of outputs, 𝐘 = (y_1, …, y_N). The model takes the form: 

    y(z)     = f_o(z) + ϵ
    ϵ   ∼𝒩(0, σ_o^2 𝐈)

where 𝐲 are the noisy observations and ϵ∈ℝ^o is modeled as independent zero-mean Gaussian white noise with variance σ_o^2. The vector-valued function f_o is assumed to follow a GP prior:

    f_o(z) ∼𝒢𝒫(m_o(z; θ_m), K_o(z, z^'; θ_k))

Here we have vector-valued mean function m_o: ℝ^d ↦ℝ^o and a matrix-valued covariance function 𝐊_o: 𝒵×𝒵↦ℝ^o × o, with vector-valued parameters θ_m and θ_k where θ_m represents parameters of the mean function and θ_k denotes kernel parameters.
Here in our formulation, the multioutput kernel K_o is assumed to be separable <cit.>. It can be expressed as:

    (𝐊_o)_dd'   = k_x(z, z') · k_f(d,d') 
    𝐊_o    = 𝐊_x ⊗𝐊_f

where ⊗ denotes the Kronecker product, 𝐊_o is the multioutput kernel in ℝ^No × No, k_f is a covariance function specifying inter-output similarities and k_x is a covariance function over inputs. In case the inputs are available on a multidimensional cartesian grid, we can further decompose the covariance matrix 𝐊_x as the Kronecker product 𝐊_x1⊗···⊗𝐊_xm. Exploiting the Kronecker product property can further improve computational efficiency. 


 §.§ Mean Function: Neural Operator

There are many available choices for the mean function in GP, each with its advantages and limitations. For instance, a linear mean function can capture linear trends in the data, while a non-linear mean function, such as a neural network, can handle more complex relationships and nonlinear patterns. By using a mean function that aligns with the characteristics of the problem, the GP can better model the underlying process as it only has to learn the discrepancy between the ground truth and the values computed using the mean function. This results in overall improved predictions. Thus, selecting a suitable mean function in GP is critical in improving the prediction capability, resulting in efficient learning of the underlying function. 
Here in this work, we have used neural operators as the mean function. The core idea of neural operators is to leverage neural networks to learn the underlying function space representation. Neural operators represent a composite function consisting of the encoder, neural network, and decoder. The encoder lifts the dimensionality of the input to a high dimensional latent space. The decoder projects the computed output back to the function space. 
Consider a fixed domain D ∈ (a,b) and x ∈ D. For a source term of a PDE f(x,t): D ×𝒯↦ℝ, with initial condition u(x,0): D ↦ℝ and the boundary condition is u(∂ D, t) : D ×𝒯↦ℝ. The solution space is u(x,t): D ×𝒯↦ℝ, where t is the time coordinate. For a Banach space with 𝒜 := 𝒞(D; ℝ^d_1) and 𝒰:= 𝒞(D; ℝ^d_2), suppose that

    𝒟: 𝒜×θ_NN↦𝒰;      θ_NN = {𝒲, 𝐛}

for arbitrarily any non-linear neural operator 𝒟, where θ_NN is the parameter space of the neural network. If we have access to N number of pairs {(a_j, u_j)}_j=1^N, then we can approximate 𝒟 using the pairs of N data points.
Although the proposed approach can accommodate any neural operator as the mean function, we restrict ourselves to the scenario where Wavelet Neural Operator (WNO) <cit.> is used as the mean function.

WNO exploits the kernel integration proposed in <cit.> in combination with wavelet transformation. WNO takes N input functions discretized on a grid and lifts them to a higher dimension using a densely connected layer. The output from the densely connected layer is then fed into the wavelet block, and its dimensions are preserved after the final operations of the wavelet block. 
Mathematically, the wavelet blocks are represented as: 

    f_j+1(x) = σ( (𝒦(a(x); ϕ∈β) * v_j(x)) + b v_j(x) ),      j = 1, …, J,
 
where σ represents a non-linear activation function, 𝒦 denotes the wavelet integral operator, v_j(x) corresponds to the transformed function from the preceding layer, and b signifies a linear transformation. The wavelet coefficients are derived through a discrete wavelet transform (DWT) and an inverse discrete wavelet transform (IDWT) of the mother wavelet. For further details, the interested reader can refer <cit.>. 


 §.§ Likelihood of the data

Consider N distinct observations with inputs, 𝐙 = (x_1, …, x_N,a_1, …, a_N)^T, and a corresponding set of outputs, 𝐘 = (y_1, …, y_N)^T as before. For brevity of representation, we consider the vectorized form of 𝐘, q = (y_11, …, y_N1,…, .
. y_1o,…, y_No), where y_ij denotes the j-th component of the output y_i corresponding to the i-th input x_i. With this setup, the likelihood can be expressed as follows:

    q|f_o,𝐙, θ_m, θ_k ∼𝒩(f_o(𝐙), σ_o^2 𝐈)

Given the fact that the mean function is parameterized using the WNO, we set θ_m = θ_NN, and the mean function computed using WNO is represented using h. In practice, we work with the Negative log-marginal likelihood (NLML),






    log p(q|𝐙,θ)     =-1/2 (q-h(𝐙; θ_NN))^T(K_o(𝐙,𝐙'; θ_k)+σ_o^2 𝐈)^-1 (q-h(𝐙;θ_NN))
          -1/2log|K_o(𝐙,𝐙';θ_k)+σ_o^2 𝐈|+ constant,

where the 𝐈 denotes an identity matrix. The unknown hyperparameters θ = {θ_NN, θ_k, σ_o} are determined by minimizing the negative log-marginalized likelihood (NLML) as given by eq. (<ref>). The optimization problem to find the parameter θ is expressed as follows:

    θ̂ = θmin  -log p(q|𝐙, θ)

Both the parameters of the mean and covariance are optimized simultaneously. Once the training is done the hyperparameter θ is obtained. More thorough and rigorous details can be found in Murphy <cit.> and Rasmussen and Williams <cit.>.
An algorithm depicting the steps involved in training the proposed model is shown in Algorithm <ref>.




 §.§ The Predictive Inference

Since we are working with Gaussian likelihood we can get the predictive distribution in its closed form. We get the following joint distribution for a new data sample z^*.

    [     q; f_o^* ]∼𝒩(
    [   h(𝐙; θ_NN); h(z^*; θ_NN) ],
    [ 𝐊_o(𝐙, 𝐙; θ_k) + σ_o^2 𝐈          𝐊_o(𝐙, z^*;θ_k);          𝐊_o( z^*,𝐙;θ_k)        𝐊_o(z^*, z^*;θ_k) ])

Here, the observed training output is denoted by q, the predictions for the test sample by f_o^*, and the mean function values for the training samples and the test sample by h(𝐙; θ_NN) and h(z^*; θ_NN).
Computing conditional distribution from the joint distribution is straightforward, refer <cit.>. Inference can be done using standard GP formulas for the mean and variance of the predictive distribution.
By conditioning the joint distribution on q the posterior predictive distribution of f_o(z^*) can be expressed as:


    p(f_o^* |𝐙, z^*, q) = 𝒩(f_o^* |m^*, cov(f_o^*)), 
       m^* = h^*(z^*; θ_NN) + 𝐊_o(z^*, 𝐙; θ_k)(𝐊_o(𝐙, 𝐙;θ_k) + σ_o^2 𝐈)^-1(q - h(𝐙; θ_NN)), 
       cov(f_o^*) = 𝐊_o(z^*, z^*; θ_k) - 𝐊_o(z^*, 𝐙; θ_k)(𝐊_o(𝐙, 𝐙; θ_k) + σ_o^2 𝐈)^-1𝐊_o(𝐙, z^*; θ_k)


where, m^* represents the predictive mean of f_o^*. The covariance of f_o^*, denoted as cov(f_o^*), is computed using the Eq. (<ref>).
An algorithm depicting the steps involved in inference is shown in Algorithm <ref>.



 §.§ Choice of the Covariance function

One key aspect of any GP model is the choice of the covariance kernel.
Out of a vast range of available kernels, the  Matérn family of kernels is quite popular as it offers a diverse variety of smoothness for modeling signals, and therefore it has been the choice for the numerical examples covered in Section <ref>. Mathematically, the family of Matérn kernels for 𝒳⊂ℝ^d, constants α > 0 and h >0 can be expressed as

    k_α, h(z, z^')=1/2^α-1Γ(α)(√(2 α)z-z^'/h)^α K_α(√(2 α)z-z^'/h),   z, z^'∈𝒳

where Γ is the gamma function, and K_α is the modified Bessel function of the second kind of order α. For specific values of α, we get different Matérn kernels namely, Matérn-1/2, Matérn-3/2, and Matérn-5/2. The following equation shows the mathematical expression for the Matérn-5/2:


    k_5 / 2, h(z, z^')=(1+√(5)z-z^'/h+5z-z^'^2/3 h^2) exp(-√(5)z-z^'/h).











§ NUMERICAL EXAMPLES

This section explores a variety of case studies ranging from fluid dynamics, gas dynamics, and phase-field modeling on which the proposed NOGaP framework is tested for its efficacy. The examples include (a) 1D Burger's equation, (b) 1D wave-advection equation, (c) 2D Darcy flow equation with a notch in the triangular domain, and (d) 2D non-homogeneous Poisson equation. The architecture of NOGaP consists of two key blocks: mean function and covariance function. The Daubechies family of mother wavelets is used for the mean function block with the GeLU which is differentiable at all points <cit.>. The information related to the architecture of the framework is summarised in <ref>.
The framework's performance is summarised in the <ref>, displaying a comparison of relative Mean Square Errors (MSE) and the comparison with WNO, GP (zero mean) and Bayesian WNO (B-WNO)  for various examples.
The results obtained from NOGaP demonstrate that the proposed framework shows improved accuracy and can learn the desired mapping between input and output space effectively. The proposed NOGaP framework also provides the uncertainty estimates associated with the predictions. Further details of individual cases are illustrated in the following subsection. 







 §.§ 1D Burger Equation

The 1D Burgers equation is a well-known partial differential equation that models one-dimensional flows in various fields, including fluid mechanics, gas dynamics, and traffic flow. One of the interesting features of this equation is that it exhibits both diffusive and advective behavior, making it a useful tool for studying a wide range of phenomena. The 1D Burger equation with periodic boundary is mathematically expressed as,


    ∂_t u(x, t) + 0.5∂_x u^2(x,t)    = ν∂_x x u(x,t),        x ∈(0,1), t ∈(0,1] 
    
    u(x = 0,t)    = u(x= 1,t),        x ∈(0,1), t ∈(0,1] 
    
    u(x,0)    = u_0(x),        x ∈(0,1).

ν here is the viscosity of the flow, which is positive and greater than zero, and u_0(x) is the initial condition. For generating the initial conditions u_0(x), a Gaussian random field is used where u_0(x) ∼𝒩(0,625(- Δ + 25I)^-2. The aim is to learn the mapping from input and output function spaces for the 1D Burger's PDE that is u(x, t = 0) ↦ u(x, t = 1). We consider ν = 0.1 and a spatial resolution of 512. The dataset is taken from Ref. <cit.>.  
Results: The predictions for three different initial conditions are presented in <ref>. The relative MSE values of the predictions are mentioned in <ref>, revealing that the proposed framework exhibits favorable performance for the 1D Burger's problem. Notably, it provides an improvement over the 
standalone WNO, GPR (zero mean) and B-WNO predictions. A visual inspection of the plot confirms the agreement between the ground truth obtained from the numerical solver and the predictions generated by the NOGaP framework. Additionally, the plot includes the corresponding 95% confidence intervals, providing a measure of uncertainty associated with the predictions. Further, a case study is carried out by considering different training sample sizes (see <ref>). The mean response for all the cases is found to have an excellent match with the ground truth. The predictive uncertainty, on the other hand, reduces with an increase in the number of training samples. This indicates a reduction in the predictive uncertainty with an increase in training samples.








 §.§ 1D Wave Advection Equation

The wave advection equation is a type of hyperbolic partial differential equation (PDE) that is commonly used in physics and engineering to describe the movement of a scalar under a known velocity field. The equation is often used to model the behavior of waves in various media, such as sound waves in air or water waves in the ocean. In this example, periodic boundary conditions are applied to simulate the behavior of the medium over time.
When the wave advection equation is considered with periodic boundary conditions, it can be expressed mathematically as:


    ∂_t u(x,t) + ν∂_x u(x,t)    = 0,        x ∈ (0,1), t∈ (0,1) 
    
        u(x - π)    = u(x+π),        x ∈ (0,1).

Here ν is positive and greater than zero,u representing the speed of the flow and the boundary condition is 2π periodic. The initial condition is chosen as,

    u(x,0) = h1_{c- ω/2, c+ω/2} + √(max(h^2 - (a(x-c))^2,0) ).

here the variable ω represents the width and h represents the height of the square wave. The wave is centred around at x = c and the values of {c,ω,h} are selected from [0.3,0.7]× [0.3,0.6] × [1,2]. For this example, the value of ν is set to 1 and we are learning the mapping from the initial condition to the final condition at t = 0.5 that is u(x,t=0) ↦ u(x,t=0.5). The spatial resolution grid is 40 and we have considered 1000 instances of training samples. 

Results: <ref> depicts the predictions obtained from the NOGaP for different initial conditions. The <ref> shows the relative MSE value from NOGaP and other frameworks. The RMSE value for NOGaP is better compared to WNO, GP (zero mean) and B-WNO. Similar to the previous example it can be observed that the NOGaP model can make predictions effectively for this example. Further, a case study is carried out by considering different training sample sizes (see <ref>). The mean response for all the cases is found to have an excellent match with the ground truth. The predictive uncertainty, on the other hand, reduces with an increase in the number of training samples. This indicates a reduction in the predictive uncertainty with an increase in training samples.










 §.§ 2D Darcy flow equation with a notch in triangular domain

The 2D Darcy flow equation is a popular tool utilized for modeling the flow of gases or liquids through a porous medium. This equation is particularly useful in fields such as environmental engineering, geology, and hydrology, where understanding fluid dynamics in porous media is essential. Without the time component, the 2D Darcy flow equation can be represented as a second-order nonlinear elliptic partial differential equation (PDE). The equation is given as,


    -∇· (a(x,y)∇ u(x,y)    = f(x,y),        x,y ∈(0,ℝ) 
    
         u(x,y)    = u_0(x,y),        x,y ∈∂(0,ℝ)

here we have u(x,y) = u_0(x,y) is the initial condition, a(x,y) is the permeability field, u(x,y) is the pressure, f(x,y) is a source function, and ∇ u(x,y) is the pressure gradient. For this example, the setup and dataset are taken from Ref. <cit.>. In this example, a complex boundary condition is considered that is a triangular domain with a notch. The dataset is generated using Gaussian Process (GP) as follows:

    u(x) ∼𝒢𝒫(0,k(x,x'), 
    
        k(x,x') = exp(-(x-x')^2/2l^2),        l=0.2,        x,x' ∈ [0,1]

A notch is introduced in the flow medium with the triangular geometry. The forcing function f(x,y) is set to -1 and the permeability field a(x,y) is set to 0.1. The objective is to learn the mapping from the input (boundary conditions) and output (pressure field) function spaces that are u(x,y)|_∂ω↦ u(x,y). The spatial resolution is taken to be 26 × 26.

Results: The <ref> shows the prediction obtained from the proposed NOGaP model. It displays predictions corresponding to one of the boundary conditions. The corresponding ground truth is obtained from the numerical solver and the mean predictions from the NOGaP. The plot suggests that the proposed framework can make accurate predictions and is able to approximate the true solution. From the <ref> we can observe, that there isn't much difference between the RMSE value for both NOGaP and WNO but since the former provides an uncertainty measure over the predictions it becomes more useful and reliable. However, the predictions are significantly improved when compared with GP (zero mean) and B-WNO. Through visual inspection, one can hardly find any discernable difference between the ground truth and the predictions.




















 §.§ 2D Non-homogeneous Poisson's Equation

The non-homogeneous Poisson equation is an elliptic partial differential equation (PDE) widely applicable in diverse physical scenarios. It serves as the fundamental mathematical model for various problems, such as the steady-state heat diffusion and the electric field generated by a specified electric charge. The Poisson equation, incorporating a source term f(x, y) and periodic boundary conditions, is mathematically expressed as:

    ∂_xxu + ∂_yyu = f(x, y),    x ∈ [-1, 1],    y ∈ [-1, 1] 
    
          u(x = -1, y) = u(x = 1, y) = u(x, y = -1) = u(x, y = 1) = 0

 where, f(x,y) is the source function, u(x,y) is the solution of the equation over the 2D domain. The objective is to learn the mapping between the input and the output space, that is f(x,y) ↦ u(x,y).
 To establish the ground truth solution, we opt for an analytical solution represented as u(x, y) = αsin(π x)(1 + cos(π y)) + βsin(2π x)(1 - cos(2π y)). For generating input data, we employ an analytical form for the source function, defined as f(x, y) = 16βπ^2 (cos(4π y) sin(4π x) + sin(4π x)(cos(4π y) - 1)) - απ^2 (cos(π y) sin(π x) + sin(π x)(cos(π y) + 1)). This expression for f(x, y) is derived by substituting the analytical expression for u(x, y) = αsin(π x)(1 + cos(π y)) + βsin(2π x)(1 - cos(2π y)) into <ref>. Consequently, the solution implicitly satisfies the boundary conditions. The training instances are generated by varying the parameters α and β such that α∼Unif(-2, 2) and β∼Unif(-2, 2). The spatial resolution of 33 × 33 is taken for this example. 

Results: The <ref> shows the predictions obtained from the proposed NOGaP model. It displays predictions corresponding to one of the source functions. The plot consists of the ground truth obtained from the numerical solver and the mean predictions from the NOGaP framework. The plot suggests that the proposed framework can make predictions effectively. From <ref>, we can observe that the Relative Mean Square Error (RMSE) value is significantly improved over WNO. Thus, we can conclude that NOGaP is an improvement over the standalone WNO, GPR (zero mean) and B-WNO. Moreover, visually, the predicted and true solutions are almost indiscernible.


















§ CONCLUSION

This paper introduced a novel framework called NOGaP for effectively learning solutions to non-linear PDEs. By harnessing the capabilities of GP, the proposed framework not only provides accurate predictions but also provides a measure of uncertainty. By incorporating mean functions to be WNO, the NOGaP framework demonstrates effective performance in learning the solutions of non-linear PDEs.
Compared to a standalone WNO, GPR (zero mean) and B-WNO, the proposed framework offers distinct advantages in terms of improved prediction capabilities.

To validate the effectiveness of the framework, four numerical examples were considered, and the predictions were compared against the ground truth obtained using numerical solvers. The results demonstrate that the NOGaP framework achieves improved accuracy, with the predictions closely matching the ground truth solutions. The following observation has been made from the results obtained: 

    
  * The predictions obtained from the NOGaP framework demonstrate a strong agreement with the ground truth, providing a good approximation of the true solutions.
    
  * The selection of a suitable mean function poses a significant challenge in constructing effective frameworks involving GP. However, by incorporating the Wavelet Neural Operator (WNO) as the mean function, the NOGaP framework achieves improved accuracy over standalone WNO, GP (zero mean) and B-WNO.
    
  * Additionally, the NOGaP framework is convenient to implement and applicable to a wide variety of examples, making it a versatile tool for solving non-linear partial differential equations.


The results suggest that the proposed NOGaP can make predictions effectively and efficiently. We have used vanilla kernels as the choice of covariance function, but it will be interesting to note the effect of other possible kernels on the proposed NOGaP framework. Working with GPs is computationally expensive, especially with large and high dimensional datasets, thus further study can be done to bring down the computational cost.



§ ACKNOWLEDGEMENTS

SK acknowledges the support received from the Ministry of Education (MoE) in the form of Research Fellowship. RN acknowledge the financial support received from the Science and Engineering Research Board (SERB), India via grant no. SRG/2022/001410. SC acknowledge the financial support received from the Ministry of Port and Shipping via letter no. ST-14011/74/MT (356529). Faculty Seed grants received from IIT Delhi are also acknowledged.



§ CODE AVAILABILITY

On acceptance, all the source codes to reproduce the results in this study will be made available to the public on GitHub by the corresponding author.









10

li2020fourier
Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar.
Fourier neural operator for parametric partial differential equations.
arXiv preprint arXiv:2010.08895, 2020.

lu2019deeponet
Lu Lu, Pengzhan Jin, and George Em Karniadakis.
Deeponet: Learning nonlinear operators for identifying differential equations based on the universal approximation theorem of operators.
arXiv preprint arXiv:1910.03193, 2019.

tripura2023wavelet
Tapas Tripura and Souvik Chakraborty.
Wavelet neural operator for solving parametric partial differential equations in computational mechanics problems.
Computer Methods in Applied Mechanics and Engineering, 404:115783, 2023.

pathak2022fourcastnet
Jaideep Pathak, Shashank Subramanian, Peter Harrington, Sanjeev Raja, Ashesh Chattopadhyay, Morteza Mardani, Thorsten Kurth, David Hall, Zongyi Li, Kamyar Azizzadenesheli, Pedram Hassanzadeh, Karthik Kashinath, and Animashree Anandkumar.
Fourcastnet: A global data-driven high-resolution weather model using adaptive fourier neural operators, 2022.

lin2023spherical
Kenghong Lin, Xutao Li, Yunming Ye, Shanshan Feng, Baoquan Zhang, Guangning Xu, and Ziye Wang.
Spherical neural operator network for global weather prediction.
IEEE Transactions on Circuits and Systems for Video Technology, 2023.

Goswami_2022
Somdatta Goswami, Minglang Yin, Yue Yu, and George Em Karniadakis.
A physics-informed variational deeponet for predicting crack path in quasi-brittle materials.
Computer Methods in Applied Mechanics and Engineering, 391:114587, March 2022.

rani2023fault
Jyoti Rani, Tapas Tripura, Hariprasad Kodamana, Souvik Chakraborty, and Prakash Kumar Tamboli.
Fault detection and isolation using probabilistic wavelet neural operator auto-encoder with application to dynamic processes.
Process Safety and Environmental Protection, 173:215–228, 2023.

garg2023vb
Shailesh Garg and Souvik Chakraborty.
Vb-deeponet: A bayesian operator learning framework for uncertainty quantification.
Engineering Applications of Artificial Intelligence, 118:105685, 2023.

Magnani2022ApproximateBN
Emilia Magnani, Nicholas Kramer, Runa Eschenhagen, Lorenzo Rosasco, and Philipp Hennig.
Approximate bayesian neural operators: Uncertainty quantification for parametric pdes.
ArXiv, abs/2208.01565, 2022.

chen2021solving
Yifan Chen, Bamdad Hosseini, Houman Owhadi, and Andrew M Stuart.
Solving and learning nonlinear pdes with gaussian processes, 2021.

batlle2023kernel
Pau Batlle, Matthieu Darcy, Bamdad Hosseini, and Houman Owhadi.
Kernel methods are competitive for operator learning, 2023.

Chen2023SparseCF
Yifan Chen, Houman Owhadi, and Florian Schafer.
Sparse cholesky factorization for solving nonlinear pdes via gaussian processes.
ArXiv, abs/2304.01294, 2023.

n2023physics
Navaneeth N, Tapas Tripura, and Souvik Chakraborty.
Physics informed wno, 2023.

lu2021learning
Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and George Em Karniadakis.
Learning nonlinear operators via deeponet based on the universal approximation theorem of operators.
Nature Machine Intelligence, 3(3):218–229, 2021.

wang2021learning
Sifan Wang, Hanwen Wang, and Paris Perdikaris.
Learning the solution operator of parametric partial differential equations with physics-informed deeponets.
Science advances, 7(40):eabi8605, 2021.

hesthaven2018non
Jan S Hesthaven and Stefano Ubbiali.
Non-intrusive reduced order modeling of nonlinear problems using neural networks.
Journal of Computational Physics, 363:55–78, 2018.

bhattacharya2021model
Kaushik Bhattacharya, Bamdad Hosseini, Nikola B Kovachki, and Andrew M Stuart.
Model reduction and neural networks for parametric pdes.
The SMAI journal of computational mathematics, 7:121–157, 2021.

li2020neural
Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar.
Neural operator: Graph kernel network for partial differential equations, 2020.

li2022fourier
Zongyi Li, Daniel Zhengyu Huang, Burigede Liu, and Anima Anandkumar.
Fourier neural operator with learned deformations for pdes on general geometries, 2022.

li2024geometry
Zongyi Li, Nikola Kovachki, Chris Choy, Boyi Li, Jean Kossaifi, Shourya Otta, Mohammad Amin Nabian, Maximilian Stadler, Christian Hundt, Kamyar Azizzadenesheli, et al.
Geometry-informed neural operator for large-scale 3d pdes.
Advances in Neural Information Processing Systems, 36, 2024.

tripura2023foundational
Tapas Tripura and Souvik Chakraborty.
A foundational neural operator that continuously learns without forgetting, 2023.

thakur2023multifidelity
Akshay Thakur, Tapas Tripura, and Souvik Chakraborty.
Multi-fidelity wavelet neural operator with application to uncertainty quantification, 2023.

rani2024generative
Jyoti Rani, Tapas Tripura, Hariprasad Kodamana, and Souvik Chakraborty.
Generative adversarial wavelet neural operator: Application to fault detection and isolation of multivariate time series data, 2024.

navaneeth2023waveformer
N Navaneeth and Souvik Chakraborty.
Waveformer for modelling dynamical systems, 2023.

langehegermann2019algorithmic
Markus Lange-Hegermann.
Algorithmic linearly constrained gaussian processes, 2019.

swiler2020survey
Laura P Swiler, Mamikon Gulian, Ari L Frankel, Cosmin Safta, and John D Jakeman.
A survey of constrained gaussian process regression: Approaches and implementation challenges.
Journal of Machine Learning for Modeling and Computing, 1(2), 2020.

gulian2022gaussian
Mamikon Gulian, Ari Frankel, and Laura Swiler.
Gaussian process regression constrained by boundary value problems.
Computer Methods in Applied Mechanics and Engineering, 388:114117, 2022.

Yang_2021
Liu Yang, Xuhui Meng, and George Em Karniadakis.
B-pinns: Bayesian physics-informed neural networks for forward and inverse pde problems with noisy data.
Journal of Computational Physics, 425:109913, January 2021.

garg2022variational
Shailesh Garg and Souvik Chakraborty.
Variational bayes deep operator network: A data-driven bayesian solver for parametric differential equations, 2022.

zou2022neuraluq
Zongren Zou, Xuhui Meng, Apostolos F Psaros, and George Em Karniadakis.
Neuraluq: A comprehensive library for uncertainty quantification in neural differential equations and operators, 2022.

pförtner2023physicsinformed
Marvin Pförtner, Ingo Steinwart, Philipp Hennig, and Jonathan Wenger.
Physics-informed gaussian process regression generalizes linear pde solvers, 2023.

chen2023sparse
Yifan Chen, Houman Owhadi, and Florian Schäfer.
Sparse cholesky factorization for solving nonlinear pdes via gaussian processes, 2023.

besginow2022constraining
Andreas Besginow and Markus Lange-Hegermann.
Constraining gaussian processes to systems of linear ordinary differential equations, 2022.

bilionis2013multi
Ilias Bilionis, Nicholas Zabaras, Bledar A Konomi, and Guang Lin.
Multi-output separable gaussian process: Towards an efficient, fully bayesian paradigm for uncertainty quantification.
Journal of Computational Physics, 241:212–239, 2013.

bonilla2007multi
Edwin V Bonilla, Kian Chai, and Christopher Williams.
Multi-task gaussian process prediction.
Advances in neural information processing systems, 20, 2007.

alvarez2012kernels
Mauricio A. Alvarez, Lorenzo Rosasco, and Neil D. Lawrence.
Kernels for vector-valued functions: a review, 2012.

murphy2022probabilistic
Kevin P Murphy.
Probabilistic machine learning: an introduction.
MIT press, 2022.

Rasmussen2004
Carl Edward Rasmussen.
Gaussian Processes in Machine Learning, pages 63–71.
Springer Berlin Heidelberg, Berlin, Heidelberg, 2004.

bishop2006pattern
C Bishop.
Pattern recognition and machine learning.
Springer google schola, 2:531–537, 2006.

hendrycks2023gaussian
Dan Hendrycks and Kevin Gimpel.
Gaussian error linear units (gelus), 2023.






