Refined Generalization Analysis of the Deep Ritz Method and Physics-Informed Neural Networks

Thank you for taking the time to review our article and for your valuable suggestions. Let us address your questions point by point.

Q1: Differences with the Weinan's section 2 (i.e., ref [1]) and CMAME section 3 (i.e., ref [2]).

A1: We have cited Weinan's work in line 238 and explained that our definition of the Barron space differs. Nevertheless, we acknowledge the lack of a detailed comparison with Weinan's work, i.e., [1]. Section 2 of [1] (and [3]) only presents approximation results under the $L^2$ norm, with a slow rate of $m^{- \frac{1}{2}}$, where $m$ is the width of a two-layer neural network. Although our definition of the Barron space differs from that in [1], our method can be naturally applied to [1], yielding approximation results under the $H^1$ norm with a faster rate, that is, $m^{-(\frac{1}{2}+\frac{1}{3d})}$. Regarding why we are able to achieve better approximation results, [1] employed the Monte Carlo method, which can only yield a slow approximation rate of $m^{- \frac{1}{2}}$. In contrast, we utilized Lemma 10. After applying Lemma 10, the challenge lies in estimating the metric entropy of the target function class in the $H^1$ norm, which we addressed using a new decoupling technique. Moreover, in addition to approximation error, we have also provided generalization error, while [1] only focuses on approximation error.

As for comparing with Section 3 of [2], [2] does not provide any theoretical analysis regarding approximation results, focusing only on the experimental part. It may be interesting to consider the theoretical analysis of the Deep Double Ritz Method presented in [2].

In summary, for the Barron space, we have provided a faster approximation rate in a higher-order Sobolev norm, and the associated constants are independent of the dimension. Moreover, our method can also be naturally extended to other Barron spaces defined differently.

Q2: Novelty of the method for deriving the fast approximation rates.

A2: [1] employed the Monte Carlo method to derive the approximation rates, hence could only obtain a slow rate of $m^{- \frac{1}{2}}$ in the $L^2$ norm. In contrast, we utilized Lemma 10 (line 2141-2150). As mentioned in A1, the difficulty lies in the estimation of the metric entropy under the $H^1$ norm of target function class, since it is not Lipschitz continuous with respect to the parameters. To address this challenge, we employed a novel decoupling technique, which involves splitting the estimation of the metric entropy in the $H^1$ norm into the independent estimation of the metric entropy in the $L^2$ norm. Specifically in the estimation of the metric entropy under the $H^1$ norm, for $(\omega _1,t _1),(\omega _2,t _2)\in \partial B _1^d(1)\times [-1,1]$, we have $|| \sigma(\omega _1 \cdot x+t _1)-\sigma(\omega _2 \cdot x+t _2)|| _{H^1(\Omega)}^2 \leq 4(|\omega _1-\omega _2| _1^2+|t _1-t _2|^2)+2\int _{\Omega} |I _{ \{\omega _1 \cdot x+t _1 \geq 0\}}- I _{ \{\omega _2 \cdot x+t _2 \geq 0\} }|^2dx .$ It is challenging to handle the first and second terms simultaneously due to the discontinuity of indicator functions, thus we turn to handle two terms separately. Note that the first term is related to the covering of $\partial B _1^d(1)\times [-1,1]$ and the second term is related to the covering of a VC-class of functions. In this way, by dealing with these two terms separately, we can address this difficulty. For more details about the proof, one can refer to Proposition 8 starting from line 928. Moreover, our method can be naturally applied to [1], because the functions defined in the Barron space in [1] have a direct integral representation, which allows Lemma 10 from our paper to be directly applied.

Q3: Comparison with [4].

A3: In [4], the authors view solving different PDEs as different tasks and consider solving a single PDE as one task. However, in our article, when solving a single PDE, we consider minimizing the internal residual and the boundary residual as different tasks, thus our criteria for classifying different tasks differ from those in [4]. Furthermore, [4] does not provide theoretical analysis, while we offer detailed error analysis.

References:

[1]: W. E, C. Ma, and L. Wu. The Barron space and the flow-induced function spaces for neural network models. Constructive Approximation, 55(1):369–406, 2022.

[2]: C. Uriarte, D. Pardo, I. Muga, J. Munoz-Matute, A deep double ritz method (d2rm) for solving partial differential equations using neural networks, Computer Methods in Applied Mechanics and Engineering 405 (2023) 115892.

[3]: Chen, Z., Lu, J., and Lu, Y. On the representation of solutions to elliptic PDEs in Barron spaces. Advances in Neural Information Processing Systems, 34, 2021

[4]: Zongren Zou and George Em Karniadakis. L-HYDRA: Multi-head physics-informed neural networks. arXiv preprint arXiv:2301.02152, 2023.

