The Challenges of the Nonlinear Regime for Physics-Informed Neural Networks
Authors: Andrea Bonfanti
Giuseppe Bruno
Cristina Cipriani
Venue: NeurIPS 2024 poster (2024)
Keywords: Physics-Informed Neural Networks
Neural Tangent Kernel
Nonlinear PDEs
Second-order optimization
Primary area: machine_learning_for_physical_sciences
TL;DR: A Neural Tangent Kernel analysis of the training dynamics of Physics-Informed Neural Networks for nonlinear PDEs.

Abstract: The Neural Tangent Kernel (NTK) viewpoint is widely employed to analyze the training dynamics of overparameterized Physics-Informed Neural Networks (PINNs). However, unlike the case of linear Partial Differential Equations (PDEs), we show how the NTK perspective falls short in the nonlinear scenario. Specifically, we establish that the NTK yields a random matrix at initialization that is not constant during training, contrary to conventional belief. Another significant difference from the linear regime is that, even in the idealistic infinite-width limit, the Hessian does not vanish and hence it cannot be disregarded during training. This motivates the adoption of second-order optimization methods. We explore the convergence guarantees of such methods in both linear and nonlinear cases, addressing challenges such as spectral bias and slow convergence. Every theoretical result is supported by numerical examples with both linear and nonlinear PDEs, and we highlight the benefits of second-order methods in benchmark test cases.

The Challenges of the Nonlinear Regime for Physics-Informed Neural Networks
Andrea Bonfanti
Giuseppe Bruno BMW AG, Digital Campus Munich Giuseppe.. GB Bruno@bmw.
BMW AG, Digital Campus Munich Basque Center for Applied Mathematics
University of the Basque Country abonfanti001 1ee00000kaleee..eeeeee ehu. eus
Cristina Cipriani
Technical University of Munich
Munich Center for Machine Learning
Munich Data Science Institute cristina. cipri an @ma. tum de
Abstract
The Neural Tangent Kernel (NTK) viewpoint is widely employed to analyze the training dynamics of overparameterized Physics-Informet Neural Networks (PINNs). However, unlike the case of linear Partial Differential Equations (PDEs), we show how the NTK perspective falls short in the nonlinear scenario. Specifically, we establish that the NTK yields a random matrix at initialization that is not constant during training, contrary conventional belief. Another significant difference from the linear regime is that, even in the idealistic infinite-width limit, the Hessian does not vanish and hence it cannot be disregarded during training. This motivates the adoption of second-order optimization methods. We explore the convergence guarantees of such methods in both linear and nonlinear cases, addressing challenges such as spectral bias and slow convergence. Every theoretical result is supported by numerical examples with both linear and nonlinear PDEs, and we highlight the benefits of second-order methods in benchmark test cases.
Introduction
PINNs have became ubiquitous in the scientific research community as a meshless and practical alternative tool for solving PDEs. The first attempts to exploit machine learning models for PDE solutions can traced back to two articles from the 90s [3] 20 while the model acquired its name and popularity through later publication (31. Due to the flexible structure of the architecture,
PINNs can used for forward and inverse problems 422 and efficiently exploited for more complex engineering practice such constrained shape and topology optimization, and surrogate modeling 35 [6]. However, the usability of PINNs for such applications is often hindered by their slow training and occasional failure to converge acceptable solutions. Due the black-box nature of PINNs, challenging to analyze their training dynamics and convergence properties mathematically [19]. Nonetheless, rapid training and reliable convergence crucial aspects of any PDE solver intended for engineering applications.
viewpoint has yielded intriguing insights, particularly in the realm of linear PDEs [40]. Although based on the assumption of overparameterized networks, this perspective has proven valuable highlighting various intrinsic pathologies in PINN training,
Related works. In this context, the NTK 15
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

such as spectral bias [39 2 29,, the complexity of the loss landscape generated by the PDE residuals [9] and the nuanced interplay among components of the loss function [38]. The salient characteristics of the NTK in the infinite-wide limit are the fact that is deterministic at initialization, constant during training, and linearizes the training dynamics due to the sparsity of the Hessian of PDE residuals
12 23
Our contributions. In this paper, we delineate the profound theoretical distinctions between the application of PINNs to linear versus nonlinear PDEs, elucidating the differences in their NTK behavior. We show that, even under the idealistic assumption of the infinite-width limit, the NTK framework fails in the nonlinear domain. Our novel contribution lies in demonstrating that the NTK stochastic at initialization, is dynamic during training, and is accompanied by non-vanishing Hessian. Given the evolution of the Hessian throughout training, we emphasize the need of employing second-order methods for nonlinear PDEs. Furthermore, we analyze their convergence guarantees, revealing that even in linear scenarios, the utilization of second-order methods proves advantageous in mitigating the issue of spectral bias. As second-order method, we employ eevenberg-Marquaae algorithm, which approximates
algorithm, stabilized version of the well-known the Hessian to make it computationally feasible even for large networks. It is important to note that goal is not propose novel training algorithm but demonstrate the benefits of using any second-order method. The reason is twofold: in the nonlinear regime, we achieve faster and better convergence, while in the linear regime, where fast convergence can be achieved by first-order methods, the advantage of second-order methods lies in their ability to alleviate spectral bias.
Our work is organized as follows: Section 2 introduces PINNs, and Section 3 covers the NTK theory, comparing its dynamics in linear and nonlinear PDEs. Section 4 examines the convergence guarantees of second-order optimization methods. Finally, Section 3 presents numerical experiments that validate our theoretical insights.
Physics-Informed Neural Networks
We address the following PDE formulated on bounded domain Ω C IRd
Ru(x) [()), x E Ω u(x) = g(z), I E an.
(I)
Here, the PDE is defined with respect to the differential operator R. while the boundary and initial conditions are collected in the function g. Notice that Ω can be either a spatial or spatio-tempore domain, depending on whether the PDE is time-dependent not. PINNs aim approximate the PDE solution u Ω + Rd with neural network ug parametrized by 8, which is a vector containing all the parameters of the network. The "Physics-Informed" nature of the neural network ug lies in the choice of the loss function employed for training
1 -0o)())(()))dddttt (o [()) = / [u)(2) - f(z) + 300
where denotes measure on the surface a2. In this work, we specifically focus scenarios where the PDE involves a nonlinear differential operator. Moreover, without loss of generality we consider the case where fx) O. Since the function f(z) does not depend on the parametrization, all of our
results hold also for the case when is nonzero. Moreover, we express (T) as
R [=]()))) = 0, Ω u(x) g(z), a2,
(2)
where D(u] gddn + doon defined as
(/u](1) = [u(z), Ox)(2), Pu(x)..... (u σ²u(z)),
(3)
denotes vector encompassing all (possibly mixed) derivatives of u until order k. while R. I& dout + IR represents a differentiable function of the components of [[u).
Remark 2.1. The importance of the function Ries in its ability to completely encode the nonlinearity of the PDE, while the term &e remains linear. Furthermore, for numerous well-known nonlinear PDEs (such as Burgers' or Anvier-Stokes equations), the function exhibits a distinctive structure as takes the form of second-order polynomial.
2

---

To illustrate this, we consider the example of the inviscid Burgers' equation, which for (T, 2) E 92 is expressed as 0,u ud u where T represents time and the space variable. It follows that
+ (T x) [u(T,a), 0,u(T, O₂ a))],
R(z1: 22, 23) = 22 + 21.3.
Neural Tangent Kernel for PINNs
We now introduce and develop the NTK for PINNs, inspired by the definition in 40 We employ fully-connected neural network featuring single hidden layer, as follows
uo(x): = W o(W'x + 60) +b²
(4)
for any x RRc Here, W Rmxd and 1 R denote the weights matrix and bias vector of the hidden layer, while W E Rdaxm and 112 E Th are the corresponding parameters of the outer layer. Additionally, hyperbolic tangent, which is a common choice for PINNs. Furthermore, we adopt the NTK rescaling ym adhere the methodology introduced in the original work [15). This is crucial for achieving a consistent asymptotic behavior of neural networks as the width of the hidden layer approaches infinity. In the following, for brevity, denote with 0 the collection of all the trainable parameters of the network, i.c. W W W', .'.
I - IR is smooth coordinate-wise activation function, such the
Remark 3.1. For the sake of brevity, we focus on the case of neural networks with a single hidden layer. However, the outcomes derived in this scenario may be extended deep networks. We leave this extension to future works and refer to (33 34) for results on finite networks with multiple hidden layers.
We consider the discrete loss on the collocation points x. E l2 and the boundary points x E En,
No
1 ( 2N, 1ro(x))2
uo(x!)) 2xx
(5)
2N
where ro(x;) = R( /b uuo) (xx)) indicates the residual term. Furthermore, N, and N₂ denote the batch size of. respectively, the collection of x" = { h and x² { which are the discrete data used for training. We now consider the minimization of (5) the gradient flow
a,0(t) V(((())
(6)
Using the following notation
un (xx) {ue (2))) (
ro(x') (racc) (x[) 1
characterize how these quantities evolve during the gradient flow, through the NTK perspective.
Lemma 3.2. Given the data and the gradient flow
6 then uo and ro satisfy the following
Oxuece @ero()(()))
(uoce) [xx))(xx) -9(x²) )
-(())
(8)
ro(e) (x")
where K(t)= = (t)) and
J(t) = Toouace) (x²)
(9)
(x)e
Proof.
The proof is presented in 40
We provide more details about the construction of J(L) in Appendix A to as Gram matrix. The analysis of Gram matrices and their behavior in the infinite-width limit [4] S yields results akin to the NTK analysis. It important to note that Lemma 3.2 applicable to any type of sufficiently regular differential operator.
The matrix also referred
3

---

