Fourier Neural Operator for Parametric Partial Differential Equations
Authors: Zongyi Li
Nikola Borislavov Kovachki
Kamyar Azizzadenesheli
Burigede liu
Kaushik Bhattacharya
Andrew Stuart
Anima Anandkumar
Venue: ICLR 2021 Poster (2021)
Keywords: Partial differential equation
Fourier transform
Neural operators

Abstract: The classical development of neural networks has primarily focused on learning mappings between finite-dimensional Euclidean spaces.  Recently, this has been generalized to neural operators that learn mappings between function spaces. For partial differential equations (PDEs), neural operators directly learn the mapping from any functional parametric dependence to the solution. Thus, they learn an entire family of PDEs, in contrast to classical methods which solve one instance of the equation. In this work, we formulate a new neural operator by parameterizing the integral kernel directly in Fourier space, allowing for an expressive and efficient architecture. We perform experiments on Burgers' equation, Darcy flow, and Navier-Stokes equation. The Fourier neural operator is the first ML-based method to successfully model turbulent flows with zero-shot super-resolution. It is up to three orders of magnitude faster compared to traditional PDE solvers. Additionally, it achieves superior accuracy compared to previous learning-based solvers under fixed resolution.

Published as a conference paper at ICLR 2021
FOURIER NEURAL OPERATOR FOR PARAMETRIC PARTIAL DIFFERENTIAL EQUATIONS
Zongyi Li zongyili caltech.edu edu
Nikola Kovachki nkovachki @ caltech.edu
Kamyar Azizzadenesheli Kamyar@purdus edu
Burigede Liu bgl ca@aalcchheedu
Kaushik Bhattacharya
Andrew Stuart astuart@ caltech.edu
Anima Anandkumar anima @ @alltech.eu
ABSTRACT
The classical development of neural networks has primarily focused on learning mappings between finite-dimensional Euclidean spaces. Recently, this has been generalized to neural operators that learn mappings between function spaces. For partial differential equations (PDEs), neural operators directly learn the mapping from any functional parametric dependence to the solution. Thus, they learn an entire family PDEs, contrast classical methods which solve one instance of the equation. In this work, we formulate new neural operator by parameterizing the integral kernel directly in Fourier space, allowing for an expressive and efñi cient architecture. We perform experiments Burgers' equation, Darcy flow, and Savier-Stokes equation. The Fourier neural operator is the first ML-based method successfully model turbulent flows with zero-shot super- resorrresolution. It is up three orders of magnitude faster compared to traditional PDE solvers. Addition- ally, it achieves superior accuracy compared to previous learning- based solvers under fixed resolution.
INTRODUCTION
Many problems science and engineering involve solving complex partial differential equation (PDE) systems repeatedly for different values of some parameters. Examples arise in molecular dynamics, micconmme and turbulent flows. Often such systems require fine discretization in order to capture the phenomenon being modeled. As consequence, traditional numerical solvers are slow and sometimes inefficient. For example, when designing materials such as airfoils, one needs to solve the associated inverse problem where thousands of evaluations of the forward model are needed. A fast method can make such problems feasible.
Conventional solvers Data-driven methods. Traditional solvers such finite element meth- ods (FEM) and finite difference methods (FDM) solve the equation by discretizing the space. There- fore, they impose trade-off on the resolution: coarse grids are fast but less accurate; fine grids are accurate but slow. Complex PDE systems, as described above, usually require very fine discretiza- tion, and therefore very challenging and time-consuming for traditional solvers. On the other hand, data-driven methods can directly learn the trajectory of the family of equations from the data. As result, the learning-basee method can be orders of magnitude faster than the conventional solvers.
Machine learning methods may hold the key revolutionizing scientific disciplines by providing fast solvers that approximate or enhance traditional ones (Raissi et al., 2019: Jiang al., 2020. Greenfeld et al. 2019, Kochkov et al. (02TT) However, classical neural networks map between finite-dimensional spaces and can therefore only learn solutions tied to a specific discretization. This is often a limitation for practical applications and therefore the development of mesh-invarian neural networks is required. We first outline two mainstream neural network- baaser--based approaches for PDEs the finit-adimeniionll operators and NeMM
Finite-dimensional operators. These approaches parameterize the solution operator as a deep convolutional neural network between fiite-dimennioonl Fuclidean spaces Guuo al. (2016): Zhu

---

Published as a conference paper at ICLR 2021
Initial Vorticity
t=15
t=20
t=25 t=30
Prediction
Zero-shot super-resolution: Navier-Stoker Equation with Reynolds number 10000; Ground truth on top prediction on bottom; trained on 64 x 64 x 20 dataset; evaluated on 256 x 256 x 80 (see Section 5.4).
Figure I: top: The architecture of the Fourier layer; bottom: Example flow from Navier-Stokes.
& Zabaras 2018); Adler & Oktem (201)); Bhatnagar al. (2019); Khoo al. [2017) Such ap- proaches are, by definition, mesh-dependent and will need modifications and tuning for different resolutions and discretizations in order to achieve consistent error (if all possible). Furthermore, these approaches are limited to the discretization size and geometry of the training data and hence, possible query solutions new points in the domain. In contrast, show, for our method, both invariance of the error to grid resolution, and the ability to transfer the solution between meshes.
Neural-FEM. The second approach directly parameterizes the solution function neural net- work & Yu 2018 Raissi et al. 20199 Bar & Sochen 2019, Smith et al. 2020 Pan Duraisamy, 2020 This approach is designed to model one specific instance of the PDE. not the solution oper- ator. It is mesh-independent and accurate, but for any given new instance of the functional parame- ter/coefficient, it requires training methods such as finite elements, replacing the linear span of a finite set of local basis functions with the space of neural networks. The Neural-FEM approach suffers from the same computational issue as classical methods: the optimization problem needs to be solved for every new instance. Furthermore, the approach limited setting which the underlying PDE is known.
neural network. The approach closely resembles classical
Neural Operators. Recently, a new line of work proposed learning mesh-free, infinite- dimensional operators with neural networks (Lu et al.. 2019 Bhattacharya al. 2020 Nelsen & Stuart 2020; Li at 2020b] a Patel et al. (021). The neural operator remedies the mesh-dependent nature of the finitedimmensons. operator methods discussed above by producing single of net- work parameters that may be used with different discretizations. It has the ability to transfer solu- tions between meshes. Furthermore, the neural operator needs be trained only once. Obtaining a solution for new instance of the parameter requires only a forward pass of the network, allevi- ating the major computational issues incurred in Neural-FEM methods. Lastly, the neural operator requires no knowledge of the underlying PDE, only data. Thus far, neural operators have not yielded efficient numerical algorithms that can parallel the success of convolutional or recurrent neural net- works in the mnite-iimensionll setting due to the cost of evaluating integral operators. Through the fast Fourier transform, our work alleviates this issue.
Fourier Transform. The Fourier transform is frequently used in spectral methods for solving differential equations, since differentiation is equivalent multiplication in the Fourier domain. Fourier transforms have also played an important role in the development of deep learning. In theory, they appear the proof of the universal approximation theorem (Hornik al. 1989) and,
empirically, they have been used to speed up convolutional neural networks (Mathieu et al. 2013).
Neural network architectures involving the Fourier transform or the use of sinusoidal activation functions have also been proposed and studied (Bengio et al.) 2007 Mingo et al.) 2004 Sitzmann ct al. 2020). Recently, some spectral methods for PDEs have been extended to neural networks (Fan et al. 2019a Kashinath al. 2020). We build these works by proposing neural operator architecture defined directly in Fourier space with quasi-linear linear time complexity and oofte--f-th-aan
approximation capabilities.
2

