On the Representation of Solutions to Elliptic PDEs in Barron Spaces
Year: 2022

Abstract:
Numerical solutions to high-dimensional partial differential equations (PDEs) based on neural networks have seen exciting developments. This paper derives complexity estimates of the solutions of d-dimensional second-order elliptic PDEs in the Barron space, that is a set of functions admitting the integral of certain parametric ridge function against a probability measure on the parameters. We prove under some appropriate assumptions that if the coefficients and the source term of the elliptic PDE lie in Barron spaces, then the solution of the PDE is -close with respect to the H norm to a Barron function. Moreover, we prove dimensionexplicit bounds for the Barron norm of this approximate solution, depending at most polynomially on the dimension d of the PDE. As a direct consequence of the complexity estimates, the solution of the PDE can be approximated on any bounded domain by a two-layer neural network with respect to the H norm with a dimension-explicit convergence rate.

heading: 1 Introduction
text: Inspired by the tremendous success of deep learning in diverse machine learning tasks including image classification, natural language processing, and artificial intelligence, there has been growing interest in exploring scientific and engineering applications of deep learning [36, 32, 34, 26, 47]. As partial differential equations (PDEs) play a fundamental role in almost all branches of sciences and engineering, numerical solutions to PDE problems based on neural networks have become an important research direction in scientific machine learning [25, 6, 23, 17, 7, 10, 22]. Among the various directions, numerical solutions to high-dimensional PDEs – the unknown function depending on many variables – are perhaps the most exciting possibility, as solving such PDEs has been a long-
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
standing challenge and breakthrough would lead to tremendous progress in fields such as many-body physics [4, 11, 18], multiple agent control [35, 17], just to name a few.
Numerical solutions to low-dimensional PDEs, such as Navier-Stokes equation in fluid dynamics, has become a standard practice after decades of work. However, the computational cost of the conventional numerical methods for PDEs grows exponentially with the dimension, as a manifestation of the curse of dimensionality (CoD). Given a target accuracy , conventional methods, such as finite element or finite difference, would need a mesh size ofO( ), and thus degree of freedom on the order of O( −d), where d is the dimension of the problem. Such complexity severely limits the numerical solutions to PDEs in high dimension, such as the many-body Schrödinger equations from quantum mechanics and the high-dimensional Hamilton-Jacobi-Bellman equations from control theory. Neural networks, in particular deep neural networks, provide a promising way to overcome the CoD in representing functions in high dimension. It is thus a natural idea to parametrize the solution ansatz to a PDE as neural networks and to employ variational search for the optimal parameters. Various neural network methods [25, 6, 7, 10, 33, 17, 41, 46, 16, 5] for PDEs have been proposed recently and some of them have demonstrated great empirical success in solving PDEs of hundreds and thousands of dimensions [7, 10, 17], much beyond the capability of conventional approaches. Question remains though on theoretical analysis of such neural-network based methods for solving high-dimensional PDEs. While there have been some recent progress on approaches including physics-informed neural networks [37, 31, 38] and the deep Ritz method [28, 27], many questions still remain open. Among them, a fundamental question is
Whether the solution of a high-dimensional PDE can be efficiently approximated by a neural network, and if so, how to quantify the complexity of the neural network representation with respect to the increasing dimension?
Our contributions The focus of the current study takes a functional-analytic approach to this question. Namely, we identify a function class suitable for neural network approximations and prove that the solutions to a class of PDEs can be well approximated by functions in this class. More specifically, the PDE we consider is a family of second-order elliptic PDEs of the form
Lu = −∇ · (A∇u) + cu = f on Rd. (1.1) We choose to work with the Barron class of functions defined in [8] (see also [1]), which is a class of functions admitting the integral of certain parametric ridge function against a probability measure on the parameters; see Definition 2.2 for a precise description. This Barron space is inspired by the pioneering work by Barron [2], where he proved that a class of functions whose Fourier transform has the first order moment can be approximated by two-layer networks without CoD. The main result of our work, stated informally, is the following; a more precise statement can be found in Section 2.3.
Main Theorem (informal version) If the coefficients A, c and the source term f of the secondorder elliptic PDE (1.1) are all Barron functions, then the solution u∗ can be approximated by another Barron function u such that ‖u − u∗‖H1 ≤ , where the Barron norm of u is upper bounded by O((d/ )C log(1/ )). Moreover, if the Barron space is defined by the cosine activation function, then the upper bound on the Barron norm can be improved to O(dC log 1/ ). We note that while the better rate is only obtained for the cosine activation function, such periodic activation function has indeed been found effective in certain PDE related tasks, see e.g., [42].
Since the Barron functions can be approximated on a finite domain Ω w.r.t. H1 norm by twolayer neural networks with a rate O(1/ √ k) where k is the network width (see Theorem 2.5), the theorem above directly implies that there exists a two-layer network uk with the number of widths k = O((d/ )C log(1/ )), or k = O(dC log 1/ ) if the activation function is cosine, such that ‖uk−u∗‖H1(Ω) ≤ . Therefore in our setting the solution can be approximated by a two-layer neural network without CoD, namely the complexity depends at most polynomially on the dimension d for fixed . Alternatively, we can rewrite the rates as O((1/ )C(log d+log 1/ )) and O((1/ )C log d) to contrast with that of conventional grid-based numerical methods for PDEs, which scales as O((1/ )d). We observe that the dependence on d is replaced with log d in the complexity bound for neural network approximations.
We emphasize that such approximation result does not follow directly from the universal approximation property of neural networks for Barron functions since it is not a priori known that the solution to
the PDE is a Barron function. In fact, directly imposing regularity or complexity assumption on the solution itself is unreasonable since the solution is unknown and its fine properties are generally inaccessible. Our main contribution is to establish the fact that the solution can be indeed approximated by a Barron function, under the assumption that coefficients and the right hand term of the PDE are Barron. From a mathematical point of view, our main theorem is in the same spirit as regularity estimates of PDEs, which are of crucial importance in the study of PDEs. While such regularity estimates are well developed in low dimension, the extension to results in high dimension is highly non-trivial and is the main focus of our work.
Related works Several theoretical work have been devoted to the above representation question. It has been established in [15, 20, 14] that deep neural networks can approximate solutions to certain class of parabolic equations and Poisson equation without CoD. The major limitation of those work lies in that the PDEs considered in those work must admit certain stochastic representation such as the Feymann-Kac formula and it seems difficult to generalize the proof techniques to broader classes of PDEs with no probabilistic interpretation. The work [28, 27] analyzed a priori generalization error of two-layer networks for solving elliptic PDEs and the Schrödinger eigenvalue problem on a bounded domain with Neumann boundary condition by assuming that the exact solutions lie in certain spectral Barron space, where the later was rigorously justified with a new regularity theory of the PDE solutions in the spectral Barron space. Similar generalization analysis was carried out in [29] for second-order PDEs and in [19] for general even-order elliptic PDEs, but without justifying the Barron assumption on the solution. Compared to those work, our work focuses on deriving complexity estimates of the solution in the integral-representation-based Barron space, which is more flexible and arguably more suitable for high-dimensional settings, see e.g., discussion in [8]. The work [9] established such estimates in the Barron space for certain specific PDEs that essentially admit explicit solution, whereas we aim to prove such estimates for general elliptic PDEs for which the analytical ansatz is not available. The work [30] is closest to ours where the authors proved that the solution of the same type of elliptic PDE with a Dirichlet boundary condition can be approximated by a (deep) neural networks with at most O(poly(d)N) parameters if the coefficients of the PDE are approximable by neural networks with at most N parameters. While our overall approach based on iterative scheme borrows idea from [30], our result differs and improves theirs in many aspects: (1) Our result shows that the solution can be well approximated without CoD by a two-layer neural network with a single activation whereas the result in [30] requires a deep network which uses a mixure of at least two activation functions; (2) Our PDE is set up on the whole space rather than a compact domain, so our setting covers some important PDEs in physics, such as the stationary Schrödinger equation; (3) The result in [30] relies on another key assumption that the source term lies within the span of finitely many eigenfunctions of the elliptic operator whereas our result completely removes such assumption. This is achieved by utilizing a novel preconditioning technique to uniformly control the condition number of the iterative scheme that underpins the proof of our main theorem.
Organization The rest of this paper will be organized as follows. In Section 2.1 we set up the PDE problem on the whole space and in Section 2.2 we introduce the definition of Barron functions and discuss their H1-approximation by two-layer networks (see Theorem 2.5). Our main theorems are stated in Section 2.3. We present the sketch proofs of the main theorems in Section 3 and defer the complete proof to Appendix. The paper is concluded with discussions on some future directions.