3.1 The difference between linear and nonlinear PDEs
In the work (40), PINNS have been thoroughly investigated using the NTK, but only in the case of linear PDEs. Additionally, [22) extensively explores the similar case of standard neural networks with linear output. particular, they show that in the unfinite-width limit, the NTK deterministic under proper random initialization and stays constant during training. Thereby, the dynamics in 8 is equivalent to kernel regression and has an analytical solution expressed in terms of the kernel. As noted in (22), the constancy of the NTK during training is equivalent to the lincarity of the model. This characteristic related to the vanishing of the (norm the) Hessian of the network's ' output in the infinite-wide convergence results for Gram matrices hold for nonlinear PDEs when using networks as in with a focus on the unexplored case where 2. The novel contribution of our paper lies in demonstrating that in this regime this phenomenon does not hold true when dealing with nonlinear PDEs, which we prove in this section. The network architecture and its associated assumptions are relatively standard, so we refer to Assumption B.2 in Appendix B However, it is essential to delineate the specific assumptions related to the nonlinear PDE.
limit. These well-known results are reported in Appendix B In 43 the same
scaling of > where > z. However, this scaling is inconsistent with the NTK model,
Assumption 3.3 (on R). The differential operator R. is nonlinear, hence the function R is nonlinear.
Moreover, the gradient R is continuous.
The first distinction with linear PDEs arises in the convergence as m - oo of the NTK at initialization.
Theorem 3,4. Consider a fully-connecte neural network given by ④ satisfying Assumption B.2 Moreover, the PDE satisfies Assumption 3.3 Then, under Gaussian random initialization ((0), holds
((0) 2 K as m 30,
where the limit is in distribution and K not deterministic, but its law can be explicitly characterized.
Proof. A detailed proof is in Appendix C However, the basic idea is to reformulate the kernel as
((0)= (R) (0) AR())
where the matrix K (0) enclose the linear components of R. hence the derivatives of the network's output, while the matrix (0) depends on the gradient of R (so its contribution is relevant just in the nonlinear case). We can establish the convergence in probability of (0) to a deterministic matrix by taking advantage of the linearity of the operator O and commuting 4 and dg (see Lemma C.2) The matrix (0) only converges in distribution, since it function of the network output and its derivatives, whose limits are Gaussian Processes at initialization by Proposition C.1
Next, we focus on the NTK behavior during training.
Proposition 3.5. Under Assumption B.2 on the network, and Assumption 3.3 on the PDE, assume additionally that R is a real analytic function. Let us be solution of the corresponding PDE and suppose that for every m I there exists tm such that
(tm) w.lle* S Em, with Em - as m - DC.
(10)
Finally, let O(t) be obtained through gradient flow as defined C6 and denote by K(t) the corre- sponding NTK. For 0(0) ~ (O, Im). the following holds:
lim sup [E[0
[[((t) [()]
0
Proof. The proof can be found in Appendix D
Remark 3.6. It is worth noticing that our result holds under the assumption that a neural network with m > do can adequately approximate the solution w² of the PDE ( and that the training process is successful in achieving this approximation. The first assumption is justified by results such as the universal approximation theorem for neural networks [... Despite this optimistic training scenario, demonstrated in Proposition 3.5 the constancy of the kernel is unattainable.

---

Linear PDES
Nonlinear PDES
NTK at initialization
Deterministic
Random (Theorem
3.4
NTK during training
Constant
Dynamic (Proposition
3.5
Hessian H,
Sparse
Not sparse (Proposition
3.7
First-order convergence bound Second-order convergence bound
Amin (K)
0 Amio [K
(t)
~
~ 0 or (Theorem 4.2)
Table I: Comparison of the theoretical results for linear and nonlinear PDEs.
In the context nonlinear PDEs, converging to a linear regime unattainable, even in the infinite- width limit, and this inability stems from the spectral norm of Hr, which is the Hessian of the residuals та with respect to the parameters 0. Indeed, in the linear scenario, the convergence of || to 0 as + is crucial for demonstrating convergence to the linear regime, established in Proposition B.3] Similar conclusions have been drawn in [22] for various deep learning architectures. However, now show that this property does not hold for nonlinear PDEs.
Proposition 3.7. Under Assumptions B.2 and 3.3 on the network and the PDE, let us further assume that R is a second-order polynomial. Then, the Hessian of the residuals H, is not sparse and
lim || 2,
m--00
where the constant è does not depend on m.
Proof. The proof can be found in Appendix together with an explicit formula for ė.
Remark 3.8. For the latter result, we additionally require that R is a second-order polynomial, which includes many classic nonlinear PDEs like Burgers' Navier-Stokes equations.
We summarize all our results and provide comparison with the linear case in Table I Motivated by the fact that the Hessian is not negligible, we shift our attention to second-ordes optimization methods and explore their convergence capabilities.
Convergence results
Before delving into second-order methods, let us revisit a convergence result for first-order ones. Traditional analyses of the gradient descent (6) often rely on the smoothness and convexity of the loss, assumptions that may not hold the context of deep learning. concentrate the infinite-width limit, particularly in connection with the NTK analysis. While we refrain from presenting a formal proof, we highlight the notable result below.
an alternative, numerous results
Theorem 4.1. Under Assumption B. on the PDE and Assumption B.2 on the network defined by 4 consider the scenario where m is sufficiently large. With high probability the random initialization, there exists a constant μ > O, depending on the eigenvalues of K. such that gradient descent, employing sufficiently small step size n, converges global minimizer of (5) with exponential convergence rate, i.e.
L(O(t)) s ((-m)))L(((0))
Proof. See
8.
22
and others.
It is noteworthy that this result is presented at the level of gradient descent, i.e. the discretization of the gradient flow (6). which explains the constant y representing its step size. Theorem 4.1 has also been extended to various types of architectures in [5). We emphasize that this convergence result is condition which, in turn, is linked
smallest eigenvalue of the tangent kernel (denoted with Amin). In the case of linear PDEs, the tangent kernel K (L) is positive definite [ for any E [0, T). leading to positive eigenvalues. The key finding in this context is that if m is sufficiently large, K(t) x K, where K is deterministic matrix, which only
rooted the applicability of the
5

---

depends on the training input and not on the network's parameters A. As a result, in the infinite- wiith regime, the dynamics (8) can be approximated by
[xug(x²) K [xo(x)) ( 9(x²)) (ira(x")
(11)
ro(x²)
In the linear case, the key steps (i.e. the fact that the NTK is deterministic and constant) of the convergence proof of Theorem 4.1 cannot be adapted to nonlinear PDEs. Indeed, the stochasticity of the matrix and its dynamic behavior during training make the reasoning of [8] or [4] inapplicable,
and it is challenging to show that the eigenvalues of (() in the nonlinear case are uniformly bounded away from zero over training time. Nevertheless, we believe this question warrants further investigation.
Another issue linked to the NTK's cigenvalues is the phenomenon recognized spectral bias by 39 2 29]. This is related to the fast decay of the NTK's cigenvalues, which characterize the rate at which the training error diminishes. The presence of small unbalanced eigenvalues leads to slow convergence, particularly for high-frequency components of the PDE solution, or even to training failure. This occurs regardless of the linearity of the PDE differential operator R. In the next section, we show that under certain assumptions, second-order methods can help mitigate both problems.
4.1 Second-Order Optimization Methods
Due to all the aforementioned reasons and Proposition 3.7 our focus turns to the investigation of second-order optimization methods. These are powerful algorithms that leverage both the gradient and the Hessian of the loss function. Within this category, Quasi-Newton methods stand out as the most natural and widely known, relying on the Newton update rule
A(t + 1) = (tt) ((00)) ~
(12)
However, the application of this update step relies on second-order derivatives, which are prohibitively expensive to compute as the number of parameters in the model increases. Indeed, the core idea behind Quasi-Newton methods involves utilizing an approximation of the Hessian as follows
L=) (t) + H.ro()) (t)
=
(13)
in the formula (12). Here, Jtt E RxXP represents the Jacobian of the loss at the training time t. and it aligns with the definition in (9). Since the Jacobian J(t) is part of the evaluation of the gradient, the approximation [3) does not necessitate the computation of higher- ooderr derivatives.
We now tackle the issues of spectral bias and slow convergence by presenting result applicable to the Gauss-Newton method. In practice, when the number of parameters p is larger than the number of samples n, the matrix J (t) is surely singular. In this case, we consider the generalized inverse (JT (t)) instead of the inverse.
Theorem 4.2. Consider the parameter (tt) obtained by the Gauss-Newton flow below
0,0(t) -J]((())))
(14)
Then, the following holds
Tieuoc (x+) -(()D((()(())) (x)
"ucce (xcc)(x()) (o(e) (x")
(15)
where (t) BRux is unitary matrix and D E Bgxx is a diagonal matrix with entries 0 or In particular, if [() is full-rank for any E [O,T], then convergence to a global minimum is attained.
Proof. The proof is presented in Appendix
This result is significant as it indicates that when utilizing second-orde methods via (1). convergence no longer depends on the eigenvalues of ((1) as (T)), but rather on the elements of the diagonal matrix (1). Consequently, the training process becomes nearly spectrally unbiased, as the nonzero eigenvalues of the controlling matrix (15) are all 1s. Let us now compare the cases of linear and nonlinear PDEs, in relation to the assumption of full-rankness of (t) and, consequently, the NTK.

---

Linear PDEs: recent research 8 has theoretically confirmed that the NTK has full-rank in this case. Hence, convergence of second-order methods is achieved with all eigenvalues equal to 1, offering a notable advantage over (11) since the training method is unaffected by the spectral bias.
Nonlinear PDEs: showing theoretically the full-rankness is a complicated task, particularly in light of Theorem 3.5 which highlights the stochastic and dynamic nature of the NTK. Similarly, verifying numerically the full-rankness of ((t) is impractical due to the matrix's 's ill-conditioning, mentioned in [0). However, even if ((t) is not full-rank, it holds that, although some singular values are zero, fast convergence for the remaining ones is attained.
Moreover, let us stress that the result in Theorem 4.2 applies to any network, including those with finite width. Thus, while the NTK model motivates the use of second-order methods, the key insights about spectral bias and convergence hold without assuming infinite width.
Remark 4.3. In practice the Gauss-Newton method becomes less computationally expensive when combined with inexact techniques such as Krylov subspace methods, conjugate gradient, BFGS, or LBFGS [27y. has been shown that BFGS and LBFGS asymptotically approach the exact Hessian under certain conditions [21]. To extend our findings to more practical inexact methods, we can lever- age these asymptotic convergence properties. However, while this approach is theoretically sound, the speed of convergence of quasi-Newton methods to the exact Newton method - specifically their matrix approximation accuracy - depends on the minimum eigenvalue of the Hessian [27]] 6]. As discussed in our paper, the Hessian in PINNs is typically very poorly conditioned. As result, quasi-Newton methods may require an impractically large number of training steps to converge to the true inverse Hessian and, thus, to begin training higher modes.
Numerical Experiments
150
mi 530
P1: 75-0000
Bi 5933
Po 10005
100601
184
2.0
2.35
1.04
1.50
0.5
1.0
117 107
15
107
10
13
5000
10000
15000
5020 12000
15003
Bandbora
Jerelann
(a)
(b)
Figure 1: (a) Mean and standard deviation of the spectral norm of K (0) as a function of the number of neurons for 10 independent experiments. Left: linear case. Right: nonlinear case. (b) Mean and standard deviation of (t == [[ K(OO]l over the network's width m, for 10 independent experiments. Left: linear case. Right: nonlinear case.
[
5.1 Empirical validation of our NTK results
First of all, we aim at numerically validate the results presented above, by comparing the NTK in case of linear and nonlinear PDEs. Our experiments are conducted the following linear equation: Lu((x)) = 16 Notably, these results exhibit consistency across various equations and experimental setups.
r sin( (x). Meanwhile, as nonlinear PDE, we consider u(x))0,uxx) = 10
= (r))
The result in Theorem 3.4 confirmed by the numerical experiments depicted in Figure I part (a): in the linear case the NTK at initialization converges to deterministic matrix when m - co, while this does not happen in the nonlinear case. The statement of Proposition 3.5 is confirmed in part (b) of Figure by showing that the constancy of the NTK during training is not attainable in the nonlinear case. Moreover, the result in Proposition 3.7 supported by part (a) of Figure [ where we compare the sparsity of the Hessian at initialization (0) both the linear and nonlinear case. Moreover, we observe that in the linear scenario || Hr, decays as m grows, contrarily to the nonlinear example. Similarly, we refer to Figure part (b) for a comparison of the eigenvalues when training with first-order or second-order methods on Burgers' equation.
7

---

10-
2.5
10-m
1.5
10m
140
40
a
(b)
Figure 2: (a) Left: in yellow the non-zero components of the Hessian matrix initialization (up in the linear case, down the nonlinear one). Center: mean and standard deviation of the spectral norm of the (0) over in the linear case (for 10 independent experiments). Right: same Center, but for nonlinear PDE. (b) Eigenvalues of K (0) for first-order optimizer and ((0) for second-order method applied Burgers' equation.
5.2 Employment of second-order methods
Among all recond-order methods, in our numerical experiments we make use of an existing variant of the Hianberr-Marquardt (LM) algorithm, as it offers further stability through the update rule
O(l 1) 0(1) (J)(()) + AId,] V(((()),
where À is a damping parameter adjusted by the algorithm. In practice, the iterative step of LM can be considered as average, weighted by X, between the Gradient Descent step and a Gauss-Newton method. This aspect of the LM algorithm represents its crucial advantage over other Quasi-Newton methods such as Gauss-Newton when the initial guess of the solution ug is close to the correct one. The update rule of LM avoids this issue by relying simil-gradient descent steps at early iteration. Moreover, the parameter À typically decreases during training, in order to converge to Quasi-Newton method when close to the optimum. Our primary aim is showcase the effectiveness of second-order methods for nonlinear PINNs, point which has been supported by findings such as those in [26]: their approach also employs a second-order method, akin to Gauss-Newton method in function spaces. For details the modified LM algorithm, along with pseudocode, we refer to Appendix G
BFGS. Indeed, Quasi-Newton methods show good performance
Details the Networks The neural network architectures adopted in the experiments standard Vanilla PINNs with hyperbolic tangent as activation function. All of the PINNs trained in our analysis are characterized by hidden layers with 20 neurons each. Every training performed for 10 independent neural networks initialized with Xavier normal distribution [TO]. All models are implemented in PyTorch (2 and trained on a single NVIDIA A10 GPU.
Test Cases
We assess our theoretical findings on the following equations:
Equation: despite being linear PDEs, they represent suitable scenario to showcase the detrimental effect of the spectral bias on the training of PINNs, due the presence of high-frequency components in the solution.
Wave/Poisson//Coneeci
Burgers' Equation: this nonlinear PDE is commonly used to test PINNs, and usually they reach valid solution even with first-order optimizer, due to the PDE's simplicity.
Navier-Stokes Equation: it poses challenges for both PINNs and classical methods, being difficult nonlinear PDEs. We test the case of the fluid flow in the wake of 2D cylinder [17]
For the sake of compactness, we refer Appendix for detailed descriptions of the mentioned PDEs, and to Appendix H for supplementary numerical experiments not included in the main text. We compare results obtained by the LM algorithm with those from commonly used optimizers for training PINNs, such as Adam [18 and L-BFGS [24]. Where not stated otherwise, Adam is trained for 10 iterations and LM for 10 iterations. Additionally, we provide comparison with other methods that are ad-hoc enhancements of PINNs, such as loss balancing [40] (also known NTK rescaling). Random Fourier Features (RFF) (39], and curriculum training (CT) [19]. Our performance
metric is the relative 2 loss the test set, detailed in Appendix H formula 28

---

10
19
UM1 377
1971
0.2
0.90
20-4
16-1
-.44
10
Addrs
2106 SAME MIRE ISARIC 1993) 2000 4382 SCED
0.6 1.5
(a)
(b)
Figure 3: (a) Poisson equation: median and standard deviation of the relative I loss for different optimizers over training iterations (repetitions over 10 independent runs). (b) Convection equation: median and standard deviation of the I loss after 1000 iterations achieved over 5 independent runs with and without CT for different values of the convection coefficient (left) and solution obtained with LM (and other enhancement) after 5000 iterations with 100 (right).
Cousaiity
Ceasily
0.30-
0052292 cerer enreren
1.1-
retorve eerer
13"1
2.20-
0.66
10.4
0.00
ml
(a)
(b)
Figure 4: (a) Burgers' equation: mean and standard deviation of the relative L loss for various optimizers over wall time (repetitions over 10 independent runs). (b) Navier-Stokes equation: mean and standard deviation of the relative L loss over the PDE time for PINNs trained with Adam and LM (10 independent runs). Both optimization methods are enhanced with causality training.
Linear PDEs affected by spectral bias In Figure 3 we demonstrate the effectiveness of second- order methods in handling equations with high spectral bias. Part (a) of Figure focuses on the Poisson equation with high-frequency components, for which is common to use RFF [39]. On the left, we show that Adam requires RFF to converge to a reasonable solution. On the right, we observe that LM not only significantly outperforms Adam combined with RFF, but also that incorporating RFF with LM leads to remarkable loss reduction from the very first iterations. In Part (b) of Figure 3 investigate the effect of high convection coefficients 3 in the convection equation discussed in [1]. where it shown that a PINN trained with Adam necessitates of curriculum training to achieve meaningful results such spectrally biased PDE. However, (b), that the LM optimizer can handle higher values 3, especially when curriculum training is introduced. Remarkably, any other enhancements, achieves high accuracy with 3 values up to 100. This level of accuracy is not feasible with Adam and curriculum training alone, which, noted in [9] manages coefficients only up to 20.
show the left Figure 3 part
the right of Figure bb), we show that PINN trained with LM, without
Nonlinear PDEs Firstly, consider the of Burgers' equation, where convergence is achiev- able even with first-order methods. To address concerns about the additional computational time required by second-order methods, in Figure 4 part (a), we display the relative L loss over wall time when training on Burgers' equation. All training methods can reach a reasonable solution, however, while the precision of PINNs trained with Adam and L-BFGS is approximately 10-3 PINNs trained with LM can consistently attain precision around 10-5 in few iterations and very short GPU time. Figure
4 also provides a qualitative estimate of the runtime of LM in comparison to Adam and L-BFGS. The intermediate performance of L-BFGS, falling between first- and second-ord -odorr methods, explained in Remark 4.. Lastly, a similar outcome can be seen part (b) of Figure where we demonstrate that employing the LM optimizer makes it possible to obtain reasonable solution even for favier-Stokes equation in terms of relative 1,2 loss over PDE time. Notice that in this case, we employ causality training [1 for both Adam and LM.

---

5.3 Limitations and possible solutions
The major limitation of our findings is related to scalability. Traditionally, second-order methods have been avoided for machine learning models due to their poor scaling with increasing number of parameters. However, one can adopt classical PDE solution approaches, such as domain decom- position, to utilize collection of smaller networks instead of single large one. Similarly, one can embrace machine sauring-basse solutions such as ensemble models [12] or mixture of experts [I]] We advocate that existing models such as [14 25 37 could already be strongly enhanced with the usage of second-order methods for training. In the scenario where these approaches are impractical, one could also resort to techniques in the field of optimization to enable the scalability of the method. For medium to large-sized networks, the challenge of storing the matrix JJ in GPU memory becomes infeasible. This can be addressed through an inexact LM method, which involves solving the equivalent system || - rel Il = 0 using a Krylov subspace iterative method (LSQR or LSMR) [27 7. These methods only require Jacobian-vector products, which can be efficiently computed through backpropagation.
6 Conclusion
In this paper, we conduct an in-depth analysis of PINNs training utilizing the NTK framework. We elucidate the distinction between linear and nonlinear cases, and reveal that even in the optimistic infinite-width limit, favorable outcomes observed with NTK in linear cases do not extend to nonlinear PDEs. Motivated by the NTK anaylsis, we emphasize the significant advantage of employing second-order methods. These seem to mitigate the spectral bias issue and to improve convergence even for challenging nonlinear PDEs. Second-order methods, such LM, consistently achieve precision comparable or even better than the state-of-the-art presented in I3 demonstrate that convergence is attainable without resorting to typical training protocols aimed at enhancing PINNs. However, combining these enhancements with second-order training methods can further improve accuracy while reducing computational time, as demonstrated in our numerical experiments. Accuracy and convergence guarantees are indeed two crucial components for the majority of real-world applications of PDE solvers. In practice, second-order methods may be preferable when the solution contains high frequencies, when the application demands high accuracy, or when the target PDE is nonlinear. A key objective of our paper is to highlight that, despite their scalability challenges, second-order methods could help bridge the gap between black-box machine learning models and PDE solutions in scientific machine learning.
Notably, our findings
References
[1] G. Cybenko. Approximation by superpositions of a sigmoidal function. Mathematics of control, signals and systems, 2(4):303-314, 1989.
[2] M. Deshpande, S. Agarwal, and A. K. Bhattacharya. Investigations on convergence behaviour of physics informed neural networks across spectral ranges and derivative orders. In 2022 IEEE Symposium Series on Computational Intelligence (SSCI), pages 1 1172-1179. IEEE, 2022.
[3] M. Dissanayake and N. Phan-Thien. keural-network-basse approximations for solving partial differential equations. communications in Numerical Methods in Engineering, 1053):19--011, 1994.
[4] S. Du and J. Lee. On the power of over-parametrizztion in neural networks with quadratic activation. In International conference on machine learning, pages 1329-1338. PMLR, 2018.
[5] S. Du, J. Lee, H. L. Wang, and X. Zhai. Gradient descent finds global minima of deep neural networks. In International conference machine learning, pages 1675-1685. PMLR, 2019.
[6] R. Fletcher. A modified marquardt subroutine for non- non-linear Ieast squares. United Kingdom Atomic Energy Authority Research Group Report, 1971.
[7] D. C.nL. and M. Saunders. Lsmr: An iterative algorithm for sparse least-squares problems. SIAM Journal on Scientific Computing, 33(5):2950-2971 2011.
10

---

[8] Y. Gao, Y. Gu, and M. Ng. Gradient descent finds the global optima of two-layer physics- informed neural networks. In International Conference on Machine Learning, pages 10676- 10707. PMLR, 2023.
[9] H. P. Gavin. The grrenberg-marquarte algorithm for nonlinear least squares curve-fitting problems. Department of civil and environmental engineering, Duke University, 19, 2019.
[10] X. Glorot, A. Bordes, and Y. Bengio. Deep sparse rectifier neural networks. In Proceedings of the fourteenth international conference on artificial intelligence and statistics, pages 315-323, JMLR Workshop and Conference Proceedings, 2011.
[11] I. C. Gormley and S. Frühwirth-Schnateee Mixture of experts models. Handbook of mixture analysis, pages 271 1--077 2019.
[12] K. Haitsiukevich and A. Ilin. Improved training of physics-informed neural networks with model ensembles. In 2023 International Joint Conference on Neural Networks (IJCNN), pages 1-8. IEEE, 2023.
[13] Z. Hao, J. Yao, C. Su, H. Su, Z. Wang. F. Lu, Z. Xia, Y. Zhang, S. Liu, L. Lu, et al. Pinnacle: A comprehensive benchmark of physics-informed neural networks for solving pdes. arXiv preprint 2rXiv:2306.08827 2023.
[14] Z. Hu, A. D. Jagtap. G. E. Karniadakis, and K. Kawaguchi. Augmented physics-informa neural networks (apinns): A gating network-based soft domain decomposition methodology. Engineering Applications of Artificial Intelligence, 126:10718 107183, 2023.
[15] A. Jacot, F. Gabriel, and C. Hongler. Neural tangent kernel: Convergence and generalization in neural networks. Advances in neural information processing systems, 31, 2018.
[16] H. Jeong, C. Batuwatta- Gamtte, J. Bai, Y. M. Xie, C. Rathnayaka, Y. Zhou, and Y. Gu. A com- plete physics- physics-informed informed neural network- based framework for structural topology optimization. Computer Methods in Applied Mechanics and Engineering, 417:116401, 2023.
[17] X. Jin, S. Cai, H. Li, and G. E. Karniadakis. NSFnets (Navier-Stokes flow nets): Physics- informed neural networks for the incompressible Navier-Stokes equations. Journal of Computa- tional Physics, 426: 109951, Feb 2021. ISSN 0021-9991. doi: 10. 10.1016/..cpp.2020..r 2099000.00000
109951. URL
http://dx. doi org/ 10. 1016/j jep. 2020. 109951
[18] D. P. Kingma and J. Ba. Adam: method for stochastic optimization, 2014. URL https /arxiiv.a org/argg////o..ooooooooo.ooo.o /412. 6980
[19] A. Krishnapriyan, A. Gholami, S. Zhe, R. Kirby, and M. W. Mahoney. Characterizing possible failure modes in physics-infeee informed neural networks. Advances in Neural Information Processing Systems, 24:26548-26560, 2021.
[20] I. E. Lagaris, A. Likas, and D. I. Fotiadis. Artificial neural networks for solving ordinary and partial differential equations. IEEE transactions on neural networks, 9(5):977-100, 1998.
[21] D. Lin, H. Ye, and Z. Zhang. Explicit convergence rates of greedy and random quasi- newtin methods. Journal of Machine Learning Research, 23(162):1-444 40, 2022.
[22] C. Liu, L. Zhu, and M. Belkin. On the linearity of large non-linear models: when and why the tangent kernel is constant. Advances in Neural Information Processing Systems, 33:15954- 15964, 2020.
[23] C. Liu, L. Zhu, and M. Belkin. Toward a theory of optimization for over-parameterized systems of non-linear equations: the lessons of deep learning. arXiv preprint arXiv:2003.00307, 7, 2020.
[24] D. C. Liu and J. Nocedal. On the limited memory BFGS method for large scale optimization. Mathematical programming, 45(1):503-528, 1989.
[25] B. Moseley, A. Markham, and T. Nissen- Finite basis physics-informed neural net- works (fbpinns): a scalable domain decomposition approach for solving differential equations. Advances in Computational Mathematics, 49(4):62, 2023.

---

[26] J. Müller and M. Zeinhofer. Achieving high accuracy with pinns via energy natural gradient descent. In International Conference on Machine Learning, pages 25471-2544 PMLR, 2023.
[27] J. Nocedal and S. J. Wright. Numerical optimization. Springer, 1999.
[28] A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein, L. Antiga, et al. Pytorch: An imperative style, high-performance deep learning library. Advances in neural information processing systems, 32, 2019.
[29] N. Rahaman, A. Baratin, D. Arpit, F. Draxler, M. Lin, F. Hamprecht, Y. Bengio, and A. Courville. On the spectral bias of neural networks. In International Conference on Machine Learning, pages 5301-5310. PMLR, 2019.
[30] M. Raissi. Deep hidden physics models: Deep learning of nonlinear partial differential equations. The Journal of Machine Learning Research, 19(1):932-95) 2018.
[31] M. Raissi, P. Perdikaris, and G. E. Karniadakis. Physics-informedd neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. Journal of Computational physics, 378:686-707, 2019.
[32] M. Raissi, A. Yazdani, and G. E. Karniadakis. Hidden fluid mechanics: Learning velocity and pressure fields from flow visualizations. Science,
2020.
[33] M. Seleznova and G. Kutyniok. Analyzing finite neural networks: Can we trust neural tangent kernel theory? In Mathematical and Scientific Machine Learning, pages 868 8955. PMLR, 2022.
[34] M. Seleznova and G. Kutyniok. Neural tangent kernel beyond the infinite-width limit: Effects of depth and initialization. In International Conference on Machine Learning, pages 19522-19560. PMLR, 2022.
[35] Y. Sun, U. Sengupta, and M. Juniper. Physics-informed informed deep learning for simultaneous surrogate modelling and pde-constrained optimization. Bulletin of the American Physical Society, 2022.
[36] M. K. Transtrum and J. P. Sethna. Improvements to the mevenaerg-marruardt algorithm for nonlinear least-squares minimization. arXiv preprint 1201v:1201.58800 2012
[37] H. Wang, R. Planas, Chandramowlishwaran, and R. Bostanabad. Mosaic flows: A transferable deep learning framework for solving pdes on unseen domains. Computer Methods in Applied Mechanics and Engineering. 389:1 14424, 2022.
[38] S. Wang, Y. Teng, and Perdikaris. Understanding and mitigating gradient flow pathologies in Physics-informed neural networks. SIAM Journal on Scientific Computing, 4315):A3055-A3081, 2021.
[39] S. Wang. H. Wang. and P. Perdikaris. On the eigenvector bias of fourier feature networks: From regression to solving multi-scale pdes with physics-informed neural networks. Computer Methods in Applied Mechanics and Engineering, 113938, 2021.
[40] S. Wang, X. Yu, and P. Perdikaris. When and why PINNs fail to train: A neural tangent kernel perspective. Journal of Computational Physics, 449:1 10768, 2022.
[41] S. Wang. S. Sankaran, and Perdikaris. Respecting causality for training physics-informed neural networks. Computer Methods in Applied Mechanics and Engineering, 421:1 16813, 2024.
[42] L. Yang, X. Meng, and G. E. Karniadakis. B-pinns: Bayesian physics-informed neural networks for forward and inverse pde problems with noisy data. Journal of Computational Physics, 425: 109913, 2021.
[43] Z. Zhou and Z. Yan. Is the neural tangent kernel of pinns deep learning general partial differential equations always convergent? Physica D: Nonlinear Phenomena, 457:133987, 2024.
12

---

Supplemental Material
This supplemental material is divided into the following eight appendices.
Appendix A: Details about the NTK Matrix
Appendix B: Standard NTK results for linear PDFs
Appendix C: Proof of Theorem 3.4
Appendix D: Proof of Proposition
Appendix E: Proof of Proposition 3.7
Appendix F: Proof of Theorem 4.2
Appendix G: Details about the Numerical Experiments
Appendix H: Further Numerical Experiments
In the following denote with || |l2 and (::-) -) the Euclidean and scalar product RR respectively, The Euclidean ball centered in æ with radius R is indicated with (xx, R). We denote with the spectral norm of matrix and with I, the identity matrix of dimension
We abbreviate with i.i.d. independently and identically distributed random variables. [[ denotes mean of the random variable E IR while Cov[X] is its covariance matrix. Convergence of
X, to X in distribution is indicated with X,, 2 X, while convergence in probability with X,, & X. GP denotes a Gaussian Process.
The operator V denotes the gradient of a function R. while f y) the partial derivative of with respect to the variable z.
A Details about the NTK Matrix
We define the following matrices
Ogue(a)) (x)) 10,, uo(x)
Oe. uz(x))
(ero(x)
[0u, ro
Oo., ro(x))
Tb, (1(u)(<))
On., @1
0, (0D(v
Bo, t (x)
don (x(u))(3))
By Ogug(x"), dere(x") and 0 we mean the same matrices as before, calculated in each xf (x, respectively) and stacked vertically, e.g.:
Oe (u0lue)
Og 0 (nol(x))
The only exception is given by:
V (((o(uo(x'))
(16)
R (x)
13

---

While the Hessians have the following structure:
0 0,0, ux(
0 0m uo(
(x) ;-
0g 0, ux(x)
om uo (x())
0g 01 ro(x)
0y nm. ro(2)
H,(x) :-
(17)
0% oe ro(x)
0g mm ro(or)
0g 0, 2, (x(
0% om (x(
Ha, (a) ==
0g 0, (x)
00m (x)
B NTK for linear PDEs
First of all, we list here all the assumptions needed on the differential operator R. and the neural network in (4).
Assumption B.1 (on R). The differential operator R is linear, which implies that R is linear.
Assumption B.2 (on the network). Given the network we assume the following properties:
(i) there exists a constant C > 0 such that all parameters of the network are uniformly bounded for & [),
sup |10()]lo SC with C. independent from m. KE [0,7]
(ii) there exists a constant
> 0 such that
(2)()) - -((مt)) dr SC,
N
Σ
c ((ua(r))<<)) dr SC.
(ii) the activation function and as well as its derivatives (8) up to power k 1 are smooth and (0x00 S for 1 where (()) denotes the i-th order derivative of a.
In order to present the results, we denote with H, the Hessian of the residuals ry with respect to the parameters 0. The Hessian plays important role in Proposition B.3 which aims to list all the prior results that can be derived by combining Theorem 4.4 of 40 Theorem 3.2 of [22]
Proposition B.3. Consider a fully-connected neural network given by (4) under the Assumption B.2 on the network and Assumption B.1 on the PDE. For the minimization of the loss function (5) through gradient flow. starting from Gaussian random initialization ((0). it holds that for any T >> O.
the randomly initialized tangent kernel K(0) converges in probability deterministic kernel K as m + oo:
the Hessian matrix H, of the residuals is sparse and
||H,|| O
hence the spectral norm converges to as m + Oo;
as a consequence, the NTK is nearly constant during training, i.e. [[((t) K (O)]12 0; 00 te[o,1]
lim sup
Proof. The proof can found in the papers mentioned above or as a special (linear) case in Appendix
B
14

---

C Proof of Proposition 3.4
First of all, we derive a result about the behavior of the vector partial derivatives $[u]. The Proposition C.I below is generalization of Theorem 4.1 in [40] for any derivative of order k. This means that there are no nonlinearities involved, since these are encoded in the function R. Moreover we study the full vector and not each component separately it is done in 40) This is needed in the following proofs.
Proposition C.1. Consider a fullly neural network of one hidden layer as in 4 under Assumption
B.2
Then, starting from (() ii.dd from (O, Id), it holds that
{(uo(o)(()) 2 GP(O, E(r,x') for any I, x, E Q
oo, where D means convergence in distribution and E is explicitly calculated.
Proof. To ease the notation, omit the initial time 0 and denote ua(0) with ug. Similarly, all the weights matrices and biases W* (O), D (0), 8 (0) indicated with W W b 89. Now according to the definition t and the fact that it is linear, we obtain that
① ((o)(())
- W Blo(W + 69) m
Σ W3 WoW
According to our assumptions, (W + b9) are i..dd. random variables. We prove below that their moments finite, hence by the multidimmesiooal Central Limit theorem (CLT) can conclude that, for every E ?
F(ua[<>) 2, N(,,
with covariance matrix:
T(x) = Cov,, [fo[ur + +
Now we compute the covariance of the limit gaussian process. In order to do first need to show that &2(uo] (x) are uniformly integrable with respect to m for every i = 1,...... k. follows from:
sup [u)))) sup E m
+ +681
sup E m
( + 62)]² = E [ [W 89)1
j=0
_
where C maxigisk | /oo and (1)) indicates the i-th order derivative of while 7
maxisisk Ey- < 0o.
Now, for any given point xr x' E we have that
(uol(()))) E(xx = lim E (x)
m-+00
lim
E
w/. W/ 2.W 6,,11, +b6,11 =
m=+00
h,1==
Σ (w 12 ( [891
lim
E
m=+00
1=1
E3 (o..)) ( v) + +)]]
15

---

Lemma C.2. Consider a conneced neural network of one hidden layer as in under
Assumption B.2 Let define
3g (o (xx)
(0)
0g 2 (ocoo) (x")
where & the collection of all the partial derivatives of u, as in (3). and ((0) ~ NN Id) iid. follows that K (0) converges in probability to a deterministic limiting kernel as m - cc.
Proof. The component Og Douooo) is linear, hence is standard as in (40, Lemma 3.1. While the rest of the matrix needs to be generalized to any derivative O; for i = 11, -...... k. For any j.j I,...,h k and every x, x' E Ω consider each entry
4m
0, (, (oP,((oro))((')²=
De t [Moco) 0g $, (='co)((²')
Am
p Σ 1m1
00. (o(o)
(x) b [0g uo(oi] (x')
where the second equality follows from Schwarz theorem (because the smoothness the derivatives of u), and the linearity of the operator 4. This sum has to be split in 4 parts, one for each possible type of 0g (in W W 50 or b²). Here we present the case when 0g = W.,. while the other cases are analogous:
Σ Φ @w? mo(o)
(x) (x)$1
(Ow; Moros] (a')=
Φ (W))(()) + Φ2 (o(W)(0))' + +f((
Eu, [ + +)] [a + +)]]
and the limit in probability in the last line comes from the law of Large Numbers.
Lemma C.3. Suppose that there exist R 0 and such that 40 E B(0(0), R) it holds
[x,,(x)) )|| <,;
|H (x')] <e Vjj = 1,...k..
Then maXoe B)) [ (t) K2 (0)] = OeeR).
Proof. Using the properties of the spectral norm, we just need to bound each block of J(t) as follows
N.
[|((t) J(O)
| 1 (uoto(()() 0g 2 (uo(o)) +
(0 - Osuoro) (x!)
< Nr max 1 (voto((()) - a B, (M0(0))))
+ No max (1000 - (ouero) (2+)
< KN, max 1.j
max (EEB(O()).R))
Ho, (z:3i] 110 Aoll
Nax
max (EA()())) [Hn(z:]]| 110 0all
$ No)
Hence:
(1) - (o)0) -0000" ()(())))))))-)---)))) - J(OII (1/(0)++ (1)()))) max(kN, No)c [|J(0))))
and the last norm is bounded on B(0(0), R)) by smoothness of the model.
16

---

Lemma C.A. Under Assumption B.1 on the PDE and Assumption B.2 on the network, then Kg. is
nearly constant during training,
lim m= rE|o,T)
sup (1) - (x(()] - 0.
Proof. The statement follows by combining Lemma C.3 and Lemma
E.1
Now we are in position to prove Theorem 3.4
Proof. (of Theorem 3.4)
By using the chain rule on the residual term, explicitly compute:
Oous(x" (xx")(x")
(0)
(ero(xe) =
(xro(x")
VR th (xo((x)) 1) () (x)) (dbu)(x)) (
(xx)
ld
(xx))(x))
[0ou0(x)) (x0
VR(4[uo|(x))
Og (oloooo(x))
V VR( (x)l)()) (T
(a(())
Ax(0)²
where we have denoted 0(0) with 0 and omitted the initial time step and ( (x")) is defined in (16). Let us first observe that the linear part, i.e. K (0), converges in probability to deterministic limit by Lemma C.2 Moreover, ["erro]] converges in distribution tion C.1 Regarding the nonlinear part denoted with (O)),, know by assumption that V R is continuous function, hence we can apply the Continuous Mapping Theorem and conclude that
gaussian process by Proposi-
v V 2 (GP((0, E(r,x'))) for x, x' E Ω.
From this, the convergence of ((0) follows by Slutsky's theorem.
D Proof of Proposition 3.5
Proof. Recall that we denote with (t) the NTK obtained with f(t), evolving according to the gradient flow (6) Similarly, K(0) is the NTK at initialization, i.e. with ((0) ~ (0, Idm). We can
rewrite the kernels in terms of their linear and nonlinear part as we did for the proof of Theorem 3.4 and obtain
lim sup [K(1) - K(O)] > lim [((mm) ([(m))(((j]] [0
lim
|An(tm) AA00 (0)
lim (4 -
m +00
- [| (K) - K4 (m)
where the last is obtained by applying the inverse triangular inequality, after summing and subtracting the needed terms. Moreover, by considering that supreJo.T] || (t) Ka (0)] as m + co by
Lemma C.4 we obtain that
lim sup 0= [c(O.T]
[[((t) K(O)]I 2 lim || (o)An
m=000
Td
ld
lim
0
V (R(F(u(tmm))
(x)())
(R(
Id
Id
0
Ko
R(
0
V
17

---

Observe that (TO) implies that (fu(tm)] + $(u*]. hence V V (2(0))))))
V (R( as + 0o
by continuity of V.R Combining this and Lemma we find
Td
Td
lim
v (R(((u(tm)]
(0)
0
V [R( O1)
m--00
Id
Id
(o(0)
v [R(D(uo(o)))
v [
Id
ld
r
Ko(0)
0
V VR(((u)))
V V [2)))))
ld
ld
(0)
v R (R((()))))))
v VR(
Finally, to prove our statement, we just need to show that the matrix above is not 0 almost surely, or at least one of its components. Let us fix a collocation point x E f and let us define the function f I I:
f(w) = V[(O(""("")) Ka (C00) VR(D(("(<)))
VR(w) (o VR(w)",
(18)
where (0)0) denotes the kernel evaluation at a fixed collocation point. The first term on the right hand side (18) is deterministic vector, so f is well defined deterministic analytic function. Moreover, if R is nonlinear. is not identically zero.
By the properties of analytic functions we can conclude that Leb( w E R R*I[(w) 0}) = where Leb denotes the Lebesgue measure. Notice that & (=x(o)](z)) ~ A(O E(x)) in the infinite-width limit as proven in Proposition C.1 and a consequence of that proof is that Σ(x) is not singular. This implies that
PPf (Pf(a()) (=)) = 0) = 0.
E Proof of Proposition 3.7
We present here some preparatory results.
Lemma F.1. For any k...k and any S2 the Hessian H4, (1) defined in
17
is such that
1
|H (2)]l = O
mm
Proof. Recall that
(H4,(()) 0g t. (x)),
where L,j 1,...,
By the linearity of the operator 4; and the smoothness of the activation function as Assumption it holds that
B.2
o T [uo] Φ; [ ug
For specific choice, e.g. first parameter is 0, = W3 and the second is 0: = WP, holds that
ohy wp (. db 1 wows)
( 69) 1 S C Vm
(19)
mm
where the last inequality follow from Assumption domain Q. Since the calculations of (9) similar for every combination of parameters W W0 50, do not report them here. Furthermore, we notice that the derivatives involving b are zeros and hence we obtain that Hl. composed by blocks (3 3 combinations of parameters). Each block is diagonal matrix, whose elements are bounded by C *mm By considering that the spectral norm of bound the spectral norm of each diagonal matrix is equal to the maximum of its components,
B.I Assumption
B.2 and the boundedness of the
18

---

block by C. Moreover the spectral norm of a matrix can be bounded by the sum of the spectral norm of its blocks, hence:
1
[H,(1)) 9C
O vm
vm
We can now prove Proposition 3.7
Proof. In the nonlinear case the Hessian of the residuals is (H,())) =bh00 ro(x) = do V R [(
=(0b, B/uol (x), V R (uo (x)) fluol(x))
V V [r)
Agg
for every collocation point J Q. The matrix Ho is defined in (77). Moreover, Lemma provides
that the spectral norm of B goes to 0 in the infinite-wide width limit. Moreover, by making use of the inverse triangular inequality, we obtaín that for any x E Q2 holds
[H) 2 lim ||| || m-+00
lim
| = lim ||]] m=000
12-200
According to the definition of spectral norm, we have that
lim ||||| lim max
||A2||2
lim
|
||
where E :
Let us now focus on the term | By using some standard
mm
incqualities and taking advantage of the fact that each entry of 2 is - obtain that
Vm
|| 2
||A=||
(Az):
Σ Ai
Vm
i,j=1
1
Oo. + uoo (0), 2
m
0o $ [ (z),
V2 R ( (x))
0o 1 (2))(())
mm j=1
Without loss of generality, we can restrict our focus to 0; = W72 and 0j W since the spectral norm of matrix is greater or equal then the norm of its submatrix, and study the term
lim
do. 2 (x) lim
m
m
(20)
-B.,)_N((.1) Bio(u +z)] =:
by the law of large numbers. In particular, deterministic. Notice that here we have considered generic f since, according to Lemma C.4 Oo, t is constant. By combining this result with the
previous one, we obtain that
lim [H+())] 2 w (= 26
where ē is a deterministic constant that does not depend on m,, but only on the value of 2 (which constant because R is second-order polynomial) and the vector defined in (20)
F Proof of Theorem 4.2
Proof. The gradient flow equation in case of Gauss-Newton methods has been defined in 14 for J(t) E Raxp where E [0, T]. It follows that
(t)) (tr
[Qquote] 0,0(t) = J()),,())= -00000) [bre(())
(o(c))
re(e)
19

---

where the last equality comes from plugging in (14 into the equation. Now, let us consider the case when p >> n, then the singular value decomposition of ((t) is as follows
J(t) U En Op-n VT
where U RR"", Σ by the square roots of the eigenvalues of the NTK. We drop the dependence on time U, and V to case the notation. Let us now study the term
R V E Rxxp and Σn E R is a diagonal matrix with elements given
(()()")))((() (²) = UY v (N
= VU
En
0p-n] )
1
[ En UT
= [En Op-n]
[In
Op
0p-n
Op-n
22
En
= [In
Op-n]
p-n
U³
Op-n
On
Op-n
= UL (E2)
UU
where D is obtained from E.,, by replacing the non-zero components with 1. In particular we can rewrite the Gauss-Newton flow as:
Opuoce) DU Orrace)
Uo(c)) re(m)
Notice that it has the same form of the gradient flow in Lemma 3.2 but the Neural Tangent Kernel is replace by a matrix with non-zeroes eigenvalues 1. This can be translated as: second-order optimizers are almost spectrally unbiased. Moreover if Jtt) stays full rank during the training, we can obtain the result of convergence regardless of the singular values of J i.c.:
[0: [ue(e)
Ocroco)
ro(e)
G Details about the Numerical Experiments
G.1
The LM Algorithm
In the following, provide more detailed deseription of the version of the Levenberg- -Marburrd algorithm along with its pseudocode and the details of the experiments whom results are shown in Section 5
The main difference between the Eevenherg-Marquarde algorithm and other Quasi-Newton method that general Quasi-Newton methods are line-search approaches, while LM is trust region approach. In practice, line search approaches determine descent direction of the loss funetion and thereinafter determine suitable step size in such direction. On the other hand, trust region method determines an area where the solution lies and computes the optimal step. If this step does not provide enough improvement in the objective function, the search area is reduced and the search is performed once more. We refer to 22 for a thorough description of trust region and line search methods.
In the following part, we drop the dependence on training time continuous function and identify f(4)) = f& for some discrete time 44- As already mentioned in Section 5 the update step e of the LM algorithm is computed follows:
Gk
[ Jv VL(O),
(21)
where Dj. is a diagonal matrix of size n. the classical LM algorithm, this matrix D₁ given by the identity matrix. Another viable alternative recommended in (6] is to use the diagonal of J Jk.
20

---

For our model, we choose Dj: to be simply the identity matrix, which appears to be more stable when J J3 is singular.
Another typical modification to the Srremberg-Marquarieee algorithm is the introduction of the geodesic acceleration 36
ak = - [J] Jx + ADx] 2y Uk
(22)
The goal of the geodesic acceleration is to introduce component which does consider all the components of the Hessian of the loss when the residuals are small and when the Hessian of the residuals is not negligible.
Moreover, every iteration, has specify a criterion C whose objective is evaluate the relative improvement of the model parameterized by 02 with respect to the update step vk- The criterion depends on the modification of the LM algorithm chosen. For our algorithm we use the same condition as (9) i.c. Cj < toll where C is defined as
C = I(0x)² L(Ok + 0g
(23)
(ok, d Da VL((%))
We provide in Algorithm [ the pseudocode of the modified LM algorithm that we chose for our numerical experiments, inspired by the implementation of 9 and modifying it by adding the component of the geodesic acceleration.
Algorithm 1 Modified Swvenherg-Maaqquaan Algorithm
Input: Maximum region area A 0, Region Radius 0 < Ao < A, Tollerance tol E ()
1) for k = 0, 1, do Compute 2/2 as in Equation 21
Compute criterion C as in Equation while Ck < tol do
23
À min(2,,
Compute va with the new value of À Compute criterion Cg as in Equation end while
23
0x+1 0g + Uk
Ax+1 = max( ( 1)
Compute dg as in Equation if 2 _ |||nk| then
Pk-1 1 + 0k
end if
end for
The main focus of the Levenberg-Marquard method is to decide the size of the trust region. In practice, at every iteration, one wants to find a better solution and afterwards reduce the size of the trust region. When this does not happen, the solution is to enlarge the trust region in order to look for a better solution. In our method we choose to include the region search as part of the inner loop, for line search approaches. This means that the iteration itself can be slower, but more accurate, which is why we include in the numerical evaluation also the computational time.
G.2 Poisson Equation
The Poisson equation that we choose for our study is monodimensional instance of the PDE defined in 39 for [0, 1] and we try find the solution Ω - R. In particular, we want solve the following equation:
8.u= f(x). x Ω u(0) = =
(24)
As in 39). the function f is constructed such way that the exact solution of Equation is given by:
1
u(x) sin (2nzz) + sin sin05.r.).
10
21

---

This approach is done to evaluate the behavior of PINNs when the target solution presents a high frequency and low frequency component. We then train the PINN model by sampling N, =03 points in Ω with latin hypercube sampling.
G.3 Wave Equation
We opt to solve the wave equation below for each (z, E === [0, 1J² and aim to find the solution u Q + IR In particular, we aim to solve the following equation:
SLu 0Lu
(z, f2
(xx 0) = sin() + sin(4T.r),
# E [0, 1],
2
(25)
u( 0) =
x E [0, 1],
(0 u
TE [0, 1].
With C' heing equal to for our case. is straightforward obtain the correct solution of this equation through Fourier transform. In particular, the exact solution of Equation (25) is given by:
1 sin (in44)))) cos(8
u( r) sin (TT) cos(2 2
We then train PINN by sampling N, = 10 training points S2 for the PDE residuals with latin hypercube sampling, and N₂- = 3 10 points for training the model against the correct solution at O..
G.4 Burgers' Equation
Burgers' equation is ID version of Navier- equations. Its solution at high times present discontinuity, which makes it challenging for spectrally biased architectures. The specific instance chosen in our numerics for Burgers' equation is the same as in [30.. In particular, we refer to the exact same data provided by the authors. particular, given (2, 7) E f [-1, x [0, 1]. we solve for u : l2 + R the following equation:
0,u + ud,u vOLu =0,
(I, E E E 11,11, TE [0. 1].
u(I; 0) - -ii((
(26)
w(-1, r) u(1, =0
with the diffusivity being equal 0.01 for this specific instance. The correct solution is provided publicly by the authors of [30].
Training is performed with N, 10 collocation points for training the PDE residuals, sampled with latin hypereube sampling, and N₂ = 3- 10 points for training the boundary and initial condition in a..
G.S Navier-Stokes Equation
The most interesting scenario taken in consideration for our experiments is that of Navier-Stokes equations. In particular, we aim to solve the fluid flow in the wake of a cylinder in 2D tackled in K7| In particular, ū 7:Ω Ω + R which is defined as u(zz y. t) = [u(z:, y, t), (ar, y. f(z,,(( are respectively the horizontal and vertical components of the fluid velocity and p is the pressure at point. Navier-Stokes equations are then expressed in vectorizer form as follows:
have (Ir, y, t) E = [2.5, 7.5] 1-2.5, 2.5] [0, 16) and we wish to find
p(z, y. t)l" In particular and
0.u + ud_u v0yu
(i u + 0.p 0. (x. T) Ω
Re
1
0,0 + u0,v + v0yv
Re (0.. + 070) + 0yp = 0, (I,8,7) E
0.u + 0,v = 0,
(2,8,7) E
(27)
u(x, = Iuo (x, v),
(x,y) E [2.5 2-2.5, 2.5],
olz,y,0 0) = gou (x, v),
x 2.5,2.5)
u(2.5, y, 7) 1,
e [0, 16).
v2..5 y, T) = O,
e x [0, 16).
22

---

30
2000
2030
3003
7300
10000
12100
15000
37503
20000
Figure 5: Mean and standard deviation the relative L loss on the test set on the Wave equation for Adam, L-BFGS and LM optimizer over iterations (repetition over 10 independent runs).
where Re represents the Reynolds' number, which an adimensional quantity defined by the problem and is set to 100 for case. The initial conditions (g 9v.) can be found in the repository published [32). as well as the correct solution. The conditions at x = 225 represents the
by the authors fluid velocity imposed at the inlet, and further conditions are given by the presence of a cylinder centered in (x. y) = (0, 0) with radius 0.25. Furthermore, an additional condition appears at the borders, namely where y +2.5, where the no-slip condition be chosen (u
0) the
correct solution can be given as boundary condition. Since the simulation provided free-flow stream, we use the correct solution the boundaries.
32 refers to a
To train our PINNs, we use N,
10 collocation points for training the PDE residuals, sampled with latin hypercube sampling, and No 2. 10 points for training the boundary and initial condition in M. Morever, at every iteration, we minimize the loss on random batches of the training data, respectively 10t points for the residuals and 5 10xh for boundary and initial condition.
H Further Numerical Experiments
In this Appendix we present some additional numerical experiments. Notice that as a performance measure we utilize the 7,2 relative loss, defined as follows
(u(a,) u(x:)) w(x:))
(28)
where u is the exact solution and the approximated one.
In Figure 5] we showcase the relative L loss obtained on the test set during training on the Wave equation with the aforementioned optimizers. While Adam and L-BFGS get stuck relatively fast in local minima, the LM algorithm is able to decrease the loss consistently, despite the complexity of the problem. The poor performance of L-BFGS can be motivated by two factors. On one hand, the Hessian computed during BFGS iterations is merely an approximation of the true Hessian; on the other hand, convergence to the true solution is heavily hindered since the initial guess is typically not close to the correct one.
In Figure 6 and Figure 7 it is possible to notice the effect of the spectral bias: the PINN trained with Adam can capture only the lower frequency components of the true solution, while the model trained with LM performs better as the spectral bias is alleviated in accordance with Theorem 4.2 It worth noticing that the same holds even when introducing the loss balancing suggested in 40 its performance is showed in Figure 7
Finally, in Figure we show that by employing the LM optimizer, is possible obtain reasonable solution even for PDE as complex as Navier-Stokes with relatively small architectures. Notice that the scale in the two plots are different.
23

---

00
000
-44
0.2
0.3
Figure 6: Experiments on the Wave equation. Left: Prediction of the parametrized solution of PINN trained with Adam (Left) and LM (Center) alongside with the true solution (Right).
1.0
Exact
LM
Acam-LLB
0.5
-0.5
Figure 7: Experiments the prediction of the solution of Poisson equation with LM and Adam (with loss balancing), both compared with the exact solution.
10c
LM
10-1
10-1
Training
20
2
10-2
10-1
10
0
500
1000
10000
20000 300n0 40000
Itarations
Itarations
Figure 8: Mean and standard deviation of the training loss over the iterations for Adam, LBFGS and LM on Navier-Stokes equation (for 10 independent runs).
NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: The Abstract and the Introduction clearly state all the claims and contributions made in the paper. This holds also for assumptions and limitations, which are shortly mentioned in the abstract and tackled more in depth the Introduction, alongside related references.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or
NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize other settings.
is fine to include aspirational goals as motivation as long is clear that these goals not attained by the paper.
24

---

2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: The practical limitations of the work are mainly connected to the scalability of the method, which is tackled in Section 5.3 Additional limitations on the theoretical analysis are clearly mentioned throughout the paper, along with related research directions and references.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create a separate "Limitations" section their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and faimesss.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?
Answer: [Yes]
Justification: All the thcoretical results are accompanied with solid proofs which are included in the appendix, for the sake of brevity, and sketched in the manuscript, in order provide an intuition to the reader. The assumptions made for each proof are also fully included (at times in the appendix).
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
25

---

4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
The pseudocode the main experimental results are given in the appendix. Moreover, the paper does mainly rely on existing algorithms and methods which are properly referenced across the paper. Furthermore, the majority of the methods referenced are also available in common Python packages.
Justification:
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how reproduce that algorithm. (b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) If the contribution is new model (e.g.. a large language model), then there should either be way to access this model for reproducing the results or a way to reproduce the model (e.g., with an open-source dataset instructions for how construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification: Making the code public is currently in discussion with the partner institutions. Due to legal reasons, it might not be possible to have released as open source. However, despite research not including any unconventional implementation, we make the code available per request to the corresponding author.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips.cd cc/
for more details.
26

---

While we encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https Codeeuomisssse asioniicc for more details.
7 //nips wtt.ccccpppppgcluieeeeee
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing much information possible in supplemental material (appended the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Details of the training and testing proceedures of all the experimental results obtained the paper are shortly provided in the paper and thoroughly discussed in the appendix
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core the paper to a level of detail that is necessary appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment Statistical Significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: All the tests which include variability (such as initialization of the networks) are obtained for several runs, and are showcased alongside the variability obtained during training.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emror bars, confi- dence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
27

---

For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out range (e.g. negative enror rates).
If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: The experimental set up used to obtain the numerical results provided in the paper is fully described.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make into the paper).
9. Code Of Ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics
Answer: [Yes]
Justification:
All the authors have reviewed the NeurlPS Code of Ethics.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., if there is special consid- eration due to laws or regulations their jurisdiction..
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification: The Conelusion delve on the potential broader impact of our work for future research direction.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation. generating fake profiles, surveillance), fairness considerations (e.g.. deployment of technologies that could make decisions that unfairly impact specific groups). privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
28

---

The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of MI.).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: [NA]
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g.. code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [Yes]
Justification: All the scientific outcome of this paper generated by the authors. Methods and algorithms fro third parties are properly referenced across the paper.
Guidelines:
The answer NA means that the paper does
existing assets. The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include a URL.
The name of the license (e.g., CC-BY 4.0) should be included for each
For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, wwpppswi thcode. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license dataset.
For existing datasets that are re- e-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
this information is not available online, the authors are encouraged to reach out to the asset's 's creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
29

---

Justification: [NA]
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create anonymized URL or include anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification: [NA]
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with
human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then much detail as possible should be included in the main paper.
According to the NeurlPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country institution) were obtained?
Answer: [NA]
Justification: [NA]
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors adhere the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
30


Reviews:
[{"review_id": "k18w9ZFB4E", "reviewer": "Reviewer_vTvV", "summary": "In this paper, the theory of the Neural Tangent Kernel(NTK) in the case of solving nonlinear partial differential equations using PINNs is investigated in detail. In particular, it is shown that typical results of the NTK framework do not hold when the simple gradient descent method is employed due to the worse behavior of the Hessian matrix compared to the linear cases. In contrast, when second-order methods are employed, it is theoretically proven that the training of the neural networks is efficient.", "questions": "In the numerical experiments, it seems that not so large neural networks are employed. Is it expected that neural networks of this size behave like the results of the theory?", "limitations": "There seems to be no problem.", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 3, "contribution": 4, "strengths": "This paper theoretically investigates the behaviors of learning dynamics for PINNs, which are known to be difficult to train. This paper provides a theoretical guarantee of the effectiveness of second-order optimization methods for training PINNs for nonlinear partial differential equations. Im my opinion, this is a significant result, which may lead to applications of PINNs to practical problems that have been inapplicable due to the training difficulties.", "weaknesses": "My concern about this paper is in the increase of the computational complexity of second-order methods; however, this concern has already been discussed by the authors in the paper.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "First and foremost, we would like to express our sincere gratitude to the reviewer for the encouraging feedback and thoughtful comments on our paper.\n\nRegarding the reviewer's question, you are correct in noting that part of our theoretical framework is developed at the NTK level, i.e., for infinitely wide neural networks. However, our result is \"negative\" in the sense that, even in this idealized scenario, a first-order method might fail to perform effectively when using PINNs to solve nonlinear PDEs. This, as you mentioned, is due to the worse behavior of the Hessian matrix. Consequently, in the more practical scenario of a network with finite width, we cannot guarantee that PINNs trained with first-order methods can accurately solve any nonlinear PDE.\n\nMoreover, we would like to emphasize that Theorem 4.2, which addresses the convergence of second-order methods, is applicable even for networks with finite width. Therefore, while the NTK model serves as motivation to employ second-order methods, the significant findings regarding spectral bias and convergence do not rely on this idealized model.", "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": "", "extra_scores": {}}, {"review_id": "Q4oS5o3wh0", "reviewer": "Reviewer_4Qrf", "summary": "The paper studies the NTK of NNs trained on non-linear PDEs, showing that they exhibit different behaviours compared to standard analysis of NTKs. The paper then discusses the issue of spectral bias that arises from first-order methods, showing that they can be alleviated by the use of second-order methods.", "questions": "- Already in PINNs, there are many works that uses NTKs in loss function scaling [1], collocation point selection [2], analysis of PINN architectures [3], and more. How would the insights in the paper be able to address the points raised in these papers, and how would they affect these proposed methods?\n- Is Theorem 4.2 general enough to be applied to regular NNs as well? How does the result compare to existing theoretical works on second-order methods in NNs or general optimisation problems?\n\n[1] Wang et al. When and why PINNs fail to train: A neural tangent kernel perspective.\n\n[2] Lau et al. PINNACLE: PINN Adaptive ColLocation and Experimental points selection.\n\n[3] Wang et. al. On the eigenvector bias of Fourier feature networks: From regression to solving multi-scale PDEs with physics-informed neural networks.", "limitations": "Limitations suggested are of the LM algorithm which is adequate.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "- The paper presents an interesting analysis of a common tool in NNs and PINNs, and presents an explanation why second-order methods (which are already used in PINNs to some extent) works better than first-order methods.\n- The paper provides both theoretical and empirical justification for the various claims, and is well-organised in that manner.", "weaknesses": "- The sections could be a bit more coherent. For example, the paper brings up the properties of the NTK in the nonlinear PDE case, but then provides less link of these properties of how it affects the convergence in terms of the spectral biases. The LM algorithm is also brought up as a second-order optimisation method, however it may warrant more description as to why it is introduced or how it differs from existing second-order methods such as LBFGS.\n- Explicit mention of LM algorithm's runtime could be mentioned for completeness.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We would like to thank the reviewer for the positive feedback and the constructive comments on our paper. In particular, we acknowledge the second weakness highlighted by the reviewer. We plan to include this information in an improved version of our paper. At present, the reviewer can get a qualitative estimate of the runtime of LM with comparison to that of Adam and L-BFGS by referring to Figure 4(a). Now, let us address the reviewer's questions in the order they were asked:\n- **Q1** We would like to emphasize that our intention is not to suggest that our method should substitute other existing approaches, such as the ones mentioned by the reviewer. On the contrary, we believe that combining second-order methods with the enhancements to PINNs mentioned above can yield highly competitive results. This potent combination is demonstrated in Figure 3(a), where we used the LM algorithm together with the method proposed in [3] (as per the reviewer's nomenclature) resulting in fast and accurate convergence on a strongly spectrally biased PDE. On the right side of the same figure, we compare and combine our method with another PINN enhancement, namely curriculum training. Thus, we believe that the LM algorithm, in conjunction with the collocation point selection presented in [2], could also result in positive improvements, especially for engineering use cases where traditional solver often relies on adaptive meshes. We appreciate the reviewer for highlighting this work. Regarding reference [1], the issue of unbalanced loss components is addressed through loss scaling. However, in our work, this issue is implicitly tackled with second-order methods according to the result presented in Theorem 4.2. Indeed, the presence of ones in the convergence matrix $D$ implies, in particular, that the various loss components are balanced during training, producing excellent results when compared to loss scaling. We already have a comparison with [1] in the Appendix, Figure 7, where we referred to the method in [1] as \"loss balancing\" instead of \"loss scaling.\"\n\n- **Q2** We believe that Theorem 4.2 can also be applied to regular NNs, as there are no specific assumptions on the loss function other than the requirement that the Jacobian $J$ has to be full-rank. However, Theorem 4.2 specifically addresses the issue of spectral bias in PINNs. This issue is a well-known problem for PINNs (even for linear PDEs), but we are not aware of any correlation between spectral bias and worse performance in regular NNs' tasks. Nevertheless, we would like to highlight that there are several papers, such as [4, 5, 6], where general second-order methods have been successfully employed to train regular NNs. Hence, it might be beneficial to extend our analysis to these cases.\n\n**[4]** Z. Yao, A. Gholami, K. Keutzer and M. W. Mahoney, \"PyHessian: Neural Networks Through the Lens of the Hessian,\" 2020 IEEE International Conference on Big Data.\n\n**[5]** Liu, G. H., Chen, T., and Theodorou, E. (2021). Second-order neural ode optimizer. Advances in Neural Information Processing Systems.\n\n**[6]** Vinyals, O., and Povey, D. (2012). Krylov subspace descent for deep learning. In Artificial intelligence and statistics, PMLR.", "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": "", "extra_scores": {}}, {"review_id": "p406cWaioO", "reviewer": "Reviewer_ocfJ", "summary": "This paper studies the training dynamics of PINNs, especially for the nonlinear PDEs. The authors find that the previous recognized NTK viewpoint is not applicable to nonlinear PDEs, although it holds for linear PDEs. Therefore, the global convergence of gradient descent on nonlinear PDEs may not be guaranteed. Moreover, the imbalance of singular values of Gram matrices, which also occurs in linear PDEs, results in slow convergence. To address this issue of spectral bias, the paper suggests using second-order methods for parameters update.  Experimental results show that LM Newton’s method can achieve lower training loss compared to Adam and L-BFGS.", "questions": "I have the following questions:  \n(1) In the previous work (e.g., https://www.sciencedirect.com/science/article/pii/S016727892300341X), it was found that NTK approach fails in solving some PDEs. Besides linear PDEs, it is still possible that NTK works for some nonlinear PDEs. However, in your paper, you exclude the special situation. It seems to be more restrictive than the published work.  \n(2) Although the second-order methods enjoy very good theoretical properties, deep learning community typically prefers first-order methods due to their computational efficiency. In practice, Newton’s method is more computational expensive even with some inexact technique (e.g., Krylob subspace, conjugate gradient, and LBFGS). However, your theorem of global convergence (Theorem 4.2) does not apply to these inexact methods. Based on your result (which I believe built upon some previous works), can you extend your findings to more practical inexact Newton’s method (e.g., replacing Hessian inverse by Krylov subspace method or Quasi-Newton’s method of BFGS or LBFGS)?  \n(3) In your experimental results, the LM method performs significantly better than LBFGS. This is surprising, as LBFGS or BFGS asymptotically approaches the exact Hessian under certain conditions. Moreover, in case where the totally number of grid points are large and the batch is small in a stochastic setting, inexact and Quasi Newtons method (e.g., Hessian averaging methods including BFGS and LBFGS) should hold the asymptotic convergence. Therefore, intuitively, I would expect LBFGS to perform at least comparably with your LM method. Am I correct? Can you explain why LBFGS performs poorly according to your results, although it can approximate Hessian and extract Hessian information?", "limitations": "Yes", "rating": 7, "confidence": 5, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "The paper is well-written and provides both theoretical and empirical analyses. The failure of NTK approach on nonlinear PDEs has not been highlighted in previous studies of the global convergence of training PINNs with gradient descent. The second-order method is essential in smooth problems and deep learning training. The authors highlight the effectiveness of Newton’s method in balancing singular values of the training dynamics.", "weaknesses": "(1) the failure of NTK is similarly investigated in some previous works (e.g., https://www.sciencedirect.com/science/article/pii/S016727892300341X), where they found that the Gram matrix does not consistently converge in some cases. Therefore, the observation seems to be not novel. They also pointed out that NTK may also work for some nonlinear PDEs.  \n(2) the second-order methods for modifying the singular values of Gram matrix of training dynamics are also not new. \n(3) the Hessian inverse is quite expensive in practice.  \n(4) regularizing the Hessian may also be intractable in practice, for high dimensional problems, e.g., training PINNs.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We would like to extend our gratitude to the reviewer for their thoughtful feedback and for bringing this reference to our attention. Now, we will proceed to address the specific questions and concerns raised by the reviewer.\n\n- **Q1** In the paper referenced by the reviewer, the convergence of the Gram matrix at initialization (Theorem 2.2) holds in the non-homogeneous case only for $s\\geq1$.  The authors state that for $s=\\frac{1}{2}$, it is impossible to guarantee universal convergence. This is why our results may seem more restrictive at first glance, but this is due to our focus on the typical NTK scaling $N^{-\\frac{1}{2}}$, which allows us to obtain meaningful characterization of the training dynamics using kernel analysis. Furthermore, we explicitly characterize the law of the NTK's limit. In addition to being applicable to any nonlinear PDE, this result enables the leveraging of random matrix theory tools to prove its positive definiteness at initialization, which is a longstanding conjecture in the research community. Hence, we would like to highlight that our results do not contradict those in the referenced paper, even if there seem to be inconsistencies that we are now going to address.\n\n    - Regarding the constancy of the Gram matrix during training (Theorem 2.3), we have noticed a problem in the proof of Lemma C.4. In equations 63-68, the authors establish a bound for a single parameter as $ O(1/N^{2s-1/2})$. However, there might be an inconsistency in equation 69, where the same bound is applied to the full vector, which has $O(N)$ components. This discrepancy implies that the bound should be $O(\\sqrt{N} \\cdot 1/N^{2s-1/2}) = O(1/N^{2s-1})$. Consequently, their theorem would hold for $ s > 1/2 $ instead of $s > 1/4$. Hence, this theorem and our Proposition 3.5 do not contradict each other, since we study $s=1/2$.\n\n    - Finally, there seem to be differences regarding the numerical experiments in the reference. However, we noticed that the equations studied in Figures 3 and 5 have a consistent linear part, while our PDE is fully nonlinear. The variation of $K(t)$  during training decreases as $N$ grows because the variation of its linear components decreases, according to our theory. However, the nonlinear components still evolve, and for even larger $N$ (we use a number of neurons several orders of magnitude greater), we observe the plateauing of the plots at a value strictly greater than 0 (note that in their Figure 5, the y-axis starts from 10).\n\n- **Q2 & Q3** To extend our findings to a more practical inexact Newton's method, we can indeed utilize the approach you suggested in Q3, specifically the asymptotic convergence properties of LBFGS and other methods. While this approach is theoretically sound, we must note the following. The theoretical guarantees on the speed of convergence of quasi-Newton methods to exact Newton methods, i.e. their matrix approximation ability, depend on the minimum eigenvalue of the Hessian ([1], Theorem 6). However, as discussed in our paper, the Hessian in PINNs is generally very poorly conditioned. Consequently, quasi-Newton methods may require a practically infinite number of training steps to converge to the true Hessian's inverse and, therefore, to start training higher modes. The aforementioned reasons would explain the intermediate performance of LBFGS, which lies between that of first-order and exact second-order methods.\n\n**[1]** Lin, D., Ye, H., and Zhang, Z. (2022). Explicit convergence rates of greedy and random quasi-Newton methods. Journal of Machine Learning Research.", "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": "", "extra_scores": {}}, {"review_id": "auSuQHbav1", "reviewer": "Reviewer_XYBU", "summary": "The authors analyze the NTK perspective for PINNs for non-linear PDEs. Previous considerations derived for linear PDEs fall short for non-linear ones and the authors attribute this difference to the non-vanishing Hessian term. Therefore, they suggest to use second order methods and show to be able to converge faster in one experiment. Second order method are shown to be useful also for linear PDEs as they alleviate the spectral bias of the NTK.", "questions": "- I would like the authors to elaborate on the applicability of the second order method. In particular, I found the explanation in lines 218-222 a bit hand-wavy. It seems like the hypothesis of Theorem 4.2 (e.g. J(t) being full-rank) are very hard to check in practice\n- No guidance or intuition is provided for practitioners on when it might be convenient to use the second order method, given that the method is expensive (the authors only mention general well-known considerations for second order methods in the limitations paragraph)", "limitations": "Limitations are discussed. However, as mentioned in the question section, I think the authors should discuss in more details the scalability of the proposed second order approach and on when it might be convenient to use it.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "- the paper is well written and flows very nicely\n- I found the analysis about the NTK in non-linear regime elegant and the experiments clearly support the findings", "weaknesses": "- the different behaviour of the NTK for non-linear PDEs is not very surprising\n- it is also not very surprising that second order method can work better, but in practice they come with significant shortcomings. Even though the method is shown to be faster in one experiment, no rate is derived so it might not be true in general.\n\nHowever, the applicability of Theorem 4.2 is just sketched (lines 218-222)", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We sincerely thank the reviewer for the kind words regarding the structure and analysis in our paper. While we agree that some of the results might be expected, most of the existing literature focuses on the linear regime and only conjectures what might happen in the nonlinear case. To our knowledge, our paper is the first to characterize the NTK dynamics in the nonlinear case. This characterization is crucial for understanding that the poor behavior of the Hessian can lead to slow and unreliable training of PINNs.\n\n- **[Q1]** We understand the reviewer's concerns regarding the applicability of Theorem 4.2. However, proving the full-rankness of $J$ in the nonlinear case is an open problem that has been widely discussed and remains unresolved. Our contribution is to highlight that the NTK is stochastic in the limit, and for this reason, if one wants to check the full-rankness of $J$, probabilistic methods need to be used. We leave this task to future works. More generally, we do not assert that $J$ is necessarily always full-rank. The remark in lines 218-222 is meant to highlight that the well-conditioning of $J$ is not crucial for a better training when using a second-order method. Indeed, one of the advantages of using second-order methods is that the diagonal of the matrix $D$, even when $J$ is ill-conditioned, only contains ones and zeros. \n    While a first-order method would yield,  in addition to zeros, values that might be extremely small (i.e. on the order of $10^{-14}$ in figure 2(b)). Regardless, this assumption is weaker than the one needed for theoretical convergence guarantees of first-order methods.\n\n - **[Q2]** We thank the reviewer for the interesting question, and we believe that adding a paragraph on this would improve the quality of our paper. In an improved version, we plan to add an intuition on when to use second-order methods. In practice, it may be convenient to use second-order methods whenever there are high frequencies in the solution, when the application requires high accuracy, and when the PDE we aim to solve is nonlinear. One of the goals of our paper is to emphasize that, despite their shortcomings, second-order method are in general a natural choice in the field of PDE solution.", "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": "", "extra_scores": {}}]

Author rebuttal:

Meta-review: The manuscript discusses the NTK regime of PINNs trained on nonlinear PDEs and shows that several results from linear PDEs cannot be shown to carry over to this setting. Specifically, the NTK is not constant during training and the Hessian does not vanish, explaining the superior performance of second-order optimization methods.

All reviewers agree that the topic is relevant and that the presentation of the paper is excellent. The experimental evidence supports the presented theory. Despite some reviewers acknowledging that NTK analyses of PINNs and insights into the performance of second-order methods for PINN training are not new, all reviewers agree that the paper should be accepted.
Decision: Accept (poster)