---

Published as a conference paper at ICLR 2021
Our Contributions. We introduce the Fourier neural operator, novel deep learning architecture able to learn mappings between dinnit-oammeniinalll spaces of functions; the integral operator is restricted convolution, and instantiated through linear transformation in the Fourier domain.
The Fourier neural operator is the first work that learns the resolution-invariant solution operator for the family of Navier-Stokes equation in the turbulent regime, where previous graph-based neural operators do not converge.
By construction, the method shares the same learned network parameters irrespective of the dis- cretization used on the input and output spaces. can do zero-shot ruper-nesolution: trained on a lower resolution directly evaluated on higher resolution, as shown in Figure
The proposed method consistently outperforms all existing deep learning methods even when fixing the resolution to be 64 x 64. It achieves error rates that are 30% lower on Burgers' Equation, 60% lower on Darcy Flow, and 30% lower on Navier Stokes (turbulent regime with Reynolds number 10000). When learning the mapping for the entire time series, the method achieves < 1% error with Reynolds number 1000 and 8% error with Reynolds number 10000.
On a 256 x 256 grid, the Fourier neural operator has an inference time of only 0.005s compared to the 2.2s of the pseudo-spectral method used solve Navier-Stokes. Despite its tremendous speed advantage, the method does not suffer from accuracy degradation when used in downstream applications such as solving the Bayesian inverse problem, as shown in Figure 6
We observed that the proposed framework can approximate complex operators raising in PDEs that are highly non-linear, with high frequency modes and slow energy decay. The power of neural op- erators comes from combining linear, global integral operators (via the Fourier transform) and linear, local actívation functions. Similar to the way standard neural networks approximate highly non-linear functions by combining linear multiplications with non-linear activations, the proposed neural operators can approximate highly non-linear operators.
2 LEARNING OPERATORS
Our methodology learns a mapping between two infinite dimensional spaces from a finite collection of observed input-output pairs. Let D C Id be a bounded, open set and A A(D; Rd)) and == Rd.) be separable Banach spaces of function taking values in &e and RRd respectively. Furthermore let G A -> be (typically) non-linear map. We study maps Gt which arise as the solution operators of parametric PDEs see Section examples. Suppose we have observations u u, ) where aj ~ is an i.i.d. sequence from the probability measure μ supported on A and uj = (az) is possibly corrupted with noise. We aim to build an approximation of Gt by constructing a parametric map
A:Ax0 or equivalently,
A 0E
for some inite-dimensional parameter space θ by choosing Øt E θ so that G(., 01) Ga Gt. This is a natural framework for learning in infinite-dimensions as one could define a cost functional C: U x x IR and seek minimizer of the problem
min ManEEu CinnEee,C,CCGG. 6EB
which directly parallels the classical fnite-dimensional setting (Vapnik 1998). Showing the exis- tence of minimizers, the infiniteddimensional setting, remains a challenging open problem. We will approach this problem in the test-train setting by using data-driven empirical approximation the cost used to determine and to the accuracy of the approximation. Because concep- tualize our methodology in the definitedimmessional setting, all unite-dimensione approximations share a common set of parameters which are consistent in infinite dimensions. table of notation is shown in Appendix 3
Learning the Operator: Approximating the operator G is a different and typically much more challenging task than finding the solution E 2 of PDE for single instance of the parame- ter a E A. Most existing methods, ranging from classical finite elements, finite differences, and finite volumes to modern machine learning approaches such physics-informed neural networks
3

---

Published as a conference paper at ICLR 2021
a(x))
Fourier layer
Fourier layer
Fourier layer T
(b)
Fourier layer
F
V(X
(a) The full architecture of neural operator: start from input 1. Lift to higher dimension channel space by neural network P. 2. Apply four layers of integral operators and activation functions. 3. Project back the target dimension by a neural network Q. Output (b) Fourier layers: Start from input e. On top: apply the Fourier transform FF; linear transform on the lower Fourier modes and filters out the higher modes; then apply the inverse Fourier transform JF On the bottom: apply a local lincar transform W.
Figure 2: top: The architecture of the neural operators; bottom: Fourier layer.
(PINNs) (Raissi et al. 2019) aim at the latter and can therefore be computationally expensive. This makes them impractical for applications where solution to the PDE is required for many different instances of the parameter. On the other hand, our approach directly approximates the operator and therefore much cheaper and faster, offering tremendous computational savings when compared to traditional solvers. For an example application Bayesian inverse problems, see Section 5.5
Discretization. Since our data (j and uj are, in general, functions, to work with them numerically, assume access only to point-wise evaluations. Let Dj {. ..p C D be n-point discretization of the domain /) and assume we have observations ajlD, Raxds, Uz lD IR axd, for a finite collection of input-output pairs indexed by j. To be discretizationiinvaiant the neural operator can produce an answer u(x) for any D, potentially æ & Dj. Such a property is highly desirable as it allows a transfer of solutions between different grid geometries and discretizations.
3 NEURAL OPERATOR
The neural operator, proposed (Li et al. 2020b is formulated as an iterative architecture up + 01 + ... vy where vy for j Rd As shown in Figure (a), the input e A is first lifted to higher dimensional representation v0(2) = P(a(x)) by the local transformation P which is usually parameterized by a shallow fully- connected neural network. Then we apply several iterations of updates e, U4+1 (defined below). The output u(x) = (v)(()) is the projection or by the local transformation Q Rdd - Rdd.. In each iteration, the update V, + Ur4 is defined the composition of non-local integral operator K. and a local, nonlinear activation function σ.
1 is sequence of functions each taking values in
Definition (Iterative updates) Define the update to the representation
0t + utt1 by
((+1(x) = a (We(()) (= + (Kia; ()).)(x))
Va D
(2)
Rd)) maps bounded linear operators on UDD Rd:) and is parameterized by o θx. W Rd. + Rd. is a linear transformation, and a RR is non-linear activation function whose action is defined component-win
where K. Ax + CUUD;
We choose ((a; Φ) to be a kernel integral transformation parameterized by neural network.
Definition 2 (Kernel integral operator K.) Define the kernel integral operator mapping in by
(Kia; ())))() = /. n(2,, U. a(x), a(v): o),()d), Vx E D
4

---