heading: 2 Problem setup and main results

heading: 2.1 Problem description
text: Notations Throughout this paper, we use ‖v‖ to denote the Euclidean norm of a vector v ∈ Rd. For a matrix A ∈ Rd×d, we denote its operator norm by ‖A‖ = supv∈Rd\{0} ‖Av‖ ‖v‖ . For R > 0, we denote by B d
R the closed ball in Rd centered at 0 with radius R, i.e., B d
R = {x ∈ Rd : ‖x‖ ≤ R}. Recall that we consider the d-dimensional second-order elliptic PDE (1.1). To guarantee the existence and uniqueness of the weak solution in H1(Rd), we make the following minimum assumptions on coefficientsA, c and right-hand side f ; this assumption will be strengthened in our main representation theorem.
Assumption 2.1. A(x) = (Aij(x))1≤i,j≤d is symmetric with ‖A(x)‖ ≤ amax <∞ and uniformly elliptic, that is for some amin > 0, it satisfies
ξ>A(x)ξ ≥ amin ‖ξ‖2 , ∀ x, ξ ∈ Rd. We also assume that 0 < cmin ≤ c(x) ≤ cmax <∞ and f ∈ L2(Rd).
Under Assumption 2.1, a standard argument using the Lax-Milgram theorem implies that there exists a unique weak solution u∗ ∈ H1(Rd), such that Lu∗ = f in H−1(Rd) which is the dual space of H1(Rd), i.e., ∫
Rd A∇u∗ · ∇vdx+ ∫ Rd cu∗vdx = ∫ Rd fvdx, ∀v ∈ H1(Rd).
Our ultimate goal is to show that the solution can be approximated by a two-layer neural network on any bounded subset of Rd with respect to the H1 norm with a rate scaling at most polynomially in the dimension. Notice that in general one cannot hope to obtain an approximation result on the whole space Rd because the asymptotic behavior of a neural network function (determined by the activation) at infinity may mismatch that of the target function u∗. On the other hand, it is well-known that the convergence rate of neural networks for approximating functions in standard Sobolev or Hölder spaces still suffers from the CoD [44, 45]. Therefore to obtain a rate without CoD for the neural networks approximation to the solution u∗, we need to argue that u∗ lies in a suitable smaller function space which has low complexity compared to Sobolev or Hölder spaces. We will work with the Barron space and show that u∗ is arbitrarily close to a Barron function which can be approximated by a two-layer neural network without CoD.

