PINNs with Learnable Quadrature
Authors: Sourav Pal
Kamyar Azizzadenesheli
Vikas Singh
Venue: NeurIPS 2025 poster ()
Keywords: PINNs
Primary area: machine_learning_for_sciences

Abstract: The growing body of work on Physics-Informed Neural Networks (PINNs) seeks to use machine learning strategies to improve methods for solution discovery of Partial Differential Equations (PDEs). While classical solvers may remain the preferred tool of choice in the short-term, PINNs can be viewed as complementary. The expectation is that in some specific use cases, they can be effective, standalone. A key step in training PINNs is selecting domain points for loss evaluation, where Monte Carlo sampling remains the dominant but often suboptimal in low dimension settings, common in physics. We leverage recent advances in asymptotic expansions of quadrature nodes and weights (for weight functions belonging to the modified Gauss-Jacobi family) together with suitable adjustments for parameterization towards a data-driven framework for learnable quadrature rules. A direct benefit is a performance improvement of PINNs, relative to existing alternatives, on a wide range of problems studied in the literature. Beyond finding a standard solution for an instance of a single PDE, our construction enables learning rules to predict solutions for a given family of PDEs via hyper-networks, a useful capability for PINNs.

PINNs with Learnable Quadrature
Souray Pal UW-Madison spal9ewisc.edi
Kamyar Azizzadenesheli NVIDIA Corporation kaazizzad@gmail. com
Vikas Singh UW-Madison
vsingh@biostat. wisc. edu
Abstract
The growing body of work on Physics-Informed Neural Networks (PINNs) seeks to use machine learning strategies to improve methods for solution discovery of Partial Differential Equations (PDEs). While classical solvers may remain the preferred tool of choice in the short-term, PINNs can be viewed as complementary. The expectation is that in some specific use cases, they can effective, standalone. A key step in training PINNs is selecting domain points for loss evaluation, where Monte Carlo sampling remains the dominant but often suboptimal low dimen- sion settings, common in physics. We leverage recent advances in asymptotic expansions of quadrature nodes and weights (for weight functions belonging to the modified Gauss-Jacobe family) together with suitable adjustments for parame- terization towards a data-driven framework for learnable quadrature rules. A direct benefit is performance improvement of PINNs, relative existing alternatives, on a wide range of problems studied in the literature. Beyond finding a standard solution for an instance of single PDE, our construction enables learning rules predict solutions for given family of PDEs via hyper-networks, useful capability for PINNs.
Introduction
Differential equations widely used across science, providing framework for modeling/analyzing diverse physical dynamics. Most real-world settings lead differential equations where analytical
solutions are not possible. but research over decades has led to mature set of numerical methods
[Ames] 2014 Trefethen and Baul [0222] Specifically, growing body of work has identified novel architectures, marrying differential equation solvers with deep learning and these formulations offer new capabilities. For example, one now has access to completely data-driven approaches [L] al.) 2020b Kovachki et al.] 2021] which use observational data to estimate the operator for a PDE. We also have a new class of differential equation solvers that exploit neural networks to encode physical laws (PINNs) [Raissi et al., 2019 Kharazmi al., 2019 without any observational data. PINNs give solutions which are mesh-independeer differentiation, and handle non-line well. As such, while classical solvers may remain the default choice in general, developments in PINN suggest that they can be an alternative in certain scenarios.
easier to implement due to automatic
Roughly speaking, the aforementioned line of work [Karniadakis al. 2021 discussed in more detail later in 7 can be broadly classified under three main threads: (a) Solving PDEs using neural networks (PINNs) [Raissi al., [019]; (b) PDE discovery (e.g., symbolic regression) [Holt et al. 2023 d' A'coolli et al. [2023] and (c) operator learning (e.g.. Fourier Neural Operator) [Azizzadenesheli 2024]. This classification is loosely based on the amount of data physics used solve/inform the Torward/inverse problem, [Boullé and Townsend 2023 Our work falls under the first category, where wish solve single (or multiple instances) of PDE with PINNs. Here, we seek to identify how learning mechanism can deliver efficiency gains solely based on the shared structure and knowledge of physics, without the use of any labeled data.
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

PDEs and Quadrature.
Consider the second-order PDE,
8u Ouu
dx +
log(2)) sin(y) + f(z y)v"
oy²
One way to find u to integrate both sides with a test function v(z, y) resulting in:
(
pu dx2 0y²
vaz yddrdy = = / (log(z) sin(y) + fzz w)y") v( ydddy
Doing so enables the use of numerical methods which build-up sums that converge to the integral's true value. This has the benefit of easing regularity conditions on u. It solves the original problem in a weighted sense. A quadrature method (Golub and Welsch 1969] chooses evaluation points (nodes) and corresponding weights to minimize approximation error. These points/nodes may be constant step/uniform or adaptive, i... either fixed or adaptive quadrature rules for estimating the integral.
Challenges in computing an integral. In many applications from fluid dynamics (turbulent flow) [Kutz] 2017] to radiation treatment planning (fluence calculation at tissue interfaces) [Lou et al.
Beckham et al., 2002] to materials science (fracture mechanics) [Aliabadi and Rooke 1991 Rice and Tracey] 1973], the associated data involves irregular behavior including singularities. Uniformly splitting the domain of integration into equal sub-domains, in many cases, is insuf- ficient owing to the singularity associated with the integrand. Further, even sophisticated par- titioning scheme, runs into difficulties in the multi-dimensional case. A common solution is to use some variant of Monte Carlo sampling. In higher dimensions, we have no choice but to sample at large and expect the estimated so- lution to converge to the true solution, given enough runtime. In lower dimensions, Monte Carlo sampling is sub-optimal [Lu et al., (22T] and several strategies to improve the speed and accuracy of the integral computation are known, the prominent ones being some variant of adap- tive quadrature scheme, essentially choosing an adaptive grid dependent on the integrand.
Generated
Integral Exiluated Quadrature
findss
1 XX
o
Epoch
1x675 X
Epoch
0
Leamable
i
Purametor:
Roy
(o000o)
in cas
Epach
Figure 1: Relevance of learnable quadratures. Given a fixed number of quadrature points, one can update a learnable module based how good/bad the numerical approximation of the integral is.
Main Idaa. We propose a learnable quadrature scheme (Fig I utilizing a rich theory of orthogonal polynomials (OP) and asymptotic expansions. Consider solving PDE with PINNS either in its strong form (PINN) [Karniadakis et al. 2021] or weakform (VPINN, hp- VPINN) [Kharazmi et al. 2019 2021]. In both cases, the goal is to either (a) determine which points in the domain evaluate the function on (this is the strong form) or (b) determine which test functions to use (for evaluation of the weak form). We tie these choices to the roots of OP w.r.t. the modified
weight
functions. Next, we use recent advances in asymptotic expansions of quadrature nodes and weights (Opsomer and Huybrechs 2023 to achieve a scheme to compute these efficiently.
For solving PDE with PINNs, we describe two separate learning modules. First is the solution function for the PDE parameterized using a neural network (standard in PINNs). The second is a parameterized weight function which induces a family of OP. Our parameterization of the weight function together with asymptotic expansions takes advantage of parallel compute on GPUs to generate a massive number (millions) of quadrature nodes and weights in near constant time. This provides an alternative to Monte Carlo sampling for low-dimensional problems. Our model can solve most PDEs from the PINN literature achieving better performance than existing adaptive/non-adee -ddaptive sampling schemes. Our parameterization of the weight function also allows learning a quadrature predictor for a family of PDEs (with shared structure), which may be of independent interest.
Contributions:
2 Preliminaries
We review some key concepts that will be useful throughout.
Orthogonal Polynomials (OP). A sequence of real-valued polynomials pa (a(2), Pi (t))) (2))),..,... where each pm (x) is a polynomial of degree are orthogonal [Olver et al. 2020 with respect to

---

continuous and son-negative weight function w(x) defined in the interval (a. b) if:
(pm Padu: = -f. [_ Pm (x) = so if m &
(3)
hn
if m n
where h, a normalization constant. In fact, if h = 1 for all n, then the family is orthonormal.
Modified Jauss-Jacobi Weight functions. We use weight functions from the modified Gauss-Jacobi
family Opsomer and Huybrechs 2023) to induce our family of OP 33) which have the form:
wo(2) (1 (1 + h(r); x x
where ax B -1. h(x) is the modifier over the standard Gauss-Jacobi weight function: w(xx
(1 (1+)) .t.::.""TTT The only restriction on h((z) is that it should be strictly positive analytic function.
Cauchy Residue Theorem. The Cauchy Residue Theorem (Stein and Shakarchi] 2010) is powerful tool to compute line integrals of analytic functions over closed curves. Let J' be a function that is holomorphic on simply connected open subset of the complex plane, except possibly at a finite set of points a1, (callee poles) and y be positively oriented simple closed curve, then we have:
f(=)d== 2Ti) Σ Res(S,ak);
Res(J, ax) = 2xd $ f(z)d=
k=1
where the quantity Res(f, (&) the complex residue the pole ag- 7k is a positively oriented simple closed curve around the pole ag not including other singularities. This is the general formula; in our case, we will only deal with simple poles, so it simplifies (Stein and Shakarchi 2010).
Res(f.ad ax) lim (z ax)f(z)
(6)
3 Strong and Weak Forms
We use u to denote the solution for given PDE. Since, u is learned (parameterized using neural networks), it is commonly called the trial function. In its most generic form, a PDE is:
Lu= f. in Ω u= g. in aΩ
where, denotes the differential operator acting on u. denotes the non-homogeneity.. Also, 9: 0.91,... denote functions corresponding the initial and/or boundary conditions needed and Ω gives the domain definition of the PDE and Ol denotes its boundary.
Solving PDEs. We examine the canonical form of a second order elliptic PDE, the Poisson's 's equation (in dimensions) as a running example. In 2D-Poisson's 's equation, the operator L is the Laplace operator, 2 For u Ω + R, where Ω C R22 is the domain of interest this is given by:
Vuu(z, v) f (a, y), (z, v) E Ω u(z y) 0, (z, 3) E an
where f is forcing function. We consider Ω to be [-1, 1| [-1, 1] and 82 as the natural boundary.
Before getting into details of our proposal, summarize two solution approaches: via the strong and weak form the PDE respectively. This will help us see how certain choices (such as quadrature rule and collocation points) will be key our parameterization and thereby, learning.
Strong Form. Solving the PDE in its strong form 7 or (8) is equivalent to asking that the equations in are satisfied exactly at several points along the domain Ω and boundary 00. This can be done by sampling a large number of points distributed uniformly (Monte Carlo). While this approach is reasonable for simple cases, in more complex settings where the solution is not smooth, uniform sampling may lead to reduced approximation quality Lu et al. (2021].
Weak Form. While the strong form enforces point-wise exactness, we may only want the property to hold in a "weighted" sense for the entire function. This yields the weak form, which involves integration with test function. For the 2D-Poisson equation, using test function o(x; y) and integrating over the domain, we have
11 8
Fu v drdy = 1 ( v drdy
(9)
3

---

Solving the PDE in its weak form moves the previous choice of points to choice of test functions. The idea is to use family of test functions based on the problem at hand and different methods emerge from this choice. These methods are called Galerkin methods. If one decides that the test functions in the weak form are Dirac-delta functions, then the Galerkin method reduces to collocation and the weak form becomes the strong form involving differential equations, which then need to be satisfied at the collocation points. When we discuss the strong form, we will interchangeably use "the choice of collocation points" and "the choice of test functions".
4 How to learn Quadrature Rules?
The above discussion underscored the importance of choice of test functions in the weak form or the choice of collocation points for solving the PDE in its strong form. One may ask: can the underlying physics inform these choices? To do so, we leverage the rich theory of OP.
Learning the weight function:
We consider weight functions to be continuous and positive functions defined in some interval I. Each such weight function w(x) induces a family of OP given by (3). By learning the weight function for a set of OP, we want to enable a learnable (or adaptive) quadrature. For this to be practical, we want to compute these efficiently. Our method exploits asymptotic expansions of quadrature nodes for efficiency which are most complete for modified Gauss-Jacobi type weight function. So, we consider the modified Gauss-Jacobi form 4 and parameterize the modifier (<x) in 4 using a neural network with parameters O. This keeps the construction simple but offers many nice benefits we will see soon. So, our learnable weight function has the form:
(x) = (1 - =1 + ho());
x E
(10)
In 10 a and B can also be parameterized/learc but in our experiments, we find that only learning (z( suffices. We consider the interval of the weight function I to be the same as Ω. the domain of the PDE. Next, we will see how this learnable weight function cleanly ties to the choice of test functions for the weak form and the collocation points in the strong form.
Remark 4.1. The method can, in principle, be extended to higher dimensions by treating each dimension independently and stacking the sampled values. Realizing such extensions in practice, however, would necessitate more elaborate constructions such as tensor- tensor-product formulations or sparse grid techniques [Gareke et al. 2006] to appropriately handle increased complexity.
Relation to Solving PDEs in Weak Form: For the weak form, consider weighted integral using our weight function Wy as the test function v, on both sides of (9).
[ʃ--v% we drdy 2
[ f ww dxdy
(11)
We must compute both sides of the integral efficiently, especially for weight functions that produce the largest error in the integral equality. Our choice of weight functions crucially helps this computation.
Use of Orthogonal Polynomial & Quadrature Rule: is well-known that a n-point Gaussian quadrature rule can be constructed to yield a very good approximation to the integral of 2n I degree polynomial, multiplied with the corresponding weight function. For example, if we use the standard (non-modified) Gauss-Jacobi weight function:
Consider a one-dimensional integral. It
w(x) = (1 - x) x²) + ) x E [-1,1] a,3> 1-3>-1
(12)
Then, the n-th order approximation to the integral given by:
[ ((0)(1=))( + x) dr = Σ wif(x;)
(13)
1=1
where 2; denotes the i-th node and w; the corresponding weight of the point Gauss-Jacobi quadrature. Here, fxx is a smooth function on []. The nodes I; are in fact the roots of the n-th degree Jacobi polynomial, which form family of OP w.r.t. the weight function in (12). Our learnable weight function not only provides a data- deta-dependent choice of test functions associated with the weak form but also induces learnable quadrature rule compute the integral for the weak form!
Relation to Solving PDFs in Strong Form: Solving a PDE in its strong form means enforcing the relationship at several collocation points. A popular method in this category is orthogonal collocation

---

where the collocation points used are roots of orthogonal polynomials. Thus, with our choice of learnable weight function, wy (which induces a family of OP (3) also provides data-dependent method to sample collocation points to solve a PDE in its strong form.
(Young/ 2019
Remark 4.2. Standard quadratures with fixed nodes and weights are easy to implement in PINNs, but their lack of adaptivity leads to suboptimal residual minimization during training. Being data independent they cannot be conveniently utilized for solving a family of PDEs via neural networks.
5 Learning Quadrature Rules Efficiently
Solving PDEs via strong form, means that must identify collocation points to evaluate the PDE. For the weak form, solving amounts to numerically approximating the integral via a quadrature rule. While $ described the relation between learnable weight function with weak and strong form solution methods, it is not obvious how this can be done efficiently. Let lay out
3(2)
Outer Region
Bulk
the task. In strong form 7 (). we must find the roots of corresponding OP. These will serve as the collocation points. In weak form, the learnable weight function(s) itself can act as test function(s). But we must find the quadrature rule (roots of OP and corresponding quadrature weights) to evaluate (TT) efficiently. Thus, while the weight function is relevant by itsllf, our main interest is in the roots of OP induced by the weight function. Next, we show how use of recent advances in asymptotic expansions of OP and their roots suggest an efficient instantiation.
R(z)
Left
Right
Figure 2: Different regions of com- plex plane for asymptotic expansion.
5.1 Instantiating Asymptotic Expansions
There is a mature literature for fast computation of quadrature nodes and weights corresponding to weight funetions of OP [Townsend 2015). Over the last few years, this race is dominated by asymptotic expansions [Bogaert] 2014] Townssed et al., 2016]. Very recently, Opsomer and Huybrechs) [2033] proposed asymptotic expansions for generalized (i.e., modified) versions of canonical weight functions including Gauss-Jacobi and Herss- type.
Division of Complex Plane. We briefly present asymptotic expansions of nodes (roots of OP)
and weights of quadrature rule for the modified Jacobi-type weight function ( [Osomer and
Huybrechs 2023 Opsomer 2018]. The details are not crucial, but useful to appreciate our choice of parameterizations. The reader can check Appendix and Opsomer and Huybrechs/ (2023) for more details on the expansions. The complex plane is divided into four regions (Fig each with a different expansion. These are: the lens covering a bulk of the interval ( (-1, 1), the two end-points referred to as left and right disks and everything else is the outer region. Let T be shorthand notation for (2n + a + B + 1) and n refers to the polynomial degree.
Left endpoint. Truncated asymptotic expansions of nodes (xx) and weights (wx) near - are: 2j²,2 2j², 02 -3a2 +1)
Wkk w(Tx)
(o)
3(T + do)
do}
(14)
where w((x) is the value of the weight function at node Fy. and co and do are expansion coefficients, described shortly. The expansion uses both (a) the zeros of Bessel functions of order / denoted as ja.k (k-th zero) and the Bessel function of order (8 1) denoted by J-31-
Right endpoint. For the right end-point, we interchange and 9 and use fu-2x) instead of h(x))
Bulk region. For expansions in the bulk region, we need find the leading order term, tg by solving:
4k + 2a +
1
log ()
π 4k + 2a 2/ +2
arccos(tx) + r
2ni Jn VE2 - 1(E (x)
(15)
Using t. the truncated asymptotic expansions of nodes and relative weights in the bulk region are:
2a²
-,,2 + (2a²2 + 2,8² - 1)
Ik 1/ +
2[ + To 2
(16)
Wkk
1 2
21
12) 2T0 +
2
w(x))
r
5

---

where To, 71 are also expansion coefficients.
Summary. In [14 (6), the values of a. correspond the used in the modified Gauss-Jacobi
weight function (4) The value of n determines the degree of the OP (x) from the family (3)
whose roots want compute. Finally, given n, k 2 .,,,,, ,.) polynomial Pa (z), which is guaranteed to exist and be unique in the interval of definition.
corresponds to the k-th root of
Expansion Coefficients. In our description above, we used several coefficients: co, do and To and Ti. While more details are in Appendix 9 a synopsis is that the coefficients cg2 dg stem from series expansion of the modulation function h(z) in (4) (or the parameterizzd version in (10) around =1. The coefficients Ti are the series coefficients resulting from the expansion of the contour integral in 115) around the leading order tg of the k-th root of the orthogonal polynomial. The above formulas involve computation contour integrals, root finding, and series expansions. Computing all these terms exactly within learnable module will be challenging. We will next perform some simplifications so that the model is amenable to learning.
5.2 Simplifications, Assumptions and Implementation
We outline key assumptions, simplifications and implementation details for an efficient instantiation.
Simple Poles. Finding the leading order of the k-th root, tg. requires solving a contour integral (15) We use the fact that roots in the bulk region are real and lie in (-1, 1), and assume the integrand has only simple poles around 4,, allowing us to compute residues via (6) as follows:
log(h(E))
((2()((2))
lim (E (x)) 2-
(17)
VI (3 (x)
V1
Assuming real roots and simple poles, the Cauchy Residue Theorem (5) simplifies
15 as
4h + 2a 3
arccos((t)) log log()((s))
(18)
4k + 2a + 2, +
r
Root finding and Implicit Function Theorem. While we avoided computing the contour integral, must solve (18) for fg- Since l, corresponds the leading order of the root of OP, it must exist in (-1,1). Thus, we use bisection method to find ty via root finding of:
4k + 2a +
-arccos(t2) + log(h(1z))
F(xx) = 4k + 2a + 23 +2
(19)
We use automatic implicit differentiation from [Bondel et al. 2022 which uses auto-diff of F(le) and the implicit function theorem to automatically differentiate through the bisection method.
Benefits of simplifications & parameterization. We use fully connected layers followed by tanh
non-lincarities to parameterize the solution function u. the modulating function / as well as to predict
expansion coefficients of the nodes and weights, cp, dg. d1, To. T1. etc. for bulk and edge regions. The sim- plifications above offer multiple benefits. First, we avoid computing contour integrals within differentiable learning framework. Second, we are able to compute all nodes and weights in parallel thereby making the process very efficient even for very large number (millions) of nodes. noting that the exact procedure to compute the nodes (beyond the leading order term) in Section 4.2 of [Opsomer and Huybrechs, 2023) has a linear time complexity due to the several re-substitutions involved to find the coefficients. Empirically, verify that the distribution of nodes and weights from the simplifications and parameterization choices does not hurt the distribution of quadrature which converges to the expected distribution.
Pn Pn+1
x Figure 3: Interlaced red and blue dots on x axis correspond to the roots of polynomial Px+1 (z) (degree n + 1) and Pr. (z) (degree r1) respectively.
worth
Interlacing of roots of orthogonal polynomial. naive application of quadrature nodes and weights for the OP is insufficient in several cases. This is because with a high degree polynomial (and so, large number of nodes), within a few epochs, jointly training the solution and quadrature functions leads to over-fitting Interestingly, this can be solved via very useful property of the family of OP, namely interlacing of roots. Given the roots of an OP of degree n + 1, the roots of the OP of degree n that belong to the same family are interlaced within the roots of Pn+1 shown in Fig. 3

---

Thus, by utilizing quadrature nodes and weights stemming from varying degree of OP (all of whom correspond to the same weight function being learned), we introduce the desired stochasticity to prevent wer-fitting.
Remark 5.1. In [Mishra and Molinaro| 2023], Theorem 2.6 decomposes the generalization error into training error ET and a quadrature term C Quaad-N -a/p, where Cquad depends on the sampling distribution. LearnQuad directly targets Cound by learning a sampling distribution that concentrates collocation points where the residual is high, without modifying the PINN architecture or PDE model.
To ensure that our learned weight function is positive, we use softplus activation the last layer of the network for (x))). Further, the presence of log in (19) can lead to
Implementation Details:
vanishing gradients, which we fixed by adding a small amount of noise (order of e -0)
). We use the
standard loss function used in PINN literature [Karniadakis et al. 2021 Cai et al. 2021 Additionally
to avoid invalid quadrature rules from numerical issues, we use:
w L uo (z(x)) + w - 2)
(20)
where the first term in 20 enforces the necessary condition that the sum of quadra- ture weights is equal to the integral of the weight function over the domain of def- inition. The second term in (20) discourages the quadrature weights from becom- ing too small. We note that the quadrature weights w; are computed numerically from the learned weight function using asymptotic expansions and not directly learned.
Instead of integrating a PDE solver into the learning framework, the PINN loss is chosen due to challenges in differentiating through ex- plicit/implicit solvers. including high memory usage, numerical instability, and computational overhead. Using PINN loss enables end-to-end learning and facilitates the development of learn- 2: for i = 1 to T do
Algorithm 1 Training for single PDE
1: Input: PDE parameter μC; #epoch: T, Learn- able model Ug, Wy; Loss L incorporates oper- ator L, inhomogeneous term, initial/boundary condition; regularizer bw-
able quadrature rules that generalize across dif-
3: Use $ to get quadrature nodes {a;}
ferent PDE instances. Our overall procedure 4: Use solution function uh on {21}
has two trainable components: one the learn- able quadrature module (I (LearnQuad) and the other is the learnable solution function for the given PDE. These can be trained jointly using the loss described above either to simultaneously decrease it or in min-max fashion where the quadrature module tries to provide hard-to-approximan function points in the domain. Empirically, we do not find a large difference in this specific choice of optimization. We provide pseudo-code for training PINN using LearnQuad in Algorithm [ As demonstrated in Appendix I0.1.8 variations of a and 33 in the modified Gauss-Jacobi weight function only lead to minor variations in performance. So, fix o, 3 to moderate values (typically in the range 1-2) for all experiments. This choice is sufficient to achieve the desired performance, without requiring extensive hyper-parameter sweeps. Our code is available at https//github.com/siiggg/ group/leam-quad
5: Loss:/ = L(u (21)) + tw (we)
6: Update Up and wy based on /
7: end for
8: Output: Learned models @ and e
Table 1: LearnQuad has lowest L2 relative error (mean standard deviation) in all PDEs.
6 Experimental Evaluations
PDE
Burgers" 2000
Allen-Cahn 1000
Wave 2000
points
We demonstrate the effectiveness of our proposed framework involving the learn- able quadrature module (LearnQuad) in solving PDEs via PINNs next. First we
Grid
0.12 1 0.04
0.88 1 0.06
0.42 + 0.09
Non-adaptive
Random
0.13 1 0.03
0.32 1 0.14
0.48 1 0.07
LHS
0.18 - 0.15
0.32 0.04
0.61 = 0.13
compare the empirical performance in solv-
Halton
0.06 = 0.02
0.18
0.46 = 0.06
Hammersley
0.07 = 0.05
0.17
0.31 = 0.09
ing single PDEs which demonstrate their effectiveness when used with PINNs. We emphasize that the use of LearnQuad
Sobol
0.08 + 0.03
0.49 + 0.09
Random-R
1.69 + 1.67
0.55
0.72 0.90
improves the performance of PINN' and helps close the gap with numerical meth-
Adaptive
RAR-G
0.12 = 0.04
0.19
0.81 - 0.11
RAD
0.02 = 0.00
0.08
+ 0.06
0.09 = 0.04
ods based solvers. Thereafter, we describe how the use of htyper-netwole can enable
RAR-D
0.03 + 0.01
0.09
1 0.03
0.29 - 0.04
LearnQuad
0.003 + 0.0020.03 + 0.0080.005
0.0006
7

---

Table 2: Ly error with 1000 collocation points. LearnQuad yields best results in 4 of 5 cases.
PDE
Convection (3
30)
Convection (3
50)
Allen Cahn
Epochs.
100k
300k
150k
300k
200k
PINN (fixed)
107.5 + 10.9%
107.5 + 10.7%
108.5 6.38%
108.7 + 6.59%
69.4 + 4.02%
PINN (dynamic)
2.81 + 1.45%
+
24.2 23.2%
56.9 + 9.08%
0.77 = 0.06%
Curr Reg
63.2 9.89%
1665+1.44%
31.5 16.6% 7.44%
CPINN (fixed)
138.8 1 11.0%
138.81
10.5% 106.5
106.5 10.5%
48.7 + 19.6%
CPINN (dynamic)
52.2 1: 43.6%
23.8
19.0
73.2133.6%
1.5 = 0.75%
RAR-G
1.
1.41%
1..71.777%
13.1428..9%
25.1 + 23.2%
RAD
1.90%
16.0+155%%
11.9%
0.78 0.05%
RAR-D
4.28%
=
5.96%
1.58%
= 51.6 0.41%
L
66.6 + 2.35%
27.9%
76.6%1004%
1.1%%
1.65 = 1.36%
R3
1.51 + 0.26%
0.78 + 0.18%
1.98 + 0.72%
2.28 + 0.76%
0.83 + 0.15%
Causal R3
2.12 + 0.67%
0.75 + 0.12%
5.99 + 5.25%
2.28 + 0.76%
0.71 + 0.007%
LearnQuad
0.78 + 0.002%
0.68 + 0.02%
0.79 + 0.02%
0.76 + 0.01%
0.87 + 0.01%
LearnQuad to efficiently solve multiple instances of a given PDE. Our extensive set of baselines include: PINN(fxee)) PINN(dynamic) [Karniadakis et al.] 202T]. Curr Reg [Krishnapriyan et al.
2021), CPINN (fixed, dynamic) Wang | [2022], RAR-G [Lu al. 2021]. RAD [Nabian et al. 202I]. RAR-D [Wu et 2023). R3 and Causal R3 Daw et al. 2023J. We performed experiments on benchmark PDEs considered in the baselines. This ensured that we able to clearly demonstrate the advantage of LearnQuad.
6.1 Solving PDEs using LearnQuad
Setup: We compare the performance of LearnQuad in -INN-based solutions of several well known PDEs. We benchmark the performance of our proposed data adaptive quadrature scheme against several other adaptive and non-adaptive algorithms [Wu al. 2023 Lu et al. 2021 Daw al.] 2023) details used in the experiment are included in Appendix [0.1 We used the exact same number of points to train the solution model in all methods for a given PDE, and the number of collocation points is chosen based on the bascline methods. The solution model in each case had the exact same number of parameters to ensure a fair comparison.
Additional experimental details including the explicit form of the PDEs, hyper-paramede
the performance metric following two different exper- imental settings from Daw et al. [2023] and Wu et al.) (2023) over 7 PDEs in Table [ and Table [] respectively. In all but one scenario, LearnQuad is able achieve the best solution function. As can be seen from the numerical results in Table [ the Le relative error for models trained using LearnQuad are better by an order of magnitude most cases and also have very little variance (results reported are an average over five runs). The performance of LearnQuad improves as the number of evaluation points increases, as presented in Table 6. Additionally, we observe that Learn can achieve similar performance to other adaptive methods with much smaller number of points in many cases. All methods have a comparable runtime and memory consumption.
Result:
We report L2 relative error (48)
As an adaptive method, LearnQuad is highly effective in PINN-style solving of PDEs, leading performance boost of PINNs in all cases. The findings highlight the usefulness of adaptive methods over non-adaptive ones specifically when the solution function is not well-behaved. LearnQuad can be used drop in replacement for any sampling strategy.
Summary:
Remark 6. 1. We include additional results on solving PDEs via LearnQuad in their strong form, weak form and also using the energy method in Appendix TO.3 These demonstrate that LearnQuad is versatile adaptive scheme which can be used to solve PDES in multiple reformulations.
Remark 6.2. We use LearnQuad in solving a 100 dimensional PDE, a Poisson equation with a very smooth solution (details in Appendix TO.2 following Yu et al. (2018) ) and achieve relative Te error of 0.085 which is similar to using naive Monte Carlo in this setting with relative T.2 error of 0.09. This illustrates the viability of LearnQuad for high dimensional PDEs. Since the solution is smooth in this particular case, there is no substantial benefit in using a data-driven adaptive method.
Remark 6.3. In Appendix 10. we show that LearnQuad remains numerically stable and improves performance by increasing the number of collocation points. We observe efficient scaling for number of collocation points runtime remains nearly flat up to large point sets, with only modest sub-
8

---

linear growth at very large scales due to hardware bandwidth limits (Appendix 11.1p We also discuss runtime comparison of introducing LearnQuad in PINNs in Appendix [1.2] While our Jax implementation is the fastest among adaptive alternatives, compared to vanilla PINNs LearQuaad is slightly more expensive due to an additional compute operations.
6.2 Comparison of LearnQuad to Classical Solvers
Table 3: L2 Relative Error Comparison: Numer- ical Solvers vs. Adaptive PINNs. Best perform- ing numerical methods are in purple and most competitive adaptive PINN methods are in blue. Learnuuad (LQ) improves performance over other adaptive methods, but there remains a performance gap with numerical solvers. Legends: LW Lax- Wendroff; BW = Beam- Warming; FE = For- ward Euler; BE Backward Euler; CN Crank-
We benchmark strong non- non-PINN solvers on multiple problems and report the results alongside adaptive PINN solvers, including LearnQuad. We used the PDEs from the Sec- tion 6. the details of which are in Appendix 10.1 The relative T-2 error is computed against the true analytical solution which is readily avail- able in these cases.
Setup:
Result: As expected, standard methods (e.g..
Nicolson; IF = Leapfrog: CFD = Centered FDM; NB Newmark-betar R3 [Daw et al., 2023] RAR-
Max-Wendroff, Crank-Nicolson), achieve the lowest relative errors as shown in [ How- ever, LearnQuad consistently improves PINN performance across all cases and notably outper- forms the strongest adaptive PINN baselines.
G [Lu al. 2021], RAD Nabian et al. 2021].
Numerical Solvers
Adaptive PINNs
PDE
Method / Le Error
Method / L2 Error
Convection (3 = 30)
Upwind
LW
BW
R3
IQQ
Summary: LearnQuad advances PINN perfor- mance over existing sampling methods, though there still exists a gap with classical solvers. We also acknowledge that PINNs in their current state are typically more expensive and less accu- rate than classical solvers for forward problems.
0,0460
0.0002
0.0002
0.0078
0.0068
Convection
Upwind
LW
BW
R3
LQ
(3 50)
0.0755
0.0003
0.0003
0.0228
0.0076
FE
BE
CN RAR-G
LQ 0.0005
Diffusion
0.0004
0.0004
0.0004
0.0009
The results do not suggest superiority over nu- merical methods but helps put the scale of gains Learn brings for PINNs context.
Wave
CFD
NB RAD
LQ 0.0050
0.0021
0.0021
0.0015
0.0900
6.3 Solving family of PDEs via LearnQuad
Given the effectiveness of LearnQuad in solving a given PDE, we now utilize our framework tackle harder problem. We consider a family of PDEs, where our end goal is solve PDE given particular choice of forcing function and/or PDE hyper- parameters and initial and/or boundary conditions. Based on 77. a family of PDEs corresponding to differential operator refers to the set of triplets, {(Si, Ii, ) where each i-th PDE satisfies Lu; = St, in Ω2 It. in a.. We only assume access to fi, g:'s and L. We note that f, denotes the forcing function and/or PDE hyper-parameter and 9i denotes the initial and/or boundary condition corresponding to the i-th PDE which is governed by operator L. As an example, if f has the following parametric form:
Setup:
[x) = (a(xT)) sin(TOr) + b(TT)) cos(TTu)); where k = {a, (a,b,0, ~
(21)
we sample p to obtain fi's (similarly for g: 's) and learn to solve for PDEs corresponding to Table 4: Absolute relative error the test set of four
C. We use two hyper-networks with parameter(s) 0 and á which provide pa- rameters of the weight function w, (x)) and solution function (x) respectively based on the input Æ ~ p. This weight function is then used to generate suit-
different PDEs. Details are in Appendix 10.4
PDE
Wave
Advection
Hcat
Burgers'
IC/BC
70
66
80
8
Error
9.9c-6 3.3e--33ee550
5| 9e-
2.1e
3..5e-
12.88
2e-1
3.3e-1
able quadrature } for the PDE corre-
sponding k. These are then used to evaluate the (x) and fx) and minimize the loss based
on the strong form as in Fig. 4 Pseudo-code and experimental details are in Appendix [0.4 We demonstrate the effectiveness of LearnQuad via: Laplace, Advection, Burger's 's, Wave and Heat equation.
We report the absolute relative error compared to the numerical solution obtained using the same number of domain points used for LearnQuad in Table As can seen, the model generalizes very well to test set PDEs. Additional visualizations are included in Appendix I0.4
Result:

---

Summary:
Learning quadratures for family of PDEs is beneficial. Once trained, given new forcing function and boundary/initial condition can generate the solution in a single forward pass.
This is highly efficient in terms of time and eliminates the need to store or process separate solutions, unlike standalone PINNs.
Ditastaiion
K
Paramettc
((a
Urpask
Gompule Guandrature
Remark 6.4 (Distinction with Op- erator Learning). While the end result of learning for family of PDES is similar operator Iearn- ing, the setting is different. While operator learning uses paired data (fi, uz), PINNs with Quar only utilize fi's and shared op- erator L without requiring any problem-solution paired data.
Expansion confficiente
Lil Pugbo)i (580),
Learmalle
(x
Weight functiors My(a)
( Fagadle
Parometon
Loarnablc
peramatera
Unpack
Solution function:
Openator
Figure 4: Solving family of PDEs governed by operator C, forcing function/external condition/PDE hyper-parameter are parameterized by known distribution p. Two hyper-networks with parameters 0 and e generate the weight and solution function.
7 Related Work
Beyond the literature described in T a large body of work focuses on discovering solutions to PDEs using neural networks. LearnQuad enables PINNs RRissi et al. 2019). , to use data-dependent sampling of collocation points. While Variational-PINN (Kharazmi et al. 2019] and hp-VPINN [Kharazmi et al.) 2021] solve PDEs in weak form, they use a careful choice of test functions, which is learnable in our case. Compared to Deep-Ritz (Yu et al. 2018). our method does not need the minimum energy principle to applicable. Our work provides a novel way to improve different variants of PINNs and is complementary to existing literature. We acknowledge recent ideas adjacent
to adaptive quadrature (Rivera al. 2022 Omella and Pardo 2024 Lau al. 2024). which either directly try to optimize node locations thereby resulting in much larger optimization problem or
fall-back to problem-specifie regularizer(s) which limit their applicability. The learnable weight function in LearnQuad induces problem-specific quadrature rules. It differs from classical adaptive methods and ML--based techniques like R3 (Daw et al.. 2023), RAR (Lu et al. [02T], RAD [Wu et al. 2023 which rely residual-based error estimators or gradient thresholding. Other works like [aag et al. 2021] introduce an annealed learning rate and a refined network architecture; [Urbán
et al. 2025
proposed modified optimizer and an adjusted loss. LearnQuad is orthogonal to these
techniques focusing on adaptive sampling and hence can be seamlessly integrated with them.
8 Conclusions
We present a data-driven approach to solve PDEs via PINNs, by exploiting new results of fast quadrature computation using asymptotic expansions and recent capabilities of implicit function dif- ferentiation. We demonstrate that incorporating our learnable quadrature scheme, LearnQuad within PINNs can lead performance improvement over exisiting adaptive and non-adaptive sampling schemes across a diverse set of PDEs. Additionally, we show that incorporation of LearnQuat is quite beneficial when solving a family of PDEs where the alternative would be to deploy Monte Carlo based scheme for each instance individually. Our hyper-network based approach generates the solution to a PDE instance from a given family in just a single forward pass without requiring any additional training.
Limitations and Broader Impact: LearnQuad improves the capabilities of PINNs on commonly studied benchmark PDEs described in the literature. However, we should acknowledge that despite active work surrounding PINNs, in their current stage, these models are not drop-in substitute for classical solvers for common real- real-world PDEs. This gap will likely get smaller with time but this context should clarify that PINNS (with or without LearnQuad) are not intended to substitute the mature body of existing work on traditional numerical methods algorithms and efficient implementa- tions deployed across many fields. Extending LearnQuad to effective higher dimensional settings will likely involve the use of tensor products or sparse grids and is ripe direction for future work. Our work describes mathematical development to improve PINNs, and as such, does not have an immediate societal impact.
10

---

Acknowledgments: Our work was motivated by interesting discussions with Daryl Nazareth (Roswell Park Cancer Institute) around the use of physics informed models [Zhou et al., 2023 for radiation therapy [Bedford 2019]. We are grateful to Daryl and Anant Gopal (also at RPCI) for their time. We thank Peter Opsomer and Daan Huybrechs for their helpful discussion during the initial phase of this project and Lopamudra Mukherjee for feedback on a first draft of this paper. S.P. and V.S. were partly supported by a contract via the DARPA Strengthen program.
References
Mohammad H Aliabadi and David P Rooke. Numerical fracture mechanics, volume 8. Springer Science & Business Media, 1991.
William F Ames, Numerical methods for partial differential equations. Academic press, 2014.
Kamyar Azizzadenesheli, Nikola Kovachki, Zongyi Li, Miguel Liu-Schiaffini, Jean Kossaifi, and Anima Anandkumar. Neural operators for accelerating scientific simulations and design. Nature Reviews Physics, pages 1-9, 2024.
WA Beckham, PJ Keall, and JV Siebers. A fluence-convolution method to calculate radiation therapy dose distributions that incorporate random set-up error. Physics in Medicine & BBooloy,, 47(19): 3465, 2002.
James L Bedford. Calculation of absorbed dose in radiotherapy by solution of the linear boltzmann transport equations. Physics in Medicine & Biology, 64(2))02TRO1, 2019.
Mathieu Blondel, Quentin Berthet, Marco Cuturi. Roy Frostig, Stephan Hoyer, Felipe Llinares- López, Fabian Pedregosa, and Jean Shu-Philippe Vert. Efficient and modular implicit differentiation. Advances in neural information processing systems, 35:5230-5242, 2022.
Ignace Bogaert. Iteration-free computation of gauss-legendre quadrature nodes and weights. SIAM Journal on Scientific Computing, 36(3):A1008-AI10) 2014.
Nicolas Boullé and Alex Townsend. A mathematical guide to operator learning. arXiv preprint arXiv:2312.14688. 2023.
Shengze Cai, Zhiping Mao, Zhicheng Wang, Minglang Yin, and George Em Karniadakis. Physics- informed neural networks (pinns) for fluid mechanics: review. Acta Mechanica Sinica, 37(12): 1727-1738, 2021.
Stéphane d' Ascoli, Sören Becker, Alexander Mathis, Philippe Schwaller, and Niki Kilbertus. Odeformer: Symbolic regression of dynamical systems with transformers. arXiv preprint arXiv:2310.05573, 2023.
Arka Daw, Jie Bu, Sifan Wang, Paris Perdikaris, and Anuj Karpatne. Mitigating propagation failures in physics-informer neural networks using retain-resample--eleass (R3) sampling. In Proceedings of the 40th International Conference on Machine Learning, 2023.
Jochen Garcke et al. Sparse grid tutorial. Mathematical Sciences Institute, Australian National University, Canberra Australia, 7, 2006.
Gene H Golub and John H Welsch. Calculation of gauss quadrature rules. Mathematics of computation, 23(106):221-230, 1969.
Samuel Holt, Zhaozhi Qian, and Mihaela der Schaar. Deep generative symbolic regression. arXiv preprint arXiv:2401.00228. 2023.
George Em Karniadakis, Ioannis G Kevrekidis, Lu Lu, Paris Perdikaris, Sifan Wang, and Liu Yang. Physics-informed machine learning. Nature Reviews Physics, 3(6))422-440, 2021.
Ehsan Kharazmi, Zhongqiang Zhang, and George Em Karniadakis. Variational physics-informed neural networks for solving partial differential equations. arXiv preprint 11X.:::112.00073.. 2019.

---

Ehsan Kharazmi, Zhongqiang Zhang, and George Em Karniadakis. hp-vpinns: Variational physics- informed neural networks with domain decomposition. Computer Methods in Applied Mechanics and Engineering, 374:1 13547, 2021.
Nikola Kovachki, Zongyi Li, Burigede Liu, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Neural operator: Learning maps between function spaces. arXiv preprint arXiv:2 08.80, 2021.
Aditi Krishnapriyan, Amir Gholami, Shandian Zhe, Robert Kirby, and Michael Mahoney. Char- acterizing possible failure modes in physics-iformedd neural networks. Advances in neural information processing systems, 34:26548-26558 2021.
J Nathan Kutz. Deep learning in fluid dynamics. Journal of Fluid Mechanics, 814:1 2017.
Gregory Kang Ruey Lau, Apivich Hemachandra, See-Kiong Ng, and Bryan Kian Hsiang Low. Pinna- cle: Pinn adaptive collocation and experimental points selection. arXiv preprint arXiv:2404.07662, 2024.
Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Neural operator: Graph kernel network for partial differential equations. arXiv preprint arXiv:2003.03485, 2020a.
Zongyi Li, Nikola Borislavov Kovachki, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew Stuart, Anima Anandkumar, et al. Fourier neural operator for parametric partial differential equations. In International Conference on Learning Representations, 2020b.
Qin Lou, Xuhui Meng, and George Em Karniadakis. Physics-informed neural networks for solving forward and inverse flow problems via the boltzmann- bagk formulation. Journal of Computational Physies, 447:1 10676, 2021.
Lu Lu, Xuhui Meng. Zhiping Mao, and George Em Karniadakis. Deepxde: deep learning library for solving differential equations. SIAM review, 63(1):208-22) 2021.
Siddhartha Mishra and Roberto Molinaro. Estimates on the generalization error of physics-informed neural networks for approximating pdes. IMA Journal of Numerical Analysis, 43(1):1-43, 2023.
Mohammad Amin Nabian, Rini Jasmine Gladstone, and Hadi Meidani. Efficient training of physics- informed neural networks via importance sampling. Computer-Aided Civil and Infrastructure Engineering, 36(8):962-97) 977, 2021.
Shechan Olver, Richard Mikaël Slevinsky, and Alex Townsend. Fast algorithms using orthogonal polynomials. Acta Numerica, 29:573-699, 2020.
Ángel J Omella and David Pardo. r-adaptive deep learning method for solving partial differential equations. Computers & Mathematics with Applications, 153:33-42, 2024.
Peter Opsomer. Asymptotics for orthogonal polynomials and fegh-frequency scattering problems. 2018.
Peter Opsomer and Daan Huybrechs. High-order asymptotic expansions of gaussian quadrature rules with classical and generalized weight functions. Journal of Computational and Applied Mathematics, 434:1 1,5377, 2023.
Maziar Raissi, Paris Perdikaris, and George Karniadakis. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. Journal of Computational physics, 378:686-707, 2019.
JR Rice and Dennis Michael Tracey. Computational fracture mechanics. In Numerical and computer methods structural mechanics, pages
Elsevier, 1973.
Jon A Rivera, Jamie M Taylor, Ángel J Omella, and David Pardo. On quadrature rules for solving partial differential equations using neural networks. Computer Methods Applied Mechanics and Engineering, 393:1 14710, 2022.
12

---

Elias M Stein and Rami Shakarchi. Complex analysis, volume 2. Princeton University Press, 2010.
Alex Townsend. The race for high order gauss-legendre quadrature. SIAM News, 48:1-3, 2015.
Alex Townsend, Thomas Trogdon, and Shechan Olver. Fast computation of gauss quadrature nodes and weights on the whole real line. IMA Journal of Numerical Analysis, 36(1):337-358, 2016
Lloyd N Trefethen and David Bau. Numerical linear algebra, volume 181. Siam, 2022.
Jorge F Urbán, Petros Stefanou, and José A Pons. Unveiling the optimization process of physics informed neural networks: How accurate and competitive can pinns be? Journal of Computational Physics, 523:1 13656, 2025.
Sifan Wang, Yujun Teng, and Paris Perdikaris. Understanding and mitigating gradient flow pathologies in physics-informed neural networks. SIAM Journal on Scientific Computing, 4355)11305--33081, 2021.
Sifan Wang, Shyam Sankaran, and Paris Perdikaris. Respecting causality is all you need for training physics-informed neural networks. arXiv preprint arXiv:2203.07404, 2022.
Chenxi Wu, Min Zhu, Qinyang Tan, Yadhu Kartha, and Lu Lu. A comprehensive study of non- adaptive and residual-based adaptive sampling for physics-informed neural networks. Computer Methods in Applied Mechanics and Engineering, 403:115671, 2023.
Larry C Young. Orthogonal collocation revisited. Computer Methods in Applied Mechanics and Engineering. 345:1033-1076, 2019.
Bing Yu al. The deep ritz method: a deep learning-basee numerical algorithm for solving variational problems. Communicutions in Mathematics and Statistics, 6(1):1-12, 2018
Jiahang Zhou, Ruiyang Li, and Tengfei Luo. Physics-informed neural networks for solving time- dependent mode-resolved phonon boltzmann transport equation. npj Computational Materials, (1):212, 2023.
13

---

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: We have presented experimental evidence to back our claims in Section 6 The abstract and introduction have well defined the scope of these experiments and set the tone.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in 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 to other settings.
is fine to include aspirational goals as motivation as long as is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification:
We discuss limitations of our work in Section g
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 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 few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low images are taken in low lighting. Or speech-to-text system might not 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 fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
the paper. The authors should use their best
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [NA]
14

---

We do not provide theoretical result in the paper. Hence, this question is not applicable here.
Justification:
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.
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]
Justification: Yes, we include all experimental details necessary to reproduce in Section and Appendix 10 Pseudo code for the main algorithms presented are in Algorithm I and
Algorithm 2
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer this question will not perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
the contribution is a 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 a specific model and empirical evaluation, may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data 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 submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how reproduce that algorithm.
(b) If the contribution is primarily
model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way access this model for reproducing the results way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reprodueibility may be tricky in some cases, in which authors are welcome to 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., to registered users), but it should be possible for other researchers have some path reproducing verifying the results.
5. Open access to data and code
15

---

Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes]
Justification: Yes, our code will be released when accepted for publication.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips.co cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (c.g., for a new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https
//nips //ipsscc//pill wttpcccpppppppppeeeeeeeeeeeeeeee httpa///giiieee/eeeSebiiisssssss /Cc////oddiiisssseeeeeee ssioiiliccyyyyy for more details.
The authors should provide instructions on data access and preparation, including how to 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 as much information as possible in supplemental material (appended the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
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 understand the results?
Question:
Answer: [Yes]
Justification:
We provide extensive details for experiments in Section and in Appendix 10
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, 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: Our performance is reported average of several runs along with mean and standard deviation as is standard in relevant literature.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, statistical significance tests, 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 with given experimental conditions).
16

---

The method for calculating the error bars should be explained (closed form formula, call a library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
It 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.
For asymmetric distributions, the authors should be careful figures symmetric error bars that would yield results that are out of range (e.g. negative rates).
show in tables or
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification:
These details are available in the Appendix 10
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 failed experiments that didn't 't make it into the paper).
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics
Answer: [Yes]
Justification: Yes, the research performed conforms, in every respect, with the NeurIPS Code of Ethics.
Guidelines:
The answer NA means that the authors have 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., there a special consid- eration due to laws or regulations in their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive socictal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification: Yes, we include this in Section
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 why the paper does not address societal impact.
17

---

Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), faireess 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 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 generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used 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.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper deseribe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification:
The paper doesn't 't invlove the release pre-traee pre-trained models or datasets.
Guidelines:
The answer NA means that the paper poses no 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 to 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 is challenging, and many papers do not 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: Yes, all sources have been cited appropriately.
Guidelines:
The answer NA means that the paper does use 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 URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website). the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms use in package should be provided. For popular datasets, paperswithcone com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
18

---

For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If 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]
Justification:
The paper does not release any new assets.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/moddl 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 to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification: The paper does not involve crowdsourcing nor research with human subjects.
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 as much detail as possible should be included in the main paper.
According the NeurIPS Code of Bthies, 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 or institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsoureing nor research with human subjects.
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 to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, not include any information that would break anonymity (if applicable), such the institution conducting the review.
19

---

16. Declaration of LLM usage
Question: Does the paper describe the usage of LLMs if it is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology. scientific rigorousness, originality of the research, declaration is not required.
Answer: [NA]
Justification: The work presented does not involve LLMs as any important, original, or non-standard components.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, or non-standard components.
Please refer to our LLM policy (https:/ /ttpss/////..p....eee. wwwptrripp..ccooofefeeeeeeee /t/CCooffeennee///22222 for what should or should be described.
20

---

Appendix
9
Asymptotic Expansion
In this section, we list the full expansion of nodes and weights used for experiments in the paper in Section [5]
For the left hard edge at x = -1 used in 14
253,8
-2j
2 -1+ (C+do)
10222 - 13a- -82 +
3(T + do)
-13.5
[16(do 3d1)03.A +3(4a² - 1)co + (12a² + 838 5)do
6(2n +
+ do)
6(4 1)d1] +... +(n
Wk ~ J/ (J dol2
8
w(Tx)
8
+ ( 1n +a + 3 do)
33a² + 62 - - 2j
2032(do 3d1 )53,5 3(402 1)c0 + (12a² + 883 5)do 6(4/2 -1143211100 1Idd]
+
O **)
3 (ig,ee 2n a 3 dols
(22)
where I' = 2n + a a 8 3+1
For the bulk region used in (16)
Dk ~ t 2x² 2/a + (2a² + 2,3 1)6k
4 t + 70 (4(a 1)xt
2 + 2
+4(0² - 1)dg + 8(a² - 33 4((a² - 5 +2(22 - 1)
+2[(2a² 28 +2 m
(23)
[[4(a² - 1)c0 4(3² 1)dg + 4(3a²
8² - 1)70 2(2 + - 1)m)(k) h.o.t.
Wk
1 - T
2 2n
12))
2To4k
+
(2n + a + B 1)2 2 + 47 1
w w(xk)
r
r
47orts + +2(7)
2T
+ 2a + 2/2 + 2 - 1) + h.o.t.
The coefficients e and de are given by:
1 log(h(E)) dg
C = 2τi Sy (E2 - 1)1/2 (E 1)
(24)
dg dx 2mi f, (E2 1)1/2 (E + 1)++1
log(h(E))
(25)
The coefficients Th and Ti are expansion coefficients in the following:
1 log h(S)dc
8
Σ 7/ ((2-tt)) tn
2Tí J V
(26)
1(5 2)
Hereafter, we refer the reader to (Opomee 2018 Opsomer and Huybrechs) 2023) for further detail
on the asymptotic expansions pertinent modified Gauss-Jacone weight functions.
In Fig. we demonstrate how different choices for hg (x: lead to different wo(x) and how these modifications are different from the standard Gauss-Jacobe weight function.
21

---

Weight Functions Comparison (s-1, p-1)
0.
We0gtt
144
00
3.5 00
30
1.00079--2 6.0
Figure 5: Comparison of weight functions. The left most plot is for the standard Gauss-Jacobi weigh function. The remaining plots are from different neural network parameterization of the modifier he(z). One can see that the weight function changes significantly when optimizing for the modified Gauss-Jacobi function. In every case, we fixed a and 3 to be 1 and note that varying these is necessary for our desired use-case.
10
Experiment Details
10.1
Solving PDEs via LearnQuad
Below describe the PDEs used the experimental results Table 2 and Table T
10.1.1
Convection Equation
For the Convection equation in Table 2 we considered the following PDE:
Ou 0t
du
0,
x E [0, 2m], t E 1]
(27)
u(r, 0) h(x))
(28)
w00, () ()
(29)
The model used in this case is fully connected neural network with hidden layers of width 50 and depth 4.
10.1.2 Diffusion Equation
We consider the following one dimensional diffusion equation:
du
8²u
e- (- (in(TT)) "" sin(Tx)), E ==(11,1), Le [0,
[0,
(30)
01
2x
u(x, 0) sin(nz))
(31)
u( u(-1,() w(1,))=
(32)
with domain [-1, 1] in space and [o, 1] in time. The exact solution to this diffusion equation given by u(x; /) sin(#x)) which smooth one and hence all methods illustrated in Table
perform reasonably well even with a small number of points. The model used in this case is a fully connected neural network with hidden layers of width 32 and depth 3.
10.1.3 Burger's Equation
consider the following Burger's equation:
du
Ou dr
82
2x [-(-1,1), LE [0,1].
(33)
0t
u(I, 0) = - sin(Tx),
(34)
u(-1,t) u(1,t) ([1,t)=0)
(35)
where v is the viscosity of the fluid and is the desired flow velocity. In our experiments, we have used 0.01/z which results in a non-smooth solution. The model used in this case is fully connected neural network with hidden layers of width 64 and depth 3.
22

---

Table 5: L2 relative enror (mean + standard deviation) of the trained solution function obtained while using different adaptive and non-adaptive methods. The lowest error for each problem is denoted in boldface. Model trained via LearnQuad achieves the lowest L2 relative error in all case.
PDE
Diffusion
# *points
30
Grid
0.004 0.001
Non-adaptive
Random
0.005 + 0.002
L.HS
0.003 + 0.002
Halton
0.002 + 0.0006
Hammersley
0.001 + 0.0007
Sobol
0.002 + 0.002
Random-R
0.12 0.06
Adaplive
RAR-G
0.0009 + 0.0008
RAD
0.0019 = 0.00097
RAR-D
0.004 + 0.0041
LearnQuad
0.0005 0.0001
10.1.4 Allen-Cahn Equation
The Allen-Cahn PDE considered in our experiments and reported in Table as follows:
du
8uu Ox2
5(u )²), --1-111], t [0.
(36)
δt
u(r, 0) = x² (T))
(37)
a(-1,t) (1,
(38)
We use a value of D = 0.001 as the diffusion coefficient in the PDE. The model used this case is a
fully connected neural network with hidden layers of width 64 and depth 3. The Allen-Cahn PDE considered in Table as follows:
du Ot
pu
0.0001
5u
5u = 0, x 111-1,11, t [0,
(39)
u(x, 0) = x² cos(TT))
(40)
u(t, 1) u(t,, 1)
(41)
Ou
Ou
(42)
a x=1
m=-1
The model used in this case is fully connected neural network with hidden layers of width 128 and depth
10.1.5 Wave Equation
We consider the following one dimensional wave equation:
du 2 2
Øt Ox.¹ x x [0, 1], E [0, 1],
(43)
u(0,t) u(1,t) = O, t [0,
(44)
ulz 0) sin(TT)) sin (TTZ)
E [0,
(45)
Ou
at (r. 0) = 0, e [0.1]
(46)
with e 2, where is the velocity of wave. The solution in this specific choice demonstrates multi-scale behavior in both space and time dimension and is as follows:
1 sin(4 cos( (oTtt
u(æ, t) sin() cos(2nt)
(47)
2
The model used in this case a fully connected neural network with hidden layers of width 100 and depth 5.
23

---

10.1.6 Other Details
The number of parameters used for the learnable weight function in the LearnQuad module was roughly 500 parameters all cases. All neural networks were implemented using fully connected layers with tanh the activation function. All experiments were performed on single NVIDIA 2080 Ti GPU. The number of epochs used for diffusion PDE was 100k while for Burger's, ss Wave and Allen-Cahn PDE they were run for 200k epochs. This was determined empirically based on convergence of the L22 relative error. We used a learning rate of 1e-3. As noted in Algorithm| [] one could either use noise sampled from the standard normal or the PDE specific parameters as input to the learnable quadrature module and results are not too different, but slightly better on using standard normal noise input. We find jointly optimizing both the LearnQuad and solution model provides very good performance without the need for a sophisticated min-max optimization scheme. The L2 relative error reported in the paper is computed the following:
lug ullp
Lzerror
(48)
||u|2
Here, ue is the learned solution function and is the "ground truth" solution. In small number of cases where the true solution is available in a closed form we use that as u or we use u to numerical solution achieved using traditional numerical scheme (finite difference). In any case, the test error is evaluated on a uniform grid of much higher density (10x) than the number of points used in the training scenario. We emphasize that the "ground truth" solution is not used in any form during the training period.
10.1.7 Performance of LearnQuad
We enumerate the performance of LearnQuad with increasing number of points in three different PDEs, (outlined previously) in Table As expected, the performance in terms of Le relative error improves increasing the number of points. Note that the solution to the diffusion equation is very smooth and hence even a very small number of points can lead to very good performance.
Table 6: L Relative Error for Different PDEs with varying number of points used by LearnQuad. Performance improves on increasing the number of points as expected.
Diffusion Equation
Allen-Cahn Equation
Wave Equation
No. of Points
L Error
No. of Points
- Error
No. of Points
L Error
20 0.0013
200
0.017
25
0.0007
200
0.0444
500
0,0076
30
0.0004
700
0.0331
1500
0.0064
35
0.0003
2500
0.0052
40
0.0002
1500
0.0280
3500
0.0044
10.1.8 Performance of LearnQuad with varying seper-par
We investigate the performance of LearQuaad with varying the hyper-parametert of a and (3 in the modified Gauss-Jacobi weight function from ((T)). We report the test performance in terms of the relative 2 relative error in Table 7 We observe minor variations in the performance of LearnQuad based on the choice of these
Table 7: L2 Relative Error for Different PDEs with varying a and 9 in the modified Gauss-Jacok weight function used by LearnQuad. We observe there are minor variations based on the choice of these
(a, 3)
Diffusion
Wave
Convection
[2, 2)
0.0004
0.0058
0.7299
(3, 3)
0.0005
0.0056
0.7163
()
0.0004
0.0065
0.7207
(2, 1)
0.0006
0.0044
0.7323
(10, 10)
0.0007
0.0062
0.6774
24

---

Figure 6: True/Predicted solution functions. 2 right-most two plots for 2 different conditions on the
ID-Laplace operator. left-most two plots for solutions to 2 settings for the operator
d²
+ d
dz
da
10.2 Solving high dimensional PDE
We consider the following high high-dimensional Poisson equation.
Au
-200,
I E (0, 1)
(49)
100
u(x)
2 I E 0(0, 1)
(50)
which is in 100 dimensional space with the true solution being u(z) Ei x². This is a Poisson
100
equation with a very smooth solution (Yu et al.] 2018]. LearnQuad achieves a relative Le error of 0.085 which is similar using naive Monte Carlo in this setting with relative L2 error of 0.09. This illustrates the viability of LearnQuad for high dimensional PDEs. Since the solution is smooth in this particular case, there is substantial benefit in using data-driven adaptive method. This experiment, demonstrates that LearnQuad is not restricted to low dimensional problems. We used fully connected neural network with hidden layers having a depth of and width of 100 as the solution model with tanh as the activation function. We used 1000 points in 100 dimensions. For this problem, our training took 18 seconds to converge in 300 epochs. After this, evaluating the trained model on any given resolution takes 0.0065 seconds. The test errors were computed with respect to the true analytical solution which is readily available in this case.
10.3 Solving PDEs in strong, weak and energy from via LearnQuad
We describe empirical evaluations of our proposed framework using LearnQuad. We show results for solving single given PDE via all three main approaches: (a) the strong form, (b) weak form and (c) minimum principle.
(A) Numerical experiments with Strong Form. We begin by deploying our learnable quadrature first in solving PDEs via the strong from described in $ We consider two operators: (a) 1D-Laplace
and (b) dx² + d. For each of these operators, we consider two different non-homogeneous conditions. As shown Fig. [6 the results of the predicted and true solution function u coincide exactly in all four cases. Both the domain and boundary loss are of the order of c-5, the same as the baseline (PINN). [Raissi al. 22119.
For the 1D-Laplace operator, we use the following two functions as the non-homogeneoos terms:
f(x) = 2 sin(2) + 60z 2((cos(x))² (sin(x))²²
(51)
f(x) = 90(x (4m²) sin(2 (4m²) cos((2T:)))
(52)
For the ID operator dæz d we the following two functions:
d²
f(x) = 3x² + 2 cos ( + + Gar + 2n cos(T=²) (ATu²) sin (ToTTT
(53)
2
f(x) 3x2 + G cos (3mx²)
+ 6x + 6t cos (3mm3) 36m² sin (3Tx²)
(54)
Remark 10.1. Using Monte Carlo based sampling to solve PDEs (as PINNs) can have undesirable outcomes when dealing with irregular boundary, hence adaptive quadrature methods have been proposed very recently [Omella and Pardo 2024]. J Our method data-driven does not suffer from such challenges since the quadratures adaptive by design.
25

---

030
Th 2
0.33
15.5155% 25 2000 22 1997777
0.00
-1.00
Our solution
Ritz Solution
Figure 7: Comparison solution curves obtained via the proposed learnable method and baseline method of Deep Ritz. Both methods perform equally well.
(B) Numerical experiments with Weak Form. We next apply our learnable quadratures to solve PDEs written in their weak form as described [3] Here, we consider the following two 1-D operators: (a) DD-Laplace and (b) dx² de
d² d Similar to the strong form, we present results for two different conditions for each operator. Again, we see from Fig. that, the predicted and true solution function in all cases coincide almost exactly. In terms of the domain and boundary loss, these are of the same order of e- as the baseline method of Hip-VPINN KKharzzmi al. 2021].
Since the weight functions can be global, in using them as test functions to solve the weak form, we can end up with a global test function. Avoiding this is possible via several schemes: one could either choose a multitude of such test functions or simply use sub-domain splitting as suggested in hp- VP-NN over VPINN Kharazmi al. 2019). Due its simplicity, choose the latter in our experiments.
For the ID-Laplace operator, the following two functions as the non-homogeneous terms:
f(x) 2 sin(x) + 60x 2(cos (xx sin (x:) f(x) = (0())) 4n² sin (2mx) x²2 cos (2Tz)
(55)
(56)
For the operator
d² d
d we use the following two functions: d
f(x)= 6z + 2T cOs(TT²)
4x. aaw.rssiinn(ta ( (3x² + 2T ( + =
(57)
f(x) 6x + GT cos(3 36x 20 sin (3Tx2) (3x² + G cos
(58)
00
-031
200
220
Figure 8: True/Predicted solution functions. 2 right-most two plots for 2 different conditions on the d² d di
ID-Laplace operator. 2 left-most two plots for solutions to 2 settings for the operator d=s
Remark 10.2. For solving PDEs in their strong and weak forms as presented above. we adopt a two stage training scheme. In the first stage, the asymptotic quadrature is learned and the second stage these learned quadratures are used to either provide orthogonal collocation points in the strong form or test function(s) for the weak form. Our overall procedure is otherwise unchanged.
(C) Energy Method. We now demonstrate the utility of learnable quadrature for solving PDE where the loss function is derived based on the minimum energy principle.
26

---

We consider the 2D-Laplace equation: Au - 100 with zero boundary conditions on a square domain: [-1,]] x [-1, 1]. In the energy form, the loss function has the form
L(u)= - (1 (1 fu) drdy + 3 / |
u dzdy
(59)
where 3 a penalty term on the second component denoting the boundary loss. The first component the loss the domain. We use learnable quadrature to approximate both integrals in (59) and compare the solution obtained with the baseline method of Deep-Ritz (Yu et al. 20188) with same number of parameters, running each for roughly 400 epochs. As can be seen from Fig. 7 our method achieves comparable performance, with approximate loss value 2000 in both case.
Remark 10.3. Since our proposed method is, in essence, data-driven way to sample points, it shows its utility in solving PDEs via all three formulations as demonstrated above, where the basic framework remains the same. In the strong form, it provides orthogonal collocation points. In the weak form, it provides test functions (which induce the quadrature rules). Finally, in the energy form it is used to directly provide a quadrature rule.
10.4 Family of PDE via LearnQuad
We specify the details of the family of PDEs which were solved using LearnQuad and the procedure
outlined in 6.3 The overall algorithm is presented in Algorithm 2
In all experiments, we used 500 parameters each for the hyper-networks predicting the weight function and solution function as outlined in Section 6..3 Specifically, we used a MI.P-based neural network with depth 5; width 100 and tanh the activation function. The number of parameters to encode the actual solution function were kept smaller than 20. Using learning rate of 0.0001, in all cases, the methods took less than 10k epochs to converge. For each family, we sampled 100 instances of the PDE and used a train/test split of 80/20. We used 600 points as standard number of points to sample from Learn Quad.
Algorithm 2 Training for a family of PDEs
1: Input: Operator L; distribution yt parametric form of forcing function F and boundary/initial conditions B. #epochs: T, Training size: n, Learnable modules 0, Φ;; PDE Loss function L and regularization loss w (20)
2: Compute: Generate Training Set, S 0
3: for = 1 to
do
4: Sample K ~
5: Get for from F, Get by from B
6: S.append fr, bx)
7: end for
8: Compute: Training Loop
9: for i 1 to do
10: for each (k, bx)) E S do
Get wx and (To, Tı, co, do, d1 from (k)
12:
Get from (xx)
Use S to get quadrature nodes {xI}x
14:
Loss:/ =
(x,)), + Le Gradient based update for and @ based on
16: end for
end for
18:
Output: Learned modules 0 and Φ
10.4.1 Family of Laplace Equation
Here, we consider the DD-Laplace operator which has the following parametric representation for the non-homogeneous function:
fx(x) = -(ax² sin( + bit ur cos(TVæ))
(60)
27

---

(a) Family corresponding to p (1,)
M
(b) Family corresponding to p U11, 4.0)
Figure 9: Test results for solving PDEs from two different family each with 1-D Laplace operator. The second family [b) has more variance than the first. In all cases, the predicted solution almost coincides with the original solution
where, k = {a, b, μu, /) belong to different distribution. In our experiments, we choose these distributions uniform, but our method can handle any distribution. In Figure we show the performance on the test set. It can seen that the predicted solution is very close to the true solution.
10.4.2 Family of Heat Equation
We consider the one dimensional heat equation and sample the heat diffusivity, c; initial distribution, f; and two boundary conditions, Ti and T,- The PDE along with initial and Dirichlet boundary conditions is given as follows:
Ou
Uu
x 1E 1], te [0,2]
(61)
OI
u( u--1,4) TI
LE [0,2]
(62)
u(1, () = Tr,
LE [
(63)
u(z, 0) = f(x), x E [-1
(64)
We perform experiments, with choices for the initial distribution:
f(z)
+
(65)
f(x) = a sin (T0x) b cos
(66)
We present a visualization the true (numerical) solution obtained using the same number of domain points as LearnQuad, the predicted solution and their relative error in Figure I0 The parameters were sampled from uniform distribution U (0.8, 2.0).
10.4.3 Family of Wave Equation
We consider the 1D wave equation and sample the wave speed, e and the initial position, f and velocity, g. The PDE along with initial conditions is given below:
pu
²²
x [-1,11 t [0,2]
(67)
ōt²
0x
u(x, 0) = f(x),
x
(68)
du (x. 0) g(x),
x E [-1,1]
(69)
Øt
We perform experiments, with the following two sets of initial conditions:
f(x) = ma. ((x) + x; x [-1,1
(70)
f(x) = m. + n, g(x) aaiannnnnsssn + cos(TT))); [-111
(71)
28

---

175
175
135
125
125
135
0.000
17
339
021
Figure 10: (Left)True solution, (Center) Predicted solution and (Right) Relative for instance from the test set of the family of heat equation using (61)- (-4 and (66)
3.90
Tree
Preskted HARLO
1.8
33
100
0000040
L2
003000
1.8
053
250
243018
1.0
-208
-2
-0.0
-026
6.00 023 558
1:00
003
-008 -2.5
226
200 -ais -8 -0.0
059
1.00
000000
Figure 11: (Left)True solution, (Center) Predicted solution and (Right) Relative error for an instance from the of the family of wave equation using (67) (69) and (70)
We present visualization of the true (numerical) solution obtained using the number of domain points as LearnQuad, the predicted solution and their relative error in Figure and Figure [2 The parameters were sampled from uniform distribution (00.5,3.00)
10.4.4 Family Advection Equation
We consider the one dimensional advection equation and sample the advection speed, c and the initial position f. The PDE along with the initial conditions is given by:
Ou
Du +
x 2E[-1,11, t [0.2]
(72)
Ot
0z
u(x,0) 0) f(x), xE[-1,1]
(73)
We conduct experiments with the two following choices for the initial displacement:
f(x) mx n, [-1,1] f(x) = a sin()0x) + boos(T x E [-1,1]
(74)
(75)
We present visualization of the true (numerical) solution obtained using the number of domain points as LearnQuad, the predicted solution and their relative error in Figure a3] Figure 14 The parameters sampled from the uniform distribution UO..5, 1.5).
10.4.5 Family Burger's Equation
consider the one dimensional viscous Burgers' equation which a non-linear PDE. We sample the diffusivity coefficient, initial velocity distribution, f: and the two boundary conditions, Ti and
29

---

175
175
472
8000000
125
125
800222
13
4:00
-0
-1.
800000
329
0:2
Figure 12: (Left)True solution, (Center) Predicted solution and (Right) Relative error for instance from the test set of the family of wave equation using (67)(69)) and (71)
179-
5.9
5.79
120
1304
1.-
1.0
6.5
1.4
0.5
025
328
800
-215
-255
0.5
8:08
-215
-035
-3225
9:0
00030
Figure 13: (Left)True solution, (Center) Predicted solution and (Right) Relative error for instance from the test set of the family of advection equation using and (74)
T,. The PDE along with the initial and boundary conditions is given by:
du du
à² 8¹
E [--1, 1]. t [0,2]
(76)
aL
Ox.2
u( 1,,t) t) Ti,
t [0, 2]
(77)
m11, 1) = Tr,
LE [0
(78)
u(x 0) f(x),
[-1,1]
(79)
We consider the following two different choices for the initial condition:
f(z)
E [-1,1]
(80)
f(x) aexp
T rE[-111]
(81)
We present visualization of the true (numerical) solution obtained using the same number of domain points as LearnQuad, the predicted solution and their relative error in Figure [5]] Figure and Figure 17 The parameters are sampled from the uniform distribution U (0.011, 0.1).
11
Additional Discussions
11.1 Parallel computation Large number nodes
Learn design, specifically our choice of parameterization avoids explicit computation of several contour integrals, allows for significant parallelization. This due to two reasons: (a) using learnable linear layers for modulating function and expansion coefficients allows bypassing compute-intensive steps such as contour integral (b) computing each quadrature node via an independent root-finding procedure. This is crucial since means that the root-finding process for all nodes and the neural
30

---

175
1.08 475
575
003100
00005
125
AZP
125
A02
1.5
13
-12
-MM
-MM
139
3.29
Figure 14: (Left)True solution, (Center) Predicted solution and (Right) Relative for instance
from the test set of the family of advection equation using 72
(3 and (75)
2.08
Predittec sgl4,
3404
1.4
120
125
125
2004
00004
3224
175
C31e
575
0.55
3,20
375
C228
035
025
0000
Figure 15: (Left)True solution, (Center) Predicted solution and (Right) Relative error an instance from the test set of the family of viscous Burgers' advection equation using (76 79 and (80)
network evaluations to get the parameters can be performed concurrently, and maps well to GPUs. checked wall-clock time by varying the number nodes on our commodity GPU (NVIDIA 2080 Ti).
Table 8: Generation time for different numbers nodes.
Number of Nodes
Generation Time (sec)
10
0.0008
100
0.0009
1000
0,0009
10000
0.0012
100000
0,0098
1000000
0.0754
We see from Table that the runtime remains roughly flat upto 10000 nodes. Beyond this we see modest sub-linear increase in runtime. This increase is due to practica/impplementaion limitations. As number of nodes reaches million,, transferring tensors associated with nodes between compute units and memory increases runtime. This memory bandwidth issue specific GPU can be mitigated (although not eliminated) by higher-end hardware and/or optimized implementation. The core algorithm, by design, remains fully parallel.
11.2 Computational Requirement for LearnQuad
We present wall-clock time per iteration of the baseline methods and LearnQuad below for the convection PDE. We used R3's Pytorch code [Daw et al.) 2023) for all methods, explaining longer runtimes our JAX version. We also tested PINN in JAX for transparency. LearnQuad's higher
31

---

175
475
3327
1.55
2204
3186
0L0044
125
1.25
2.5
C373
3186
E122
3150
13
011
3122
3.9
035
021
Figure 16: (Left)True solution, (Center) Predicted solution and (Right) Relative for instance from the test set of the family of viscous Burgers' equation using (76) (7) and (81)
Traa vou.d
Iraditted apt4,
Sellative STTC
2.0
100
854
0403
144
6.9
159
330
84
+ 1.0
100
1.00
48322
631
0023
124
150
0.59
030
00004
800
-000
+070
-206
-005
040
825
050
675
L00
031
0.0
-380
-275
-003
025
650
035
1.0
813
0:0
15.00
-00
-925
000
005
10k
500
Figure 17: (Left)True solution, (Center) Predicted solution and (Right) Relative error for an instance from the of the family of viscous Burgers' equation using 766 (9 and (81)
runtime over efficient JAX implementation PINN is likely due to the extra learnable model that is serialized before the PINN model. We find the trade-off fair given better performance. Our model converges faster than R3 and other baselines.
11.3 PINNS and Classical Solvers
Our proposed method is not designed to compete with classical solvers. While for simple problems, classical methods are indeed effective, one motivation behind the sizable PINN literature is an alternative which is advantageous in many scenarios: (a) PINN based solutions are mesh-independentt (b) they rely on conoomati-difffreniien linearity effectively given the universal function approximation properties of neural networks, We note that PINN based methods never use "classical solution" as the ground truth in the training procedure at all. It is only used to evaluate test time performance metric. This is needed in cases where the PDE solution is not given in closed form, which true for most scenarios.
which are easier implement; and (c) can handle non-
11.4 PINN loss over PDE solver
Our object interest in this work PINN. PINNs provide a mesh independent solution, are more amenable non-linearities and are easier to scale and implement. Hence, PINNs offer many benefits in several cases and for this reason, are being studied extensively. Next, we justify our choice of PINN loss instead of PDE solver.
While the learnable quadrature rule amenable for classical solvers, there are several issues. Suppose that we use classical solver instead of PINN loss. This means that for each update of parameters learnable weight function we (which induces the quadrature) we will need (i) generate quadrature points using current wg. (ii) solve the system of equations (either implicitly using solver
32

---

Table Average time per iteration for different methods.
Method
Avg. Time per Iteration (sec)
LearnQuad
0.00461 (lax)
R3
0.01172
PINN
0,00274 (Jax)
Causal R3
0.01444
PINN (fixed)
0.01291
PINN (uniform)
0.01261
Curriculum Regularization
0.01278
Causal PINN
0.01317
RAR-G
0,01305
RAR-D
0.01316
RAD
0.01342
L-inf
0.01337
iteratively), (iii) compute some loss/quality and (iv) update to improve this metric. This poses several challenges. Explicit (iterative) solvers are memory-intensive when unrolling across time steps, sensitive numerical instabilities, thereby requiring fine time steps and increased computational cost. When diffcrentiating through numerically unstable solver, the gradients can become inaccurate or blow up. Implicit solvers demand solving linear or nonlinear systems. Computing Jacobians for implicit differentiation requires significant computational resources. Furthermore, matrix inversion or solving linear systems as part of implicit differentiation introduces high computational overhead. Hence, we can agree that integrating PDE solvers into neural network modules presents challenges for both explicit and implicit solvers due to the above mentioned issues in computing gradients which are necessary to update the models via back-propagation. Therefore, in order to make learnable quadrature feasible the main goal of this work we leverage the PINN loss which is more suited for the end-to-end learning framework.
Another aspect worth mentioning regarding the setup for family PDEs. Without scheme to learn the common structure shared between different instances of the PDE, it would require solving each instance separately at each desired resolution. To conclude, the choice of PINN loss is important not only for solving individual PDEs, but also for permitting the learning of quadrature rules that can then be used across multiple problems/solution schemes.
11.5 Contrast with other adaptive methods for PINNS
Our main contribution is not just solving PDEs, but learning how to optimally sample points based on the PDE's structure. We emphasize that advantages stem directly from our core theoretical contribution: the learnable weight function that induces problem-specine
quadrature rules. This is fundamentally different from both classical adaptive methods and other existing ML approaches like R3Daww et al. 22033] RARLu et al. [222T) RADWu et al. [2023]; all of whom invariably rely on computing error estimates through residual-based estimators or gradient thresh-holding which are problem-specific, need to be chosen carefully, and sometimes may need to solve additional local problems. Instead we adaptively learn where refinement may needed in an end-to-end fashion in conjunction with the PINN loss and no additional explicit enror estimation is required.
33


Reviews:
[{"review_id": "ADjtrHxlhp", "reviewer": "Reviewer_oyRH", "summary": "This paper combines the numerical approach with the machine learning approach to learn solutions of PDEs. It proposes the method LearnQuad that enhances the physics informed neural networks (PINNs) with learnable quadrature rules, that are data-driven learned and parametrized by neural networks. By using asymptotic expansions, the method efficiently obtains quadrature nodes and weights from the learned weight function. This method can be applied in both the strong form and weak form of solving PDEs. Additionally, compared to operator learning, this method is able to generalize to solving a family of PDEs, and it does not need solution and PDE parameter pairs. For solving a single PDE, this method yields improved performance over multiple baselines including PINN, CPINN, RAD, etc. For solving a family of PDEs, this method is well generalized to the test set.", "questions": "1. Is it possible for the authors to show the total training time for all methods? Since LearnQuad introduces additional components, it is important to quantify the computational cost compared to PINNs or other baselines. I am also curious how the number of collocation points would scale the training time.\n\n2. It is unclear how the \\alpha and \\beta in equation (10) is chosen. How much effort do the authors put on for tuning the two parameters? How do different choices of \\alpha and \\beta affect the results? An ablation study and analysis would be appreciated.\n\n3. How does the method work in the high dimensional PDE problem in terms of the accuracy and computational cost? Is there any optimization challenge or training instability issue for the learnable quadrature rules?\n\n4. For learning a family of PDEs, the proposed method can be used without labelled data, which is needed for neural operators. I am curious how neural operators perform in the test set. I think the authors can train neural operators either with the solutions learned by the proposed model or with the ground truth solutions.", "limitations": "Yes.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n\n1. For choosing collocation points and their weights, several prior works have been proposed. This method proposes a data-driven method to dynamically learn collocation points and weights. Unlike fixed quadrature or random sampling, the learned quadrature rules is tailored to the structure of the PDEs, enabling a more effective use of limited collocation points. \n\n2. This paper bridges the gap between PINNs and neural operators for solving a family of neural operators. However, the proposed method doesn’t need the use of labelled pairs of solutions and PDE parameters.\n\nWeaknesses:\n1. The biggest weakness is that experiments use 1000 collocation points. Though this paper claims the method is efficient in obtaining a large number of collocation points and weights by leveraging the asymptotic expansions, it is unclear the performance gap would still exist by gradually increasing the number of collocation points. It is also unclear if the asymptotic expansions would be numerically stable at larger scales. \n\n2. Second, the PDEs solutions used in the experiments are smooth. It concerns me that the method based on orthogonal polynomials might be less effective for problems with PDEs solutions that are more complex, e.g., solutions that have multi-scale and are more chaotic. Therefore, it is unclear if the proposed method would keep its accuracy advantage in those settings. It is interesting to see if with 1000 collocation points or by gradually increasing the number of points, how would the proposed method perform compared to baselines.", "was_revised": true, "final_justification": "Previously, I had concerns about the number of collocation points and the solutions being too smooth. Those concerns have been solved. I increased my score to borderline accept.", "rebuttal": "We thank the reviewer for the thoughtful and constructive feedback. We address each point in detail below and will incorporate clarifications and additional results into the final version of the paper.\n\n**Q1. Performance scaling with number of collocation points: Is the proposed method still advantageous when more points are used? Is asymptotic expansion numerically stable at large scale?**\n\nOur experimental setup (Tables 1 and 2) uses $1000–2000$ collocation points only to maintain consistency with baseline works cited in line 291. However, LearnQuad, by design, scales efficiently to a much larger number of points, enabled by the use of asymptotic expansions. To empirically validate this, we reported runtime and scaling behavior in Table 7 (Appendix Section 11.1) and request the reviewer to take a look. To summarize, the wall-clock time to generate quadrature nodes and weights remains nearly flat up to $10,000$ points. Beyond that, we observe modest sub-linear growth. This can be attributed to GPU memory bandwidth bottlenecks, not to any instability or complexity of the core method. As the number of nodes reaches millions, transferring tensors associated with nodes between compute units and memory increases runtime. This memory bandwidth issue specific to our GPU can be mitigated (although not eliminated) by higher-end hardware and/or optimized implementation. The core algorithm, by design, remains fully parallel.\nIn Table 5 and Section 10.1.7 in the Appendix we include numerical results obtained on increasing the number of collocation points. As demonstrated the results improve on using more points in general. We are happy to expand on this section and include larger scale experiments in the main paper.\n\n**Q2. Smoothness of PDE solutions: Does LearnQuad retain its advantage for more complex or chaotic PDEs?**\n\nWe appreciate the reviewer’s concern regarding generalization to more complex PDEs. Our current benchmarks (Table 1 and 2) do include non-trivial settings like Burgers' and Allen-Cahn equations. As mentioned in Appendix 10.1.3, the value of $\\nu=0.01/\\pi$ leads results in a non-smooth solution with steep gradients. Furthermore, the solution to the wave equation used in the experiments as outlined in Appendix 10.1.5 demonstrates multi-scale behavior in space and time. We mention in l291 how existing baseline papers informed our choice of the PDEs. That said, we agree that applying PINNs (with or without LearnQuad) on much more complex PDEs remains open.\nWe note that the adaptivity of LearnQuad driven by a learnable weight function allows collocation points to cluster dynamically/adaptively in regions of high solution complexity. As the number of collocation points increases, the method redistributes these adaptively, we show this in Appendix Table 5.\nWe would be happy to include additional PDEs with more challenging behavior – we have covered a mostly exhaustive set covered in adaptive PINN papers (see lines 291) but welcome specific suggestions from the reviewer.\n\n**Q3. Computational overhead of LearnQuad: How does the training time compare to baselines? How does it scale with the number of points?**\n\nThanks for bringing up this point. We compare wall-clock training time across all methods in Table 8 (Appendix §11.2). We find that our Jax implementation of LearnQuad is significantly faster than other adaptive methods. As reported, compared to vanilla PINNs LearnQuad is bit more expensive due to an additional compute operations but converges faster apart from improving the accuracy of the solution.\nAdditionally, we analyze the scaling behavior of LearnQuad with increasing collocation points in Table 7. The time to generate quadrature nodes and weights remains nearly flat up to $10,000$ points. Beyond that, we observe modest sub-linear growth. Quadrature generation and solution training are fully parallelized and GPU-compatible.\n\n\n\n**Q4. Hyper-parameter tuning for $\\alpha$ and $\\beta$ in Eq. $(10)$: How sensitive is performance to these parameters?**\n\nWe reported an ablation study in Table 6 (Appendix §10.1.8), which explores variations of $\\alpha$ and $\\beta$ in the modified Gauss-Jacobi weight function. We observe only minor variations in performance. Based on this, we fixed $\\alpha$ and $\\beta$ to moderate values (typically in the range 1–2) for all experiments in the main paper. This choice turned out to be sufficient to achieve the desired performance across tasks, without requiring extensive hyper-parameter sweeps. We are happy to expand this section in the final version of the paper.\n\n\n**Q5. Applicability to high-dimensional PDEs: Does LearnQuad remain accurate and efficient? Any training challenges?**\n\nThis is a very good point. An extension of LearnQuad to higher dimensions is indeed possible and incurs some increased memory (which would be the case for PINNs in general) due to the requirement of sampling more points in multiple dimensions to reach a desired error. In higher dimensions techniques of tensor-products or sparse grids can help with reducing the computational complexity. The time complexity of LearnQuad does not explode with increasing dimensions due to parallel compute. However, as mentioned earlier (Table 7) one can expect a a moderate increase with sampling an extremely large number of collocation points. When we performed a simple experiment on a 100 dimensional Poisson equation with a closed form ground truth solution: $-\\Delta u = -200,\\quad x \\in (0,1)^{100};\\quad u(x) = \\sum_{i=1}^{100} x_i^2,\\quad x \\in \\partial(0,1)^{100}$. \n\nWe observe LearnQuad achieves a relative $L_2$ error of 0.085 which is similar to using naive Monte Carlo (with the same number of collocation points) in this setting with relative $L_2$ error of 0.09. This illustrates the viability of LearnQuad for high dimensional PDEs although this is a preliminary finding. Since the solution is smooth in this particular case, there is no substantial benefit in using a data-driven adaptive method. For this problem, our training took 18 seconds to converge in 300 epochs. After this, evaluating the trained model on any given resolution takes 0.0065 seconds. During the testing phase errors were computed with respect to the true analytical solution which is readily available in this case.\n\n**Q6. Contrast with Neural Operators: How does its performance compare?**\n\nWe appreciate the reviewer’s interest in this connection. The reviewer is right that with the help of LearnQuad PINNs can be extended to provide the functionality of neural operators albeit without the use of any labelled data. We believe that this is a step towards bridging the gap between the two modelling paradigms. In the problems considered in the paper, neural operators perform equally well but as the reviewer would agree we needed to provide it with labelled data. We are happy to include these results in the final version.\n\nWe once again thank the reviewer for their detailed evaluation and helpful suggestions. We believe these clarifications and additions strengthen the contributions of LearnQuad and demonstrate its applicability across a broad range of PDE-solving tasks."}, {"review_id": "QXxYep3RAP", "reviewer": "Reviewer_D9vu", "summary": "The authors propose \"LearnQuad\", an approach to sample quadrature points for physics-informed neural network training (PINN training).\nThey demonstrate the performance of the approach in comparison to other sampling strategies in a simple setting (Diffusion equation), and separately demonstrate that quadrature points can also be found for families of PDEs.", "questions": "Also see \"weaknesses\".\nQ1) how would the method perform when using the best PINN optimizers instead of Adam?\nQ2) how accurate are classical methods for the given PDEs, with the same number of sample points?\nQ3) how do the other sampling approaches perform when the number of points increases, or when the PDE is different?", "limitations": "The limitations section does not describe a limitation of LearnQuad, but of PINNs in general. No limitations of the present work are discussed.", "rating": 4, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths: \n1) Several numerical experiments for low-dimensional PDEs, and comparison to other sampling techniques.\n2) Appendix with description of all PDEs and settings.\n\nWeaknesses:\nThe paper introduces a method to sample quadrature points to improve PINN training. In general, improving methods that are still sub-par compared to the state of the art (here: compared to all PDE solvers) is fine, but the comparison should still make it clear how far away the current method still is from the SOTA. This has been discussed extensively in the PINN literature at this point, see for example:\n * McGreivy, Nick, and Ammar Hakim. \"Weak Baselines and Reporting Biases Lead to Overoptimism in Machine Learning for Fluid-Related Partial Differential Equations.\" Nature Machine Intelligence 6, no. 10 (September 25, 2024): 1256–69. https://doi.org/10.1038/s42256-024-00897-5.\n * Wang, Sifan, Yujun Teng, and Paris Perdikaris. \"Understanding and Mitigating Gradient Flow Pathologies in Physics-Informed Neural Networks.\" SIAM Journal on Scientific Computing 43, no. 5 (January 2021): A3055–81. https://doi.org/10.1137/20M1318043.\n * Jorge F. Urbán. \"Unveiling the Optimization Process of Physics Informed Neural Networks: How Accurate and Competitive Can PINNs Be?\" Journal of Computational Physics 523 (February 15, 2025): 113656. https://doi.org/10.1016/j.jcp.2024.113656.\n\nThe first citation includes important suggestions to developing new methods for PINNs, the second and third are important contributions to understanding and improving PINNs. Especially the third demonstrates that the accuracy can be drastically improved by choosing a different optimizer. These developments have to be taken into account when developing new methods for PINNs, otherwise it is impossible to understand how well the new method performs in practice.\n\n1) Weak baselines.\n2) Comparison to classical solvers / solvers in latest PINN literature is missing; points in appendix 11.3 are not enough to address this appropriately.\n3) No theoretical analysis of the approach. Stating that \"no additional explicit error estimation is required\" is not enough to convince the reader why the proposed method is providing \"better\" samples than monte carlo approximation.\n4) A lot of formatting concerns (see below).\n5) The only comparison to other methods is done on a very simple benchmark (diffusion equation 1d) with only 30 sample points. It is impossible to judge if the method performs better than others (especially classical methods) with just this experiment.", "was_revised": true, "final_justification": "The authors now compare to classical solvers, which all outperform the presented method, as expected. The new method (LearnQuad) however outperforms all related PINN training point sampling approaches for these low-dimensional PDEs. My points in the rebuttal were all addressed, but it is hard to judge the new version of the PDF without looking at it - thus I only update to 4.", "rebuttal": "We thank the reviewer for their thoughtful comments. We begin by summarizing our understanding of the core concern and then offer detailed responses.\n\n**Core Concern:** We appreciate the pointer to McGreivy/Hakim’s important critique of weak baselines in ML-for-PDEs. We agree with the call for transparent comparisons and support benchmarking against classical solvers at equal accuracy or runtime. We already align with many of these principles—comparing all sampling/PINN methods under identical setups (same network, optimizer, and compute budget). We also clearly acknowledge that PINNs are typically more expensive and less accurate than classical solvers for forward problems. What we had not done is explicitly quantify this gap. Since this is a key factor in your assessment, we now address it directly. We have benchmarked strong non-PINN solvers on multiple problems and report these results prominently. Rather than undermining our work, we believe these comparisons reinforce our contribution: LearnQuad advances PINN performance over existing sampling methods, though we do not claim it closes the gap with classical solvers.\n\n**Q1. Other Reference:**\n\nWe find the developments in the second and third papers to be complementary to us. The second introduces (a) an annealed learning rate and (b) a refined network architecture; the third proposes (a) a modified optimizer and (b) an adjusted loss. LearnQuad is orthogonal to these techniques—our focus is adaptive sampling via a principled quadrature framework using orthogonal polynomials, which precedes the solution network. In contrast, the cited methods enhance optimization and modeling. Thus, LearnQuad can be seamlessly integrated with them. To validate this, we ran quick experiments combining LearnQuad with proposed setups and observed consistently improved or comparable results. We will include these findings in the final version. This reinforces our message: LearnQuad is a modular, drop-in enhancement that boosts PINN performance across architectures and optimizers.\n\n**Q2. Comparison (Choice of PDE and baselines):**\n\nWe urge the reviewer to look at Table 1 and Table 2 in the main paper which has the following PDEs: Burger’s, Allen-Cahn, Wave, Convection with 1000-200 points. The Diffusion equation with 30 points is indeed a toy setting and hence included in the Appendix. The multiple baselines are chosen following major works which introduce adaptive sampling methods as mentioned in l291. \n\n**Q3. Comparison to classical solvers:**\n\nWe fully agree with the reviewer that PINNs, in their current form, often fall short of classical numerical solvers in accuracy, efficiency, and robustness—a well-known limitation in the field. Our work does not claim that this is not the case. In fact, we explicitly acknowledge that classical solvers (e.g., finite difference, finite element, spectral methods) generally outperform PINNs on well-posed PDEs.\n\nHowever, our goal is not to argue for PINN superiority, but to contribute to ongoing efforts to improve PINN performance and understanding. The gap between PINNs and traditional solvers is real, and our work aims to help narrow it—not by surpassing classical methods, but by enhancing generalization and accuracy through adaptive, learnable quadrature schemes.\nImportantly, our contribution is orthogonal to the broader debate about whether PINNs can eventually match classical solvers. We focus on one key aspect: sampling efficiency, which directly impacts the quadrature error term in generalization bounds (Mishra & Molinaro, 2023). By improving training-time sampling distributions, LearnQuad yields significant empirical gains. In short, LearnQuad is not meant to outperform standard numerical methods like finite elements or finite differences—but to improve over standard PINNs. We believe this is a necessary step if PINNs are to become viable PDE solvers.\n\n**$L_2$ Relative Error Comparison: Numerical Solvers vs. Adaptive PINNs**\n\n| **PDE**            |            **Numerical Solvers: (Relative $L_2$ Error)**                            |             **Adaptive PINNs: (Relative $L_2$ Error)**                    |\n|--------------------|------------------------------------------------------------------------|---------------------------------------------------|\n| Conv. (β = 30)     |  Upwind: 0.0460 &#124;   **L-W: 0.00015** &#124;   **BW: 0.00015**         |   R3: 0.0078 &#124;   **LQ: 0.0068**              |\n| Conv. (β = 50)     |  Upwind: 0.0755 &#124;   **L-W: 0.00025** &#124;   **BW: 0.00025**         |   R3: 0.0228 &#124;   **LQ: 0.0076**              |\n| Diffusion          |  FE: 0.00043 &#124;   **BE: 0.00044** &#124;   **CN: 0.00044**             |   RAR-G: 0.0009 &#124;   **LQ: 0.0005**           |\n| Wave               |  LF: 0.00208 &#124;   CFD: 0.00208 &#124;   **NB: 0.00151**            |   RAD: 0.0900 &#124;   **LQ: 0.0050**             |\n\n**Legend:**  \nL-W = Lax-Wendroff; BW = Beam-Warming; FE = Forward Euler; BE = Backward Euler; CN = Crank-Nicolson; LF = Leapfrog; CFD = Centered FDM; NB = Newmark-beta;  R3, RAR-G, RAD = existing adaptive PINNs; LQ = LearnQuad (ours)\n\nIn response to the reviewer’s suggestion, we benchmarked LearnQuad against classical solvers on three PDEs with known analytic solutions (details in Appendix §10). As expected, standard methods (e.g., Lax-Wendroff, Crank-Nicolson), still achieve the lowest relative $L_2$ errors. However, LearnQuad consistently improves PINN performance across all cases and notably outperforms the strongest adaptive PINN baselines. These results are not to suggest superiority over classical solvers, but to contextualize the scale of gains LearnQuad brings within the PINN framework. We hope this clarifies the scope of our contribution. We will highlight these results in the final version.\n\n**Q4. Theoretical Analysis:**\n\nWe thank the reviewer for pushing the need for theoretical justification. While our main contribution is algorithmic, LearnQuad is grounded in recent theory on PINN generalization—most notably Mishra & Molinaro (2023). Their Theorem 2.6 decomposes the generalization error into training error $E_T$ and a quadrature term $C_{\\text{pde}} C_{\\text{quad}}^{1/p} N^{-\\alpha/p}$, where $C_{\\text{quad}}$ depends on the sampling distribution. This shows that sampling is a central factor in solution quality—not a secondary detail.\n\nLearnQuad directly targets the quadrature term by learning a sampling distribution that concentrates collocation points where the residual is high. This lowers $C_{\\text{quad}}$ without modifying the PINN architecture or PDE model. It offers a practical means to realize the theoretical insights of (Mishra & Molinaro). We will make this connection explicit in the final version. The link is direct—Theorem 2.6 gives an interpretable bound that LearnQuad implements in a scalable, theory-aligned way. LearnQuad is a principled, theoretically grounded adaptive sampling method that minimizes one of the two main sources of PINN error.\n\n*Mishra, Siddhartha, and Roberto Molinaro. \"Estimates on the generalization error of physics-informed neural networks for approximating PDEs.\" IMA Journal of Numerical Analysis 43.1 (2023):1-43*\n\n**Q5. Why are the sampled points better than Monte Carlo?**\n\nThe case for adaptive sampling outperforming uniform Monte Carlo in PINNs has been well established by prior works cited as baselines. These show that allocating points in regions of high residual or solution complexity improves generalization. Our results align with this: LearnQuad consistently outperforms vanilla Monte Carlo PINNs and other adaptive methods across benchmarks. We introduce a principled approach to adaptive sampling via a learnable quadrature scheme. During training, the sampling distribution is not fixed—it evolves by learning a weight function informed by the PDE residual. This function induces a probability distribution that biases collocation points toward informative regions. A key challenge in adaptive sampling is that drawing from such tailored distributions is often costly or complex. We address this by using asymptotic expansions to generate quadrature points quickly and stably from the learned weights, making our method both scalable and easy to integrate with existing PINN models.\n\nWe believe the reviewer will agree that adaptive sampling from problem-specific distributions generally outperforms uniform sampling, especially when residuals or features are spatially non-uniform. By formulating sampling via asymptotic expansions and learnable weights, we enhance performance and offer a general-purpose tool that integrates smoothly into broader learning pipelines. We hope this is seen not just as an empirical gain, but as a step toward scalable, theory-informed adaptive sampling in physics-informed learning.\n\n\n**Q6. Limitations:**\n\nThe reviewer is correct in pointing this out, since PINNs are a vessel for LearnQuad we had outlined the limitations when used in this setting, despite LearnQuad improving the performance of PINNs are demonstrated empirically. We are happy to expand the section to include a discussion of challenges in extension to higher dimensions which we believe is a ripe direction for future research.\n\n\n**Q7. Paper formatting Concern:**\n\nWe have fixed typesetting errors and inserted commas after “where”.\nL101 we meant it is being parameterized using neural networks, the weights of which are learnt.\nWe adjusted the abstraction level, made style consistent throughout the paper.\nWe are happy to remove the color.\nWe have updated the tables in the Appendix to a more standard format.\n\nWe thank the reviewer for their thoughtful critique, which helped refine our scope and positioning. We hope our clarifications—on benchmarking, theoretical grounding, and compatibility with recent PINN advances—highlight that LearnQuad is a principled, modular improvement that strengthens PINNs and moves them closer to bridging the gap with classical solvers."}, {"review_id": "9WWERzoy7u", "reviewer": "Reviewer_yhsF", "summary": "This paper proposes a principled method for learning adaptive quadrature points and weights, which improves the accuracy of PINN. The quadrature point is parameterized as the roots of a parameterized Jacobi polynomial induced from a modified Gauss-Jacobi type weight function, where a learnable multiplicative factor $h_{\\theta}(x)$ is introduced. To compute the roots efficiently, the authors utilize the recent advances in asymptotic expansion of OP described in section 5.1.  Since the quadrature points are implicit functions of the parameters $\\theta$ in $h_{\\theta}(x)$, the author uses implicit differentiation over the simplified root equation (eq. 17-19). Learning is done by backpropagating the PINN loss over the parameters of the solution function and the parameters of the quadrature point. The author also proposed an additional loss for the quadrature point in equation 20. The method was tested on various common PDEs and demonstrated superior accuracy (measured in L2 relative error) against other methods. I believe that this work made a great contribution in improving PINN and making it closer to being comparable to the existing PDE solvers.", "questions": "1. $\\alpha$ and $\\beta$ in the weight function are not learned. How are they chosen? Same question for the expansion coefficients. What kind of information about the equation do we need to derive these?\n2. To derive equation 18, you assume the roots are real and there are only simple poles. When does this condition hold, and could this assumption be violated in some PDEs?\n3. It appears that the first term in the quadrature loss in equation 20 is a relaxation of a hard constraint. How much is this error? And what effects does this error have on the solution?\n4. You mentioned that the quadrature point could be learned adversarially, but in practice, no significant difference was observed. What could be the reason?\n5. In line 230, it was mentioned that \"we use linear layers...\", but from the appendix, it seems that the author is actually using a multi-layer perception with 3~4 layers. I fail to understand this paragraph.", "limitations": "Yes", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 2, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strength:\n1. The method is very principled and demonstrates strong empirical performance\n2. The writing was very clear, and the paper is easy to follow despite the heavy mathematical content\n3. To the best of my knowledge, the method is novel, as previous methods for learning quadrature points do not generate the quadrature points as roots of OP.\n\nWeaknesses:\n1. There should be a comprehensive description of the parameterization. Also, Algorithm 1 provides a sketch of the training process; the bulk of the work in this paper is in how the quadrature points are computed. In this current form, I had to gather various pieces from the paper to understand how exactly the quadrature points are computed. Also, I couldn't find the network architecture description\n2. The method's accuracy should also be benchmarked against non-PINN PDE solvers. It is good to show that the method outperforms other PINN variants, but it is also good to see how much gap there still is.\n3. The method is only tested on low-dimensional standard PDEs. Applicability to real-world problems could be challenging, especially for irregular domains", "was_revised": true, "final_justification": "I am maintaining my initial assessment of the paper and hence keeping my initial score of 5. The learned quadrature method clearly improves upon previous works in PINNs, but the limitation is still apparent. For example, the method cannot be applied to high-dimensional problems yet, and there is still a significant gap in accuracy compared to traditional solvers.", "rebuttal": "We thank the reviewer for the thoughtful and encouraging evaluation. We deeply appreciate the recognition of the novelty and principled nature of our approach, as well as the constructive suggestions. Below, we address each of the reviewer’s comments and questions in detail.\n\n**Q1. Clarification on parameterization and network architecture**\n\nWe appreciate the reviewer’s request for a more cohesive description of the parameterization strategy and network architecture. The architecture details are currently included in Appendix §10.1.6. To summarize: for both the solution function and the weight function (i.e., the LearnQuad module), we use fully connected multilayer perceptrons with depth $3–5$ and width ranging from $50$ to $200$, depending on the PDE. We use tanh as the activation function throughout. We acknowledge the wording in line 230 (“linear layers”) is misleading and will revise this to explicitly state the use of MLPs. Additionally, we will expand Algorithm 1 in the final version to provide a more comprehensive description of how quadrature points and weights are computed and how they are integrated into the end-to-end training loop.\n\n**Q2. Benchmarking against non-PINN PDE solvers**\n\nThank you for this suggestion. In response, we conducted a quick experiment comparing the performance of existing numerical solvers and PINN with LearnQuad on three representative PDEs with known analytic solutions (the details of which are in Appendix 10). We find that LearnQuad consistently closes the gap between existing (best adaptive) PINNs and classical solvers in terms of $L_2$ relative error, while maintaining computational efficiency. These results will be included in the final version of the paper, and we welcome further recommendations for such comparisons.\n\n**$L_2$ Relative Error Comparison: Numerical Solvers vs. Adaptive PINNs**\n\n| **PDE**            |            **Numerical Solvers: (Relative $L_2$ Error)**                            |             **Adaptive PINNs: (Relative $L_2$ Error)**                    |\n|--------------------|------------------------------------------------------------------------|---------------------------------------------------|\n| Conv. (β = $30$)     |  Upwind: $0.0460$ &#124;   **L-W: $0.00015$** &#124;   **BW: $0.00015$**         |   R3: $0.0078$ &#124;   **LQ: $0.0068$**              |\n| Conv. (β = $50$)     |  Upwind: $0.0755$ &#124;   **L-W: $0.00025$** &#124;   **BW: $0.00025$**         |   R3: $0.0228$ &#124;   **LQ: $0.0076$**              |\n| Diffusion          |  FE: $0.00043$ &#124;   **BE: $0.00044$** &#124;   **CN: $0.00044$**             |   RAR-G: $0.0009$ &#124;   **LQ: $0.0005$**           |\n| Wave               |  LF: $0.00208$ &#124;   CFD: $0.00208$ &#124;   **NB: $0.00151$**            |   RAD: $0.0900$ &#124;   **LQ: $0.0050$**             |\n\n**Legend:**  \nL-W = Lax-Wendroff; BW = Beam-Warming; FE = Forward Euler; BE = Backward Euler; CN = Crank-Nicolson; LF = Leapfrog; CFD = Centered FDM; NB = Newmark-beta;  R3, RAR-G, RAD = existing adaptive PINNs; LQ = LearnQuad (ours)\n\n**Q3. Applicability beyond low-dimensional problems and standard PDEs**\n\nWe agree with the reviewer that evaluating LearnQuad on more complex and real-world PDEs including irregular domains and higher dimensions is a valuable direction ripe for future work. As noted in line 291, our current PDE choices follow established baselines to enable direct empirical comparison. We also mention that extension to higher dimensions can involve the use of tensor-products or sparse grids to make it computationally efficient. If there is a specific PDE the reviewer recommends we are happy to include it in the final version.\n\n\n**Q4. Choice of $\\alpha$, $\\beta$, and expansion coefficients:**\n\nThank you for raising this point. In Appendix §10.1.8 (Table 6), we performed an ablation study on various choices of $\\alpha$ and $\\beta$ in the modified Gauss-Jacobi weight function. We observed that LearnQuad's performance shows minor variations to these hyperparameters. Based on this, we chose $\\alpha, \\beta \\in \\[1,2\\]$ for all experiments reported in the paper. We emphasize that the choice was not arbitrary, but rather empirically sufficient to achieve strong performance across PDEs, without requiring exhaustive hyper-parameter sweeps. We will revise the text to reflect this clearly.\nAs for the expansion coefficients (e.g., $c_0, d_0, \\tau_0, \\tau_1$​), these are not fixed but are learned from data as part of the weight function network, as described in lines 232–236. We will make this explicit in the main text to improve clarity.\n\n\n**Q5. Assumption of real roots in asymptotic expansion (Eq. 18)**\n\nThis is an important technical point. Our derivation assumes that the domain of the PDE is real-valued, which ensures that the roots of the orthogonal polynomials are real and simple. This assumption holds for all PDEs considered in the paper and for most physical PDEs encountered in practice. However, as the reviewer rightly points out, this assumption could break down in special cases where PDEs are defined on the imaginary plane. We will include this clarification and limitation explicitly in the final version.\n\n**Q6. Impact of quadrature loss term in Equation (20)**\n\nThe first term in Equation (20) enforces that the sum of quadrature weights approximates the integral of the learned weight function. This is indeed a relaxation of the exact constraint. Empirically, we observe that this term typically converges to a value of order $e^{-3}$ or lower. If we don’t include this term, the final performance is worse and also, its inclusion enforces that the quadrature weights are valid.\n\n\n**Q7. Adversarial training of the quadrature module**\n\nWe explored the idea of learning quadrature points adversarially — where one module proposes points and the other learns under this challenging distribution. However, in our experiments, training both modules jointly (end-to-end) or in an alternating min-max fashion resulted in comparable performance. Adversarial training is known to be more unstable and harder to tune. Without careful balancing (e.g., learning rates, update schedules), it may converge to a similar fixed point as joint training, effectively reducing to the same solution. This could have been an issue specific to our implementation but we were unable to get it working. If there are specific papers we can use as guidance, we welcome suggestions. \n\n**Q8. Clarification on “linear layers” in line 230**\n\nWe apologize for the inconsistent wording. As detailed in the appendix, we use multi-layer perceptrons (MLPs) with $3–4$ hidden layers for both the solution and quadrature modules. Each hidden layer consists of fully connected (linear) operations followed by tanh non-linearities. We will correct the text in line 230 to reflect this architecture accurately. Thank you for pointing this out.\n\n\nWe sincerely thank the reviewer for their thoughtful engagement and strong endorsement. We believe that the clarifications and additions outlined above further strengthen the paper and reinforce the rigor and applicability of LearnQuad."}, {"review_id": "3R4zI2iEcU", "reviewer": "Reviewer_q5zG", "summary": "This paper proposes an enhancement to physics-informed neural networks (PINNs) by introducing learnable quadrature rules. Specifically, the weight function in the quadrature rule is parameterized by a neural network. The loss function includes regularization to discourage small weights and enforces that the learned quadrature weights sum to the integral of the weight function. The method shows improved empirical performance on several 1D benchmark problem. The authors also extend their methods to train neural operators without using numerical solutions.", "questions": "- How does the learned weight function compare to classical choices such as Gauss-Jacobi weight functions? Please provide plots or examples of the learned weight functions for illustration.\n\n- Equation (30), which defines $l_w$, is unclear. It is not evident which parameters are trainable. Are the quadrature weights $w_i$ themselves learnable? It seems that for any given weight function $w(x)$, the weights $w_i$ and nodes $x_i$ are computed using a numerical procedure described in Section 5.\n\n- The paragraph on the interlacing of roots of orthogonal polynomials is confusing. What exactly does “overfitting” mean in this context? A plausible failure mode would be the learning of a trivial weight function (e.g., $w(x) = 0$), which seems to motivate the penalty term in Equation (20). However, the paper states: “By utilizing quadrature nodes and weights stemming from varying degrees of orthogonal polynomials (all of which correspond to the same weight function being learned), we introduce the desired stochasticity to prevent over-fitting.” This seems to be a key component of the algorithm, but it is unclear how it is implemented or integrated into the training process.", "limitations": "Yes.", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 2, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths\n- The method is compared with lots of existing methods and shows empirical success.\n- The approach is well motivated and is principled in adaptive quadrature.\n\nWeaknesses\n- All examples are 1D in space.\n- The presentation lacks clarity in several aspects. The use of $\\theta$ and $\\phi$ is confusing. For example, $w_\\theta$ in Equation (2) and (10), $w_\\phi$, $u_\\theta$ in Algorithm 1, and $\\phi$ goes to the solution function in Figure 4. In addition, in Algorithm 1, why is $\\phi(\\mu)$.", "was_revised": true, "final_justification": "The rebuttal addressed major issue on clarity, including inconsistent notation and the explanation of the interlacing of roots. The planned revision, especially visualization of the learned weight function, will help readers better understand the methods. The method remains evaluated only on 1D problems. The approach is technically sound, well motivated, and shows consistent empirical improvement over baselines.", "rebuttal": "We thank the reviewer for the constructive and insightful feedback. We address each point below and will incorporate the relevant clarifications and improvements into the final version of the paper.\n\n**Q1. Scope of Experiments: Why are all examples 1D in space?**\n\nThanks for bringing this up. Our primary intent was to show the effectiveness of LearnQuad on benchmark problems used in prior PINN literature. As mentioned in line 291, we evaluated our method on standard PDEs from several baselines to enable direct comparison. That said, our formulation (Remark 4.1) naturally can be extended to higher-dimensional settings via tensor-product or sparse-grid constructions. We have checked that in principle the implementation works for simple cases, but these need to be developed in a separate follow-up work. We are happy to include additional results in the final version if the reviewer recommends a specific benchmark PDE.\n\n**Q2. Notation clarity: Confusion around $\\theta, \\phi$, and their roles in Algorithm 1 and Figure 4.**\n\nMany thanks for pointing this out. Yes there was some notational inconsistency in the original draft. We have revised the text and figures to consistently denote:\n\n$\\theta$: parameters of the weight function $w_\\theta​$\n\n$\\phi$: parameters of the solution function $u_\\phi​$\n\n\nIn Algorithm 1, the line $w_\\phi$ from $\\phi(\\mu)$ has a typo and was intended to show that the weight function may optionally depend on PDE parameters $\\mu$. However, this caused confusion since the dependence was not consistently carried forward. We have now removed this line. Figure 4 has been revised to indicate more clearly that $\\theta$ governs the quadrature module and $\\phi$ the solution function.\n\n**Q3. How does the learned weight function compare to classical ones like Gauss-Jacobi? Can you show plots or examples?**\n\nThis is a nice suggestion. We observe that the learned weight functions differ quite a bit from classical choices. While Gauss-Jacobi weights have fixed singularities at the endpoints, our learned modifier $h_\\phi(x)$ dynamically adapts the weight distribution based on the PDE solution’s characteristics. Empirically, we find that the collocation points tend to cluster in high density regions and re-distribute them as training proceeds, just as desired. We will include representative plots comparing learned versus classical weights and their induced node distributions in the final version to illustrate these differences clearly.\n\n**Q4. Equation (20) and clarity about trainable parameters: Are the quadrature weights $w_i$ themselves learnable?**\n\nWe agree that Equation (20) could be better explained. Yes your interpretation is correct, only the weight function $w_\\theta(x) $(parameterized via a neural network) is learnable. The quadrature weights $w_i$ and nodes $x_i$​ are computed numerically from the learned weight function using asymptotic expansions (Section $5$), and not directly learned. We have revised the text to clearly emphasize that the model learns the modulator $h_\\phi(x)$, from which quadrature rules are derived numerically.\n\n**Q5. Clarification on “interlacing of roots” and how it mitigates overfitting. What does “overfitting” mean here?**\n\nWe appreciate the question in this regard and make an effort to clear the confusion. Given a weight function, it induces a family of orthogonal polynomials each of a different degree. LearnQuad utilizes the roots of these polynomials by computing them using asymptotic expansions. While the weight function is being learnt, we observed that only sampling the roots of a fixed degree lead to local minima and the solution function didn’t generalize well to unknown points in the domain. To mitigate this, we leverage a classical property of orthogonal polynomials: interlacing of roots. The roots of the degree $n+1$ polynomial interleave those of degree $n$, as shown in Figure 3. Hence, by using different degrees to sample roots in different epochs leads to a drastically better generalization of the learnt solution. We thank the reviewer for pointing this out and have updated the manuscript to clearly describe this mechanism, its theoretical motivation, and its empirical benefit to generalization.\n\n**To Summarize:**\n\nWe will improve clarity throughout the manuscript, particularly in Section 5 and Algorithm 1.\nWe will revise Equation (20) to explain more precisely how it prevents degenerate quadrature weights and enforces normalization.\nAs suggested, we will add visualizations of the learned weights and node distributions.\nCode will be released to enhance reproducibility and clarity.\n\nWe thank the reviewer again for their thoughtful evaluation and helpful suggestions, which we believe will greatly improve the final version of the paper."}]

Author rebuttal:

Meta-review: (1) Summary of this work: This paper proposes LearnQuad to enhance PINNs with learnable quadrature rules. With asymptotic expansions, LearnQuad can efficiently obtain quadrature nodes and weights from the learned weight function. In addition to experimenting with strong and weak forms of PDEs, LearnQuad can also be extended to a family of PDEs, surpassing previous methods.

(2) Strengths and weaknesses: This paper proposes a novel and reasonable method to choose quadrature nodes and weights, which is different from previous residue-based sampling methods and demonstrates strong performance. In the rebuttal, the authors also compare LearnQuad with numerical methods, where their work also performs well. There are no significant weaknesses in this work, but there is also no evidence to demonstrate that the proposed method is principal or "optimal" to the sampling problem. Besides, as mentioned by some of the reviewers, the writing of this paper needs further improvement.

(3) Summary of rebuttal: The authors made a solid rebuttal and provided extensive experiments, clarifications and discussions to the reviewers. Especially, the experiment compared to numerical methods is significant and surprising. All the reviewers thank the authors for their rebuttal and there are no serious issues remaining.

**Final decision:** This paper presents a novel and reasonable approach to enhance PINN's sampling process, which can fit the structure of targeted PDEs better and show solid performance. After the rebuttal, all the reviewers acknowledged the authors' clarification. Thus, I would like to recommend acceptance. The authors should include the newly added experiments in the final version.
Decision: Accept (poster)