Published as a conference paper at ICLR 2021
where Kg (2(()))) Rddx xdu is neural network parameterized by ø Okc.
Here Ki plays the role of kernel function which we learn from data. Together definitions 1 and 2 constitute a generalization of neural networks to infinit--dimensional spaces as first proposed in Li et al. (020b). Notice even the integral operator is linear, the neural operator can learn highly non-linear operators by composing linear integral operators with non-linear activation functions, analogous to standard neural networks.
If we remove the dependence on the function and impose Kg (2, v) (a y), we obtain that (3) a convolution operator, which is a natural choice from the perspective of fundamental solutions. We exploit this fact in the following section by parameterizing Ks directly in Fourier space and using the Fast Fourier Transform (FFT) to efficiently compute This leads to a fast architecture that
obtains state-of of-the-art results for PDE problems.
4 FOURIER NEURAL OPERATOR
We propose replacing the kernel integral operator in 3 by a convolution operator defined in Fourier space. Let F denote the Fourier transform of function f D Rd and F- its inverse then
(k) for 1, ,,..,, dy where = is is the imaginary unit. By letting Kg (Z, (0 a(r)) a(y)) kx(x in 3 and applying the convolution theorem, we find that ())4 (x) = F (F(r)) ((),
(F)) J. f3 fz()) (2,
- dx,
(F (1,00)
Vx E D.
We, therefore, propose to directly parameterize Ki in Fourier space.
Definition 3 (Fourier integral operator X.) Define the Fourier integral operator
(K(()),)xx)= = =F RR (Fu))
Vo E
(4)
where Rs is the Fourier transform of a periodic function 1s D - Rd. xd. parameterized by @ E Ox. An illustration is given in Figure (b).
For frequency mode i D, have (F) E Cd. and (k) E Cd.xd. Notice that since we assume k is periodic, admits a Fourier series expansion, so we may work with the discrete modes k E dd. We pick a dinite-dimensional parameterization by truncating the Fourier series at -
a maximal number of modes kimax = 12 = lik E 2d |kil I kmax.j, for l,......). d)|. We thus parameterize R. directly as complex-valuee
(Kmax x du x ) comprising a collection of truncated Fourier modes and therefore drop á from our notation. Since k is real-valued, we impose conjugate symmetry. We note that the set 7km.. is not the canonical choice for the low frequency modes of uy. Indeed, the low frequency modes are usually defined by placing an upper-bound on the fn-nrrm of 1: E Zd. We choose Zkm as above since it allows for an efficient implementation.
The discrete case and the FFT. Assuming the domain D is discretized with n E N points, we have that v. E Raxd, and F(v,) E Cnxd, Since we convolve v, with a function which only has kmax Fourier modes, we may simply truncate the higher modes to obtain F(vz) E Ch xd. Multiplication by the weight tensor R E Xdi
xd. xd= is then
de
FF.FF)) k.
Σ 1=1
k 1...., Kmax, j j= 1,...,,d de.
(5)
When the discretization is uniform with resolution 81 x... x X = n, F can be replaced by the Fast Fourier Transform. For f E Raxdu, (kI...., ka) E Z., x.. x...x22 Zn, and x (œ1 Ia) E D, the FFT ff and its inverse f are defined
#1-1 Σ
dd 2g fi Ia)e
(F))(K))
=
===0
1
25 Edh 2./2
(f-1))(()) = Σ
k: 0
0==0
5

---