heading: 2.2 Barron spaces
text: The definition of Barron space is strongly motivated by the two-layer neural networks. Recall that a two-layer neural network with k hidden neurons is a function of the form
uk(x) = 1
k k∑ i=1 aiσ(w > i x+ bi), x ∈ Rd. (2.1)
Here σ : R → R is some activation function and (ai, wi, bi) ∈ R × Rd × R, i = 1, 2, . . . , k are the network parameters. If the parameters are randomly chosen accordingly to some probability distribution, then in the infinite width limit the averaged sum in (2.1) formally converges to the following probability integral
uρ(x) :=
∫ aσ(w>x+ b)ρ(da, dw, db), x ∈ Rd, (2.2)
where ρ is a probability measure on the parameter space R× Rd × R. Observe that (2.1) is a special instance of (2.2) if we take ρ(a,w, b) = 1k ∑k i=1 δ(a− ai, w − wi, b− bi).
The Barron norms and Barron spaces are then defined as follows, where we require the marginal measure in w to have compact support. This is because that the (formal) first-order and second-order partial derivatives of uρ(x) would involve with components of w by chain rule. By adding some uniform bounds on w, we can to control the Barron norms after taking derivatives. In the subsequent discussion, we may also need to restrict our attention on functions defined on a bounded set. Therefore we present below the formal definition of a Barron function defined any domain Ω ⊂ Rd. Definition 2.2. Fix Ω ⊂ Rd and R ∈ [0,+∞]. For a function g = uρ with some probability measure ρ, we define the Barron norm of g on Ω with index p ∈ [1,+∞] and support radius R by
‖g‖BpR(Ω) = infρ
{(∫ |a|pρ(da, dw, db) )1/p : g = ∫ aσ(w>x+ b)ρ(da, dw, db) on Ω,
ρ is supported on R×BdR × R } ,
where B d
R = {x ∈ Rd : ‖x‖ ≤ R}. The corresponding Barron space is then defined as BpR(Ω) = { g : ‖g‖BpR(Ω) <∞ } .
It is worth making some comments on the definition above. Our definition of Barron space adapts a similar definition in [8] (see also [1]) with several important modifications for the purpose of PDE analysis. First we require that the w-marginal of the probability measure ρ has compact support in order to control the derivatives of a Barron function defined in (2.2); in fact differentiating the integral of (2.2) leads to an integral of the product of the ridge function with w (or its powers) and enforcing ρ has a compact w-marginal thus controls the Barron norm of the derivatives of uρ. In addition, our definition of Barron norm only involves the p-th moment of ρ with respect to a parameter whereas the Barron norm in [8] takes the moments in all parameters into account. This is because [8] uses the unbounded ReLU activation function, which requires the moment condition in all parameters to make the integral in (2.2) well-defined; whereas we will only consider bounded σ (see Assumption 2.3) and the integral is guaranteed to be finite under such assumption.
Both our notion of Barron space and the one in [8] are motivated by the seminal work of Barron [2] where he proved that if the Fourier transform F(f) of a function f satisfies that∫
Rd |F(f)(ξ)||ξ|dξ <∞,
then there exists a two-layer network uk with k hidden neurons such that ‖f − uk‖L2(Ω) ≤ Ck− 1 2 . Since Barron’s original function class is defined via the Fourier transform, we call such function class the spectral Barron space to distinguish it from our Barron space based on the probability integral. We refer to [24, 3, 39, 40, 28] for recent developments on the spectral Barron space.
As we investigate the solution theory of the second-order PDE in the Barron space, we expect to differentiate the integral representation (2.2) up to the second order. Therefore, we assume that the activation function σ as well as its first-order and second-order derivatives are all bounded in R. Assumption 2.3. σ : R → R is smooth with C0 := supy∈R |σ(y)| < ∞, C1 := supy∈R |σ′(y)| < ∞, and supy∈R |σ′′(y)| <∞.
Thanks to the Hölder inequality, it is clear that BpR(Ω) ⊂ B q R(Ω) when p ≤ q. The following useful proposition (see also [8, Proposition 1]) shows that the reverse is also true and that the Barron norms and the Barron spaces are in fact independent of p. Proposition 2.4. For any function g ∈ B1R(Ω), it holds that ‖g‖B∞R (Ω) = ‖g‖BpR(Ω) = ‖g‖B1R(Ω) for any 1 ≤ p ≤ ∞. As a consequence, B∞R (Ω) = B p R(Ω) = B1R(Ω) for 1 ≤ p ≤ ∞.
The proof of Proposition 2.4 can be found in Appendix B.
The most important property that makes Barron functions distinct from Sobolev or Hölder functions is that they can be approximated by two-layer neural networks with a dimension-independent approximation rate in H1 norm as shown in Theorem 2.5. Theorem 2.5 (Approximation theorem inH1 norm). Suppose that Assumption 2.3 holds and that g ∈ B1R(Ω). Then for any open bounded subset Ω0 ⊂ Ω and any k ∈ N+, there exists {(ai, wi, bi)}ki=1 satisfying ∥∥∥∥∥1k k∑ i=1 aiσ(w > i x+ bi)− g(x) ∥∥∥∥∥ 2
H1(Ω0)
≤ 2(C20 +R 2C21 )m(Ω0) ‖g‖ 2 B1R(Ω)
k , (2.3)
where C0 and C1 are the constants in Assumption 2.3, and m(Ω0) is the Lebesgue measure of Ω0.
Theorem 2.5 provides an H1-approximation rate for Barron functions defined by the integral representation (2.2). The proof is deferred to Appendix B. Similar approximation results in the sense of L2 for Barron functions (including formulations based on spectrum and integral representation) have been proved in [2, 24, 3, 39, 8]. H1-approximation results for spectral Barron functions were previously obtained in [40] and [28].

