Nonlinear Reconstruction for Operator Learning of PDEs with Discontinuities
Authors: Samuel Lanthaler
Roberto Molinaro
Patrik Hadorn
Siddhartha Mishra
Venue: ICLR 2023 notable top 25% (2023)

Abstract: Discontinuous solutions arise in a large class of hyperbolic and advection-dominated PDEs. This paper investigates, both theoretically and empirically, the operator learning of PDEs with discontinuous solutions. We rigorously prove, in terms of lower approximation bounds, that methods which entail a linear reconstruction step (e.g. DeepONets or PCA-Nets) fail to efficiently approximate the solution operator of such PDEs. In contrast, we show that certain methods employing a non-linear reconstruction mechanism can overcome these fundamental lower bounds and approximate the underlying operator efficiently. The latter class includes Fourier Neural Operators and a novel extension of DeepONets termed shift-DeepONets. Our theoretical findings are confirmed by empirical results for advection equations, inviscid Burgers’ equation and the compressible Euler equations of gas dynamics.

Published as a conference paper at ICLR 2023
NONLINEAR RECONSTRUCTION FOR OPERATOR LEARNING OF PDES WITH DISCONTINUITIES
Samuel Lanthaler
Computing and Mathematical Science
California Institute of Technology Pasadena, CA, USA
olanth@calttch..e. edu
Roberto Molinaro, Patrik Hadorn & Siddhartha Mishra
Seminar for Applied Mathematics
ETH Zurich
Zurich, Switzerland
froberto mo inaro, siddhartha sidhhartha.missrra@eetnn
ABSTRACT
A large class of hyperbolic and advection-dominate PDEs can have solutions with discontinuities. This paper investigates, both theoretically and empirically, the operator learning of PDEs with discontinuous solutions. We rigorously prove, in terms of lower approximation bounds, that methods which entail linear re- construction step (e.g. DeepONet or PCA-Net) fail to efficiently approximate the solution operator of such PDEs. contrast, we show that certain methods employing nonlinear reconstruction mechanism lower bounds and approximate the underlying operator efficiently. The latter class includes Fourier Neural Operators and novel extension of DeepONet termed shift- Dettt- Our thcoretical findings are confirmed by empirical results for advection equation, inviscid Burgers' equation and compressible Buler equations aerodynamics.
overcome these fundamental
INTRODUCTION
Many interesting phenomena in physics and engineering are described by partial differential equations (PDEs) with discontinuus solutions. The most common types of such PDEs are nonlinear hyperbolic systems of conservation laws (Dafermos) 2005). such as the Euler equations of aerodynamics, the shallow-wate water equations of occanography and MHD equations of plasma physics. It is well-known that solutions these PDEs develop finite-time discontinuities such as shock waves, even when the initial and boundary data smooth. Other examples include the propagation of waves with jumps in linear transport and wave equations, crack and fracture propagation in materials (Sun & Jin 2012). moving interfaces in multiphase flows (Drew Passman 1998) and motion of very sharp gradients as propagating fronts and traveling wave solutions for reaction-diffusion equations (Smoller 2012). Approximating such (propagating) discontinuities in PDEs is considered to be extremely challenging for traditional numerical methods (Hesthaven) 2018) as resolving them could require very small grid sizes. Although bespoke numerical methods such as high-resolution finite- volume methods, discontinuous Galerkin finite-element and spectral viscosity methods (Hesthaven) 2018) have successfully been used in this context, their very high computational cost prohibits their extensive use, particularly for many-query problems such UQ optimal control and (Bayesian) inverse problems (Lye et al. 2020), necessitating the design of fast machine learning- basedd surrogates.
As the task at hand in this context is to learn the underlying solution operator that maps input functions (initial and boundary data) to output functions (solution at given time), recently developed operator learning methods be employed in this infinite-dimension
setting (Higgins) 2021). These methods include operator networks (Chen Chen 1995) and their deep version, DeepONer (Lu et al. 2019 2021) where sets of neural networks (branch and trunk nets) combined in

---