Published as a conference paper at ICLR 2021
for 1 = 1, .... de. In this case, the set of truncated modes becomes
Kknm {Kk1..... ................................ .ua) E Ze,
Znd k; I kmax.j or 8j - k2 S kmax.j, for j l..... .d..
When implemented, R is treated as a (81 X ... X sd x de d,) -tensor and the above definition of Zy corresponds to the "corners" of R. which allows for straight-forward parallel implementation of via matrix-vector multiplication. In practice, we have found that choosing kaux.j = 12 which yields kmax = 12d parameters per channel to be sufficient for all the tasks that we consider.
Parameterizations of R. In general, R can be defined depend on (Fa) parallel (3). Indeed, can define Rs 2d x Rd, + Xd. as parametric function that maps (k, (F))(k) to the values of the appropriate Fourier modes. We have experimented with linear as well as neural network parameterizations of Rgg. We find that the linear parameterization has similar performance to the previously described direct parameterization, while neural networks have worse performance. This likely due to the discrete structure of the space Z Our experiments in this work focus the direct parameterization presented above.
Invariance to discretization. The Fourier layers are discretization-innvriinn because they can lear from and evaluate functions which are discretized in an arbitrary way. Since parameters are learned directly in Fourier space, resolving the functions in physical space simply amounts to pro- well-defined everywhere
jecting the basis 2 which zero-shot super-resolution as shown in Section 5.4] Furthermore, our architecture has a consistent error at any resolution of the inputs and outputs. On the other hand, notice that, in Figure 3 the standard CNN methods we compare against have an error that grows with the resolution.
Rd This allows
achieve
Quasi-linear complexity. The weight tensor contains kimax < n modes, so the inner multipli- cation has complexity O(kimux). Therefore, the majority of the computational cost lies in computing the Fourier transform F(ve) and its inverse. General Fourier transforms have complexity (nn²). however, since truncate the series the complexity is in fact O(nkmax), while the FFT has com- plexity On log n). Generally, we have found using FFTs to very efficient. However a uniform discretization is required.
5 NUMERICAL EXPERIMENTS
In this section, we compare the proposed Fourier neural operator with multiple architectures as well as operator-based approximation methods on the 1-d Burgers' equation, the 2-d Darcy Flow problem, and 2-d Navier-Stokes equation. The data generation processes are discussed and A.3.3 respectively. We not compare against traditional solvers (FEM/FDM) neural FEMM type methods since goal is to produce efficient operator approx- imation that can be used for downstream applications. We demonstrate one such application to the
Appendices I A.3.2
Bayesian inverse problem in Section 5.5
We construct our Fourier neural operator by stacking four Fourier integral operator layers as spec- ified in (2) and (4) with the ReLU activation well as batch normalization. Unless otherwise specified, we use N = 1000 training instances and 200 testing instances. We use Adam optimizer to train for 500 epochs with initial learning rate of 0.001 that is halved every 100 epochs. We set kmax. 16, d, = 64 for the 1-d problem and kmax.. = 12, d, 32 for the 2-d problems. Lower resolution data are downsampled from higher resolution. All the computation carried on a single Nvidia V 100 GPU with 16GB memory.
Remark on Resolution. Traditional PDE solvers such as FEM and FDM approximate a single function and therefore their error to the continuum decreases as the resolution is increased. On the other hand, operator approximation independent the ways data discretized as long as all relevant information is resolved. Resolution-invarian operators have consistent error rates among different resolutions as shown Figure Further, resolution-invarti operators can do zero-shot super-resolution, shown in Section 5.4
Benchmarks for time-independent problems (Burgers and Darcy): NN: a simple point-wise wise feedforward neural network. RBM: the classical Reduced Basis Method (using POD basis) (De-

---

Published as a conference paper at ICLR 2021
(a) Burger's 's Equation
(b) Darey Flow
(c) Navier-Stokes
10
ECN
KBM
FNO.3D INO
GUN PCANN
10
ECN
GND
U-Net
Relative entr
TNO
Relative citor
TNO MGNO
Relative crto
TT-Na
MGNO ENU
10
10
10
256
2044
4096
8192
85
141
211
421
100
200
300
400
500
Resolation
Resolution
Epochs
Left: benchmarks on Burgers equation; Mid: benchmarks on Darey Flow for different resolutions; Right: the learning curves Navier-Stokes 1e-3 with different benchmarks, Train and test on the same resolution.
For acronyms, see Section 5
details in Tables I
Figure 3: Benchmark on Burger's equation, Darey Flow, and Navier-Stokes
Vore 2014). FCN: a the-state-on of- ahe neural network architecture based on Fully Convolution
Networks (Zhu & Zabaras, 2018). PCANN: an operator method using PCA as an autoencoder on both the input and output data and interpolating the latent spaces with a neural network (Bhat- tacharya et al., [020). GNO: the original graph neural operator (Li al.] 2020b). MGNO: the multipole graph neural operator (Li et al. 2020a). LNO: neural operator method based on the
low-rank decomposition of the kernel K(Z, y) == Ej= (1)م, similar to the unstacked Deep-
ONet proposed in (Lu et al. 2019). FNO: the newly purposed Fourier neural operator.
Benchmarks for time-dependent problems (Navier-Stokes): ResNet: 18 layers of 2-d convo- lution with residual connections (He et al., 2016). U-Net: A popular choice for image-to to- to-image regression tasks consisting of four blocks with 2-d convolutions and deconvolutions (Ronneberger eet al. 2015). TF-Net: A network designed for learning turbulent flows based on a combination of spatial and temporal convolutions (Wang et al.) 2020). FNO-2d: 2-d Fourier neural operator with RNN structure in time. FNO-3d: 3-d Fourier neural operator that directly convolves in space-time.
5.1 BURGERS' EQUATION
The 1-d Burgers' equation is non-linear PDE with various applications including modeling the one dimensional flow of viscous fluid. It takes the form
0,u(x, t) + 3. (u² (x. 2) = vd uxx,t), 1(z 0) = un(
x E (0, e (0,
(6)
x (0,
with periodic boundary conditions where uh E I (0, 1); is the initial condition and v E RR the viscosity coefficient. We aim to learn the operator mapping the initial condition to the solution at time one, G ((0, 1); R) + ((,, 1); R) defined by ug + u(.,1) for any r >> 0.
The results of our experiments are shown in Figure 3 (a) and Table 3 (Appendix A.3. Our pro- posed method obtains the lowest relative error compared to any of the benchmarks. Further, the error is invariant with the resolution, while the error of convolution neural network based methods (FCN) grows with the resolution. Compared to other neural operator methods such as GNO and MGNO that use Nyström sampling in physical space, the Fourier neural operator is both more accurate and more computationally efficient.
5.2 DARCY FLOW
We consider the steady-state of the 2-d Darcy Flow equation on the unit box which is the second order, linear, elliptic PDE
V (a(z)Vu(z)) f(x)
E (0, 1)2
u(x)
E a(0, 1)²
7

---

Published as a conference paper at ICLR 2021
with Dirichlet boundary where ((0,1))) is the diffusion coefficient and f 1() is the foreing function. This PDE has numerous applications including modeling the pressure of subsurface flow, the deformation of linearly elastic materials, and the electric po- tential in conductive materials. We are interested in learning the operator mapping the diffusion coefficient to the solution, G I 1)2; IR+) ((, 1)2; R2) defined by + u. Note that
although the PDE is linear, the operator G' is not.
The results of our experiments are shown in Figure (b) and Table 4 (Appendix A.3.2 The pro- posed Fourier neural operator obtains nearly one order of magnitude lower relative error compared to any benchmarks. We again observe the invariance of the error with respect to the resolution.
5.3 NAVIER-STOKES EQUATION
We consider the 2-d Navier-Stokes equation for viscous, incompressible fluid in vorticity form on the unit torus:
0, t) + u(x,t) Vw(a, t) = vAw(e, f(z),
= (0,1) (0,7]
V u(x t) = O,
x E (0,1) E [0, x x (0,
(8)
w(x, 0) wo(x),
where u E C((0, 7]; H ((0, 1)); for any r > 0 is the velocity field, = V x u is the vor- 14 is the viscosity coefficient, and f e
ticity, wo E Lier (0, 1)²; R) is the initial vorticity, ((0, 1)2;R) is the forcing function. We are interested in learning the operator mapping the ticity up to time 10 to the vorticity up to some later time > 10, G C((O, 10]; ((,, 1)?;))) C((10, T); ((0, IR): defined by w| (011))xx + w| (0.11)x easy to derive the velocity. While vorticity harder model compared to velocity, provides more information. By formulating the problem on vorticity, the neural network models mimic the pseudo- spectral method. We experiment with the viscosities v 1c 1c-3, 3, Ic 4, 1,-5, decreasing the final time 7 as the dynamic becomes chaotic. Since the baseline methods are not resolution-invariant, we fix the resolution to be 64 64 for both training and testing.
Given the vorticity it is
Table 1: Benchmarks on Navier Stokes (fixing resolution 64 64 for both training and testing)
Parameters
Time per epoch
1c-3
1e-4
== 1e-4
1e-5
Config
T 50
T = 30
T 30
T 20
N 1000
N 1000
N 10000
N 1000
FNO-3D
6, 558, 537
38.99s
0.0086
0.1918
0.0820
0.1893
FNO-2D
414,517 24,950, 491
127.80s
0.0128
0.1559
0.0834
0.1556
U-Net
48.67s
0.0245
0.2051
0.1190
0.1982
TF-Net
7,451, 724
47.21s
0.0225
0.2253
0.1168
0.2268
ResNet
266, 641
78.47s
0.0701
0.2871
0.2311
0.2753
As shown in Table 1 the FNO-3D has the best performance when there is sufficient data (v 1e-3, N 1000 and y 1e = 10000). For the configurations where the amount of data insufficient
1e-4, . 1000 and y 1e-5, N 1000), all methods have 15% error with FNO-2D achieving the lowest. Note that we only present results for spatial resolution 64 x 64 since all benchmarks we compare against are designed for this resolution. Increasing it degrades their performance while FNO achieves the same errors.
2D and 3D Convolutions. FNO-2D, U-Net, TF-Net, and ResNet all do 2D-convolution in the spatial domain and recurrently propagate in the time domain (2D+RNN). The operator maps the solution at the previous 10 time steps to the next time step (2D functions to 2D functions). On the other hand, FNO-3D performs convolution in space-time. It maps the initial time steps directly to the full trajectory (3D functions to 3D functions). The 2D+RNN structure can propagate the solution to any arbitrary time T in increments of fixed interval length At, while the Conv3D structure is fixed to the interval [0, T] but can transfer the solution to an arbitrary time-discretizatin. We find the 3-d method to be more expressive and easier train compared to its RNN-structured counterpart.
8

---

Published as a conference paper at ICLR 2021
5.4 ZERO- SHOT SUPER- REEEERRRSOLUION.
The neural operator is invariiva so it can be trained on lower resolution and evaluated at higher resolution, without seeing any higher resolution data (zero-shot super-resolution). Figure shows an example where we train the FNO-3D model on 64 x 64 x 20 resolution data in the setting above with (v 1e-4, = 10000) and transfer to 256 x 256 80 resolution, demonstrating super-resolution in space-time. Fourier neural operator is the only model among the benchmarks (FNO-2D, U-Net, TF-Net, and ResNet) that can do zero-shot super-resolution. And surprisingly. can do super-resolution not only in the spatial domain but also in the temporal domain.
5.5 BAYESIAN INVERSE PROBLEM
In this experiment, we use a function space Markov chain Monte Carlo (MCMC) method (Cotter et al. 2013) to draw samples from the posterior distribution of the initial vorticity in Navier-Stokes given sparse, noisy observations at time T 50. We compare the Fourier neural operator acting as a surrogate model with the traditional solvers used to generate our train-test data (both run on GPU). We generate 25,000 samples from the posterior (with a 5,000 sample burn-in period), requiring 30,000 evaluations of the forward operator.
As shown in Figure C (Appendix A.5). FNO and the traditional solver recover almost the same posterior mean which, when pushed forward, recovers well the late-time dynamic of Navier Stokes. In sharp contrast, FNO takes 0.005s to evaluate a single instance while the traditional solver, after being optimized to use the largest possible internal time-step which does not lead to blow-up, takes 2.2s. This amounts to 2.5 minutes for the MCMC using FNO and over 18 hours for the traditional solver. Even if we account for data generation and training time (offline steps) which take 12 hours, using FNO is still faster! Once trained, FNO can be used to quickly perform multiple MCMC runs for different initial conditions and observations, while the traditional solver will take 18 hours for every instance. Furthermore, since FNO is differentiable, it can easily be applied to PDE-constrained optimization problems without the need for the adjoint method.
Spectral analysis. Due the way we parameterize Rg, the function output by 4 has at most kmax.. Fourier modes per channel. This, however, does not mean that the Fourier neural operator can only approximate functions up to kmax. modes. Indeed, the activation functions which oceur between integral operators and the final decoder network Q recover the high frequency modes. As an example, consider a solution to the Navier-Stokes equation with viscosity v = 1e-3. Truncating this function at 20 Fourier modes yields an error around 2% while our Fourier neural operator learns the parametric dependence and produces approximations to an error of < 1% with only kmax. = 12 parameterized modes.
Non-periodis boundary condition. Traditional Fourier methods work only with periodic bound- ary conditions. However, the Fourier neural operator does not have this limitation. This is due to the linear transform W (the bias term) which keeps the track of non-periodic boundary. As example, the Darcy Flow and the time domain of Navier-Stoker have non-periodic boundary conditions, and the Fourier neural operator still learns the solution operator with excellent accuracy.
6 DISCUSSION AND CONCLUSION
Requirements on Data. Data-driven methods rely on the quality and quantity of data. To learn Navier-Stokes equation with Reynolds number Re 1e+4, we need to generate N 10000 training pairs {ay, 4,} with the numerical solver. However, for more challenging PDEs, generating a few training samples can be already very expensive. future direction is to combine neural operators with numerical solvers to levitate the requirements on data. Recurrent structure. The neural operator has an iterative structure that can naturally be formulated as a recurrent network where all layers share the same parameters without sacrificing performance. (We did not impose this restriction in the experiments.) Computer vision. Operator learning is not restricted to PDEs. Images can naturally be viewed as real-valued functions on 2-d domains and videos simply add a temporal structure. Our approach is therefore a natural choice for problems in computer vision where invariance discretization crucial is important (Chi et al. 2020

---

Published as a conference paper at ICLR 2021
ACKNOWLEDGEMENTS
The authors want to thank Ray Wang and Rose Yu for meaningful discussions. Z. Li gratefully acknowledges the financial support from the Kortschak Scholars Program. Anandkumar is sup- ported in part by Bren endowed chair, LwLL grants, Beyond Limits, Raytheon, Microsoft, Google, Adobe faculty fellowships, and DE Logi grant. K. Bhattacharya, N. B. Kovachki, Liu, and A. M. Stuart gratefully acknowledge the financial support of the Army Research Laboratory through Cooperative Agreement Number W91 1011NN-12-0022 Research was sponsored by the Army Re- search Laboratory and accomplished under Cooperative Agreement Number W91 12-2- 0022. The views and conclusions contained in this document are those of the authors and should not be interpreted as representing the official policies, either expressed or implied, of the Army Research Laboratory or the U.S. Government. The U.S. Government is authorized to reproduce and distribute reprints for Government purposes notwithstanding any copyright notation herein.
REFERENCES
Jonas Adler and Ozan Oktem. Solving ill-posed inverse problems using iterative deep neural networks. Inverse Problems, nov 2017. doi: 11.188//1611-6422/aa9888 //db..oog/1.1108822/7. 200.
URL https:
Leah Bar and Nir Sochen. Unsupervised deep learning algorithm for pde-based forward and inverse problems. arXiv preprint arXiv:1904.0544 2019.
Yoshua Bengio, Yann LeCun, et al. Scaling learning algorithms towards ai. Large-scale kernel
machines, 34(5):1-41 1-11) 2007.
Saakaar Bhatnagar, Yaser Afshar, Shaowu Pan, Karthik Duraisamy, and Shailendra Kaushik. Predic- tion of aerodynamic flow fields using convolutional neural networks. Computational Mechanics, pp. 1-21, 2019.
Kaushik Bhattacharya, Nikola Kovachki, and Andrew M. Stuart. Model reduction and neural networks for parametric pde(s). preprint, 2020.
Lu Chi, Borui Jiang, and Yadong Mu. Fast fourier convolution. Advances in Neural Information Processing Systems, 33, 2020.
S. L. Cotter, G. O. Roberts, A. M. Stuart, and D. White. Meme methods for functions: Modifying old algorithms to make them faster. Staristical Science, 28(3):424 446, Aug 2013. ISSN 0883-4237. doi: 11t211411-- URL http: doi wtt.:/...ooo..1111111111111-1-10
Ronald A. DeVore. Chapter 3: The Theoretical Foundation of Reduced Basis Methods. 2014. doi:
10.113771.9781611197488
[RL/pttps://ppubss sianuussssia...
hess//oi/aaa /bt///b///ii orgg/ei/abs///aaaeeeeor 1137/
1. 1.9781611974823.chm
Weinan E and Bing Yu. The deep ritz method: A deep learning-based numerical algorithm for solving variational problems. Communications in Mathematics and Statistics, 3 2018. ISSN 2194-6701. doi: 10.10077/443300001---107--
Yuwei Fan, Cindy Orozco Bohorquez, and Lexing Ying. Ber-net: A neural network based on the nonstandard wavelet form. Journal of Computational Physics, 384:1-15, 2019a.
Yuwei Fan, Lin Lin, Lexing Ying, and Leonardo Zepeda-Núnez. A multiscale neural network based hierarchical matrices. Multiscale Modeling & Simulation, 17(4): 18931 2019b.
Daniel Greenfeld, Meirav Galun, Ronen Basri, Irad Yavneh, and Ron Kimmel. Learning to optimize multigrid pde solvers. In International Conference Machine Learning, pp. 2019.
PMLR,
Xiaoxiao Guo, Wei Li, and Francesco Iorio. Convolutional neural networks for steady flow ap- proximation. In Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, 2016.
10

---

Published as a conference paper at ICLR 2021
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep residual learning for image recog- nition. In Proceedings of the IEEE conference on computer vision and pattern recognition, pp- 770-778, 2016.
Kurt Homik, Maxwell Stinchcombe, Halbert White, et al. Multilayer feedforward networks are universal approximators. Neural networky, 2(5):359-366, 1989.
Chiyu Max Jiang, Soheil Esmaeilzadeh, Kamyar Azizzadenesheli, Karthik Kashinath, Mustafa Mustafa, Hamdi A Tchelepi, Philip Marcus, Anima Anandkumar, et al. Meshfreeflownet: A physics-constrained deep continuous space-time super-resolution framework. arXiv preprini arXiv:2005.01463, 2020.
Karthik Kashinath, Philip Marcus, et al. Enforcing physical constraints in cnns through differen- tiable pde layer. In ICLR 2020 Workshop on Integration of Deep Neural Models and Differential Equations, 2020.
Yuehaw Khoo, Jianfeng Lu, and Lexing Ying. Solving parametric PDE problems with artificial neural networks. arXiv preprint arXiv: 2017.
Dmitrii Kochkov, Jamie A Smith, Ayya Alieva, Qing Wang, Michael P Brenner, and Stephan Hoyer. Machine learning accelerated computational fluid dynamics. arXiv preprint arXiv:2102.01010, 2021.
Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, An- drew Stuart, and Anima Anandkumar. Multipole graph neural operator for parametric partial differential equations, 2020a.
Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, An- drew Stuart, and Anima Anandkumar. Neural operator: Graph kernel network for partial differ- ential equations. arXiv preprint arXiv:2003.0344 03485, 2020b.
Lu Lu, Pengzhan Jin, and George Em Karniadakis. Deeponet: Learning nonlinear operators for iden- tifying differential equations based on the universal approximation theorem operators. arXiv preprint arXiv: 1XXv:1900.03193, 2019.
Michael Mathieu, Mikael Henaff, and Yann LeCun, Fast training of convolutional networks through ffts, 2013.
Luis Mingo, Levon Aslanyan, Juan Castellanos, Miguel Diaz, and Vladimir Riazanov. Fourier neural networks: An approach with sinusoidal activation functions. 2004.
NH Nelsen and AM Stuart. The random feature model for input-output maps between banach spaces. arXiv preprint vrXiv:2005.10224 2020.
Shaowu Pan and Karthik Duraisamy. Physics-informee probabilistic learning of linear embeddings of nonlinear dynamics with guaranteed stability. SIAM Journal on Applied Dynamical Systems, 19(1):480-509 2020.
Ravi G Patel, Nathaniel A Trask, Mitchell A Wood, and Eric C Cyr. A physics-informed operator gression framework for extracting data-driven continuum models. Computer Methods in Applied Mechanics and Engineering, 373:1 13500, 2021.
Maziar Raissi, Paris Perdikaris, and George Karniadakis. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. Journal of Computational Physics, 378:686-707, 2019.
Olaf Ronneberger, Philipp Fischer, and Thomas Brox. U-net: Convolutional networks for biomedi- cal image segmentation. In International Conference on Medical image computing and computer- assisted intervention, pp. 234-241. Springer. 2015.
Vincent Sitzmann, Julien NP Martel, Alexander W Bergman, David B Lindell, and Gordon Wetzstein. Implicit neural representations with periodic activation functions. arXiv preprint arXiv:2006.09665 2020.
11

---

Published as a conference paper at ICLR 2021
Jonathan D Smith, Kamyar Azizzadenesheli, and Zachary E Ross. Eikonet: Solving the eikonal equation with deep neural networks. arXiv preprint arXiv:2004.00361, 2020.
Vladimir N. Vapnik. Statistical Learning Theory. Wiley-Interscience,, 1998.
Rui Wang, Karthik Kashinath, Mustafa Mustafa, Adrian Albert, and Rose Yu. Towards physics- informed deep learning for turbulent flow prediction. In Proceedings of the 26th ACM SIGKDD International Conference Knowledge Discovery & Data Mining, pp. 1457-1466, 2020.
Yinhao Zhu and Nicholas Zabaras. Bayesian deep convolutional encoder-decoder networks for surrogate modeling and uncertainty quantification. Journal of Computational Physics, 2018. ISSN 0021-9991. doi: http://kko.oog////.0001ooooo2220 URL http:://ww. sciencedirect. htt.e//ww....a..atatoooooooo///a
12

---

Published as a conference paper at ICLR 2021
A APPENDIX
A.1 TABLE OF NOTATIONS
A table of notations is given in Table 2
Table 2: table of notations
Notation
Meaning
Operator learning
D C Rd
The spatial domain for the PDE Points in the the spatial domain The input coefficient functions The target solution functions The discretization of (dj, uÿ)
x D a 8 (D:; Rd.))
uEU- (D:Rd..
Dj
G A +
The operator mapping the coefficients to the solutions
1
A probability measure where a, sampled from.
Neural operator
v(x) E Rd.
The neural network representation of u(x)
da
Dimension of the input a(x).
du
Dimension of the output u(z).
d,
The dimension of the representation e(z)
2(d+1)
1dd xd.
The kereel maps (z, y, a(x), a(y)) to a d, d, matrix
The parameters of the kernel network K
O,....7 T
The time steps (layers) The activation function
0
Fourier operator
F,F
Fourier transformation and its inverse.
The linear transformation applied on the lower Fourier modes.
W
The linear transformation (bias term) applied on the spatial domain. Fourier modes wave numbers. The max Fourier modes used in the Fourier layer.
k
kmaz
Hyperparameters
N
The number of training pairs. The size of the discretization. The resolution the discretization (sd nt). The viscosity.
T
The time interval [0, 2" for time-dependent equation.
A.2 SPECTRAL ANALYSIS
The spectral decay of the Navier Stokes equation data is shown in Figure 4 The spectrum decay has slope k 5/33 matching the energy spectrum in the turbulence region. And we notice the energy spectrum does not decay along with time.
We also present the spectral decay in of the truncation kmaz defined in aas shown in Figure5 We note all equations (Burgers, Darcy, and Navier-Stokes with S 1e-4 ) exhibit high frequency modes. Even we truncate at kmax 12 in the Fourier layer, the Fourier neural operator can recover the high frequency modes.
A.3 DATA GENERATION
In this section, we provide the details of data generator for the three equation used in Section 5
13

---

Published as a conference paper at ICLR 2021
10
1-50
1690
101
10
12-
33"
16
11
60
16
12'
12
The spectral decay of the Navier-stokes equation data we used in section 5.3 The y-axis is the spectrum: the x-axis is the wavenumber [k] kn
lrz.
Figure 4: Spectral Decay of Navier- Stokes equations
10
Error of truncation
Darcy
-011
NS 1c 3
1e-d
NS le:5
0.2
:0
:0:4
10
15
20
25
Fourler modes per dimension
The ertor of truncation in one single Fourier layer without applying the linear transform R. The y-axis is the normalized truncation error; the x-axis is the truncation mode kmaz-
Figure 5: Spectral Decay in term of kmaz
A.3.1 BURGERS EQUATION
Recall the 1-d Burger's equation on the unit torus:
del(z, t) (x, (2)
vd ((,
(0,1), E (0, x E (0, 1).
u(x 0) = (z))
The initial condition uo(x) is generated according un ~ μ where μ NO, 625 251) with periodic boundary conditions. We set the viscosity to = 0.1 and solve the equation using split step method where the heat equation part solved exactly in Fourier space then the non-linear part is advanced, again in Fourier space, using very fine forward Euler method. We solve on spatial mesh with resolution 2 8192 and use this dataset subsample other resolutions.
A.3.2 DARCY FLOW
The 2-d Darcy Flow is a second-order linear elliptic equation of the form
V (a(x)) = f(x) u(x) 0
(0, 1)² 000,
The coefficients a(z) are generated according to where μ ( (-A 97) with zero Neumann boundary conditions on the Laplacian. The mapping v R takes the value 12 on the positive part of the real line and 3 on the negative and the push-forward is defined pointwise. The forcing is kept fixed f(x) = 1. Such constructions are prototypical models for many physical systems such as permeability in subsurface flows and material microstructures in elasticity. Solutions u are obtained by using a recond-order finite difference scheme on a 421 421 grid. Different resolutions downsampled from this dataset.
14

---

Published as a conference paper at ICLR 2021
A.3.3 NAVIER-STOKES EQUATION
Recall the 2-d Navier-Stokes equation for viscous, incompressible fluid in vorticity form on the unit torus:
O, u(x, Vu(r vAw(e, + I(z), u(2 t) 0, w(2, 0) = wo(r),
x x (0,1)), (0,
x x (0, 1), [0.T]
# E (0,1)²
The initial condition wo(z) is generated according to wn μ where μ N/ (--
with periodic boundary conditions. The forcing is kept fixed f(x) = 0. (21(si((27(1) + x2)) + cos(2 + x2))). The equation is solved using the stream-function formulation with pseu- dospectral method. First a Poisson equation is solved in Fourier space to find the velocity field. Then the vorticity is differentiated and the non-lincar term is computed is physical space after which is dealiased. Time advanced with -ricol-Ni update where the non-linear term does not enter the implicit part. All data are generated We use a time-step of le-4 for the Crank-Nicolson scheme in the data-generated process where we record the solution every L = 1 time units. The step increased 2e-2 when used in MCMC for the Bayesian inverse problem.
256 x 256 grid and downsampled to 64 64.
A.4 RESULTS OF BURGERS' EQUATION AND DARCY FLOW
The details error rate on Burgers' equation and Darcy Flow are listed in Table 3 and Table
Table 3: Benchmarks on 1-d Burgers' equation
Networks
256
512
1024
8 2048
4096
8192
NN
0.4714
0.4561
0.4803
0.4645
0.4779
0.4452
GCN
0.3999
0.4138
0.4176
0.4157
0.4191
0.4198
FCN
0.0958
0.1407
0.1877
0.2313
0.2855
0.3238
PCANN
0.0398
0.0395
0.0391
0.0383
0.0392
0.0393
GNO
0.0555
0.0594
0.0651
0.0663
0.0666
0.0699
LNO
0.0212
0.0221
0.0217
0.0219
0.0200
0.0189
MGNO
0.0243
0.0355
0.0374
0.0360
0.0364
0.0364
FNO
0.0149
0.0158
0.0160
0.0146
0.0142
0.0139
Table 4: Benchmarks on 2-d Darcy Flow
Networks
85
141
211
421
NN
0.1716
0.1716
0.1716
0.1716
FCN
0.0253
0.0493
0.0727
0.1097
PCANN
0.0299
0.0298
0.0298
0.0299
RBM
0.0244
0.0251
0.0255
0.0259
GNO
0.0346
0.0332
0.0342
0.0369
LNO
0.0520
0.0461
0.0445
MGNO
0.0416
0.0428
0.0428
0.0420
FNO
0.0108
0.0109
0.0109
0.0098
A.5 BAYESIAN INVERSE PROBLEM
Results of the Bayesian inverse problem for the Navier-Stokes equation shown in Figure 6 can be seen that the result using Fourier neural operator as a surrogate is as good as the result of the traditional solver.
15

---

Published as a conference paper at ICLR 2021
F.4
The top left panel shows the true initial vorticity while bottom left panel shows the true observed vorticity at T 50 with black dots indicating the locations of the observation points placed 7 x grid. The top middle panel shows the posterior mean of the initial vorticity given the noisy observations estimated with MCMC using the traditional solver, while the top right panel shows the same thing but using FNO surrogate model. The bottom middle and right panels show the vorticity at T = 50 when the respective approximate posterior means are used as initial conditions.
Figure 6: Results of the Bayesian inverse problem for the Havier-Stoker equation.
16


Reviews:
[{"review_id": "wyQ9GjXPhcR", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This work attacks the important problem of learning PDEs with neural networks. To this extent, it proposes a neural network architecture where (in essence), the linearity is not applied in the Fourier space but in the state space. This allows the network to implicitly learn the partial differentials of the equation, and results in a parametrization which is invariant to spatial discretization. This approach yields good results when regressing to data derived from numerical solutions of complex differential equations. \n\n\nStrengths:\n\nOverall, the paper is well written. The theoretical and experimental sections are, for the most part, clear and concise, although some important details remain unclear/lacking. \n\nThe method seems to be quite generic and can be applied to a large range of PDEs. \n\nThe inverse problem experiment is interesting and highlights a potentially useful application of the proposed method.\n\n\nWeaknesses:\n\nPerhaps most importantly, the data generation is not always clear to me. The distinction between training data and test data is not clearly specified. Generalization to different useful initial conditions is obviously of utmost importance in order to properly evaluate the quality of data-driven methods for dynamical systems. \n\nThe loss function and training method is not clearly specified: are observations acquired every t=k \\in N steps for NS? Do you use teacher forcing or compute an error on the full sequence? What about the other experiments?\n\n\nComments and questions:\n\nIt seems that in the proposed model, the Fourier coefficients R do  not vary with the layers. This seems to be quite constraining, could you comment on this?\n\nIt would have been interesting to analyze the expressiveness of the architecture. \n\nYou mention that you are not restricted to periodic boundary conditions due to your final decoder (which is never clearly defined), however, does the finite support of the data not lead to artifacts in the Fourier transform, causing you to have to learn a high number of Fourier modes? Could you comment on this? \n\nIn the NS experiment, why have you used N=10 000 when viscosity is set to 1e-4, and N=1000 for the others? It seems as if you have selected the number of samples where your method outperforms the baselines.\n\nMoreover, for the NS experiments, test data seems to be future data (T>10), however this seems to not be the case for the Burgers’ equation, even though it is also time dependent.\n\nAs one of your motivations behind this work is to learn the operator, it could have been interesting to test your approach using different sample points as in the training data. \n", "extra_scores": {}, "extra_text": {}}, {"review_id": "KBqilna8RHb", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The manuscript proposes a neural network for solving partial differential equations. This is performed by a supervised learning of a mapping between the PDE coefficients (a) and the PDE solution (u), where a is drawn from a certain distribution. This work follows  a work by Li et al. 2020b.  Inspired by the properties of Green’s function, in Li et al. the authors suggest to learn a kernel function in an iterative manner.  In this work, the general kernel is replaced by a convolution kernel which is represented in the Fourier space. \n\n\n\n\nPros. \n\n1.\tThe authors address an important and practical problem \n2.\tI like the idea of learning a mapping for a class of PDEs, rather than optimizing per instance\n3.\tThe numerical results are impressive\n\nCons. \n1.\tNovelty. I feel that the novelty is incremental with respect to Li et al. The motivation for replacing the kernel function (3) by a convolution operator is not entirely clear. I guess that efficiency is one of the reasons. If so, I would suggest to discuss the complexity of the proposed scheme vs. Li et al. complexity.  Moreover, how does it affect the expressivity of the network? \n2.\tClarity of the paper.  I didn’t like the idea that on the one hand entire sentences are taken form Li et al. but on the other hand, the paper itself is not self-contained. In practice, the reader should follow the intuition and the motivation given in Li et al in order to grasp the idea of neural operator.  In addition, a concise detailed description of the network architecture is missing. \n3.\tMissing references and comparison.  There are two works, aiming at solving PDEs, which are very relevant in the context of the proposed approach. The first copes with the structured case, Greenfeld at al. “Learning to optimize multigrid PDE solvers” and the second handles the unstructured case, Luz et al. “Learning algebraic multigrid using graph neural networks”. It resembles the proposed approach in the sense that the training is performed for a class of problems, yielding a mapping from the coefficients to the solution. In these works, the learning is unsupervised, while generalizing over boundary conditions and domains.\nTo summarize, I believe that this work should be published. I hope that the authors are able to address my concerns. \n", "extra_scores": {}, "extra_text": {}}, {"review_id": "HGDtY8le-kE", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Paper Summary: The authors proposed a novel neural Fourier operator that generalizes between different function discretization schemes, and achieves superior performance in terms of speed and accuracy compared to learned baselines.\n\nI have to admit that my understanding of this paper is rather limited and I have lots of conceptual questions about the implementation.\n- The authors claimed to novel contributions: one is that the method is fast, one is that the method is discretization invariant. These two seem to be in conflict. The authors leveraged FFT for performing fast convolutions differentiably. However FFT assumes the data to be present on a uniform grid (of 2^n grid cells in each dimension). I understand that since the learned content is a frequency \"filter\" that can be applied to any discretization, but this can be said for just about any CNN (in that case a learned 3x3 filter in physical space can also just be projected to the spectral space and work with \"any discretization\"). Furthermore (correct me if I'm wrong), it doesn't seem that the authors experimentally illustrated a case of learning on non-discretized data, and performing inference on discretized versions of it, or vice versa.\n- From Fig. 1, it seems that the actual parameters that are being learned is the convolutional kernel R (represented as spectral coefficients), and W. That basically amounts to a single convolution plus bias operation. The expressivity of a single such operation should be pretty limited, since this operation itself is a linear operator, while the underlying equations that the authors demonstrated on are mostly nonlinear. \n\nMore insights into the two questions above are certainly appreciated.\n\nAdding a reference that seems related to this work. The paper below uses a spectral solver step as a differentiable layer within the neural network for enforcing hard linear constraints in CNNs, also taking the FFT -> spectral operator -> IFFT route.\n\nReferences\n---------------\nKashinath, Karthik, and Philip Marcus. \"Enforcing Physical Constraints in CNNs through Differentiable PDE Layer.\" ICLR 2020 Workshop on Integration of Deep Neural Models and Differential Equations. 2020.\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "GMj3B9nCFzC", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Paper summary:\n\nBuilding on previous work on neural operators, the paper introduces the Fourier neural operator, which uses a convolution operator defined in Fourier space in place of the usual kernel integral operator. Each step of the neural operator then amounts to applying a Fourier transform to a vector (or rather, a set of vectors on a mesh), performing a linear transform (learnt parameters in this model) on the transformed vector, before performing an inverse Fourier transform on the result, recombining it with a linear map of the original vector, and passing the total result through a non-linearity. The Fourier neural operator is by construction (like all neural operators) a map between function spaces, and invariance to discretization follows immediately from the nature of a Fourier transform (just project onto the usual basis). If the underlying domain has a uniform discretization, the fast Fourier transformation (FFT) can be used, allowing for an O(nlogn) evaluation of the aforementioned convolution operator, where n is the number of points in the discretization. Experiments demonstrate that the Fourier neural operator significantly outperforms other neural operators and other deep learning methods on Burgers’ equation, Darcy Flow, and Navier Stokes, and that that it is also significantly faster than traditional PDE solvers.\n\n------------------------------------------\nStrengths and weaknesses:\n\nMuch of the theoretical legwork for this paper, namely, neural operators, was already carried out in previous papers (Li et al.). The remaining theoretical work, namely writing down the Fourier integral operator and analysing the discrete case, was succinctly explained. The subsequent experimentation was extremely thorough (e.g. demonstrating that activation functions help in recovering high frequency modes) and, of course, the results were very impressive. I liked the paper a lot, and it’s definitely a big step-forward in neural operators. I’m assigning a score of 8 (a very good conference paper), and I think that the paper is more or less ready for publication as is. I’ve included a few questions below (to help my own understanding), as well as some typos I spotted whilst reading the paper.\n\n------------------------------------------\nQuestions and clarification requests:\n\n1)\tSection 4, The Discrete Case and the FFT – could you explain the definition of bounds in the definition of Z_{k_{max}}?\n2)\tSection 4, Parametrizations of R, sentence 2 – could you explain the definition R_{\\phi}? At present I can’t see how the function signature of R matches the definition given.\n------------------------------------------\nTypos and minor edits:\n- Page 3, bullet point 3 – “solving Bayesian inference problem” -> “solving Bayesian inference problems”\n- Section 1, final paragraph, sentence 2 - “approximate function with any boundary conditions” -> “approximate functions with any boundary conditions”\n- Section 4, The discrete case and the FFT, final paragraph, last sentence - “all the task that we consider” -> “all the tasks that we consider”\n- Section 4, Parametrizations of R, last sentence - “while neural networks have the worse performance” -> “while neural networks have the worst performance”\n- Section 4, final sentence – “Generally, we have found using FFTs to be very efficient, however a uniform discretization if required.” -> “Generally, we have found using FFTs to be very efficient. However, a uniform discretization is required.”\n- Section 5, final paragraph, sentence 2 – “FNO takes 0.005s to evaluate a single instances while the traditional solver” -> “FNO takes 0.005s to evaluate a single instance while the traditional solver”\n- Section 6, final sentence – “Traditional Fourier methods work only with periodic boundary conditions, however, our Fourier neural operator does not have this limitation.” -> “Traditional Fourier methods work only with periodic boundary conditions. However, our Fourier neural operator does not have this limitation.”", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: Pros:
- Provides a practical technique which can dramatically speed up PDE solving -- this is an important and widely applicable contribution.
- Paper is simultaneously clearly written and mathematically sophisticated.
- The experimental results as impressive.

Cons:
- There were concerns that the paper lacks novelty compared to Li et al 2020b, where the underlying theoretical framework was developed. The primary novelty would seem to be:
- - using Fourier transforms as the specific neural operator
- - the strength of the experimental results

Overall, I recommend acceptance. I believe the techniques in this paper will be practically useful for future research.
Decision: Accept (Poster)