heading: 2.3 Main theorems
text: To state our main theorems, we need to make some additional complexity assumption on the coefficients A, c and the source term f of the PDE (1.1), which is reasonable as otherwise there is no hope that the solution would lie in a smaller function class.
Assumption 2.6. For some RA, Rc, Rf ∈ (0,+∞), we have `A := max1≤i,j≤d ‖Aij‖B1RA (Rd) < ∞, `c := ‖c‖B1Rc (Rd) <∞, and `f := ‖f‖B1Rf (Rd) <∞.
We remark that Assumption 2.6 is compatible with our earlier Assumption 2.1 on the coefficients A, c and the source f . In fact, it is easy to see that constant coefficients A, c satisfy both assumptions if im(σ) 6= {0}, i.e., σ is not constantly zero. As for f , we provide in Proposition A.1 of Appendix A a concrete class of f that satisfies both assumptions.
We also need two additional technical assumptions on the activation function.
Assumption 2.7. The function h : R2 → R, (y1, y2) 7→ σ(y1)σ(y2) satisfies that `m := ‖h‖B1Rm (R2) <∞, for some Rm ∈ (0,+∞).
Assumption 2.8. It holds that `d,1 := ‖σ′‖B1Rd,1 (R) <∞ and `d,2 := ‖σ′′‖B1Rd,2 (R) <∞, for some Rd,1, Rd,2 ∈ (0,+∞).
Assumption 2.7 and Assumption 2.8 guarantee that Barron spaces are closed under multiplication and differentiations (up to the second order) respectively; see Lemma 3.3 (iii)-(iv) for a precise statement. These operations and the associated closeness will be useful for constructing approximation to the exact solution u∗ of the PDE (1.1) in Barron spaces. Proposition A.2 shows that Assumption 2.7 and Assumption 2.8 hold for a relatively large class of activation functions including cosine.
With the preparations above, we are ready to state our main theorems below. The first main theorem concerns the complexity estimate of the exact solution u∗ in the Barron space.
Theorem 2.9. Suppose that Assumption 2.1, 2.3, 2.6, 2.7, and 2.8 hold. For any ∈ (0, 1/2), there exists u ∈ B1R(Rd) withR ≤ γ1
( 1 )γ2 and ‖u‖B1R(Rd) ≤ β1 (d )β2|ln |, such that ‖u− u∗‖H1(Rd) ≤ . Here γ1, γ2, β1, and β2 only depend on ‖f‖H−1(Rd) and constants in Assumptions 2.1, 2.6, 2.7, and 2.8.
Furthermore, if σ = cos, then ‖u− u∗‖H1(Rd) ≤ can be achieved with R ≤ γ′1|ln | and ‖u‖B1R(Rd) ≤ β ′ 1d β′2|ln |, where γ′1, β ′ 1, and β ′ 2 only depend on ‖f‖H−1(Rd) and constants in Assumption 2.1 and 2.6.
Theorem 2.9 shows that the exact solution u∗ is -close (in the sense of H1) to a Barron function u ∈ B1R(Rd). In addition, the Barron norm of u grows at most polynomially in d, indicating that the complexity of u dose not suffer from the CoD. Also the complexity estimate gets substantially improved when the activation function is cosine. In fact, advantages of periodic activation functions have been empirically observed in some earlier works, see e.g., [42]. It remains an open question whether results similar to Theorem 2.9 can be established for activation functions that do not satisfy Assumption 2.7 and Assumption 2.8. This will be investigated in future works.
Thanks to Theorem 2.5 and Theorem 2.9, it is easy to conclude that the PDE solution u∗ can be approximated on any bounded subset Ω ⊂ Rd using two-layer neural networks with the number of hidden neurons k scaling at most polynomially in d.
Theorem 2.10. Under the same assumptions as in Theorem 2.9, given any ∈ (0, 1/2) and any open bounded subset Ω ⊂ Rd, there exists a two-layer neural network uk(x) with k ≤ γm(Ω)
( d )β|ln | such that ‖uk − u∗‖H1(Ω) ≤ , where γ and β only depend on ‖f‖H−1(Rd) and constants in Assumptions 2.1, 2.3, 2.6, 2.7, and 2.8.
Furthermore, if σ = cos, then ‖uk − u∗‖H1(Ω) ≤ can be achieved with k ≤ γ′m(Ω)dβ ′|ln |, where γ′ and β′ only depend on ‖f‖H−1(Rd) and constants in Assumptions 2.1, 2.3, and 2.6.