Published as a conference paper at ICLR 2023
linear reconstruction procedure to obtain an siiinite-ddmessinnaal output. DeepONets have been very
successfully used for different PDEs (Lu et al.) 2021; Mao et al., 2020b Cai et al. 2021 Lin et al.
202I) An alternative framework is provided by neural operators (Kovachki et al.
202 wherein
the affine functions within DNN hidden layers are generalized to infiite-dimensions by replacing
them with kernel integral operators as in (Li et al. 2020a Kovachki el. 2021a [L et al.
2020b). A
computationally efficient form of neural operators is the Fourier Neural Operator (FNO) (Li et al.
2021a), where a translation invariant kernel evaluated in Fourier space, leading many successful applications for PDEs (Li et al. 202 b
Pathak et al. 2022).
Currently available theoretical results for operator learning (e.g. Lanthaler et al. (2022) ; Kovachki et al.
(202Ia]b); De Ryck & Mishra (2022b);(eam) et al.) (2022)) leverage the regularity (or smoothness) of solutions of the PDE prove that frameworks such as DeepONet, FNO and their variants approximate the underlying operator efficiently. Although such regularity holds for many elliptic and parabolic PDEs, it obviously destroyed when discontinuities appear in the solutions of the PDEs such as in the hyperbolic PDEs mentioned above. Thus, a priori, it is unclear if existing operator learning frameworks can efficiently approximate PDEs with discontinuous solutions. This explains the paucity of theoretical and (to a lesser extent) empirical work operator learning of PDEs with discontinuous solutions and provides the rationale for the current paper where,
using lower bound, we rigorously prove approximation error estimates to show that operator learning architectures such as DeepONet (Lu et al. [02I] and PCA-Net (Bhattacharya et al. 202I), which entail linear reconstruction step, fail to efficiently approximate solution operators of prototypical PDEs with discontinuities. In particular, the approximation error only decays, at best, linearly in network size.
We rigorously prove that using nonlinear reconstruction procedure within an operator learning architecture can lead to the efficient approximation of prototypical PDEs with discontinuities. In particular, the approximation error can decay exponentially in network size, even after discontinuity formation. This result shown for two types of architectures with nonlinear reconstruction, namely the widely used Fourier Neural Operator (FNO) of Li et al. 2021a) and for a novel variant DeepONet that we term as Sritt-DeepONet.
We supplement the theoretical results with extensive experiments where FNO and shift- DeepONet are shown to consistently outperform DeepONet and other baselines for PDEs with discontinuous solutions such linear advection, inviscid Burgers' equation, and both the one- and two-dimensional versions of the compressible Euler equations of gas dynamics.
2 METHODS
Setting. Given compact domains D C Rd U C Rd we consider the approximation of operators G X + Y, where X C LDD and y' C (U) are the input and output function spaces. In the following, we will focus on the case, where ú + (ūū) maps initial data ú to the solution at some time t > 0, of an underlying time-dependent PDE. We assume the input to be sampled from probability measure μ E Prob(X').
DeepONet. DeepONet (Lu et al. 202T) will be our prototype for operator learning frameworks with linear reconstruction. To define them, let æ := (#1...., Fmm E D be a fixed set of sensor points. Given an imput function ū E X, we encode by the point values E(ū) (ū(z)).... ū(*m)) R
DeepONet is formulated in terms of two neural networks: The first is the branch-net B, which maps the point values (() to coefficients ((E(T)) = (E)))))), ... ... (,(((() resulting in mapping
3: R" + R', E(ū) (.((E())......8(.))))) (,(E))(
...
(2.1)
The second neural network is the so-called trunk-net ( = (v), ..... Ty()), which is used to define a mapping
7: T y (T) 7p()))
(2.2)
While the branch net provides the coefficients, the trunk net provides the "basis" functions in an expansion of the output function of the form
DON (D)(( = Σ (o) âEX, y yEu U,
(2.3)
2

---

Published as a conference paper at ICLR 2023
with {x(a)= (E(E) The resulting mapping NDDON :'' + a (D)N is DeepONet.
Although DeepONet were shown to be universal in the class of measurable operators (Lanthaler et al. 2022). the followingfunndamentn lower bound on the approximation error also established,
Proposition 2.1 (Lanthaler et al. (2022) Thm. 3.4)). Let X' be separable Banach space, ) a separable Hilbert space, and let μa be probability measure on X. Let G X + y be Borel measurable operator with Ba ()) [3] < oo. Then the following lower approximation bound holds for any
DeepONet DDN with trunk-/branch-net dimension p:
@(VDON) = En-y [[LADOON (ū) - -(m)] 1/2
2 8opt ==
Σ j>p
X
(2.4)
V
where the optimal error Eopt is written in terms of the eigenvalues A1 2 A2 ... of the covariance operator I ugu == Eu-Gmallu @ u)] of the push-forward measure 9gH..
We refer to SM for relevant background on the underlying principal component analysis (PCA) and covariance operators, well example illustrating the connection between sharpness of gradients and the decay of the PCA eigenvalues X; (S....)) The same lower bound (.4 in fact holds for any operator approximation of the form (ū) 2 B& (ū)TTk where 8g X R are arbitrary functionals. In particular, this bound continues to hold for e.g. the PCA-Net architecture Hesthaven & Ubbiali (2018); [Bhattacharya al. (202T). We will refer any operator learning architecture of this form as method with "linear reconstruction", since the output function ((ū)) is restricted the linear p-dimensional space spanned by the 71..... Tp Y. In particular, DeepONer are based on linear reconstruction. To overcome the lower bound (2.4) which sets fundamental limitations on DeepONets, the basis therefore needs to additionally depend on the input u.
shift- Deet-eDeeeOet. The lower bound 2.4 shows that there are fundamental barriers to the expressive power of operator learning methods based on linear reconstruction. This is of particular relevance for problems in which the optimal lower bound &opt. in (2.4) exhibits a slow decay in terms of the number of basis functions p. due to the slow decay of the eigenvalues X of the covariance operator. It is well-known that even linear advection- or transport-dominated problems can suffer from such a slow decay of the eigenvalues (Ohlberger & Rave 2013, Dahmen et al. 2014 Taddei et al. 2013 [Peherstorfer] 2020), which could hinder the application of ineaarreconstruution based operator learning methods to this very important class of problems. In view of these observations, it is thus desirable develop nonlinear variant of DeepONet which can overcome such lower bound in the context of transport-dominated problems. We propose such an extension below.
A shift-DeepONet //=DONNeN//D M/DON, X y is an operator of the form
MOD (z) Σ Ba((ū)T (a(a))+7((z))
(2.5)
where the input function ú is encoded by evaluation at the sensor points E(u) E R. We retain the DeepONet branch- and trank-nets B, T defined in (2.)) (2.2). respectively, and we have introduced a scale-net A (Ax) consisting of matrix-valued functions
A& : RM - Red xd'
E(u) + As(u) : (c()))),
and a shift-net 7 = (/ with
7k : 1RT
dd
E(τ) H ( = 70 (E()),,
All components of a DeepOet are represented by deep neural networks, potentially with different activation functions.
Remark 2.2. The form of shift- Seitt-DeeepNett (2.5) is very natural from a theoretical perspective. Practical experimentation indicates that an extended architecture based on trunk-net T" RRd' xp - R". depending jointly on all values A1 Az))) + 71 (ū) ( Ap + (1) and defining mapping
DON (z) = Σ (u(T)T [((T)y + 7(0)
(A()))+>(1))
(2.6)
3

---

Published as a conference paper at ICLR 2023
Ap(u)y + Tp 7p(ū)) achieves better with concatenated input A(ā)y + +ū( := (ū) A.). accuracy. Our numerical results will reported for (2.6). We emphasize that all theoretical results in this work apply to both architectures, 2.. and (.6))
Since Shift-DeepONets reduce to DeepONets for the particular choice A = and Y = 0. the universality of DeepONets (Theorem 3.1 of Lanthaler et al. (2022) is clearly inherited by shift- DeepONets. However, Shift-DeepONets do not use a linear reconstruction (the trunk nets in (25) depend on the imput through the scale and shift nets), the lower bound (2.4) does not directly apply, providing possible space for Shift-DeepONet to efficiently approximate transport-dominated problems, especially in the presence of discontinuities.
Fourier neural operators (FNO). FNO ANNO (Li al. 2021a) a composition
NPNO NHy: A = Q 6 ....o o R.
(2.7)
consisting of "lifting operator" u(x) + R(ū(x),)) x), where represented by a (shallow) neural network R I x R + fd with d, the number of components of the input function, d the dimension of the domain and dh, the "lifting dimension" (a hyperparameter), followed by L hidden layers Le: vx() p(+1(x) of the form
(x) (We v'(x) be(x) (Kee") (z))
with W2 E Rd, weight matrix (residual connection), I + be(x)) E Rd. a bias function and with a convolution operator Kev (x) JTu Ke(æ - o dy, expressed in terms of a (learnabee) integral kernel x + (x) E Rd. xd. The output function is finally obtained by linear projection layer (x) + ANNO (z) Q (xx)
The convolution operators Ke add the indispensable non-local dependence the output on the input function. Given values on equidistant Cartesian grid, the evaluation of Kevº can be efficiently carried out in Fourier space based on the discrete Fourier transform (DFT), leading to a representation
Kev FN (Pe(k) FNu
where Fruf(K) (k) denotes the Fourier coefficients of the DFT of of (z). computed based on the given N grid values in each direction, Pe(k) E Cd xds is complex Fourier multiplication matrix indexed by l zd, and FN' denotes the inverse DFT. In practice, only finite number of Fourier modes can be computed, and hence coefficients of (x()) as well as the Fourier multipliers, (k) = O and Pz(k) = O, vanish whenever be efficiently computed in O(((2kmax + 1) operations (i.e. linear in the total number of grid points). The output space of FNO 2.7 is manifestly nonlinear it is not spanned by fixed number of basis functions. Hence, FNO constitute a nonlinear reconstruction method.
introduce a hyperparameter kmax E N. such that the Fourier
1k > kmax- In particular, with fixed kmax the DFT and its inverse
3 THEORETICAL RESULTS.
Context. Our aim in this section is to rigorously prove that the nonlinear reconstruction methods (shift-DeepOeet, FNO) efficiently approximate operators stemming from discontinuous solutions of PDEs whereas linear reconstruction methods (DeepONet, PCA-Net) fail to do so. To this end, we follow standard practice numerical analysis of PDEs (Hesthaven 2018) and choose two prototypical PDEs that are widely used to analyze numerical methods for cennspoor--mmiittte PDEs. These are the linear transport or advection equation and the nonlinear inviscid Burgers' equation, which the prototypical example for hyperbolic conservation laws. The exact operators and the corresponding approximation results with both linear and nonlinear reconstruction methods are described below. The computational complexity of the models is expressed in terms of hyperparameters such as the model size, which are described in detail in SMB
Linear Advection Equation.
We consider the one-dimensional linear advection equation
Oyu ad 0, u(., 0) ū
(3.1)
4

---

Published as a conference paper at ICLR 2023
on 2x -x-periodic domain D T. with constant speed E R. The underlying operator is Gady + (T)
L' (TT) n with initial data i up to any final time = T. We note that X' = L (T) L ( measure μa E Prob(X), we consider random input functions ū ~ lr given by the square (box) wave of height h, width and centered at &,
ū + Gadv == u( . T), obtained by solving the PDE (3.1)
(T) (T). As input
a(x) /2/ (r E),
(3.2)
where h E [h, hl. Jw, t E [0, 2m] are independent and uniformly identically distributed. The constants 0 < h S h, < w I w are fixed.
DeepONet fails at approximating Gudv efficiently. Our first rigorous result is the following lower bound on the error incurred by DeepONet 2.33 p in approximating Gadv.
Theorem 3.1. Let N. There exists constant C O, independent of m, p and T, such that for
any D DON 2.3 with supa-a IINDON (ū)ll I M < oc, we have the lower bound
C
E Eazy [llgadv (ū) ADDON [w)lla>} min(m,p)
Consequently, to achieve 6 (ADDN) S e with DeepONet, we need p, m Z E trunk and branch net basis functions and sensor points, respectively, entailing that Sizo(NDON) 2 pm 2 (cp. BM
The detailed proof is presented SM D.2 relies on two facts. First, following Lanthaler et al. (2022), one observes that translation invariance of the problem implies that the Fourier basis is optimal for spanning the output space. As the underlying functions are discontinuous, the corresponding eigenvalues of the covariance operator for the push-forward measure decay, most, quadratically in p. Consequently, the lower bound (2.4) leads to linear decay of error in terms of the number of trunk net basis functions. Second, roughly speaking, the linear decay of error in terms of sensor points is a consequence of the fact that one needs sufficient number of sensor points to resolve the underlying discontinuous inputs.
Shift-DeepONet approximates Gado efficiently. Next and in contrast to the previous result on DeepONet, we have following efficient approximation result for shift- Siftt-DeeONNet (2.5).
Theorem 3.2. There exists constant C >> O, independent of T, such that for any 0 there exists shift-DeepONet. A/:DNN (2.5) such that
Eazu [aly (ū) -NDDON (a))la) Se,
(3.3)
with uniformly bounded p < C, and with the number of sensor points m S Ce-1, Furthermore, we have
width(A S C, depth(A/DDON) S C' C'log(E-1)2 size (D S C1
The detailed proof, presented in SM D.3 is based on the fact that for each input, the exact solution can be completely determined in terms of three variables, i.e., the height h. width w and shift & of the box wave (3.2). Given an input u, we explicitly construct neural networks for inferring each of these variables with high accuracy. These neural networks are then combined together to yield hhift-DeepONet that approximates Gado, with the desired complexity. The nonlinear dependence of the trunk net Seitt-DeepONet (2.5) on the input the key encode the shift in the box-wave (.. and this demonstrates the necessity of nonlinear reconstruction in this context.
FNO approximates Gedv efficiently. Finally, we state an efficient approximation result for Gads with FNO (2.7) below, where the constant C > 0 is again independent of the final time T:
Theorem 3.3. There exists
0, such that for any e > 0, there exists an FNO FNNO 2.7 with
Ea-n [l] A,FFNO (a)l(( Se.
with grid size N < Ce and with Fourier cut-off kmax. lifting dimension d., depth and size:
kinax = 1, d, < C, depth(A/FNO) < C
size(AENO) < C'og(e-1))
5

---

Published as a conference paper at ICLR 2023
A priori, one recognizes that Gade can be represented by Fourier multipliers (see SM D.4. quently, a single linear FNO layer would in principle suffice in approximating Godu- However, the size of this FNO would be exponentially larger than the bound in Theorem 3.3 To obtain a more efficient approximation, one needs to leverage the nonlinear reconstruction within FNO layers. This is provided in the proof, presented in DMD where the underlying height, wave and shift of the box-wave inputs (3.2) are approximated with high accuracy by FNO layers. These are then combined together with a novel representation formula for the solution yield the desired FNO.
Conse-
Comparison. Observing the complexity bounds in Theorems 3.1 3.2 we note that the DeepONet size scales least quadratically, size 2 in terms of the error in approximating Gude. whereas for shift-DeepONet and FNO, this scaling is only linear and logarithmic, respectively. Thus, we rigor- ously prove that for this problem, the nonlinear reconstruction methods (FNO and shiftDDeepONet) can be more efficient than DeepONet and other methods based on linear reconstruction. Moreover, FNO is shown to have a smaller approximation error than even Shift-DeepONet for similar model size. We provide two remarks on extensions of these results to the approximation of the time-cvolution and to higher dimensions in SMC
Inviscid Burgers' equation. Next, we consider the inviscid Burgers' equation in one-space di- mension, which is considered the prototypical example of nonlinear hyperbolic conservation laws (Dafermos 2005):
0.u + Ox (2²²== =0, u(.,t= 0) = ū,
(3.4)
on the 2 br--peiioiic domain D T. It is well-known that discontinuities in the form of shock
waves can appear in finite-time even for smooth ū. Consequently, solutions of (3.4) are interpreted in the sense of distributions and entropy conditions are imposed to ensure uniqueness (Dafermos) 2005). Thus, the underlying solution operator is g (T) n 1 - (T) n (T).
Burg (7) == u(). with being the entropy solution of (3.4) at final time T. Given E ~ (n0 2T)), we define the random field
ū(x) == : sin(s - E),
(3.5)
and define the input measure μ E Prob T.- (T) 1, (T)) as the law of ū. We emphasize that the difficulty in approximating the underlying operator Burg arises even though the input functions are smooth, in fact analytic. This is in contrast the linear advection equation.
First, we recall the following result, which follows directly from Lanthaler et al. (2022) (Theorem 4.1 and the lower bound (2.4)
DeepONet fails at approximating Grung efficiently.
Theorem 3.4. Assume that Burg = ul . T), for T> x and is the entropy solution of (3. with initial data ū ~ μ.. There exists a constant C > 0, such that the - for any DeepONet DON with trunk-/branch net output functions is
by
Eang [[] (ū) - -DDON (T)]l C
((NDON)
(3.6)
Consequently, achieving an error & (NDDNN) S requires at least size (NDON) 2E
Chift-DeepONet approximate Borr efficiently. In contrast to DeepONet, we have the following result for efficient approximation of Borg with Nett-e
Theorem 3.5. Assume that > т. There is a constant C > 0, such that for any c > 0, there exists a
shift-DeepONet MaeeNNet.ttll MeDON such that
((A-DON)- Ba-n Burg (ū) - DON (a)lle>} Se, with a uniformly bounded number - C' of trunk/branch net functions, the number of sensor points can be chosen 3, and we have width (DNNNN < C, depth(A..Do)) < C'log(e") 1)) )2, size(A < C'log(e 1
(3.7)
The proof, presented in SM D.5 relies on an explicit representation formula for Buurg obtained using the method of characteristics (even after shock formation). Then, we leverage the analyticity of the underlying solutions away from the shock and use the nonlinear shift map in (2.5) to encode shock locations. Careful inspection of the proof implies that the constant C' is independent of T >

---

Published as a conference paper at ICLR 2023
ENO approximates SBurg efficiently Finally prove (in S.6 the following theorem,
Theorem 3.6. Assume that T иτ then there exists constant C (again independent of T), such that for any e > 0 and grid size N 2 3, there exists an FNO NFNO (27 such that
8(N = Ea-n Burg (ū) - AFNO ④ [ S,
and with Fourier cut-off kimax: lifting dimension d, and depth satisfying,
kmax 1, du C, depth(Av < C' (og(()))))
size(A)NN)) S C'log(E+1),,
Comparison. A perusal of the bounds in Theorems 3.4 3.5 and revells that after shock formation, the accuracy e of the DeepONet approximation of gBung scales at best as ~ n= in terms of the total number of degrees of freedom n = size(AVOON) of the DeepONet. In contrast, shift- DeepONet and FNO based
nonlinear reconstruction can achieve exponential convergence rate 5 exp(-eri 1/2) in the total number of degrees of freedom size(AN*DoN) size(A\ even after the formation of shocks. This again highlights the expressive power of nonlinear reconstruction methods in approximating operators of PDEs with discontinuities.
ResNet
FCNN
DecpONet
Shift DeepONet
FNO
Advection
14.8%
11.6%
7.95%
2.76%
0.71%
Equation
Burgers' Equation
20.16%
23.23%
28.5%
7.83%
1.57%
Shocktube Problem
4.47%
8.83%
4.22%
2.76%
1.56%
Riemann Problem
2.6%
0.19%
0.89%
0.12%
0.12%
Table Relative median-L error computed over 128 testing samples for different benchmarks with different models.
EXPERIMENTS
In this section, illustrate how different operator learning frameworks approximate solution operators of PDEs with discontinuitics. To this end, we will compare DeepONet (2.3) (a prototypical operator learning method with linear reconstruction) with (the extended) shift-DDeppONet (2.6) and FNO (2.7) (as nonlinear reconstruction methods). Moreover, two additional baselines (described in detail SM EEI) are also used, namely the well-known ResNet architecture of al. (2016) and fully convolutional neural network (FCNN) of Long et al.| Below, we present results for the best performing hyperparameter configuration, obtained after grid search, for each model while postponing the description of details for the training procedures, hyperparameter configurations and model parameters SM E
Linear Advection. We start with the linear advection equation 3.1 in the domain [0, 1] with wave speed 0.5 and periodic boundary conditions. The initial data is given by (3.2) corresponding square waves, with initial heights uniformly distributed between h 0.2 and h 0.8, widths between w = 0.05 and w = 0.3 and shifts between 0 and 0.5. We seek to approximate the solution operator Gade final time T 0.25. The training and test samples are generated by sampling the initial data and the underlying exact solution, given by translating the initial data by 0.125, sampled on very high-resolution grid 2048 points (to keep the discontinuities sharp), see SM Figure TO for examples of the input and output of Gadv- The relative median test error for all the models are shown in Table [ We observe from this table that DeepONet performs relatively poorly with a high test error of approximately 8%, although its outperforms the ResNet and FCNN baselines handily. As suggested by the theoretical results of the previous section, shift- Dett-DeeeeNtttt is significantly more accurate than DeepONet (and the baselines), with at least two-fold gain in accuracy. Moreover, as predicted by the theory, FNO significantly outperforms even shift- Seftt-DeepONet on this problem, with almost five-fold gain in test accuracy and thirteen-fold gain vis vis DeepONet.
7

---

Published as a conference paper at ICLR 2023
Inviscid Burgers' Equation. Next, consider the inviscid Burgers' equation (3.4) in the domain D = [0, 1] and with periodic boundary conditions. The initial data is sampled from a Gaussian Random field i.e., Gaussian measure corresponding to the (periodization of) frequently used covariance kernel,
(- x'12
k( 2') exp
28
with correlation length & = 0.06. The solution operator Burg corresponds to evaluating the entropy solution at time = 0.1. We generate the output data with a high-resolution finite volume scheme, implemented within the ALSVINN code Lye (2020), at spatial mesh resolution of 1024 points. Examples of input and output functions, shown in SM Figure [T] illustrate how the smooth yet oscillatory initial datum evolves into many discontinuities in the form of shock waves, separated by Lipschitz continuous rarefactions. Given this complex structure of the entropy solution, the underlying solution operator is hard to learn. The relative median test error for all the models is presented in Table and shows that DeepOnet (and the baselines Resnet and FCNN) have an unacceptably high error between 20 and 30%. In fact, DeepONet performs worse than the two baselines. However, consistent with the theory of the previous section, this error is reduced more than three-fold with the nonlinear Shift-DeepON outperforms DeepOnet by factor of almost 20 and learns the very complicated solution operator with an error of only 1.5%
The error is reduced even further by FNO and in this case, FNO
Compressible Euler Equations. The motion of an inviscid gas is described by the Euler equations of aerodynamics. For definiteness, the Euler equations two space dimensions are,
pu , F(U)= = pu² +
pv
U.+F(UU,++(U)) 0, U
pu pv
G(U)
pwv p + p
puw (E + p)u
F
E + p)u
(4.1)
with P. u, v and p denoting the fluid density, velocities along æn-and y-axis and pressure. E represents the total energy per unit volume
E = p0( + + P
7-1
where 7 = cp/cu, Cv is the gas constant which equals 1.4 for diatomic gas considered here.
Shock Tube. We start by restricting the Euler equations (4.1 to the one-dimensional domain
D D=[55,5] by setting 0 in (4.1) The initial data corresponds a shock tube of the form,
-f ur,
Pr. po(x) {pn
x S 20 x > zo
S zo x
x S x0
uo(x) (
1 unn
po(x)=
(4.2)
> xo
PR
x > z0
parameterized by the left and right states (pz, U4, UL.),,, (pR, un, PR), and the location of the initial discontinuity As proposed in Lye et al. (2020). these parameters are, in turn, drawn from the measure; Pr. 0.75 + 0.45G(=)), pR = 0.4 0.3G us 0.5 + 0.5G( (2)),
3 un 0,Pr. = 2.5 1.6G(22)2 PR = 0.375 0.325G((2)) 20 Oz5G(ze). with = [21, 22. z6] U ( 1] and G(z)= = 2- 1. We seek approximate the operator G [po. Polio Ep] E(1.5). The training (and test) output are generated with ALSVINN code Lye (2020). using a finite volume scheme, with a spatial mesh resolution of 2048 points and examples of input-output pairs, presented in SM Figure [2] show that the initial jump discontinuities in density, velocity and pressure evolve into a complex pattern of (continuous) rarefactions. contact discontinuities and shock waves. The (relative) median test errors, presented in Table I reveal that shift- hhett-DeepONet and FNO significantly outperform DeepONet (and the other two baselines). FNO is also better than shift- neitt-DeepONet and approximates this complicated solution operator with median error of R 1.5%.
Four- Quadrant Riemann Problem. For the final numerical experiment, we consider the two-dimensional Euler equations (4.1) with initial data, corresponding to a well-known four-quadrant
Riemann problem (Mishra & Tadmor 2011 with Uo(z. y) U Usw if x. < 0, Uo(z, (2, y) = U se-
æ < 0,y > 0, Uo(2, y) = U nnw if r
0. < O and Un(2, y) U Une- if 2, > 0,
with states given by Pa,ne Po,sw
Po.ne
Po,se 1.1, Pa.nw
Po.se 0.5065.
8

---

Published as a conference paper at ICLR 2023
Po,nw se 0.35 ne www ne, No,se] O.35 G(2).. ((z)), G(z)) and lus ww V. sw,
0.8939 + (.35 G(ze), ( G(Zz). G(z18), with 2 = [21, 22, za] U (0, and G(z) = 2z - 1. We seek to approximate the opera- tor 9 [po: Poyto: Pollo Egl H E(1.5). The training (and test) output are generated with the ALSVINN code, on spatial mesh resolution of 256 points and examples of input-output pairs, presented in SM Figure 13 show that the initial planar discontinuities
the state variable evolve into very complex structure of the total energy final time, with a mixture of curved and planar discontinuities, separated by smooth regions. The (relative) median test errors are presented in Table I. We observe from this table that the errors with all models are significantly lower in this test case, possibly on account of the lower initial variance and coarser mesh resolution at which the reference solution is sampled. However, the same trend, vis a vis model performance, is observed i.e. DeepONet is significantly (more than seven-fold) worse than both shift-DeepONet and FNO, On the other hand, these two models approximate the underlying solution operator with a very low error of approximately 0.1%.
5 DISCUSSION
Related Work. Although the learning of operators arising from PDEs has attracted great interest in recent literature, there are very few attempts to extend the proposed architectures to PDEs with discontinuous solutions. Empirical results for some examples of discontinuities or sharp gradients were presented in Mao et (2020b) (compressible Navier-Stokes equations with DeepONets) Kissas equations with an attention based framework) and in Seidman et al. (2022) (free-surface waves based on a "operator learning manifold hypothesis"", " However, with the notable exception of Lanthaler et al. (2022) where the approximation of scalar conservation laws with DeepONets is analyzed, theoretical results for the operator approximation of PDEs are not available. Hence, this paper operators arising in PDEs with discontinuous solutions has been presented, particularly for FNOs. On the other hand, there is considerably more work the neural network approximation of parametric nonlinear hyperbolic PDEs such as the theoretical results of De Ryck & Mishra (2022a) and empirical results of Lye al. 2020 202I). Also related are results with physies informed neural networks or PINNs for nonlinear hyperbolic conservation laws such as De Ryck et al.] (2022); Dagtap et al. (2022) ; Mao et al. (020a)). However, in this setting. the input measure assumed to be supported on a finite- dimensioilii subset of the underlying dimensoni input function space, making them
et al. (2022) (Shallow-War)
be considered to be the first where rigorous analysis of approximating
too restrictive for operator learning as described in this paper.
Conclusions.
A priori, it could be difficult approximate operators that arise in PDEs with discontinuities. Given this context, we have proved a rigorous lower bound to show that any operator learning architecture, based linear reconstruction, may fail at approximating the underlying operator efficiently. In particular, this result holds for the popular DeepONet architecture.
On the other hand, we rigorously prove that the incorporation of nonlinear reconstruction mechanisms can break this lower bound and pave the way for efficient learning of operators arising from PDEs with discontinuities. We prove this result for existing widely used architecture i.e., FNO, and a novel variant of DeepONet that we term as shift-DecpONet. For instance, we show that while the approximation error for DeepONets can decay, at best, linearly in terms of model size, the corresponding approximation errors for shift-DeepONet and FNO decays exponentially in terms of model size, even in the presence or spontaneous formation of discontinuiries.
These theoretical results are backed by experimental results where we show that FNO and shift- DeepONet consistently beat DeepONet and other ML baselines by wide margin, for variety of PDES with discontinuities.
Moreover, we also find theoretically (compare Theorems 3.2 and 3.3) that FNO is more efficient than even shift-DeepONet. This fact is also empirically confimmed in our experiments. The non-local as well as nonlinear structure of NNO is instrumental in ensuring its excellent performance this and and SM
context, see Theorem 3.3 for FNOs.
B.2 for further demonstration of the role of nonlinear reconstruction

---

Published as a conference paper at ICLR 2023
REFERENCES
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.
Shengze Cai, Zhicheng Wang, Lu Lu, Tamer A Zaki, and George Em Karniadakis, Deepm&mnet: Inferring the electroconvection multiphysics fields based operator approximation by neural networks. Journal of Computational Physics, 436:1 10296, 2021.
Tianping Chen and Hong Chen. Universal approximation to nonlinear operators by neural networks with arbitrary activation functions and its application to dynamical systems. IEEE Transactions on Neural Networks, 6(4):911-917, 1995.
C. M. Dafermos. Hyperbolic Conservation Laws in Continuum Physics (2nd Ed.). Springer Verlag, 2005.
Wolfgang Dahmen, Christian Plesken, and Gerrit Welper. Double greedy algorithms: Reduced basis methods for transport dominated problems. ESAIM: Mathematical Modelling and Numerical Analysis, 48(3):623-665 2014.
De Ryck and Mishra. Error analysis for deep neural network approximations of parametric hyperbolic conservation laws. arXiv preprint arXiv:2207.07366 2022a.
De Ryck, Mishra, and R. Molinaro. weinns:eeak physics informed neural networks for approximating entropy solutions of hyperbolic conservation laws. arXiv preprint arXiv:2207.08883, 2022.
Tim De Ryck and Siddhartha Mishra. Generic bounds on the approximation error for physics- informed (and) operator learning. In Advances in Neural Information Processing Systems (NeurlPS), 2022b.
Beichuan Deng, Yeonjong Shin, Lu Lu, Zhongqiang Zhang, and George Em Karniadakis. Ap- proximation rates of deeponets for learning operators arising from -ifffiinn equa- tions. Neural Networks, 153:411 1-66,, 2022. ISSN 0893-6080. doi: http://ddo..rg//0.1006// neumet.2022.06.015 URL https::://// hattps://ww.. wtw.e//eneeeeeee..toooooeeeeeeee
D. A. Drew and L. Passman. Theory of Multicomppnent Fluids. Springer Verlag, New York, 1998.
Dennis Elbrächter, Dmytro Perekrestenko, Philipp Grohs, and Helmut Böleskei. Deep neural network approximation theory. IEEE Transactions on Information Theory, 6755)2258,, 2021.
Ian Goodfellow, Yoshua Bengio, and Aaron Courville. Deep learning. MIT press, 2016.
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep residual learning for image recognition. In Proceedings of the IEEE conference on computer vision and pattern recognition, pp. 770-778, 2016.
S. Hesthaven. Numerical methods for conservation laws: From analysis to algorithms. SIAM, 2018.
Jan Hesthaven and Stefano Ubbiali. Non-intrusive reduced order modeling of nonlinear problems using neural networks. Journal of Computational Physics, 363:55-78, 2018
Irina Higgins. Generalizing universal function approximators. Nature Machine Intelligence, 3()) 192-193, 2021.
Ameya D. Jagtap, Zhiping Mao, Nikolaus Adams, and George Em Karniadakis. Physics-informed neural networks for inverse problems in supersonic flows. Journal of Computational Physics, 466: 14402, October 2022. ISSN 0021-9991. doi: 10. 10.10116//eppp222221111111 URL https
7/ w/www.c. sci enceddieieiinccccoo rect com ott.aa/..........aaa..aaaaaaaaaa
10

---

Published as a conference paper at ICLR 2023
Georgios Kissas, Jacob H Seidman, Leonardo Ferreira Guilhoto, Victor M Preciado, George Pappas. and Paris Perdikaris. Learning operators with coupled attention. Journal of Machine Learning Research, 23(215)1163,, 2022.
N. Kovachki, Z. Li, B. Liu, Azizzadenshelii Bhattacharya, A. Stuart, and A. Anandkumar. Neural operator: Learning maps between function spaces. arXiv preprint arXiv:2108.8888/33, 2021a.
Nikola Kovachki, Samuel Lanthaler, and Siddhartha Mishra. On universal approximation and error bounds for fourier neural operators. Journal of Machine Learning Research, 22:Art-No,, 2021b.
Samuel Lanthaler, Siddhartha Mishra, and George E Karniadakis. Error estimates for DeepONets: deep learning framework in infinite dimensions. Transactions of Mathematics and Its Applications. 6(1):taac001, 2022. URL https: 1/tnac001/425 0ataacc0014278554444 2t127785554//ac011......
oup. com/ imatrm/ wttp://aaanmmaat.ccleepaaa////// article-pdf///
Zongyi Li, Nikola B Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattachayya, Andrew M Stuart, and Anima Anandkumar. Neural operator: Graph kernel network for partial differential equations. CoRR, abs/2003.03485. 2020a.
Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Andrew M Stuart, Kaushik Bhattacharya, and Anima Anandkumar. Multipole graph neural operator for parametric partial differential equations. In H. Larochelle, M. Ranzato, R. Hadsell, M. Balcan, and H. Lin (eds.), Advances in Neural Information Processing Systems (NeurIPS), volume 33, pp. 6755-6766. Curran Associates, Inc., 2020b.
Zongyi Li, Nikola Borislavov Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhat- tacharya, Andrew Stuart, and Anima Anandkumar. Fourier neural operator for parametric partial differential equations. In International Conference on Learning Representations, 2021a. URL
Sttpp::/opperevieee net/ Forumrid=ccee PtN.iidd 9NOVImnO
Zongyi Li, Hongkai Zheng, Nikola Kovachki, David Jin, Haoxuan Chen, Burigede Liu, Kamyar Azizzadenesheli, and Anima Anandkumar. Physiee-informed differential equations. arXiv preprint arXiv:2111.03794 2021b.
informed neural operator for learning partial
Chensen Lin, Zhen Li, Lu Lu, Shengze Cai, Martin Maxey, and George Em Karniadakis. Operator learning for predicting multiscale bubble growth dynamics. The Journal of Chemical Physics, 154 (10):104118, 2021.
Jonathan Long, Evan Shelhamer, and Trevor Darrell. Fully convolutional networks for semantic segmentation. In Proceedings of the IEEE conference on computer vision and pattern recognition, pp. 3431-3440, 2015.
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 1eXiv:1910.03//1 2019.
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.
K. 0. Lye. Computation of Statistical Solutions of hyperbolic systems of conservation laws. ETH Dissertation 26728, 2020.
Kjetil O Lye, Siddhartha Mishra, and Deep Ray. Deep learning observables in computational fluid dynamics. Journal of Computational Physics, pp. 109339, 2020.
Kjetil O Lye, Siddhartha Mishra, Deep Ray, and Praveen Chandrashekar. Iterative surrogate model optimization (ISMO): An active learning algorithm for pde constrained optimization with deep neural networks. Computer Methods in Applied Mechanics and Engineering, 374:1 2021.
Z. Mao, A. D. Jagtap. and G. E. Karniadakis. Physics informed neural networks for high-speed flows. Computer Methods in Applied Mechanics and Engineering, 360:1 12789, 2020a.
11

---

Published as a conference paper at ICLR 2023
Z. Mao, L. Lu, O. Marxen, T. Zaki, and G. E. Kamiadakis. DeepMandMMet for hypersonics: Predicting the coupled flow and finite-rate chemistry behind a normal shock using neural-network approximation of operators. Preprint, available from 1.XX::0011.03349v 2020b.
Siddhartha Mishra and Eitan Tadmor. Constraint preserving schemes using potential based fluxes- ii: Genuinely multi-dimensional schemes for systems of conservation laws,. SIAM Journal on Numerical Analysis, 49:1023-100 1045, 2011.
Mario Ohlberger and Stephan Rave. Nonlinear reduced basis approximation of parameterized evolution equations via the method of freezing. Comptes Rendus Mathematigue, 351(23):901- 906, 2013. ISSN 1631-073X. doi: http::/ddbi.or///10001......m 2013. 10.028. URL https
www. www.ssccieceiiree.. rect. com cttp:///ccccee/aarr.aol////////. 30028477
J. Pathak, S. Subramanian, P. Harrington, S. Raja, A. Chattopadhyay, M. Mardani, T. Kurth, D. Hall, Z. Li, K. Azizzadenesheli, p. Hassanzadeh, K. Kashinath, and A. Anandkumar. Fourcastnet: A global data-driven high-resolution weather model using adaptive fourier neural operators. arXiv preprint arXiv:2202. 11214, 2022.
Benjamin Peherstorfer. Model reduction for transport-dominate problems via online adaptive bases and adaptive sampling. SIAM Journal on Scientific Computing, 4280-4328,,,,,836, 2020.
Jacob H Seidman, Georgios Kissas, Paris Perdikaris, and George Pappas. NOMAD: Nonlinear manifold decoders for operator learning. arXiv preprint 22XX:v:206.033550 2022.
J. Smoller. Shock waves and reaction-diffusion equations. Springer, 2012.
C.T. Sun and Z-H. Jin, Fracture Mechanics. Elsevier, 2012.
Tommaso Taddei, Simona Perotto, and ALFIO Quarteroni. Reduced basis techniques for nonlinear conservation laws. ESAIM: Mathematical Modelling and Numerical Analysis, 49(3):787-814 2015.
Qingcan Wang et al. Exponential convergence of the deep neural network approximation for analytic functions. Science China Mathematics, 61(10):1733-1744 1740, 2018.
Dmitry Yarotsky. Error bounds for approximations with deep ReLU networks. Neural Networks, 94: 103-114, 2017. Publisher: Elsevier.
12

---

Published as a conference paper at ICLR 2023
Supplementary Material for:
Nonlinear Reconstruction for operator learning of PDEs with discontinuities.
A PRINCIPAL COMPONENT ANALYSIS
Principal component analysis (PCA) provides a complete answer to the following problem (see e.g. Bhattacharya al. (2021); Lanthaler al. (2022) and references therein for relevant results in the infinitedimmnsionall context);
Given a probability measure Prob(3)) on Hilbert space y and given p N, would like to characterize the optimal linear subspace Vp C y of dimension p, which minimizes the average projection error
Erww [Iw n w 13]
min kpp
Ewse [|| Hyy w11/)
(A.1)
where Il, denotes the orthogonal projection onto Vp and the minimum is taken over all p- dimensional linear subspaces V X.
Remark A.1. A characterization of the minimum in A.I is of relevance to the present work, since the outputs of DeepONet, and other operator learning frameworks based on linear reconstruction (u) = Σ ] Bg (u) Tk. are restricted to the linear subspace V := span (T....., Tp}. From this, it follows that (Lanthaler et al. (2022)):
Eusu [l[[(u) (N (V))l Enzy [1)(()) (())l Euwaan [|| 11, w1b
lower bounded by the minimizer in A.I with G f the push-forward measure of under G.
To characterize minimizers of AA.1 introduces the covariance operator r, y V, by
r, := Ewoe [w @ w), where @ denotes the tensor product. By definition, r, satisfies the following relation
(v',rue)y = Ewne [(v,, w) (w, e)yl, VE, v' E e
It is well-known that r, possesses a complete set of orthonormal eigenfunctions Φ1, 02. with corresponding eigenvalues A 2 X 2 ... 0. We then have the following result (see e.g. Lanthaler
et al. 2022 Thm. 3.8)):
Theorem A.2. A subspace V. C can be written as the span of the first p eigenfunctions of an orthonormal eigenbasis Φ1, 02,... of the covariance operator Ty. with decreasing eigenvalues A1 2 X2 2 ..... Furthermore, the minimum in A.II is given by
minimizer of [A.I] I if, and only if, V, = span op
min Bwnn [|]l HI,, w 11:3 =B dim V
op
j>P
in terms of the decay of the eigenvalues of
A.1 ILLUSTRATIVE EXAMPLE
To illustrate the close connection between the decay of PCA eigenvalues and the sharpness of gradients in that distribution, we consider the solution operator G 1.2(T)) (.² of the advection equation 0.u + aQ_u= = 0 with 1, mapping the initial data + ( u(t = 1) the solution at time t = 1. We consider the input probability measure ()) E P(L?) which is defined as the law of random indicator functions ú = 1|=w/2.w/2[ <) on the periodic torus T = [0, 2n), where E [4, 3T/44 and & [0, 2TT are drawn uniformly random and independently. For we define "smoothened" probability measure (()) on input functions, whose law is obtained by mollifying random draws from (0) against Gaussian mollifier gx() with Fourier coefficients To (k) =
# = laww (ū ga(y) ~ a w0)))
13

---

Published as a conference paper at ICLR 2023
1.0
6=0.00
100
6=0.00
6=0.05
0=0.05
0.8
0
8=0 20
0=0.50
10
0m0.50
0.6
k-2
0.4
0.2
10
0.0
-2
0
109
101
10a
10a
k
(a) Representative sample from (())
(h) Decay PCA eigenvalues x
Figure 1: Sharper gradients (at length scale d) cause slow decay of the corresponding PCA
eigenvalues (
We denote by v()) ( the push-forward measure under the solution operator G. Observing that v(d) is a translation-invariant probability measure, and following (Lanthaler et al. 2022 Proof of Lemma 4. 44), we note that the PCA eigenbasis is given by the standard Fourier basis, and the eigenvalue associated with the eigenfunction eike is given by
(
2T
3x/4
1 (1:)2 dw,,
( 11:-/2,1//1 * (x()),
(A.2)
Aw
x/4
v/ (x)e- dr denotes the k-th Fourier coefficient of vi and Aw = x/2 a normalizing factor. We approximate X1e)) numerically by a (trapezoidal) quadrature in with No 4001 guadrature points. For each value of w, approximate the Fourier coefficients (6 (k) via the fast Fourier transform of wil (x) computed on a fine grid of 10 equidistant points. The approximate PCA eigenvalues are finally obtained by sorting the computed eigenvalues in decreasing order and retaining only the first 1000, which yields decreasing sequence X > 1
(k) = (2T) 50th where (1)
X
Figure D compares and contrasts the decay of the PCA eigenvalues for the probability measures v( with gradients different length scales 8 {0, 0.05, 0.2, 0.5). We note that the theoretically predicted asymptotic A ~ Ck decay is recovered for = 0.
B MEASURES OF COMPLEXITY FOR (HHFFT)DDEPPNEE AND FNO
As pointed out in the main text, there several hyperparameters which determine the complexity of SEANNN//hhft-DeepOON and FNO, respectively. Table 2 summarizes quantities of major importance for (shift-) Deeptt))eepONett and their (rough) analogues for FNO. These quantities are directly relevant to the expressive power and trainability these operator learning architectures, and are described in further detail below.
(shift-) Deeift)DeepONet
FNO
spatial resolution
Na kid kiaas size(N)
"intrinsic" function space dim.
p size(AN) depth(N)
# trainable parameters depth
depth(AV) Aed
width
(idth(AV)
dy
Figure 2: Approximate correspondence between measures of the complexity (d=dimension of the domain D C Rd of input/output functions).
14

---

Published as a conference paper at ICLR 2023
(shift-DDeepONet:: Quantities of interest include the number of sensor points m, the number of trunk-/branch-nee functions p and the width, depth and size of the operator network. We first recall the definition of the width and depth for DeepONet,
width(ADDNN == width(33) width(T),
depth(ADON) == max {depth(3)) depth(r))
where the width and depth of the conventional neural networks on the right-hand side are defined in terms of the maximum hidden layer width (number of neurons) and the number of hidden layers, respectively. To ensure a fair comparison between DeepONet, thift-DeepOn and FNO, we define the size of a DeepONet assuming a fully connected (non-sparse) architecture, as size(ADoN) : (m + w)dth(((DD))))) +
where the second term measures the complexity of the hidden layers, and the first term takes into account the input and output layers. Furthermore, all architectures consider have width which scales at least as width(A (DD)))
2 min(p, m), implying the following natural lower size bound,
size(ADON) 2 (m p) min(p. m) + min(p, m) (neddethAADDOO
(B.1)
We also introduce the analogous notions for shift-DeepONet:
width(A>sDON) :- width + width(()) width + width(7),
depth(Assbon, = max depphi(())
depth(r))
sizo(AADON) == (m + (mm+p)wwiiii) p)width(W/DO)) + Srddth(//DoNN/deeptttttttt
FNO: Quantities of interest for FNO include the number of grid points in each direction N (for total of O(N)) grid points), the Fourier cut-off kmax (we retain a total of O Fourier coefficients in the convolution operator and bias), and the lifting dimension d,, We recall that the lifting dimension determines the number of components of the input/output functions of the hidden layers, and hence the "intrinsic" dimensionality of the corresponding function space in the hidden layers is proportional dy. The essential informational content for each of these dy components encoded in their Fourier modes with wave numbers [k| S kmax (a total of O(keae) Fourier modes per component), and hence the total intrinsic function space dimension of the hidden layers is arguably of order ~ ked dna The width of an FNO layer is defined in analogy with conventional neural networks as the maximal width of the weight matrices and Fourier multiplier matrices, which of order ~ kd kmax d,, The depth is defined as the number of hidden layers L. Finally, the size is by definition the total number of tunable parameters in the architecture. By definition, the Fourier modes of the bias function by (z) restricted to wavenumbers |k| I kmax (giving total number of OKku du)) parameters), and the Fourier multiplier matrix is restricted to wave numbers [6] S kmax (giving O (kmlaxdd) parameters). Apriori, it is easily seen that if the lifting dimension d, is larger than the number of components of the input/output functions, then (Kovachki et 2001bb)
size(A/NNO) S (dl. + de + d,, depth(AVFN))
where the first term in parentheses corresponds to sin) = &e,, the second term accounts for size(Pt) Odd kmax)) and the third term counts the degrees of freedom of the bias, size(be(zz)) Nd), The additional factor depth(AVFNO) takes into account that there are L (eeth((NNNO) such layers.
If the bias by is constrained to have Fourier coefficients (k() = 0 for [k] > Kmax (as we assumed in the main text), then the representation of by only requires size (ize (k)) Od max) degrees of freedom. This is of relevance in the regime kimax < N. reducing the total FNO size from ( L) to
size(AFNNO S Me dephh(AT)
(B.2)
Practically, this amounts to adding the bias in the hidden layers in Fourier k-space, rather than physical r-space.
C EXTENSIONS OF THEORETICAL RESULTS
In this section, we remark two forward extensions of our theoretical results in Section 3
15

---

Published as a conference paper at ICLR 2023
C.I TIME-EVOLUTION
We first consider the approximation of the time-evolution + u(+,() () for solutions of the linear advection equation
0,u + adgu = 0, u(.,t= 0) i,
(C.1)
with ü drawn from the probability measure μ specified in (3.2 in the main text.
One way to apply neural operators to this time-evolution setting is by recursive application, where one fixes At > and learns an upproximation ((ū) = Sar (ū) for the given time-step, with Sar(i) given by the cata-to-solution mapping u + u( At) of the PDE (C.T). Given ū, approximation of
the time-evolution u( ,4,) at the discrete time-steps 1j jAi for j = 1,..., Ny with Ny T, is then obtained by iterative evaluations
u(,,t,) (.,1))nNoa)) NN(a) = ( (Wo...oN(())
I times
In this context, Theorems 3.1 3.2 and 3.3
can be extended to show that
For DeepONets, we have lower bound
8 = sup Eu-p (a) DON (a)la.) (o)l min j=
with C > 0 independent of T, At. Consequently to achieve an enror & c requires p. m. 2 6 trunk and branch net basis functions and sensor points, entailing lower size bound size(A/DON)
Z mp 2 -
For shift-DeepONets, there exists constant C > 0, independent of Al and Ny, such that for any there exists shift-DeepONet. A such that
Eawa [[]S,, (ū) NNDDNN ( (m)lla) sc,
sup J=I...,,Ny Nr
with uniformly bounded p S C, with number of sensor points m S CNr/c and size(AssDO)) I CNN/c.
For ENO, there exists C > 0 independent of At and Ny, such that for any c > 0, there exists FNO AFNO with
sup 8 Ny
=
Easy [||S (ū) - -N |1 Se
with grid size N < CNT /e, Fourier cut-off kmax 1, lifting dimension d, < C and size(AFFOO) & CloogNNr/e).
The additional factors of Ny in these estimates stem from the fact that the errors iterative time-steps accumulate, requiring an accuracy of order e/Ny per time-step in order to achieve an cumulative error of most e.
The above results imply that for any fixed choice of At and T = At, shift-DeepONets and FNO can approximate the time-evolution of (C.T) more efficiently than DeepONets. As an avenue for future work, it would be interesting to consider the approximation of the solution operator S : ū + where the output function u u(2, t) depends on both position æ, as well as on time f. But this outside of the scope of the present work.
C.2 EXTENSION TO HIGHER DIMENSIONS
While the analysis becomes considerably more cumbersome in higher dimensions, the main insights of this work also apply to problems on higher-dimensional domains, as we indicate for the linear advection example in the following: We consider the PDE
Oeu + Σ aarr, u =0, 0, u( = 0) = ū.
(C.2)
==1
16

---

Published as a conference paper at ICLR 2023
where ū = ((11, ....,....... Id) is function defined on the d-dimensional torus T which assume to be given by a random box wave of the form
2(
rd) = 1 wy /,//2,0,/2 -5))
(C.3)
Here E [L, Mi,, wy [e, w. and E E [0, 2n]] are independent, uniform random variables. As the 1 1-dimensional case, we consider the input probability measure defined as the law of this random box-wave in d-dimensions. Fixing T > O, and by a slight abuse notation, we will write the solution operator of C22) as Godv (ū) TT).
DeepONet: Following an analysis analogous to the one-dimensional case, it can be shown that the d-dimensional Fourier basis provides optimal PCA basis for the push-forward measure Gadv. # fs. Furthermore, an argument based on the decay of the Fourier coefficients of the box wave (C.3) implies that the PCA eigenvalues A satisfy the lower bound A 2 k- 2 In particular, this provides a similar lower bound on the number of required basis functions p for the DeepONet approximation also the d-dimensional case. Furthermore, extension of the argument in the proof of Proposition D.10 also provides a similar lower bound in terms of m. i.e. if ess sup u- ||9adv (4)]l S M, then
(()adv()))l
C
E = Eanu [lICadv (α) - - DDD (+)ll. minn(m,))
for a constant C C(M) > 0 that is independent of p and T. Again this shows that a large number of basis functions is necessary to approximation Gady by a DeepONet. We note that establish lower bound 2 1/m rather than ž 1/m-/d, as one might have expected to appear from the "curse of dimensionality" associated with
only
problems.
shift- Deitt-DeepONet: In contrast to the case of DeepONet, for neftt-DeepONetii it can be shown that there exists a constant C = C(d)= > 0, such that for any > 0, there exists a shift-DeepONet Me A/DDN such that
8 = Euw [ll (z) - _A=DON (a)]l c,
and with bounded number of basis functions p S C, number of sensor points m S Ce 1, and such that size(A\ssoo.) < C.
Sketch of proof: The idea underlying this construction is to first reduce the multi-d problem to d one- dimensional problems, to apply the known result in the one-dimensional case (along each coordinate) and to finally reconstruct the output wave-form from the solutions of the mieddimeniional problems. To reduce to the one-dimensional case, we fix a coarse uniform grid 2: {= -....," Neaasse on = [0, 2T) with step size < wt (smaller than the smallest possible box width in (C.3) ), as well as fine uniform grid 2.fine, L t=1 1,..,, Nine. on T = [0, 2T) with step size e. In terms of these one-dimensional grids, and for any given coordinate direction j E {I,...,,d] d]. we fix sensor points
2 Ei
(2 (2oaas
coarse Aun coarse ....., coRrse) T rT,
where E 1,..... Nines and where the other Eg (k # j) run over the coarse index set l 1...., Nooarse. Since Ncoarse is independent of while Nane e1 the total number of sensor points
{
jE (l......d.... df. ey E fl..... NRnee, E (L.....,eoaare for + 1}.
can be bounded by m S Cc- where C(d, w)) independent of Importantly for the reduction to the -me-dieeniional case:
the height h of the d-dimensional box can be obtained by taking the maximum over Lg) with f,...., La belonging to a coarse subset of indices (cp. Step in the proof in SM D.3 this requires neural network of size O ())-)
the sensor values a(x(() (7)
for any box wave of the form (C.). and for any direction j = 1,... d, a summation over
the coarse indices (i.e. sum over Es for k + j) of the encoded values ū ((..a)....a)) combined
with a ReLU truncation, allows us to construct DNN mapping
d
ū = h II
1 ( (-wx/2 wwt/2/w/////// ( - Ek) + 1 1 1/w5//2w,//2(zz; - E2),
k=
17

---

Published as a conference paper at ICLR 2023
where the output function is encoded by evaluation at the fine grid points a:ine,
Given this reduction to the one-dimensional case, we can then apply our one mie-dimenioonal results to find suitable approximations of h, Wj, Ey. (j = -I..... d) as in the one-dimensional case. Based on this, we construct an approximation of the box wave solution u(- T) h IT, -1 12,,/ /2) S ()T) by defining the trunk net as a suitable approximation of the dddimensiional unit box T(x) = I j== (1,,, defining the scale-net to scale each coordinate direction by A(ū) = diag(1/ W1... ................................ 1/wa), setting the shift-net to be shift by 7()) ~ (Ex + a,T"), and finally defining the branch net to be equal to the box height ((a) = h, so that
d
..)) = h II 1 1//5////////333
E a,T) B(ū)T (A(a) x (())
provides the desired approximation of the box wave solution. We will not provide the precise details and required estimates here.
Similarly, for FNO it can be shown that there exists constant C = ((d) 0, such that for any > 0, there exists an FNO APNO such that
FNO:
= En-p [[lSnav (ū) - (NNO Se
and with bounded truncation parameter kimxx = 1, lifting dimension d, S C, a number of grid points N S Ce-d, and such that size(A)FNO) C' ('looe--)))
Sketch of proof: The idea is very similar to the case of shift-DeepONet, and follows by reduction to the one-dimensional case. Note that in this case, fine grid needs be chosen in each coordinate direction (requiring many grid points, of order ~ e ed but crucially the number of weights and biases of the FNO architecture independent of this discretization parameter. Hence this does not affect the overall size. In the case of FNOs, we furthermore note that the box wave
ū(x) h (2 - E,),
j=1
be uniquely reconstructed (by nonlinear reconstruction procedure) from knowledge of it's Fourier coefficients
Fu(k,,..., ka), for bi.....,ky-. kitt....................... ka 0, (-1,0,1), and j E {l....... dd.
following the same approach detailed in the one-dimensional case in SM below.
Therefore, given the Fourier coefficients Fū(k) Fū(kl.... ka)) for k = (ki,..., k.)) with |Kl < we can first compute the corresponding Fourier coefficients the solution u(., 7) ū( aT).
by simple phase shift Fiu(-,T1 (k) Fu(k)e a)k and then reconstruct ū( aT) from
knowledge of these Fourier coefficients. Thus kmax = 1 also suffices in this case (corresponding to
3d Fourier coefficients which are retained). The lifting dimension d, of the detailed construction outlined above scales linearly in the dimension of the domain d, but is independent of
D MATHEMATICAL DETAILS
In this section, we provide detailed proof's of the Theorems in Section 3 preliminary results below,
We start with some
D.1 RELU DNN BUILDING BLOCKS
In the present section we collect several basic constructions for ReLU neural networks, which will be used as building blocks in the following analysis. For the first result, we note that for any 8 > 0, the
Similar to the one-dimensional case, additional shifted copies of these basis functions are needed in practice account for the peepiiiiy of the output function.
18

---

Published as a conference paper at ICLR 2023
following approximate step-function
0,
x < 0,
(x(x) ==
0 5
6
can be represented by a neural network:
<5(x) = 18
(~ -m--0
where (xx) max(x, 0) denotes the ReLU activation function. Introducing additional shift multiplying the output by and choosing > sufficiently small, we obtain the following result:
Proposition D.1 (Step function). Fix an interval [a, 6] C. I. & E [a, 6). h R. Let h Hir be step function of height h. For any and p [1, 00), there exist a ReLU neural network &. + R. such that
depth (.) = 1, width( = 2,
and
p - h [alee)ll S.
The following proposition is an immediate consequence of the previous one, by considering the linear combination 4 a) 2s(x b) with suitable choice of 5 >
Proposition D.2 (Indicator function). Fix an interval [a, b] C R. Let (a, (x) be the indicator function of [a, b]. For any > 0 and p E [1, oo), there exist ReLU neural network o R - R, such that
depth( = width(#.. 4,
and
2 (a.ollLL((a..))) Sc.
1.0
A1
A2
0.6
As
0.4
A
0.2
As
0.0
0.00
0.25
0.50
0.75
1.00
X
Figure 3: Illustration of partition of unity network for / 5, [a, 6] = [0,
useful mathematical technique glue together local approximations of given function rests the use of a "partition of unity". In the following proposition we recall that partitions of unity can be constructed with ReLU neural networks (this construction has previously been used by Yarotsky (2017); cp. Figure 3):
Proposition D.3 (Partition of unity). Fix an interval [a, b] C R. For N. let Ar := (b (b- a)/J and let x, jAr, j 0,..., J be equidistant grid on (a, b). Then for any (0, Ar/2), there exists ReLU neural network A R Rdd + (A₁ ... Az(x), such that
width(A) 4J, depth(A)
19

---

Published as a conference paper at ICLR 2023
each A, is piecewise linear, satisfies
0, (zS -c, 1, t -),
Ay(0)
0, (x
and interpolates linearly between the values 0 and the intervals [xz-1 6, + c) and (ry 5 + e). In particular, this implies that
supp(A,) c [zz-1 e, + e). for all I,...,,,
Az(x) 2 0 for all x E IR,
The {Ay)j=1.....J form a partition of unity, i.e.
E 1,)=) Vo [a,
We also recall the well-known fact that the multiplication operator (ar, y) + ay can be efficiently approximated by ReLU neural networks (cp. Yarotsky (2017):
Proposition D.4 (Multiplication, (Yarotsky] 2017 Prop. 3)), There exists a constant 0, such that for any E (0, 1. M 2 2, there exists neural network Xe,M |-M, M] x [-M, M] + R, such that
ardth(&,,M)) S C, depth( R..Mr S C'log(Me=)), size(8e.M) S C Colg(MMe
and
sup 2.VE[-M, M.M [&,, xzl Se
We next state general approximation result for the approximation of analytic functions by ReLU neural networks. To this end, we first recall
Definition D.5 (Analytic function and extension). function F (a, B) R is analytic, for any 3o E (a,88 (3) there exists radius > 0, and a sequence (aR)REN, such that E |ax|r" < and
F(x) - Σ cala ro)*, V|r- zol <T.
If f (a, b] + R is function, then we will say that f has an analytic extension, if there exists F: (a, + R. with (a, b) (a, (3), with F(x) f(x) for all x (a, and such that analytic.
We then have the following approximation bound, which extends the main result of Wang et al. (2018). In contrast to Wang al. (2018), the following theorem applies analytic functions without globally convergent series expansion.
Theorem D.6. Assume that f [a, bj - R has an analytic extension. Then there exist constants
C, 0, depending only on J. such that for any L E N, there exists ReLU neural network
Φ1= R. with
sup [f(x) 2x(2)) < C Cexp(-71\13 xE(a,b]
and such that
depth(+r.) CL, width(0 C.
20

---

Published as a conference paper at ICLR 2023
Proof. Since f has an analytic extension, for any [a. there exists radius Tu O, and an analytic function Fg : [x T2, 2 +z] + RR which extends f locally. By the main result of (Wang et al. 2018 Thm. 6), there constants Cz: Yi > 0 depending only on E (a, b), such that for any
E N. there exists ReLU neural network + F. IR + R, such that
[F(E) 22. (2) C2 exp(-
sup 1s-2/sro
and depth((..) < C. (x) I Cr. For J E N. set Ar = (b a)/J and consider the equidistant partition 2j == + jAr of [a, Since the compact interval [a can be covered by finitely many of the intervals (x Tz: I + rz), then by choosing Az sufficiently small, we can find x0f) E [a, bj. such that [g)] C (x r ; x!1) + +ran)) for each j I,............... J. J.
By construction, this implics that for C == max Czo., and 7 == minj= Yoor, we have that for any L N, there exist neural networks 3 (())))) :RR IR, such that
[f(x) (x)) < C exp( exp(-)
sup * =
and such that depth((9,.)) z CL, width( SC.
Let now A : R + R the partition of unity network from Proposition D.3 and define
(2 = Σ M, (Az(z), (z,(()))
where X., denotes the multiplication network from Proposition D.4 with == 1 Ceep(--)) supxe(a.0) [x]]] and e == J - (eceeppp<>L)) Then we have
depth(#r)) s depth( Mr,e) + depth(A) max (epth((),1.) =
< Logle +1+CL
L)
where the constant C > 0 on the last line depends on supze(a.M [f()), 7 and on Ò, but is independent of Similarly, we find that
width( S width(A) + max width s.,.) S 4J + C,
j= j=I.....
is bounded independently of L. After potentially enlarging the constant C' > 0, we can thus ensure that
(epth((rr.) S CL, width( S C,
with a constant C >> O that depends only on f. but is independent of L. Finally, we note that
| - f(x) <
XM. - -
<
M.c (,L((,)) (,(z)
A5(2) (z2(()) f(x))
By construction of X M,e. and since 1 (x)) M, A, M. the first can be bounded by Je Č' (exp(->T Furthermore, each term in the second sum is bounded by Az(x) 1/11)))) and hence
Σ A,(0)) Cexp(-7L) 1/3
1=
21

---

Published as a conference paper at ICLR 2023
for all x E (a, bl. We conclude that su [x] f(x) - 2C 2C exp((>)) with constants
C, 7 > independent L. Setting y : : and after potentially enlarging the constant > 0 further, we thus conclude: there exist C, > O, such that for any L E N, there exists neural network Φr + R with depth (1)) - CL. width( S C. such that
sup
|
[(=) S C Ceep(=>X''
xE(a,b]
This conclude the proof of Theorem
D.6
By combining suitable ReLU neural network approximation of division + 1/a based on Theorem D.6 (division is an analytic function away from 0), and the approximation of multiplication (cp. Proposition D.4 above), we can also state the following result:
by Yarotsky 2017)
Proposition D.7 (Division). Let < < b be given. Then there exists C = C(a, b) > 0, such that for any c C (0, 2). there exists ReLU network ta. : IR x IR + R, with
depth(Fa.b.e) S C'log 15= width(Fa.dd S SC, size( a, & log 1
satisfying
Fa. (x:
x. Sc
sup a.yE [.
y
We end this section with the following result.
Lemma D.8. There exists constant C > 0, such that for any > 0, there exists a neural network
E, TR such that
sup <E (E(
3. (cos()) S
with (cos sin(E) [0, 2n] for all & E [0, 2x], and such that
dcpth(E.) s C' C'log(e 1)2 width(E.) S C.
Sketch of proof. We can divide up the unit circle (ino)))))) [0, 2x]} into subsets, where
{ E n/4],
- x xZ 1/v2, y20,
{ E (T/4,
= y 1/v2,
E [3T/4, /4, 5 57/4]. /4],
- I IS -1/V2,
E (5T/ 7744),
= y < -1/V2,
{ E [7x/4, 2m),
- x S -1/V2 < O,
On each of these subsets, one of the mappings
1
+ R. x = cos(E) +
2
y E
V2 - R y == sin(E) + E.
V2
is well-defined and possesses analytic, invertible extension to the open interval (-1, 1) (with analytic inverse). By Theorem D.6 follows that for any > 0, we can find neural networks an open set containing the corresponding
i,,-.., Ss. such that (cos(()()) (E) domain, and
depth( s C' log(e 1))
width (4)) C.
By a trraight-forwann partition unity argument based on Proposition D.3] we can combine these mappings to global map
2 which is represented by neural network =, lg [0, 2m), such that
sup
[ (cos(f))) sin(E)) El <
6E(0.2x-e)
this step is where the e-gap at the right boundary 2x is needed, as the points angles & and 2T identical the circle.
22

---

Published as a conference paper at ICLR 2023
and such that
depth(E.) C'log(e"1), width(E.)) S C.
D.2 PROOF OF THEOREM 3.1
The proof of this theorem follows from the following two propositions,
Proposition D.9 (Lower bound in p). Consider the solution operator Gady L n (T) - L' n (T) the linear advection equation, with input measure μ given as the random law of box functions of height
[h, li), width w [ue, w) and shift & E [0, 2 Let M > There exists and M, with the following property: If (ū = Ek=1 (ū)Ta The any operator approximation with linear reconstruction dimension such that suppu- [ I M co, then
constant C C(M, ()) > 0, depending only
Ea=n [I - Gadv(u)llz) C
Sketch of proof. The argument is an almost exact repetition of the lower bound derived in Lanthaler ct al. (2022) and therefore we will only outline the main steps of the argument, here: Since the measure translation-invarian,, it can be shown that the optimal PCA eigenbasis with respect to the L no (cp. A) is the Fourier basis. Consider the (complex) Fourier basis kkzz, and {etk>)
denote the corresponding eigenvalues (Ax JREz The k-th eigenvalue d of the covariance-operator r₁ Eg (u @ u)] satisfies
Fa (e (e- tke)
dre -izz
short calculation, in (Lanthaler al. 2022 Proof of Lemma 4.14), then shows that
Th
1 11 /
107
Aw Ah h² [
dw dh
dw (k))
Aw
where (k) denotes the k-thh Fourier coefficient of (xr == /2, 2 (I). Since ( has jump discontinuity of size for any O, it follows from basic Fourier analysis, that the asymptotic decay of (w C/|k/,, and hence, there exists constant C C(h, wt w)) > such that
C
Re-ordering these eigenvalues X in descending order (and renaming). À 2 A2 2 it follows that for some constant C' C(h, (h, w. w) > we have
Cp-
j>p
1>8
By Theorem 2.1 this implies that (independently of the choice of the functionals (e(ū)), in the Hilbert space [TT), we have
ΣトZ C
Eavu [I[((0) - ((M))a)
j>p
for a constant C that depends only on μ. but independent of p. To obtain a corresponding estimate with respect to the L -norm, we simply observe that the above lower bound on the L - together with the a priori bound supa~a 1 (ū)ll S on the underlying operator and the assumed -bound supa~u [I I M, imply
C
< E Tay
[[) [(
p
< B
[
-
< (M + B (T) (())l)))
This immediately implies the claimed lower bound.
23

---

Published as a conference paper at ICLR 2023
Proposition D.10 (Lower bound in m.) Consider the solution operator Gadv L' LT) n (T) L' (T))) of the linear advection equation, with input measure 1l given as the law of random box functions of height E [h, hi. width w [w, and shift ç E [0, 2m.. There exists an absolute constant C with the following property: If ADDON is DeepONet approximation with m sensor points, then
C
E B (1))(a) Gadv (u)l m
Proof. We recall that the initial data ú is randomly shifted box function, of the form
ū( [11--/22+m//]] (x E),
where ç E [0, 2T], h E 14, h] and w variables.
w u) are independent, uniformly distributed random
Let #1...... Im E (0, 2m] be an arbitrary choice of m sensor points. In the following, we denote a(X) := (ū(r1),.... .... ((z))) e R"" Let now (x, ((X)) (2;; ū(X)) be any mapping, such that x + ((x; (X)) E L (T) for all possible random choices of à (i.c. for all p). Then we claim that
((X))))) C
E Ea-pl (ewpullluual Fudv (1)
((:;;(())
(D.1)
m'
for a constant C(h, w) 0. holds for all m N. Clearly, the lower bound D.I immediately implies the statement of Proposition D.10 upon making the particular choice
(xx) = No))(<)) = ( (Txm))
To prove D.I,, we first recall that ū(x)) =[(;; h, w, {) depends on three parameters h, w and &. and the expectation over ū p in [.11) amounts to averaging over h E [h. 7.. E [w w| and
E E [0. 2x). In the following, we fix and h, and only consider the average over Ç. Suppressing the dependence on the fixed parameters, we will prove that
2n
C
[[(x;;) E) ((x; ((x 2 m'
(D.2)
2T
with constant that only depends h. w. This clearly implies D..] To prove D.2), we first introduce two mappings & + I(E) and - J(E), by
W E (2,,
+ w
I(E) = - = i
J(Z) = j =
E
2
2
where we make the natural identifications on the periodic torus (e.g. #m+1 is identified with 01 and Ç + w/2 is evaluated modulo 2x). We observe that both mappings
T(E), J({) cycle exactly once through the entire index set (1, (1.....,) m} as & varies from 0 to 2v. Next, we introduce
Ang = ( E [0, 2T) ((E) = i, (E) Vi,j E mL....,m..
Clearly, each & [0, 2T) belongs to only one of these sets Atj. Since E + I{E) and E++ J({) have m jumps on [0, 2m). it follows that the mapping & + (T(E), J()) can have at most 2m. jumps. In particular, this implies that there are at most 2m non-empty sets As * (these are sets of the form A [0, 2m)), i.e.
#(A, #01ij E (1,... m)) 2m.
(D.3)
Since ū(2; {) h 1[- 2 / (x E) one readily sees that when varies in the interior of As, then all sensor point values af x: remain constant, i.e. the mapping
interior(A.)) + R", E ū(X:E) = const.
We also note that A, = [x; + w/2, 2i+1 + w/2) n [x, - w/2, w/2) is in fact an interval.
Fix i, j such that Aij * for the moment. We write Aij (a A a + A) for some a, A E T,
24

---

Published as a conference paper at ICLR 2023
and there exists constant U E R" such that U ((X;E) for Ja A, A). It follows from the triangle inequality that
rtt
[ū(z;
((;; (xx de
|Í([;; {) - ((;;;;;; 12 dg
Ay
a-2
1ū(2ça E) 0(x; +1x dg'
[ ||ū(z; + xx:;; U) de'
[ +5) ū(x; l' de'
Since ū();{) h (= E), we have, by simple change of variables
|ū(2; 5¹) uiz;a E' = /11 - 1///2, + 29') dr.
The last expression is of order E'. provided that 5' is small enough to avoid overlap with a periodic shift (recall that we extension). To avoid such issues related to periodicity, we first note that g S A S and then choose a (large) constant Cp = Cu(w), such that for any g' S T/Co and w S w. we have
working on the torus, and /1=//2,w/2 (ar) is identified with its periodic
1. /1=1
w0/2 [
/2 ( 2(') =
1 +
de 48'.
-1/22-2
2//2-2
From the above, we can now estimate
S1p [[ū(1;))) (a;a((;;00 d
lū(x; + --(x; a E') dg'
A/Co
[ū(;; ū(z;a 5') de'
A/Co
2 h
4g' de'
A
2h 2 AlA
C
where C(h, w a constant only depending on the fixed parameters h, w. Summing over all Aij * 0. we obtain the lower bound
2x
|ū(x; { (x X: 1l)l de
Σ
/
[a(x;
((;;)( dg
A17701
Agg
2 C |Aty1
A,
We observe that 2x)) = As is a disjoint union, and hence Σ EAnyo vy AAil = 2v. Furthermore, as observed above, there are at most 2m non-zero summands |Ayy * 0. To finish the proof, we claim that the functional > |ag/2 is minimized among all a.,..... 0.2m satisfying the constraint 22m |ass = 2n if, and only if, |ai| = Jazml = π. /m. Given this fact, it then immediately
follows from the above estimate that
20T²"
2x Jo
()) - & (X: 2
Σ A
m
Ag70
where C' = Chh, w) > 0 is independent of the values w E [w, w| and h E [h, FJ. This suffices to conclude the claim of Proposition D.10
25

---

Published as a conference paper at ICLR 2023
It remains prove the claim:
We argue by contradiction. Let di...., 4m2 be minimizer of Σx aagl under the constraint Σg [ax = 2n. Clearly, we can wlog assume that 0 0 01 ....
a.mm are non-negative numbers. If the claim does not hold, then there exists minimizer, such that on < 022m Given > 0 to determined below, we define 3x by
B = 01 + 6, 32m 22m
and Ork, for all other indices. Then, by simple computation, observe that
Σ Σ
oh 20(a.m
- 8).
K
Choosing 5 > 0 sufficiently small, we can ensure that the last quantity is > O, while keeping 3p. 2 0 for all In particular, it follows that Σn 13
aggl
2n, but
Σ aak
in contradiction to the assumption that 01,,... 02m minimize the last expression. Hence, any minimizer must satisfy Jon| =... |02ml = क t/m.
D.3 PROOF OF THEOREM 3.2
Proof. We choose equidistant grid points .i,..., 2m for the construction of shift-DeepONet approximation to Gadv. We may wlog assume that the grid distance Ax < w. as the
statement asymptotic m - oo. We note the following points:
Step 1: We show that h can be efficiently determined by max-pooling. First, observe that for any two numbers a, b, the mapping
mmax (0, a - -b)) max (0, b)
(
+ max(0, a 6) max(0, b) max (0, -b) = max(a,b) 6),
max (0, bb)
exactly represented by ReLU neural network max(a, b) of width 3, with a single hidden layer. Given & inputs .1..... ag. we can parallelize O(k/2) copies of max, to obtain a ReLU network of width S 3k and with a single hidden layer, which maps
a2
max(aj, 02)
ak-1
max max((x+11,0)
0k
Concatenation of Ollog2 (k)) such ReLU layers with decreasing input sizes k. [k/2], [k/4],... provides ReLU representation of max-pooling
@1
max(a1, a2)
max(ai, 02, 03, a4)
max(a ax).
max mnax(ax-1,40 (x)
max mnax(ak-3. (ak-3, Gk ak 1 (x))))))))
This concatenated ReLU network maxpool I E has width < 3k, depth (log(k))) and size Ok log(k)).
Our goal is to apply the above network maxpool the shift-Deep input u(. u() to determine the height h. To this end, we first choose .,..... l E {1....,m., m}. such that 2t -It, < w, with & minimal. Note that is uniformly bounded, with a bound that only depends on (not on m)) Applying the maxpool construction above the u(zr)), ...., ...,,,,)) we obtain a mapping
u(x1)
(u(eee))
+
+ maxpool(u(22,).... ,((r,)) = h.
u(Fm))
(((ree))
26

---

Published as a conference paper at ICLR 2023
This mapping can be represented by ((log(k)) ReLU layers, with width I 3k and total (fully connected) size O(k² log(k)). In particular, since & only depends on w, we conclude that there exists C' = C(w) > O and neural network h mh - R with
depth(h) S C, width(h) S C, size(h) SC,
(D.4)
such that
M(ū(X)) = h,
(D.5)
for any initial data the form ū(x) = h1 h1/-w/2w/////t E). where E [h, M., w E [w, w.. and ç E [0, 2x]
Step 2: To determine the width we can consider linear layer (of size m.) followed by an approximation division, {a; b) xs a/b (cp. Proposition D.7
((rz) +
Ar
a(2,);
(((x):((((())))
Denote this by (((ū(X)), Then
w
1
2n
u ((x)
dx Σ ū(rig)
Ar
h
Ar ū(x;)) - hhhe
Ax
ū(x,)) h
2
+ E.
m
And we have (c) S (clo((c1))) 1), width(w) S C, sizo(w) S C. (m log(c-1)3), by the complexity estimate of Proposition D.7
Step 3: To determine the shift & [0, 2T). we note that
m
Ax
ū (z
ū(x)e dx +0
m
2sin(w/2)) -i + +0 ( )
m
Using the result of Lemma and the observation that w E (w, π implies that sin (w/2) 2 sin(w/2) > is uniformly bounded from below for all [wN. w.. it follows that for all (0, z). there exists neural network
D.8 combined with the approximation of division of Proposition D.7
: Rm > [0, 2x.. of the form
]
u u(xx)e ) ia;; XM.. (h., 2sin (, u8 (1/22))
{(u(X)) = E.
=r Ax Σ
such that
1: ((=(((()) SE
for all & [0, 2n ee. and
depth(()= < C' Clog(e-1)) width(E) SC size(E) & C' (m + log((
Step 4: Combining the above three ingredients (Steps 1-3). and given the fixed advection velocity
E RR and fixed time i, define shift-DeepONet with 6, scale-net A₁ 1, and shift-net
27

---

Published as a conference paper at ICLR 2023
with output 7/ (a) + at, as follows:
NaDoN(()) = Σ (a(7))((-~<((1))
Σ hl (0,0)) ( f at + 1/2 00/22/2+111 112++211 27j)
x- (,,000) (x
at 8/2 + 2TJ
(0,00)
where ñ = (())))),
ã- ((((())
x(((X)) )) and & = (()), and where 1
0.00) is a sufficiently accurate
I appapoooxinttoon of the indicator function 1 10,00) (a)) (ep. Proposition
D.2 To estimate the approxi-
mation emorr we denote 1 w//.22//2 (x)
= [ x) (or + 7/2) 1 (0))00) (x w/2) and identify it with
it's periodic extension to T. so that we can more simply write
ANDON (z) (1
at).
We also recall that the solution u(z, 1) of the linear advection equation 0,u + ad = 0. with initial data u(2; 0) u(xr) is given by u(x, l) ū(x al), where al fixed constant, independent of the input ū. Thus, we have
Gadv (x) = ū(x - = le - 1100/200/222 - al).
We now write (()))()(()) -NDDON (z)()) [ - - at) -
at)
We next recall that by the construction of Step 1, we have h(a) h for all inputs u. Furthermore, upon integration over z, we can clearly get rid of the constant shift at by change of variables. Hence, we can estimate
1| (u) ANDON (i) S h $ - 5) -
-allao------aaaaaaa -if-w/220/2(1 3) de.
(D.6)
Using the straight-forwan bound
1 1// E) - (1 < 1 w {) 1 w|/ 5)
(x-
1-w/2, { 1- 0/2,0/2 E
1/ Q) if / //2. (I
one readily checks that, by Step 3, the integral over the first term is bounded by /
I10) sc [. -
2m
K El de
+
d S (C' + 2)e.
2
where C. C(w, w) 0. By Step 2, the integral over the second term can be bounded by (I1) < Clw / < C' (1/m +e).
Finally, by Proposition D.I by choosing e sufficiently small (recall also that the size of İ is independent of e), we can ensure that
(IIIT))) = 1 //2,//2 - 1/=/22 S
holds uniformly for any W. Hence, the right-hand side of (D.6) obeys an upper bound of the form
dw f de (Z) - -VDDON (ū)llra
Einu [1]Gadv - DDN (m)l
dh
J
h
dw
dg (l(T))lz + II(77) + [[(IT)))
2T
1
C
m
28

---

C N Published as a conference paper at ICLR 2023
for constant C' (w,, w. h)) O. We also recall that by construction,
depth(A>soo)) v C Clog(e )
width(A<DON) C, size(AaDON)
S C (m + log(E ()))
Replacing c by c/2C and choosing ct, we obtain
Ban [advv ANDON (ū)ll>) cc
with
depth(A>sDON) S C log(( 1)2 width(AsDOON) < C, size(\>DON) S SCE
where C' depends only on μ., and is independent of c. This implies the claim of Theorem 3.2
D.4 PROOF OF THEOREM
For the proof of Theorem 3.3 we will need few intermediate results:
Lemma D.11. Let ū = h 1 11-/22,w/0) /2 Б) and fix a constant at E IR There exists a constant C > such that given N grid points, there exists an FNO with
kinax 1, de S C, depth S C, size S C,
such that
C
sup N h
[NNNO(T)(a) (z)( sin(w//2 cos(ar E - - S N
Proof.
We first note that there a ReLU neural network & consisting of two hidden layers, such that @ ((((x)) min(1,
((:)) (=))) min h (x))
||= (2))
for all h E [h. Rl. Clearly, & be represented by FNO layers where the convolution operator K₂ = O.
Next, note that the Fourier coefficient of ū 1/=//2w//2 is given by
1
FN
w / - E)etie
N
2m
11-//2w/21122 - Eie de + O(N-1)
0
sin w 2/))eFFi
O(N-1),
where the O(N error is bounded uniformly in g [0, 2T]] and w [wu w. It follows that the FNO AENO defined by
PFx((() a (F
where implements a projection onto modes | and multiplication by etior (the complex exponential introduces a phase-shift by at), satisfies
FNO (z)(x) sin(w/2) cos(x (o-ss==>) - E- s-at)]
sup
zET
where C' independent of N, w. h and &.
Lemma D.12. Fix 0 < h < h and 0 < w There exists constant C(h, w. W) > 0 with the following property: For any input function ū(x)) h 1/-//2w//2{(t - {) with [1. hi and
[w, w). and given N grid points, there exists FNO with constant output function, such that
C
sup LNFNO (z)(x) w S h,wag
and with uniformly bounded size,
kimax 0, d, 5 C, depth C, size & C.
29

---

Published as a conference paper at ICLR 2023
Proof. We can define FNO mapping
2T
[ 1 (2,) w O(N
where we observe that the first mapping is just a H> max(h (11(h----p ((x), 1), which is easily represented by an ordinary ReLU NN of bounded size. The second mapping above is just projection onto the O-th Fourier mode under the discrete Fourier transform. In particular, both these mappings can be represented exactly by FNO with Amax = 0 and uniformly bounded d,, depth and size. To conclude the argument, we observe that the error O(N-!) depends only on the grid size and is independent of w [w, w].
Lemma D.13. Fix w < W. There exists a constant = C(w, w) 0, such that for any >
there exists an FNO such that for any constant input function ū(x)) = [w, W. we have
NPNO (z) 2
sin(w) s.
Vx [0, 2m],
and
kmax = O, dg CC depth I C'log(e 1)2
size S C'log(e-1).
Proof. It follows e.g. from (Elbrächter al. 2021 Thm. III.9) (or also Theorem D.6 above) that there exists a constant C = C(w,, w >> O, such that for any > O, there exists a ReLU neural network S. with size(S.) S Cloo(ee 1)2 depth(S.) $ C'log(e (12 and width(S.)) S C, such that
S. sin(w) S.
sup
w fuu T0
To finish the proof, simply note that this ReLU neural network S. can be easily represented by an FNO S. with kmax = 0, d, S C. (S.th((.)) I C Coog(e" ))2 and size(S,) I C' C'log(e )): it suffices copy the weight matrices We of Se, set the entries of the Fourier multiplier matrices Pz(k) and choose constant bias functions (e(x)= = const. (with values given by the corresponding biases in the hidden layers of S.).
Lemma D.14. Let ū = h h w//1//////2 (x E). Assume that /N/N= S w. For any e > 0, there exists an FNO with constant output function, such that
sup (VFNO(O)(a) Se. h,w.g
and
kmax 0, d, S, depth < Clog(e" 1)2
Proof. The proof follows along similar lines as the proofs of the previous lemmas. In this case, we can define an FNO mapping
h 1 (=1 12-0/2,002 (I)
1 11-/2
h
1 /2 (Zz), Σ
j=1 11-//2 (o,))
jml
The estimate on kmax de, depth follow from the construction of _ in Proposition D.7
Proof of Theorem 3.3
We first note that (the 2 Serpeiiodization of) 1 1 (x - -- at)
and only if
cos(æ - E at) 2 cos(w/2)
> sin(w/2) cos(x - - at) 2 sin(w:).
(D.7)
The strategy of proof is as follows: Given the input function ū(x) h 1-1
{ with
unknown [w, w). E [0, 2n] and [4, hj. and for given a, t E (these are fixed for this problem), we first construct an FNO which approximates the sequence of mappings
h
T +
w
+
2
sin sin(w) ww)
+
sin(w//2 (onww2) cos(x - - - at) - sin(w)
sin(w//2)oom E)
(u /2) cos( - at)
30

---

Published as a conference paper at ICLR 2023
Then, according to (D7.p the identity
approximately reconstruct 1 (/ ( by approximating
1 sin(w))
1[ / - --a)) 10, sin(ww/2 - 5 at)
2
where oe) the indicator function of [0, oo). Finally, we obtain Gadv (u) h 1 w/2 - at) by approximately multiplying this output by h. We fill in the details of this construction below. Step 1: The first step is to construct approximations
the mappings above. We note that we can
choose a (common) constant Co = (o, h, w. 3) > 0. depending only on the parameters h, h, w and W. such that for any grid size N E N all of the following hold:
There exists an FNO HN with constant output (cp. Lemma D.14), such that for ū(x) h 1 /2 E).
1
sup wch
HN ( 11| S N
(D.8)
and with
kmax $ 1, dg S Co.
depth < C. log(N) size S Co log(()))
2. Combining Lemma D.12 and D.13 output, such that for x = h 1-w//2 2/ (x {), we have
conclude that there exists FNO Sy with constant
Co N
sup
(N) sin(w) S
(D.9)
wc [ p 1,74+11
and with
kmax O, dg S Co. depth < Co log( size S Co logg(
3. There exists an FNO Cy (cp. Lemma D.IT), such that for h 1 11-w/2,w/2)
sup [Cnx (x) sin(w/2) cos(ar - - Co x,E,w
(D.10)
N
where the supremum is over vr & E [0, 2m] and w E jw, w). and such that kynax 1, d, S Co, depth S Co. size S Co.
4.
There exists a ReLU neural network I(0.>0) (cp. Proposition
D. such that
0, Mooo)lle 1,, IN (2) =
(x 0). (x 2 N )
N
(D.11)
1,
with
width (,,00) Co. depth((i I Co-
5. there exists a ReLU neural network x (cp. Proposition D.4 such that
1 (a, abl s N'
sup
(D.12)
a.b
where the supremum is over all |a]. 1b < h 1, and
width(Xxx) S Co. depth( < Co log(N).
Based on the above FNO constructions, we define
N FNO (ū) : xN (HN( TM (Cx(O) - SN (a))
(D.13)
31

---

Published as a conference paper at ICLR 2023
Taking into account the size estimates from points 1-5 above, well as the general FNO size estimate follows that AFNN can be represented by an FNO with
B.2
kima =1
dy S C,
depth S Clog(NN², size I Cloog(N
(D.14)
To finish the proof of Theorem 3.3 it suffices to show that AFNO satisfies an estimate
C
sup IAFNO Gadv (u)ll 2 a-
N
with C > independent of N.
Step 2: We claim that if a [0, 2m] is such that
2Cu + 1
sin(w/2) cos(x - -) sin(w) 2
N
with Co the constant of Step I, then
, (C)( -Sx(t)) 1- - -)..
To see this, we first assume that
2Co
sin(w/2) cos(x - -E) sin(w) (w)
>
N
Then
1
Cx(T)(x) Sx(u)= 2 iinww/2 cos(x - E) sin(w)
(C)( - sin(w/2 cos(1 - E) - SN Sw((a) sw)
2Co +1
Co
Co
1
0.
N
N
N
N
Hence, it follows from (.111 that
IN T(0.,0) (Cx((())
Sx(u))
(sin(w/2) cos(x -E) 2 sin(w)
10.000)
-1 (I)).
The other case,
2C + E) sin(w)= -
sin(w/2) cos(æ
2
N
is shown similarly.
Step We note that there exists C(w, w. w > O. such that for any 8 O, the Lebesgue measure
meas x E [0, 2x] sinn(w/2 cos(x -) - sin(w)) < s Co.
Step 4: Given the previous steps, we now write
Gadv - FNOO ())
h 111 xx (H.N. I(),0e) (CN - -Sx(a))
NN
h 1111/22w/2 hNN No. (CN- -SN)
+ (h - H )1/ 1,000 (CN - SNN)
+ H
TN 0,xx)
CN - Sx((i) - xx (HHN, IN To (cN - Sx(())
(7) + (IT) + (IIT).
32

---

Published as a conference paper at ICLR 2023
The second (II) and third (IIT) terms uniformly bounded by NN1 by the construction of X and H.. By Steps and 3 (with 8 = (2Ca + 1)/N), we can estimate the L -norm of the first term as
cos(æ æ - E) ---- 2-' sin(w)) < 8) S C/N, I s 2h ( (in(w//)
where the constant C' is independent of N, and only depends on the parameters h. h, w and W. Hence, A/FNO satisfies
C
sup JA i Jadv (a)ll S N'
Roops
for constant C' 0 independent of N, and where recall (cp. (.14 above):
6max de C,
depth(AFNO) I Cllog(V)?, size(AVNNO) I C' log((
The claimed error and complexity bounds of Theorem N N.l..
immediate upon choosing
D.5 PROOF OF THEOREM 3,5
To motivate the proof, first consider the Burgers' equation with the particular initial data ū(z) si), with periodic boundary conditions on the interval x E [0, 2T]. The solution for this initial datum can be constructed via the well-known method of characteristics; u(x, t) with initial data ū(x) is smooth for time E [0, 1), develops a shock discontinuity at x = (and x = 2T) for & 1, but remains otherwise smooth on the interval x E (0, 2T) for all times. In fact, fixing a time 2 0, the solution u(x, L) can be written down explicitly in terms of the bijective mapping (cp. Figure
observe that the solution
V, [Ke, 2T xz] - [0, 2m),
Ve(2r)) xp -ssin(r))
where
Ft 0,
for S 1, O is the unique solution of xa = ssnn(z;), for > 1.
(D.15)
We note that for given zo, the curve r (2ro) traces the characteristic curve for the Burgers' equation, starting at zo (and until it collides with the shock). Following the method of characteristics, the solution u(x, t) is then given in terms of V by
u(x, t) = sin (V, (x)), for æ [0, 2m..
(D.16)
We are ultimately interested in solutions for more general periodic initial data of the form ū(x) -sin(x E); these easily be obtained from the particular solution (D.16) via shift. We summarize this observation in the following lemma:
Lemma D.15. Let & E [0, 2T) ) be given, fix a time 2 O. Consider the initial data u(x) sin Et) Then the entropy solution u(x, 1) of the Burgers' equations with initial data á is given by
sin( (a E 2x)), sin( V (I 1 E)), (x 2 E),
(x +
t(x,t)
(D.17)
for x. E [0,2m], 1 0.
Lemma D.16. Let and define U [0, 2m R by U2x) Gin (x)). There exists
A. > 0 (depending on t), such that x + U(x) can be extended to an analytic function Ú Ar, 2n A.) + R, U(x); i.e., such that Ü (x:) U(x) for all [0, 2m).
33

---

Published as a conference paper at ICLR 2023
2
2m
2nT
3m Si
3n
FNO
77
3m
2n
377
277
17
377
2
2
X
X0
XD
(a) & + (To)
(h) V. (2o) at t 0..8
(c) 1/ (e at 1.4
Figure Illustration of V (zo): (a) characteristics traced out by t (x())) (until collision with shock), (b) (x())) before shock formation, (c) (2)) after shock formation, including the
interval [Ty 2T - _t (red limits) and the larger domain (blue limits) allowing for bijective analytic continuation, (A., 2T At).
Corollary D.17. Let U: [0, 2T] + R defined as in Lemma D.16 There exists a constant C > such that for any > 0, there exists ReLU neural network O R - R. such that
(2) - U(x) e, =
sup
and
depth(#.) S C'log(e-)) width( SCC
D.16 U can be extended analytic function Ü (-Ar, 2n A.) R. Thus, by Theorem D.6 there exist constants C, > 0, such that for any L E N, there exists neural network
Proof. By Lemma
1 R
such that
sup (U(x) - xx()) & C Cexp(->111 =E[0.2m]
depth S CL, width( S C.
Given e O, we can thus choose L L -2 log(E -1)2, 2 obtain a neural network 4, Lr, such that (0,2x [Uxx) (x)) S c, and
depth( (2.) S C' Cogg(--1) width($.) C,
for constant C independent of c.
Lemma D.18. Let 1, and let U [ 2n, [-2T, 2 + R be given by
U(x) = sin( v (a 2x)),
(x < <), (x 2 0).
sin
(,- (x)),
There exists constant C = C(t) > depending only f, such that for any (0, }]. there exists a neural network 4. - R, such that
depth(@.) s C C'lo(e- width( SC
and such that
2 - -U S.
Proof. By Corollary D.17 there exists a constant C 0, such that for any e >> 0, there exist neural networks th
R + R. such that
(cpth(#=) s C'log(c-133 width( SC
and
c (x) - -()))(l)=(222000)1 S6, | (x) - Sc.
34

---

Published as a conference paper at ICLR 2023
This implies that
U (2 $ + 10.22=0+11
(1-2x S..
By Proposition D.2 (approximation of indicator functions), there exist neural networks x : R - I with uniformly bounded width and depth, such that
X 11-28,0jll.- S, |xx= - [0.220) S.
and (x))ll>) S 1. Combining this with Proposition D.4 (approximation of multiplication),
follows that there exists neural network
B R +R
$.(x)= ( xe(xx x*) xelxe x)),
such that
4.-
[1-22,00 x 10.2000**1 (1((22224))
xx, - x
(112227277)
+ XM,cc (8) 1 (2(((22.28))
By construction of X have
1xM.((x.. ,
110,20 Φ+
S
x x, x+) - x ++
11
1.1
+ ||| (1++x(10000 - (
111
< 4Te + l1xt 1 ||2 +
$ (4T + c.
And similarly for the other term. Thus, it follows that
| [1]-27,0] + 10,2x; ++]
2(4T + 2)e,
((((222,55)
and finally.
IIU Dellza((-2 2x, S 4T ||U - [1-20,,0) 2- + 10,27 x+] 00 (1))
+ ||[ x- + x+1 D, (1((-22
< 4TE + 2(4 + 2)E (12m + +)e.
for a neural network P.: RR + R size:
depth( s C' log( width( S C.
Replacing e with c/(12m + 4) yields the claimed estimate for Pa..
Based on the above results, can now prove the claimed error and complexity estimate for the shift-DeepONet approximation of the Burgers' equation, Theorem 3.5
Fix >>> 1. Let U : [-2m, 2x] + IR be the function from Lemma D.18 By Lemma D.15 the exact solution of the Burgers' equation with initial data u(x) sin {) at time i, is given by
Proof of Theorem 3.5
u(x,t) Uxx E), Vz E [0, 2T).
From Lemma D.18 (note that - E E [-2 2n ). follows that there exists constant C > such that for any e O. there exists a neural network 4. IR + R. such that
[lu(-,t) - T.( - 3) Sc.
and
depth((.) & log(( -1)2 width( CC
35

---

Published as a conference paper at ICLR 2023
We finally observe that for equidistant sensor points 20, 21. E [0, 2, 2 there exists
matrix A E TR2x3, which for any function of the form g(x) oo + on sin(z) + or cos(x) maps
((1), s(x2)l² + A g(11) = a2)
Clearly, the considered initial data w(x) sin(x {) is of this form, for any & E [0, 2x), or more precisely, we have
u (z) = -in(x - {) = cos(4) sin(E) cos(2))
so that
A (a(2o), ū ((x1), a(x2))" [Cos(E), sin
As next step, recall that there exists C 0, such that for any there exists a neural network
E, : IR - [0, 2n] (cp. Lemma D.8). such that
sup & (cos(()), sin(E)) çe(o,2m-c)
such that (cc) E [0, [m] for all ç. and
dopth(E.) v C' C'log(c-1) width(E,)) S C.
Based on this network E., now define shift- Shift-DeepONet approximation of of ASDON (ū) x
9urgg(ü of size:
depth(AssDO)) s C'log(e width(A>DON) S SC
by the composition
M (T)(x) - EA ((X))
(D.18)
where ū(X)= == [ū( (0). ū(x1). ((22) and note that (denoting E. == Ee(A (())), we have for & [0, 2 - c|:
Burg (7) AN (ū) = |IU( 5) - -.( - =
$ |IU( - -(- -E.--
+ [IU[ - -Ee) -0.(x - -
< < - +e
(C++
where only depends on U, and is independent of > On the other hand, for & 2 we have
1 (1) - DON (ū)ll
IIU( 5) a =e
< |IT 1) + 1. - =
2n (IIU(- =) + +| =
< 6TT,
uniformly bounded. It follows that
24-c
Easu [lrg (a) - -DDON (T)lla
| (γ) ADDON
(T)ll
S 2m (C + 1)e 6TE,,
with a constant C > 0, independent of > 0. Replacing b C' for sufficiently large constant C' (depending only on the constants the last estimate above), one readily sees that there exists shift- Deli-DeepNeeteeeeeee ADDON, such that
8 ((N=DON) = En-n ( (ū) - -AaDON(N (T)]l:. Sc
and such that
width(AADOON < C, depth(AssDON) < Clog(e -1)2
and size(A>ooo))) & C'log(e-1)) for a constant > 0, independent e > 0. This concludes our []
proof.
36

---

Published as a conference paper at ICLR 2023
D.6 PROOF OF THEOREM 3.6
Proof. Step I: Assume that the grid size is N 2 3. Then there exists an FNO Wi, such that
Ni( [cos(E), sin(E)).
and depth(AV1) S C. d. I C, kiax = 1.
To see this, we note that for any ç E [0. 2 the input function ū(æ) -si((r E) -cos(E) sin(x) - sin(E) cos(æ:) can be written in terms of sine/cosine expansion with coeffi- cients cos(g) and Sin(E). For 2 grid points, these coefficients can be retrieved exactly by discrete Fourier transform. Therefore, combining a suitable lifting to dy = 2 with Fourier multiplier matrix P. we can exactly represent the mapping
[cos(E) sin(ar) sin(() cos(1r)
á 7 F
u)) (P-
by a linear Fourier layer. Adding a suitable bias function b(x) = [sin (æ), cos(x))) and composing with an additional nonlinear layer, it is then straightforward to check that there exists (ReLU-)FNO, such that
cOs(E) sin (r) + sin(r)
(coo(E) + 1) sin (z) (sin(E) + cos (z)
sin(E) cos (z) co)
=
cos(E) + 1 sin [x)
sin(E) + 1|| cos (z)]
cos(E) + Σ .v ==
sin (2)
sin({) +
cos(2;))
cos(E) + 1 1-11-11
sin(E) + - = coo((E)
sin(E)
Step 2: Given this construction of M1. the remainder of the proof follows essentially the same Theorem argument as in the proof Theorem 3.5 We again observe that the solution u(æ, t) with initial data u(x) sin(æ 5) well approximated by the composition
NFNO == (x (N))
such that (by verbatim repetition of the calculations after (D.18 for shift-DeepONets)
(NNPNO) Ba-n [ (ū) - AFNO (ū) [ S
and where & : a ReLU neural network of width
depth(#.) C' Clooge-1)), width (b) 5 SC
and E. R + [0] is an ReLU network with
depth(E.) S C' C'log(=-)),
(e.th(E.)
Being the composition of an FNO A1 satisfying kmax 1, d, C, depth(Ai) S C' with the two ordinary neural networks t and E., it follows that that AFNO can itself be represented by an FNO with kmax d, I and depth(A/FNO) > C' C'log(--)22 By the general complexity estimate (B.2).
(iNo(A/FNO) S d. denaa
also obtain the claimed an upper complexity bound size(A (NNNO) I (looge--1)
LI
E DETAILS NUMERICAL EXPERIMENTS AND FURTHER EXPERIMENTAL RESULTS.
E. TRAINING AND ARCHITECTURE DETAILS
Below, details concerning the model architectures and training discussed.
37

---

Published as a conference paper at ICLR 2023
E.1.1 FEED FORWARD DENSE NEURAL NETWORKS
Given an input y E _ feedforward neural network (also termed multi-layer perceptron), transforms it to an output, through layer of units (neurons) which compose of either affine- maps between units (in successive layers) or scalar nonlinear activation functions within units Goodfellow et al. (2016), resulting in the representation.
uo(y) = Ct. a C1-1.... o C2 o Ci(v).
(E.1)
Here, o refers to the composition of functions and o is scalar (nonlinear) activation function. For any < I., we define
C+ = Weze + be, for W E Rede xdz E &R=, by E Rdeit.,
(E.2)
and denote,
(We,
(E.3)
to be the concatenated set of (tumable) weights for the network. Thus in the terminology machine learning, feed forward neural network (E.I) consists of an input layer, output layer, and L hidden layers with de neurons, < <. L.. In all numerical experiments, we consider a uniform number of neurons across all the layer of the network de d d, < l < L. The number of layers I., neurons d and the activation function a are chosen though cross-validation.
E.1.2 RESNET
A residual neural network consists of residual blocks which use skip or shortcut connections to facilitate the training procedure of deep networks He et al. 2016). residual block spanning le layers is defined as follows,
r(ze, =xxx = (W,2 + be) + ze-k-
(E.4)
In all numerical experiments we set k = 2.
The residual network takes input sample function (ā(z1).... ū(zm)) R" and outputs the output sample ( E y encoded
encoded at the Cartesian grid points
(zi...., ), ,. {(ū) the same of points, (G((i) = (((()((1)) (())(am)) IR. For the compressible Euler equations the encoded input is defined as
E(ū) = (o(
Po uo ( (m Eo (o(x1) ..... Eo (2m))
E(U) - ... (p
po u0 (Im²),
po (o(21)00(21) (nfx))) Eo (Wmi)) RR 1m4m
(E.5)
for the and 2d problem, respectively.
E.1.3 FULLY CONVOLUTIONAL NEURAL NETWORK
Fully convolutional neural networks are special class of convolutional networks which are inde- pendent of the input resolution. The networks consist of an encoder and decoder, both defined by composition linear and nonlincar transformations:
Eg_ (v) Ci, a o Ch.11
C (i()) o C oo (2),
Do, (z) Ct oa C
(E.6)
ug (w) Do Eg. (v).
The affine transformation C2 commonly corresponds to convolution operation in the encoder, and transposed convolution (also know deconvolution), in the decoder. The latter can also be performed with simple linear (or bilinear) upsampling and a convolution operation, similar to the encoder.
The (de)connolttion is performed with kerne We R (for Id-problems. and We R for 2d ad-problems), stride and padding p. takes as input a tensor 2; E " " (for 1d- Id-problems, and 24 IR wy xhexee for 2d-problems), with e being the number of input channels, and computes
38

---

Published as a conference paper at ICLR 2023
2051
R04
16
N
26
2715
H
(
Figure 5: Fully convolutional neural network architecture for the linear advection equation and shock tube benchmarks. (2²) = BN o a o C* (2), with BN denoting a batch normalization and C a convolution defined by the tuple (3, 2, 1, Cin c), with Cin for the advection equation and Cin = shocktube benchmarks. F(z) BN o a o C1 o BN o o C3 o BN c a C5 BN oos o (2) with Cf, C2 C3 C₁ identified with (3, 2, Cins 20n)) (1, 1,0, 0, 20m, 20n)). ) (1,1, 0, 20m 20m)). (3, 2, 1, 2ci,, 20m.. (z) BN a a C4 BN oa 03 BN oa O2 o BN oa (z), with CT
Cid. Cg Oj transposed convolutions defined by (3,2, 1, Ciny Cn), ) (1, 1, Ciny Cin). (1, 1,1,0, O, Cins Cin)
(3, 2, Cin, 0.3cmn). G is a transposed convolution defined by the tuple (3, 2, 1,1,,,
1024
H
((2))
Figure 6: Fully convolutional neural network architecture for the Burgers' equation benchmark: B, F. G follow the same definition the caption of figure 5 H (z) o BN of Poo BN o o
UPoa BN o C2 oPo BN a U (z), with Cf. C5. Cg. O2 being standard convolutions
defied by the tuples (3. 2, 1, Can: Cin), (1, 1, 0, Cin: Cin), (1, 1, 0, Cany Cin), (3, 2, 1, Ciny 0.5cin). and
U P denoting up-sampling operation with scaling factor 2.
255
120
π7
D
(Tz)
Figure 7: Fully convolutional neural network architecture for the 2-dimensional Riemann problem.
B²2) = a BN C5 a BN o Cf, with Cf and Cj being standard convolutions defined by the tuples (3, 1, 0, Cm. 2cm)). (3, 1, 0, 20 2cmn). D(z) B(z) M P(z), with MP being max pool with kernel size 2. U(Z)= = B(2) UP),,, where U P denotes an up-sampling with scale factor 2. G is convolution defined by 1,,,,
2t-1 1wr+1xee Xce+1 (for Id-problems, and 2 _ Xcr- for 2d-problems). Therefore, (de)convolutional affine transformation can uniquely identified with the tuple (ke, s, Cts ct+1).
The main difference between the encoder's and decoder's transformation is that, for the encoder het1 < he, Wi++1 < We, ce+1 > ce and for the decoder h4+1 > fgg Wi+1 > wig ce+1 < cc.
For the linear advection equation and the Burgers' equation we employ the same variable encoding of the input and output samples as ResNet. On the other hand, for the compressible Euler equations, each input variable is embedded in an individual channel. More precisely, E(T) IR m for the shock-tube problem, and (()) E R mm x4 for the 2d Riemann problem. The architectures used in the benchmarks examples shown in figures 6
In the experiments, the number of channel c (see figures 5 6
for an explanation of its meaning)
and the activation function selected with cross-validation.
39

---

Published as a conference paper at ICLR 2023
E.1.4 DEEPONET AND SHIFT- DEEPONET
The architectures of branch and trunk are chosen according to the benchmark addressed. In particular, for the first two numerical experiments, we employ standard feed-forward neural networks for both branch and trunk-net, with skip connection between the first and the last hidden layer in the branch.
On the other hand, for the compressible Euler equation composition L blocks, each defined as:
convolutional network obtained as
B BN co <<t<
(E.7)
with BN denoting batch normalization. The convolution operation instead defined by h 8 = 2. = 1, Cg C++11 = d, for all < L < L - 1. The output is then flattened and forwarded through a multi layer perceptron with 2 layer with 256 neurons and activation function ø.
For the shift and scale-nets of beitt-DeepNNet we use the same architecture as the branch.
Differently from the rest of the models, the training samples for DeepONet and shift-DeepONet are encoded at m and n uniformly distributed random points, respectively. Specifically, the encoding points represent randomly picked subset of the grid points used for the other models. The number of encoding points m and n. together with the number of layers L, units d and activation function of trunk and branch-nets are chosen through validaaali....n..
E.1.5 FOURIER NEURAL OPERATOR
We use the implementation of the FNO model provided by the authors of et al. (2021a). Specifically, the lifting defined by a linear transformation from Id xm to Rd,.m, where dy is the number of inputs, and the projection Q to the target space performed by neural network with single hidden layer with 128 neurons and GeL.U activation function. The same activation function used for all the Fourier layers, as well. Moreover, the weight matrix We used in the residual connection derives from a convolutional layer defined by (ke = L, = 1, 0, Ce dy. Ct+11 = d,),, for all 1 < l < L-1. 1. We use the same samples encoding employed for the fully convolutional models. The lifting dimension du the number of Fourier layers L and kimzs defined in 2 are the only objectives of cross- vassdvaaiit
E. TRAINING DETAILS
For all the benchmarks, a training set with 1024 samples, and validation and test set each with 128 samples, are used. The training is performed with the ADAM optimizer, with learning rate 5 10-4 for 10000 epochs and minimizing the ll-so function. We use the learning rate schedulers defined in table We train the models in mini-batcher of size 10. A weight decay of 10 6 is used for ResNet (all numerical experiments). DON and sDON (linear advection equation, Burgers' equation, and shock-tube problem). On the other hand, no weight decay is employed for remaining experiments and models. At every epoch the relative L Lerror is computed on the validation set, and the of trainable parameters resulting in the lowest erтor during the entire process saved for testing. Therefore, no early stopping is used. The models hyperparameters are selected by running grid searches over range of hyperparameter values and selecting the configuration realizing the lowest relative L eeror on the validation set. For instance, the model size (minimum and maximum number of trainable parameters) that are covered in this grid search are reported in Table
The results of the grid search i.e., the best performing hyperparameter configurations for each model and each benchmark, are reported in tables
and
E.2 FURTHER EXPERIMENTAL RESULTS.
In this section, we present some further experimental results which supplement the results presented in To the main text. We start by presenting more statistical information about the median emrors shown in Table I To this end, in Table we show the errors, for each model on each benchmark, corresponding the 0.25 and 0.75 quantiles, within the test set. We observe from this table that the same trend, as seen in Table [ also holds for the statistical spread. In particular, FNO and shift- Deftt-eeeOONe outperform DeepONet and the other baselines every experiment, Similarly, FNO
40

---

Published as a conference paper at ICLR 2023
ResNet
FCNN
DeepONet
Shift DeepONet
FNO
Step-wise decay 100 steps,
Step-wise Decay 50 steps, 0.99
Step-wise decay 100 steps,
Advection Equation
Exponential decay 7 0.999
None
0.999
7 0.999
Step-wise decay 100 steps,
Step-wise Decay
Step-wise decay 100 steps,
Burgers' Equation
Exponential decay
50 steps,
None
= 0.999
7 0.999
? 0.99
y 0.999
Step-wise decay
Step-wise Decay steps. 0.99
Step-wise decay 100 steps.
Shocktube Problem
Exponential decay 0.999
steps.
None
0.999
7 = 0.999
Step-wise decay 100 steps.
Step-wise Decay
2D Riemann Problem
Exponential decay
Exponential decay
50 steps.
None
y 0.999
0.999
0.999
= 0.99
Table 2: Learning rate scheduler for different benchmarks and different models: 7 denotes the learning rate decay factor
ResNet
FCNN
DeepONet
Shift DecpONet
FNO
Advection
576,768
1,156,449
519,781
1,018,825
22,945
Equation
1,515,008
1,8240,545
892,561
1,835,297
352,961
Burgers' 313,600
1,155,025
519,781
1,018,825
22,945
Equation 989,696
18,219,489
892,561
1,835,297
352,961
Shocktube
1,101,056
1,156,497
1,344,357
3,492,553
23,009
Problem 2,563,584
18,240,737
3,190,673
8,729.633
353,089
2D Riemann
42,059,008
442,985
361,157 2,573,513
821,961 7,082,785
268,833 13,132,737
Problem
84,416,000
7,066,529
Table 3: Minimum (Top sub-row) and maximum (Bottom sub-row) number of trainable parameters among the grid-search hyperparameters configurations.
Advection
Burgers Equation
Shocktube Problem
Riemann
Equation
Problem
Leaky ReLU
Leaky ReLU
Leaky ReLU
Leaky ReLU
L
4
8
8
d
128
256
256
256
Trainable
576,768
989,696
2,563,584
84,416,000
Params
Table 4: ResNet best performing hyperparameters configuration for different benchmark problems.
Advection
Burgers' Equation
Shocktube Problem
2D Riemann Problem
Equation
Leaky ReLU
Leaky ReLU 16
Leaky ReLU 16
ReLU 16
8
Trainable
1.156.449
4,576,033
4,581,537
1,768,401
Params
Table 5: Fully convolutional neural network best performing hyperparameters configuration for different benchmark problems.
41

---

Published as a conference paper at ICLR 2023
Advection Equation
Burgers' Equation
Shocktube Problem
2D Riemann Problem
256
512
512
64
n
512
256
256
256
p
200
50
200
100
Liranch
3
3
Ltrunk
6
4
6
doranch
256
256
256
32
derunk
256
256
256
256
Obannhh
Leaky ReLU
Leaky ReLU
Leaky Le
SoftSign
Trrunk
Leaky ReLU
Leaky ReLU
Leaky ReLU
Leaky ReLU
Trainable
761,233
618,085
3,190,673
607,433
Params
Table 6: DeepONet best performing hyperparameters configuration for different benchmark problems.
Advection Equation
Burgers' Equation
Shocktube Problem
Riemann
Problem
WA
$12
512
512
256
512
512
512
128²
P
50
200
100
Loranch
3
4
4
Lteunk
6
6
doranch
256
256
256
32
derunk
256
256
256
256
Chranch
Leaky ReLU Leaky ReLU
Leaky ReLU Leaky ReLU
Leaky ReLU Leaky ReLU
Leaky ReLU Leaky ReLU
Orrunk
Trainable
1,445,321
1,835,297
6,047,633
6,851,785
Params
Table 7: Shift DeepONet best performing hyperparameters configuration for different benchmark problems.
Advection Equation
Burgers' Equation
Shocktube Problem
2D Riemann Problem
kmaz
15
19
7
15
d,
64
32
32
64
2
4
4
Trainable Params
148,033
90,593
41,505
8,414,145
Table 8: Fourier neural operator best performing hyperparameters configuration for different bench- mark problems.
42

---

Published as a conference paper at ICLR 2023
ResNet
FCNN
DeepONet
Shift DeepONet
FNO
Advection Equation
10.1% 22.8%
8.2%
17.3%
4.9% 13.8%
1.4% 5.4%
0.35% 1.25%
Burgers' Equation
18.8% 22%
20.3%
25.7%
25.4%
32.4%
6.7% 9.6%
1.3% 1.9%
Shocktube Problem
3.6%
5.6%
7.0%
10.25%
3.4%
5.4%
2.0%
3.75%
1.2% 2.1%
2D Riemann Problem
2.4% - 2.9%
0.17%
0.21%
0.77%
1.1%
0.10%
0.15%
0.10% 0.14%
Table 9: 0.25 and 0.75 quantile of the relative L error computed over 128 testing samples for different benchmarks with different models.
Test Samples
ResNet
FCNN
DeepONet
Shift DeepONet
FNO
128
10.1, 14.8, 22.8
8.2, 11.6, 17.3
4.9,8.0, 13.8
1.4, 5.4
0.25, 0.71, 1.25
256
11.1, 15.5, 22.8
8.0, 11.6, 16.8
5.4,9.2, 15.1
1.5, 25,,28,65 5.1
0.36, 0.63, 1.3
512
11.3, 15.9, 23.6
8.3, 11.6, 16.8
9.6,9.6 15.2
1.4, 6.5
0.36, 0.65, 1.36
1024
11.0, 15.25, 23.1
8.0, 11.35, 16.0
5.4,9.1,555 15.4
1.4,288 5.9
0.35, 0.62, 1.29
2048
11.31, 15.7, 7, 24.0
8.2, 11.5, 16.8
5.4, 9.3, 15.8
1.5, 2.9, 6.6
0.36, 0.67, 1.37
Table 10: 0.25, 0.5 and 0.75 quantile of the relative I. error computed over 128, 256, 512, 1024 and 2048 testing samples for the linear advection equation with different models.
FNO Median Testing L' -error
FFT Median L'-error
kmar
dg
Trainable Params
192
247,169
2.30%
164.2 %
160
252,097
1.21 %
137.9
128
263,169
0.48 %
63.3 %
241,113
0.54 %
38.9
Table 11: Testing error obtained by training FNO with different number of modes for the linear advection problem and corresponding error obtained with linear Fourier projection.
FNO Median Training L -error
FNO Median Testing -e-error
FFT Median L--error
kmaz
L
d,
Trainable Params
0
3
64
33,409
2.2%
2.38% %
164.2 %
1
3
64
45,697
0.70%
0.80 %
137,9 %
3
3
64
70,273
0.54%
0.55 %
63.3 %
7
3
64
119,425
0.42%
0.46 %
38.9 %
Table 12: Training and testing error obtained by training FNO with different number of modes for the linear advection problem and corresponding error obtained with linear Fourier projection.
43

---

Published as a conference paper at ICLR 2023
Trainable Parameters
DON Median Training L -error
SDON Median Training L -error
DON Median
SDON Median
p
DON / SDON
Testing -error
Testing L -error
658947/ 1054213
100.1%
4.1%
100.1%
13.1%
2
659461 / 1055497
100.0%
3,5%
100.0%
5.0%
4
660489/ / 1058065
40.9%
2.5%
37.5%
5.0%
8
662545 / 1063201
26.6%
1.6%
29.7%
4.1 %
16
666657/ / 1073473
16.5%
3.1%
15.6%
4.1%
32
674881 1094017
10.9%
3.0%
9.8%
4.6%
64
691329 1135105
10.8%
2.2%
11.1%
3.9%
Table 13: Training and testing error obtained by training DON and SDON for the linear advection problem with different number of basis functions.
10a
Training Error DON
Tasting Emor DON
Training Error
Tasting Emor SDON
10a
10a
10a
Figure 8: Testing error obtained by training DON and SDON with different number of basis function p for the linear advection problem.
outperforms shift-DeepOOi handily on each experiment, except the four-quadrant Riemann problem associated with the Euler equations of gas dynamics.
The results presented in Table 9 further demonstrate that FNO is the best performing model on all the benchmarks. In order to further understand the superior performance of FNOs, we consider the linear advection equation. As stated in the main text, given the linear nature of the underlying operator, a single FNO Fourier layer suffices to represent this operator. However, the layer width will grow linearly with decreasing error, requiring kmax use the nonlinear reconstruction, as seen in the proof of Theorem 3.3 to obtain good performance. To empirically illustrate this, compare FNO with different kmax (number of Fourier modes) with corresponding error obtained by simply projecting the outputs of the operator into the linear space spanned by the corresponding number of Fourier modes. This projection onto Fourier space amounts to the action of a linear version of FNO. The errors, presented in Table TT (keeping the total number of degrees of freedom roughly constant) and Table d, = 64 fixed), clearly show that as predicted by the theory, very few Fourier modes with a kmax = 1 are enough obtain an error of approximately 1%. On the other hand, the corresponding error with just the linear projection is two orders of magnitude higher. In fact, one needs to project onto approximately 500 Fourier modes to obtain an error of approximately 1% with this linear method. This stark contrast is further illustrated in Figure 9 which shows that FNOs relying on nonlinear reconstruction mechanism can outperform a competing method based on linear Fourier reconstruction (Fourier linear projection) by one two orders of magnitude in terms of approximation accuracy. This experiment clearly brings out the role of the nonlincarities in FNO in enhancing its expressive power advection-dominall problems.
Fourier modes. Hence, it is imperative to
the lifting dimension
44

---

Published as a conference paper at ICLR 2023
10a
Projection
10
Iraining Enro HNO Testing Error FNO
10a
10a
10a
Figure 9: Testing error obtained by training FNO with different number of modes for the linear advection problem and corresponding error obtained with linear Fourier projection.
Besides our analysis of FNOs, we have shown by means of approximation theoretical arguments that shift-DeepONets, which are based on nonlinear reconstruction, overcome the fundamental limitations faced by DeepONets. Clearly, the approximation theoretic arguments employed in this work neglect some aspects of the practical training of neural operators. most notably possible errors due to finite number of training samples and errors due to numerical optimization procedure (e.g. convergence to sub-optimal local minima, etc.). In order to further understand the superior performance of shift-DeepONee over DeepONets, carry out scan in the number of basis functions p.
For the linear advection problem, first recall that the fundamental limitation of DeepONets is due to the lower bound of Theorem 3.1] which shows that DeepONets can achieve at- ar-best linear decay error ~ 1/p in terms of the number of DeepONet basis functions p. To compare this theoretical prediction with empirical observation, we carry out a scan over different choices of = 1,2, 4..... 64 in Table 13 The table clearly shows that the DeepONet training/testing errors both decrease with increasing p. Furthermore, inspection of the corresponding an approximately linear decay, which is only slightly worse than the (optimal) lower bound 1/p. consistent with approximation theory.
reveals that this decrease follows
In contrast, our approximation theoretic result shows that shift-DeepONets can overcome the lower bound, and achieve high accuracy even with bounded number of basis functions S C' (cp. Theorem 3.2). This reflected by the almost immediate saturation of the error as function of in Figure 8 (cf. also Table 13). Consistent with the theoretical insight, our numerical experiment demonstrates that (a) the error for Dest (based on nonlinear reconstruction)
considerably smaller than the corresponding error for DeepONets (based on linear reconstruction), and (b) increasing the number basis functions beyond very moderate values of p ~ 4, 8 does not improve the accuracy of shift-DeepONets, consistent with the fact that uniformly bounded number p S C is sufficient.
45

---

Published as a conference paper at ICLR 2023
M0,025)
0.0
0.0
0.0
0.
644
0.6
1.0
0.
0.4
0.6
0.0
a) Sample
1.0
0.4
08
0.0
0.2
0.4
0.6
1.5
a) Sample 2
Figure 10: Illustration of two input (blue) and output (orange) samples for the advection equation.
2.0
at. 0.1)
13
0.5
G
-18
-18
-15
-2.8
00
9.
1.0
60
02
0.4
0.
08
a) Sample
G
-05
-15
-15
a) Sample 2
Figure 11: Illustration of two input (blue) and output (orange) samples the Burgers' equation.
46

---

Published as a conference paper at ICLR 2023
63-
0.4
03
02
6.1
0.0
E(.15)
3.00
+3.00
0.4
033
-04
(a) Sample 1
0.7
2.0
05
02
5.0
11.15
8.15
-3.05
-0.10
0.2
0.3
(b) Sample
Figure 12: Illustration of two input (bluc) and output (orange) samples for the shock-tube problem.
47

---

Published as a conference paper at ICLR 2023
1.80
A
Pally
1.00
E11.0)
1
1236
0.50
0.5
0.20c
0.500
0.12
07100
1.844
6.55
-035
0.500
0.666
0.50
7.00
0332
0.2
6.5
-0.75
0.264
0.000
11.00
0.504
0.352
5.79
1.009
0.75
0.882
030
0.756
025
0.830
3.82
0.00
0.33
2.56
2.24
-0 50
2.43
-0.50
0 136
0.252
1.82
3.02
-0.5
Sample
0.35
2.004
8.50
5.976
0.75
0.914
5311
6.50
0.00
0.544
0.50
6.744
-0.00
6.524
0354
2454
23356
0.75
1.80
A/e
Eu
0.10
4.72
0.75
6,755
3.435
2.000
-2.25
0.50
0.2
8.555
8,700
21122
325
6.510
2.56
0.25
0330
3.855
6:100
5
1.40
434
-63
02
03
1.0
0335
-45
5.125
23.00.0 -0.75
-0.53
-025
0.00
0.25
0.50
1.00
1.22
(b) Sample 2
Figure 13: Illustration of two input (left) and output (right) samples for the 2-dimensional Riemann problem.
48


Reviews:
[{"review_id": "pJqxjNme9xx", "reviewer": "Reviewer_vjUn", "summary": "", "questions": "", "limitations": "", "rating": 8, "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": "This paper provides rigorous analysis of operator learning of PDEs with discontinuities. The main contributions are two-fold. First, author prove that frameworks using linear reconstruction (DeepONet or PCA-Net) fail to efficiently capture the discontinuities by proving the lower bound of the approximation error decay merely linearly in terms of size of DeepONet. Second, the authors prove that two nonlinear reconstruction architectures, Shift-DeepONet and FNO, outperform the linear reconstruction architectures in the sense that much smaller model size is required. And the theoretical results are firmly supported by extensive numerical experiments.", "summary_of_the_review": "This paper is the first paper rigorously investigate the error estimate of operator learning of PDEs with discontinuities. The authors carefully illustrate why linear reconstruction approaches fail to efficiently handle discontinuities, and how nonlinear reconstruction fix this theoretically and numerically. Additional numerical experiments that support the quantitative result (linear rate in DeepONet and exponential rate in FNO and Shift-DeepONet) are expected.", "strength_and_weaknesses": "Strengths:\n1.This paper is the first theoretical paper investigating the operators learning of PDEs with discontinuities.\n2.The authors provide quantitative theoretical results and rigorous proof.\n3.The paper is well-organized and well-written.\n4.Extensive and supportive numerical experiments are provided.\n\nWeakness:\nNumerical examples do not illustrate the quantitative result, for example, the linear decay of model size of DeepONet, exponential decay rate of FNO and Shift-DeepONet.\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works.", "empirical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works."}, "extra_text": {}}, {"review_id": "ZhlK4sNs_f", "reviewer": "Reviewer_koS7", "summary": "", "questions": "", "limitations": "", "rating": 3, "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": "This paper tries to investigate the operator learning of PDEs with discontinuous solutions. They prove that the linear method cannot approximate efficiently the solution operator. Furthermore, they propose a nonlinear reconstruction mechanism to overcome these fundamental difficulties.   And they propose new algorithms such as Fourier Neural Operators and a novel extension of DeepONet termed shift-DeepONet. Finally, they provide a rigorous theory for the empirical results. ", "summary_of_the_review": "Using the method of deep learning to solve PDE is a very interesting topic.   The basic numerical method for the nonlinear PDE is\n also unknown. Moreover, the motivation, the dynamic process and any kind of approximation still need to be presented clearly.  ", "strength_and_weaknesses": "Strength: The topic of learning the solution of PDE, especially the discontinuous solution, is very interesting.\n\nWeaknesses: (1) The results in all theorems, from Theorem 3.1 to Theorem 3.3,  look to learn the solution to the PDE at some time T. They \n                            do not tell the readers what kind of numerical method to get G_adv = u(\\cdot, T).\n\n                       (2) The result is about the solution at some T. Does the approximation depend on the time T? \n\n                       (3)  How does the dynamic approximation estimation evolve with time?\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "2: Several of the paper’s claims are incorrect or not well-supported.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "1: The contributions are neither significant nor novel."}, "extra_text": {}}, {"review_id": "WMZiuz9pP6", "reviewer": "Reviewer_ZZzU", "summary": "", "questions": "", "limitations": "", "rating": 8, "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": "The authors consider the problem of learning the solution operator of PDEs with discontinuous solutions (such as hyperbolic conservation laws) from initial data/solution pairs with initial data draw from a given measure.\nThe authors consider the DeepONet approach for operator representation that uses a sort-of finite sum over basis functions leveraging a \"branch-net\" for learning the coefficients and a \"trunk-net\" for learning the \"basis-functions\".\nBased on an existing lower bound on the expected error of DeepONet given the eigenvalues of the covariance of the push-forward measure obtained from the solution operator the authors motivate a variation of DeepONet called shift-DeepONet that through a shift-net enables learning shift relations between initial data and solution that are typical for hyperbolic conservation laws where shock discontinuities in the solution can shift as a result of perturbing the initial data.\nThe authors also consider Fourier neural operators (FNOs) as an alternative for better learning solution operators for discontinuous solutions because these also do not suffer from a similar dependence on the covariance spectrum of the push-forward.\nIn the theory section of the paper the authors prove that for a simple linear transport with initial data measure representing box-type initial data DeepONet size has to scale quadratically as error decays while their shift-DeepONet size only has to scale linearly and an FNO only has to scale logarithmically.\nFurthermore, for a simple nonlinear scalar conservation law (Burgers' equations) with sinus initial data and solution after shock-formation the authors prove that DeepONet has to scale linearly, while shift DeepONet and FNO only have to scale logarithmically.\nIn the experiments section the authors again demostrate that both shift-DeepOnet and FNO have better relative approximation error than DeepONet and some other alternatives for a linear advection, a Burgers, an 1D Euler (shock tube) and a 2D Euler (Rieman) problem.", "summary_of_the_review": "I recommend to accept the paper and would ask the authors to try to make it more accessible by providing a better intuition on the issue of the eigenspectrum decay for the covariance of the push-forward with discontinuous solutions. I like the combination of proposing the new shift-parameterization of the model which seems well motivated, providing a theoretical analysis and a comprehensive set of experiments. I'm not sure whether the ICLR audience is the right audience for this sort of problem but the fourier neural operator is an ICLR paper so it could be okay. I also feel like although the theoretical analysis is high quality and I like the insight it provides it is not all that general since it still only pertains a set of examples and does treat any more general class of operators.", "strength_and_weaknesses": "Strength:\n- Clear structure and motivation, original contribution of the shift idea to the DeepONet operator learning approach\n- High quality mathematical analysis of the given example problems\n- A comprehensive set of example PDEs / initial value instances for the experiments section\n\nWeakness:\n- Intuition of push forward covariance eigenvalue distribution does not become clear to someone who is not already familiar with what that looks like. Some explanation of what eigenvalues would look like\nfor a simple example of a solution operator for a smooth solution vs discontinuous solution would be helpful.\n- In the experiments sections it's not entirely clear why only 128 test examples were used to compute relative errors (that would still allow for relatively high error variance). A training set of 1024 example also seems quite small. Also adding standard defiation bars around the errors would be helpful.\n\nOther remarks:\n- Typo: \"In particular, this bound continuous to hold\" -> should be \"continues\" instead?\n- Typo: \"described in detailed\" -> \"described in detail\"\n- Typo: \"Riemman\" -> Riemann\n- I'm not sure if e.g. Theorem 3.2 should be a theorem or just a proposition as you're considering an example initial data measure and perhaps the generality is not as great as to warrant calling it a theorem. There are perhaps possibilities for achieving more generality by treating general piecewise constant solutions of scalar conservation laws (although I'm not sure having more generality is necessary for making the point that the paper makes).\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works.", "empirical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works."}, "extra_text": {}}, {"review_id": "BEnua9oQ15", "reviewer": "Reviewer_sk7n", "summary": "", "questions": "", "limitations": "", "rating": 8, "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": "The main premise of the paper is that for hyperbolic conservation laws PDEs, operator learning frameworks like DeepOnets are inefficient as compared to counterparts like Neural Fourier Operators. This is because DeepOnets (and PCA-Nets) are linear reconstruction methodologies (i.e, the trunk and branch net outputs are linearly combined) which effects their representational capabilities, whereas FNO based methodology have a nonlinear reconstructions step. \n\nThe authors therefore extend the DeepOnet framework to a shift-deepOnet Framework with an extra shift and scale terms that are learned. \n\nThe proofs show that for one dimensional Linear Advection equation and inviscid Burgers equation, the size of the linear reconstruction methods scale at least quadratically in the reconstruction error whereas for the nonlinear reconstruction method size scales linearly only.", "summary_of_the_review": "This paper theoretically and empirically shows that conditions under which DeepOnets are less efficient than methodologies like FNOs and also provides an extension of DeepONets that can make up for the difference.", "strength_and_weaknesses": "The paper shows a useful results which delineates DeepONets are not as efficient as FNOs for approximating PDEs with discontinuous solution. \nThey introduce shift-DeepONet, a variant of deepOnet which performs as well as FNOs on pdes with discontinuous solutions. \nThe derive a lower bound (based similar proof techniques used in prior works) for the approximation error of DeepONet for advection equation and burgers equation (in one dimension) and show that FNO and shift DeepONet can upper bound the solution more efficiently. While I didn't go into the details of all the proofs, the overall sketch and techniques seem correct. \n\nThe authors back their theory with empirical results where they show that shift-DeepONet and FNO perform better than normal DeepOnet and resnets on PDEs like linear advection, burgers equation etc. \n\nOne small question i have is how would the upper bound change with increasing the dimension? A trivial way would be to treat each dimension independently, but I wonder if the authors have any other (more efficient way) that they have in mind?", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}]

Author rebuttal:

Decision: Accept: notable-top-25%