Please check if the authors' response addresses your concerns.

Thank you for taking the time to review our article and for your insightful comments and suggestions. Let us address your questions point by point.

Q1: Comparison with [2] on the Deep Ritz Method aspect.

A1: First, the results of [2] are unrealistic, meaning there is still a gap between theory and practice. Specifically, the strong convexity of the DRM objective in Proposition 2.1 of [2] shows that for any $u\in H _0^1(\Omega)$, $\mathcal{E}^{DRM}(u)-\mathcal{E}^{DRM}(u^{* })\lesssim ||u-u^{* }|| _{H^1}^2 \lesssim \mathcal{E}^{DRM}(u)-\mathcal{E}^{DRM}(u^{* }),$ where $u^{* }\in H _0^1(\Omega)$ is the true solution and $\mathcal{E}^{DRM}$ is the Ritz functional. However, this only holds for the functions in $H _0^1(\Omega)$ and the neural network functions fail to belong to $H_0^1(\Omega)$. If we force the neural network functions to belong to $H _0^1(\Omega)$, for example, by multiplying them with a function that vanishes on the boundary, then the approximation results in [2] no longer hold. In this article, we consider the zero Neumann boundary conditions instead of the zero Dirichlet boundary conditions in [2]. The assumption that $\int _{\Omega}fdx=0$ is made to ensure the existence of the solution under the zero Neumann boundary condition. In fact, when dealing with zero Dirichlet boundary conditions, we would prefer to add a boundary penalty term, that is, to consider the following loss function: $\mathcal{E}^{DRM}(u)=\int _{\Omega} ||\nabla u|| _2^2+V|u|^2dx-2\int _{\Omega}fudx+\lambda \int _{\partial \Omega} u^2ds.$ Then combining this with our Section D.2 (line 2845-2871), which discusses the non-zero Neumann boundary condition and the results in [3], we may achieve a better generalization bound that is more aligned with reality under the zero Dirichlet boundary condition.

Second, the method in [2] for deriving fast rates does not work for the setting of Poisson equation, since the expectation of empirical loss is not equal to the variational formulation. This fact also limits the use of some popular methods, such as the local Rademacher complexity. As for why the method in [2] and the method in the original paper of local Rademacher complexity (i.e., [5]) are not feasible, for example, the core lemma in [2], i.e., Lemma B.4, requires the condition (from [5]) that $R _n(\{f\in \mathcal{F}| \mathbb{E}[f]\leq r\}) \leq \phi(r),$ where $R _n$ is the empirical Rademacher complexity and $\phi$ is a sub-root function. In the subsequent proof of [2], an appropriate class of functions class $\mathcal{F}$ can be chosen such that $\mathbb{E}[f] = \mathcal{E}^{DRM}(u)-\mathcal{E}^{DRM}(u^{* })$ for some $u$. By doing so, the strong convexity of the Ritz functional can be utilized. However, in the setting of Poisson equation, the Ritz function cannot be expressed in the form of an expectation, i.e., there does not exist a function class $\mathcal{F}$ such that $\mathbb{E}[f] = \mathcal{E}^{DRM}(u)-\mathcal{E}^{DRM}(u^{* })$. Therefore, in this work, we use a new error decomposition method and new peeling method to derive the fast rates.

Q2: Contribution of this paper to the DRM

A2: We apologize for the error in Remark 2. In fact, it should be "the convergence rate associated with the DRM improves from $n^{-\frac{2k-2}{d+4k-4}}$ to $n^{-\frac{2k-2}{d+2k-2}}$." We have corrected this in the revised version, and thank you for pointing out this error. Regarding the contributions of this paper to DRM, in Answer 1 we have already explained the difficulties in dealing with the Poisson equation and our new method. Although this method also works for the static Schrödinger equation, it is too complicated. On the other hand, the case of the static Schrödinger equation represents a large class of important problems with strongly convex properties like equation (6), such as the $L^2$ regression problem under bounded noise in [1]. In fact, the methods used for the Poisson equation in this article, as well as those in [2] and [1], are applicable to the class of problems mentioned above. However, these methods are too complicated. In this article, by combining a new error decomposition with the use of local Rademacher complexity, we provide a simpler method for solving this class of problems. Finally, we also discuss other boundary conditions for the DRM (see Section D.2, line 2845-2871) and show that our methods can be extended to a broader class of PDEs. Moreover, in the Discussion section starting on line 2873, we investigate the complexity of over-parameterized two-layer neural networks when approximating functions in Barron space and demonstrate meaningful generalization errors in the setting of over-parameterization.

Please check if the authors' response addresses your concerns.