heading: 3 Proofs of the main results
text: We sketch the proof ideas in this section and present the full details in the Appendix.

heading: 3.1 Preconditioned functional iterative scheme
text: The key ingredient of our proof of Theorem 2.9 is a functional iterative scheme for solving the elliptic PDE, which can be viewed as an infinite dimensional analog of the preconditioned steepest descent algorithm to solve linear algebra equations. Recall when solving the linear equation Ax = b with A ∈ Rn×n and x, b ∈ Rn, the preconditioned steepest descent algorithm [13] runs the iteration
xt+1 = xt − αP (Axt − b),
where P is a preconditioning matrix, α is the step size, and t = 0, 1, 2, · · · indicates the iteration index. The purpose of the preconditioned iteration is to reduce the condition number of the iteration κ(PA) by choosing a suitable P and hence accelerate the convergence of the iterative algorithm.
In the case of solving the elliptic PDE (1.1), we generalize the preconditioned steepest descent iteration to the functional setting by considering the following iteration scheme in H1(Rd):
ut+1 = ut − α(I −∆)−1(Lut − f), (3.1)
where the inverse operator (I − ∆)−1 plays the role of preconditioner. As a matter of fact, we will show that the condition number of (I −∆)−1L is bounded and this directly implies that the iterative scheme (3.1) converges exponentially to the exact solution u∗. Indeed, we have the following contraction estimate for the iteration (3.1), whose proof can be found in Appendix C.
Proposition 3.1. Recall the constants amin, amax, cmin, cmax defined in Assumption 2.1. For any α > 0 and any u ∈ H1(Rd),∥∥(I − α(I −∆)−1L)u∥∥
H1(Rd) ≤ Λ(α)‖u‖H1(Rd), (3.2)
where the contraction factor Λ(α) = supλ∈[λmin,λmax] |1− αλ| with λmin = min{amin, cmin} and λmax = max{amax, cmax}.
In particular, minimizing Λ(α) with respect to the step size α yields an optimal choice of step size
α∗ := 2
λmin + λmax .
With α = α∗ in (3.2), we obtain that∥∥∥∥(I − 2λmin + λmax (I −∆)−1L ) u ∥∥∥∥ H1(Rd) ≤ λmax − λmin λmax + λmin ‖u‖H1(Rd). (3.3)
As a direct consequence, we obtain the following estimate for the number of iterations required to achieve a given error tolerance.
Corollary 3.2. Let u∗ be the exact solution of the PDE (1.1). Under Assumption 2.1, consider the iteration scheme (3.1) with α = α∗ = 2λmin+λmax . Then for any
T ≥ (
ln λmax + λmin λmax − λmin
)−1 ln ‖u0 − u∗‖H1(Rd) ,
the iterate uT satisfies ‖uT − u∗‖H1(Rn) ≤ .
Let us remark that the idea of using iterative scheme to establish neural network representation results of solutions to PDEs is not new, see e.g., [23, 30], similar ideas have been also used to construct neural network architectures inspired from iterative schemes, see e.g., [43, 12]. Closely related to our setting, the work [30] uses a steepest descent iteration with the right hand side of the equation assumed to be in the span of first several eigenfunctions of the elliptic operator, while [23] considered general right hand side, but only after discretization which also effectively truncates the problem onto a finite dimensional subspace. These restrictions were made to limit the condition number of the iteration. Unlike those works using standard steepest descent iterations, by using the preconditioning technique, we can deal with general right hand side without restricting to a finite-dimensional subspace.

heading: 3.2 Algebra of Barron functions and representation of the solution
text: Corollary 3.2 in the previous subsection shows that we can obtain an approximate solution by running the iteration (3.1). To complete the proof of Theorem 2.9, we show in this subsection that the iteration (3.1) can be carried out in the Barron space B1R(Rd), i.e. each iteration ut ∈ B1R(Rd) (with the support radius R potentially depending on t). To this end, we first need to establish the closeness of Barron space under function operations involved in the iteration. In fact, by decomposing each of the iteration step in (3.1) into two steps, we can write{
vt = Lut − f = − ∑ i,j(∂iAij∂jut +Aij∂ijut) + cut − f, ut+1 = ut − α(I −∆)−1vt. (3.4)
Thus, to show that the iterate ut remains in Barron space, it suffices to establish that addition, scalar multiplication, product, differentiation, and action of (I − ∆)−1 are closed in the Barron space. The closedness of Barron functions under those operations are not only useful for proving our main results, but also of its own interest. The next two lemmas summarize the algebras and the stability estimate of the inverse (I −∆)−1 in the Barron space. Their proofs can be found in Appendix D. Lemma 3.3 (Algebras in Barron spaces). The followings hold:
(i) (Addition) Suppose that ‖gi‖B1Ri (Rd) <∞, i = 1, 2, . . . , k. Then ‖g1 + · · ·+ gk‖B1R(Rd) ≤∑
1≤i≤k ‖gi‖B1Ri (Rd) , where R = max1≤i≤k Ri.
(ii) (Scalar multiplication) Suppose that ‖g‖B1R(Rd) <∞ and that λ ∈ R. Then ‖λg‖B1R(Rd) = |λ| ‖g‖B1R(Rd).
(iii) (Product) Suppose that Assumption 2.3 and Assumption 2.7 hold and that ‖gi‖B1Ri (Rd) <∞
for i = 1, 2. Then ‖g1g2‖B1R(Rd) ≤ `m ‖g‖B1R1 (Rd) ‖g‖B1R2 (Rd) , where R = Rm(R1 +R2) with Rm and `m being constants in Assumption 2.7.
(iv) (Derivatives) Suppose that Assumption 2.3 and Assumption 2.8 hold and that ‖g‖B1R(Rd) < ∞ with R < ∞. Then ‖∂ig‖B1Rd,1R(Rd) ≤ `d,1R ‖g‖B1R(Rd) and ‖∂ijg‖B1Rd,2R(Rd) ≤
`d,2R 2 ‖g‖B1R(Rd) for any i, j ∈ {1, 2, . . . , d}, whereRd,1,Rd,2, `d,1, and `d,2 are constants
in Assumption 2.8. Lemma 3.4 (Applying (I − ∆)−1 on Barron functions). Suppose that ‖g‖B1R(Rd) < ∞. Then∥∥(I −∆)−1g∥∥B1R(Rd) ≤ ‖g‖B1R(Rd). we include a proof for Lemma 3.4 in Appendix D that uses similar arguments from [9], though the analysis in [9] is for d ≥ 3. The lemmas above lead to the following recursive estimate on the Barron norm of ut. Lemma 3.5. Suppose that Assumption 2.3, Assumption 2.7, and Assumption 2.8 hold. If ‖u‖B1Ru,t < ∞ with Ru,t <∞, then ut+1 defined in (3.1) or (3.4) satisfies that ‖ut+1‖B1Ru,t+1 (Rd) ≤ ( α`m`A(` 2 d,1RARu,t + `d,2R 2 u,t)d 2 + α`m`c + 1 ) ‖ut‖B1Ru,t (Rd) + α`f , (3.5) for any
Ru,t+1 ≥ max{RmRd,1(Ru,t +RA), Rm(Rd,2Ru,t +RA), Rm(Ru,t +Rc), Ru,t, Rf}. (3.6)
The proof of Lemma 3.5 is deferred to Appendix D. One observation is that the amplification factor of the Barron norm in Lemma 3.5 increases as the support radius R increases. The reason is that differentiating the function would introduce components of w and hence the amplification depends on how large ‖w‖ can be and thus the support of the measure. One possible direction to improve the estimate is to realize that the preconditioner (I −∆)−1 can counteract the action of taking derivatives. It is indeed possible to to remove the R dependence from the amplification factor, at least for some specific activation functions, through a more careful analysis. In particular, we have the following lemma for the cosine activation function, the proof of which can also be found in Appendix D.
Lemma 3.6. Suppose that Assumption 2.6 holds. If σ = cos and ‖u‖B1Ru,t (Rd) <∞ with Ru,t <∞, then ut+1 defined in (3.1) or (3.4) satisfies
‖ut+1‖B1Rt+1 (Rd) ≤ ( 6α`A max{R2A, 1}d2 + α`c + 1 ) ‖ut‖B1Ru,t (Rd) + α`f , (3.7)
for any Ru,t+1 ≥ Ru,t + max{RA, Rc, Rf}. (3.8)
Lemma 3.5 and Lemma 3.6 estimate the amplification of the Barron norm in each iteration of (3.1). Combining them with the control of number of iterations, Corollary 3.2, we are ready to finish the proof of Theorem 2.9.
Proof of Theorem 2.9. Fix u0 = 0 and α = 2λmin+λmax . According to Corollary 3.2, it holds that ‖uT − u∗‖H1(Rn) ≤ for any
T ≥ (
ln λmax + λmin λmax − λmin
)−1 ln ‖u∗‖H1(Rn) .
Moreover, thanks to the estimate λmin ‖u∗‖2H1(Rd) ≤ ∫ A∇u∗ · ∇u∗dx+ ∫ c|u∗|2dx = ∫ fu∗dx ≤ ‖f‖H−1(Rd) ‖u ∗‖H1(Rd) ,
we have ‖u∗‖H1(Rd) ≤ 1 λmin ‖f‖H−1(Rd). Therefore, it suffices to take
T = ⌈( ln λmax + λmin λmax − λmin )−1 ln 1 + ( ln λmax + λmin λmax − λmin )−1 ln ‖f‖H−1(Rd) λmin ⌉ .
Set Ru,0 = max{RA, Rc, Rf , 1} and Ru,t+1 = max{2RmRd,1, 2RmRd,2, 2Rm, 1} ·Ru,t ≥ Ru,t. Then (3.6) is satisfied for any t. Let us define a sequence {Xt}t≥0 via X0 = 1 and Xt+1 =( α`m`A(` 2 d,1 + `d,2) + α(`m`c+`f )+1 d2 ) R2u,td
2 ·Xt. By (3.5), we have ‖ut‖B1Ru,t (Rd) ≤ Xt for any t. Therefore, it holds that
Ru,T = max{RA, Rc, Rf , 1} ·max{2RmRd,1, 2RmRd,2, 2Rm, 1}T , and that
‖uT ‖B1Ru,T (Rd) ≤ XT
= ( α`m`A(` 2 d,1 + `d,2) + α(`m`c + `f ) + 1
d2
)T d2T (Ru,0 · · ·Ru,T−1)2
≤ ( α`m`A(` 2 d,1 + `d,2) + α(`m`c + `f ) + 1
d2
)T d2T
· (max{RA, Rc, Rf , 1})T ·max{2RmRd,1, 2RmRd,2, 2Rm, 1}T 2 .
The first part of Theorem 2.9 is established by setting u = uT and R = Ru,T .
If σ = cos, (3.8) is satisfied by setting Ru,t = max{RA, Rc, Rf} · t.
Define Y0 = 0 and Yt+1 = ( 6α`A max{R2A, 1}d2 + α`c + 1 ) Yt + α`f . By (3.7), we obtain that ‖ut‖B1Ru,t (Rd) ≤ Yt for any t, and in particular that
‖uT ‖B1Ru,T (Rd) ≤ YT =
α`f (( 6α`A max{R2A, 1}d2 + α`c + 1 )T − 1) 6α`A max{R2A, 1}d2 + α`c ,
which finishes the proof by setting u = uT and R = Ru,T .
Theorem 2.10 is then a corollary of Theorem 2.9 and Theorem 2.5 (the approximation theorem).
Proof of Theorem 2.10. Theorem 2.10 follows directly from applying Theorem 2.5 with error tolerance /2 and applying Theorem 2.9 with error tolerance /2.

heading: 4 Conclusion
text: In this work, we establish the approximation rate for the solution of a second-order elliptic PDE by a Barron function and by a two-layer neural network. Under the assumption that the coefficients and the source of the PDE are all in the Barron spaces with some compact support property on the underlying probability measure, the approximation rate is shown to depend at most polynomially on the dimension. Therefore, our results indicate that even a neural network as simple as a two-layer network with a single activation function can have adequate representation ability to encode the solution of an elliptic PDE, without incurring the CoD. Our result provides theoretical guarantee for numerical methods for solving high-dimensional PDEs using neural networks.
For future directions, it is of interest to extend the functional analysis framework to more general activation functions (such as unbounded ones) and more general neural network architectures. One interesting direction is to establish depth separation result for representing PDE solutions. Our analysis also indicates some potential benefit of using periodic activation function such as cosine in terms of approximation, further studies and understanding of the choice of activation function and architecture are crucial. Moreover, while we focus on approximation error, generalization error and analysis of training should also be considered in future works.
It is possible to extend the approximation results to a wider range of high-dimensional PDEs such as parabolic PDEs, PDE eigenvalue problems, and nonlinear equations such as those arise from control theory. The analysis tools and characterization of Barron space we establish in this work would be useful for these future studies.
Acknowledgments and Disclosure of Funding
The work of Z.C. and J.L. is supported in part by the National Science Foundation via grants DMS2012286 and CCF-1934964. Y.L. thanks the National Science foundation for its support through the award DMS-2107934.

Review summaries:
The objective of this paper is to prove that the solution $u$ of a certain class of PDEs, of the form  \begin{equation}\label{eq:PDE} - \nabla\cdot (A \cdot \nabla u) + cu = f \mbox{ on } \mathbb{R}^d   \end{equation} for some matrices $A$, coefficient $c$ and right-hand side $f$ can be well-approximated by functions belonging to some Barron spaces, provided that $A$, $c$ and $f$ belong to some Barron spaces as well. The precise approximation results are stated  in Theorems~2.9 and~2.10 of the manuscript.   Such kinds of approximation result are very interesting since they can give some insight on how well the solutions of some high-dimensional PDEs can be approximated by neural networks-based numerical methods. Indeed,  it is proved in earlier works that elements of these Barron spaces can be approximated in the $H^1$ norm by two-layer neural networks.    The article is organized as follows: notation and assumptions on the data of the problem are introduced in Section~2.1. The definition of Barron spaces is recalled in Section~2.2, together with approximation results in the $H^1$ norm of elements of Barron spaces  by two-layer neural networks (Theorem~2.5). The proofs of these results is postponed in Appendix~B. The main results of the paper concerning the approximation of the solution of PDE (\ref{eq:PDE}) by elements of a Barron space are stated in Theorem~2.9 and Theorem~2.10.  The main idea of the proofs of these results is sketched in Section~3, together with the statement of useful auxiliary lemmas. The detailed proofs of these results are given in the Appendix. 

This paper studies the neural network approximation of high-dimensional elliptic PDEs. Specifically, it shows that when the source term and the coefficient functions lie in a Barron space, the solution can be approximated by two-layer neural networks without the curse of dimensionality.  The Barron norm characterizes the approximability of functions by two-layer neural networks.  The primary technique is to represent the PDE solution as the fixed point of the preconditioned steepest descent iteration.  Then,  it shows that the iteration satisfies two properties: (1) The iteration converges to the fixed point exponentially fast; (2) Each iteration only amplifies the Barron norm by a small fact.  To guarantee the first property, the PDE must be strongly elliptic. To guarantee the second property, the Barron space must be closed under multiplication, which excludes the most commonly-used activation functions. 

This paper studies the representation of solutions of second order elliptic PDEs using two layer neural networks, and gives sufficient conditions so that this can be done without curse of dimension (the number of units grows at  most polynomially in the dimension).  Precise theorems to this effect are stated and proven.  The mathematical framework includes the definition of a Barron space of functions (essentially, finite norm of the final layer weights), many structural results about this function space, and an iterative solution scheme with its analysis.