The reviewer would like to thank the authors for their detailed response. However, it seems to the reviewer that current organization of the whole manuscript still seems to be a bit disorganized and disconnected, especially for the "Related Works" section. This makes it hard for the readers to tell the main contribution and focus of this paper. Therefore, the reviewer is inclined to maintain the score of 5 for now.

Thank you for taking the time to review our article. Let us address your questions point by point.

Q1: What is the intention behind studying the two methods using different mathematical tools for different PDEs?

A1: We consider these two methods, DRM and PINNs, because they are both very popular for solving PDEs using neural networks. The reason we used different approaches for the generalization bounds of these methods is due to the nature of each method. Section 2 deals with DRM, which involves only one task and has strong convexity, while Section 3 deals with PINNs, which involves two tasks and lacks strong convexity. Of course, the methods for handling DRM and PINNs can also be combined, as demonstrated in Section D.2, starting on line 2844, where we considered non-zero Neumann conditions.

Q2: The assumption that the solutions of the PDEs lie within the Barron space is too strong.

A2: In this article, we consider not only the case where the solutions of the PDEs lie in the Barron space but also the case where the solutions belong to the general Sobolev space. The consideration of the Barron space scenario is aimed at demonstrating that the Deep Ritz Method (DRM) and Physics-Informed Neural Networks (PINNs) do not suffer from the curse of dimensionality in this context. That is, the constants in the generalization bounds are at most polynomially related to the dimension, and the convergence rate is almost independent of the dimension. In fact, higher-order Sobolev spaces can be embedded in Barron spaces, for example, $H^{s+\frac{d}{2}+\epsilon}(\Omega) \subset \mathcal{B}^{s}(\Omega)$ for any $s>0, \epsilon>0$ and for any bounded domain $\Omega\subset \mathbb{R}^d$. Additionally, [1] and [2] have considered the solution theory of certain PDEs in Barron spaces, demonstrating when the solutions of certain PDEs are in the Barron space.

Q3: Assumptions are restrictive when deriving the generalization bounds under certain Sobolev norms.

A3: For the DRM, we only require boundedness and strong convexity. The assumption of boundedness is important in the field of machine learning for deriving better generalization bounds. Regarding the assumption of strong convexity, the conditions required by the traditional Ritz method can lead to strongly convex structures similar to equations (4) and (6). Specifically, let $H$ be a real Hilbert space with norm denoted by $||\cdot|| _{H}$ and let $\langle \cdot, \cdot \rangle$ denote the inner product. Let $a(u, v)$ be a real bilinear form on $H\times H$, $f$ be an element of $H$ and consider the following problem: find $u\in H$ such that $a(u, v) = \langle f, v\rangle, \ \forall v \in H.$ When the conditions of the Lax-Milgram theorem are met, that is $|a(u,v)|\leq M||u|| _{H}||v|| _{H}$ for any $u,v \in H$ and $a(v,v)\geq \alpha ||v|| _{H}^2$, and $a$ is symmetric, we have that the solution $u^{* }$ of the above variational problem is also the only element of $H$ that minimizes the following quadratic functional (also called energy functional) on $H$: $J(u)= a(u,u)-2\langle f, u\rangle.$ At this point, we naturally have for any $u\in H$, $||u-u^{* }|| _{H}^2 \lesssim J(u)-J(u^{* }) \lesssim ||u-u^{* }|| _{H}^2,$ since $J(u)-J(u^{* })=a(u-u^{* }, u-u^{* })$, this is exactly the strongly convex structure we need. For the PINNs, equation (30) only requires the strict elliptic condition and smoothness of the boundary. In fact, these assumptions are not restrictive.

Q4: Limitation of RePU activation functions.

A4: Thank you for highlighting the limitations of using the RePU activation function in our theoretical analysis. We acknowledge that in deeper networks, the power of RePU increases with depth, which can lead to floating-point precision issues. We will point out this limitation in the paper.

Q5: This article considers PDEs that are too simple and traditional methods can solve them faster and more accurately.

A5: Classical methods may suffer from the curse of dimensionality when solving high-dimensional PDEs. In contrast, part of our results are meaningful in high dimensions as well. As we mentioned in answer 3, the results regarding DRM in the article can be generalized to a broader class of PDEs. The strong convexity required can naturally be derived from the conditions of the Lax-Milgram theorem, which are also needed by the traditional Ritz method. Similarly, for PINNs, our conclusions do not require strong conditions, so they can be naturally extended to other types of PDEs.

Q6: The Poincaré constant is dimension-dependent.

A6: The Poincaré constant of the domain $[0, 1]^d$ is dimension-independent, since $[0,1]^d$ is a space of tensor-product form. Please see Theorem 4.2 in [3] or the content in line 1223-1232 of this paper for more details.

I appreciate the authors' efforts in addressing the concerns I raised. However, I am still unclear about the main focus and contributions of the paper. The discussions on DRM, PINNs, and the expressivity of neural networks seem somewhat disconnected, making it difficult to discern the central message of the paper.

While mathematical analysis of deep learning models is important, I would like to remind you that the majority of readers of the ICLR proceedings may not be mathematically inclined. With this in mind, I remain uncertain about the key takeaway of this study for the broader ICLR audience. DRM and PINNs have been explored in previous works, so it would be helpful to clarify what new insights and contributions your work offers over existing research, particularly for the broader deep learning community. I am also unclear about how your mathematical contributions provide a useful tool compared to existing research.

You also mention that the model addresses the curse of dimensionality. However, I believe this claim may be premature. Deep learning models involve complex interdependencies between network expressiveness, optimization, and generalization. It is not entirely clear that generalization bounds alone can be used to definitively argue that your model resolves the curse of dimensionality. For instance, the convergence of gradient descent in PINNs is known to be exponentially dependent on the dimension (C.H.Song at el., How does PDE order affect the convergence of PINNs?, NeurIPS, 2024.)

I believe it would be helpful to clearly position your work within the broader context of deep learning-based PDE solvers and to clarify its implications for general advancements in deep learning research. As of now, I still have some reservations regarding the paper's contribution to the field and its potential impact on the broader community.

Please check if the authors' response addresses your concerns.

Thank you for taking the time to review our article and for your recognition of it. Let us address your questions point by point.

Q1: Extension to other PDEs for the Deep Ritz Method (DRM).

A1: We first explain why we consider these two simple PDEs, and then illustrate how to generalize our results to more general PDEs. When considering the Poisson equation, the situation is different from most cases in the field of machine learning, because the expectation of empirical loss is not equal to the Ritz functional. Therefore, we need to employ new peeling techniques to obtain better generalization bounds. As for the static Schrödinger equation, the scenario represents a large class of problems with a similar strongly convex structure like equation (6) (line 199), and for such problems, we provide a simple and effective theoretical framework to achieve better generalization bounds.

Then we will demonstrate how our theory can be extended to a broader range of PDEs. In the discussion section starting on line 2844, we have illustrated how to extend our results to non-zero Neumann boundary condition. And for other boundary conditions, our results can also be similarly extended. As for extending to other PDEs, note that our method relies on two essential conditions, one is the boundedness condition and another is the strongly convex structure (see line 188 or line 199). For strongly convex structure, it is not only the Poisson equation and the static Schrödinger equation that meet the conditions. In fact, even the traditional Ritz method has similar requirements, which stem from the need to satisfy the Lax-Milgram theorem. Specifically, let $H$ be a real Hilbert space with norm denoted by $||\cdot|| _{H}$ and let $\langle \cdot, \cdot \rangle$ denote the inner product. Let $a(u, v)$ be a real bilinear form on $H\times H$, $f$ be a element of $H$ and consider the following problem: find $u\in H$ such that $a(u, v) = \langle f, v\rangle, \ \forall v \in H.$ When the conditions of the Lax-Milgram theorem are met, that is $|a(u,v)|\leq M||u|| _{H}||v|| _{H}$ for any $u,v \in H$ and $a(v,v)\geq \alpha ||v|| _{H}^2$, and $a$ is symmetric, we have that the solution $u^{* }$ of the above variational problem is also the only element of $H$ that minimizes the following quadratic functional (also called energy functional) on $H$: $J(u)= a(u,u)-2\langle f, u\rangle.$ Hence, for any $u\in H$, we naturally have $||u-u^{* }|| _{H}^2 \lesssim J(u)-J(u^{* }) \lesssim ||u-u^{* }|| _{H}^2,$ since $J(u)-J(u^{* })=a(u-u^{* }, u-u^{* })$, and this is exactly the strongly convex structure we need. As a simple example, we can replace the Laplacian operator in the static Schrödinger equation with operator $A$, satisfying $Au(x)= \sum\limits _{i,j=1}^{d}\frac{\partial}{\partial x _j} (a _{ij}\frac{\partial u}{\partial x _i}(x) ),$ where $a _{i,j}\in L^{\infty}(\Omega)$, $a _{ij}=a _{ji}$, the matrix ($a _{i,j}$) is uniformly (a.e.) positive definite.

In summary, the assumption of strong convexity is not overly restrictive, as even the traditional Ritz method requires it.

Q2: No numerical experiments.

A2: Apologies, this article primarily focuses on the theoretical aspects; in the future, we plan to design efficient algorithms to validate our conclusions.

Please check if the authors' response addresses your concerns.