pith. sign in

arxiv: 2607.00560 · v1 · pith:Q6PQQFFGnew · submitted 2026-07-01 · 🧮 math.NA · cs.NA

A multilevel stochastic-gradient neural solver for boundary integral equations

Bing-Ze Lu , Richard Tsai This is my paper

Pith reviewed 2026-07-02 08:09 UTC · model grok-4.3

classification 🧮 math.NA cs.NA
keywords boundary integral equationsneural networksmultilevel methodsstochastic gradient descentNyström discretizationneural tangent kernelspectral bias
0
0 comments X

The pith

A neural network trained level-by-level on refining quadrature grids solves second-kind boundary integral equations with total work proportional to the finest grid.

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The paper establishes a multilevel stochastic-gradient solver in which a multilayer perceptron represents the unknown density and minimizes the Nyström-discretized residual on a sequence of successively finer quadrature grids. Each level is warm-started from the parameters obtained on the previous grid, so that the spectral bias of the network, which normally stalls training on a fixed grid, instead acts as a smoother by deferring high-frequency modes until finer grids make them visible to the optimizer. The empirical neural tangent kernel sampled on the current grid governs the contraction rate. Under a uniform regularity bound preserved across levels, the cumulative work remains a fixed multiple of the work on the finest grid, and the uniform conditioning of the discrete second-kind operator converts the final training residual into an a posteriori error bound on the solution density.

Core claim

Successive quadrature refinement brings additional modes of the continuum neural tangent kernel within reach of the optimizer, converting the frequency principle into the mechanism of a multigrid-type iteration whose total cost is a constant multiple of the finest-grid cost and whose residual directly bounds the error because the second-kind operator remains uniformly well-conditioned.

What carries the argument

Multilevel warm-started stochastic-gradient training on a ladder of refining Nyström quadrature grids, where the empirical neural tangent kernel sampled on each grid controls residual contraction.

Load-bearing premise

The neural network obeys a uniform regularity bound that remains valid when the quadrature grid is refined.

What would settle it

An experiment that records the total number of network evaluations across all levels and finds the count grows faster than linearly with the number of points on the finest grid, or that measures the actual solution error and finds it exceeds a fixed multiple of the training residual.

Figures

Figures reproduced from arXiv: 2607.00560 by Bing-Ze Lu, Richard Tsai.

Figure 1
Figure 1. Figure 1: Initial NTK spectral distributions at Levels 1–4, aggregated over 20 independent runs. Top: NTK [PITH_FULL_IMAGE:figures/full_fig_p021_1.png] view at source ↗
Figure 2
Figure 2. Figure 2: NTK spectral distributions upon termination at Levels 1–4, aggregated over 20 independent runs, [PITH_FULL_IMAGE:figures/full_fig_p021_2.png] view at source ↗
Figure 3
Figure 3. Figure 3: Single-level training on the finest grid: NTK spectra at iterations 1, 1000, 2000, and the final step. [PITH_FULL_IMAGE:figures/full_fig_p022_3.png] view at source ↗
Figure 4
Figure 4. Figure 4: Box plots of epoch counts at each training level comparing the progressive and cyclic schedules [PITH_FULL_IMAGE:figures/full_fig_p022_4.png] view at source ↗
Figure 5
Figure 5. Figure 5: Stanford Bunny example. The first two panels show the surface density obtained by the neural [PITH_FULL_IMAGE:figures/full_fig_p023_5.png] view at source ↗
Figure 6
Figure 6. Figure 6: Total wave field |u + u inc|reconstructed from the MLSG density. (left) n = 256 points per dimension, κ = 4. (right) n = 512 points per dimension, κ = 8. 4.3.1 Scattering by multiple spheres The multi-sphere configuration is the canonical three-dimensional Helmholtz benchmark and the natural setting for our wall-clock comparison with GMRES. We report two experiments: how the MLSG runtime scales with the nu… view at source ↗
Figure 7
Figure 7. Figure 7: shows the computed density and the recovered total field at κ = 8 with N = 1024 (system size N2 = 1,048,576); the MLSG solve took 4,621.6 s of wall-clock time. -4 -2 0 2 -5 4 0 5 -5 0 5 10 20 30 40 50 (a) Computed density |ρ param θ |. 10 -5 5 -5 0 0 0 5 5 -5 0 0.5 1 1.5 2 2.5 (b) Total field |u + u inc| on a constant-z plane [PITH_FULL_IMAGE:figures/full_fig_p026_7.png] view at source ↗
Figure 8
Figure 8. Figure 8: Linked-tori Helmholtz problem (κ = 4) with a zero Neumann boundary condition and an incident plane wave in the x-direction. Figure 9a shows the 3-genus computational domain, where the incident plane wave propagates in the x￾direction. The domain is discretized with n = 128 grid points per dimension (h ≈ 0.2197), and the Helmholtz equation is solved with wavenumber κ = 3. The neural network consists of 10 h… view at source ↗
Figure 9
Figure 9. Figure 9: Three-genus Helmholtz problem (κ = 3) with a zero Neumann boundary condition and an incident plane wave in the x-direction [PITH_FULL_IMAGE:figures/full_fig_p028_9.png] view at source ↗
read the original abstract

We develop a multilevel stochastic-gradient neural solver for boundary integral equations of the second kind. The unknown density is represented by a multilayer perceptron, trained by minimizing the Nystr\"om-discretized residual on a ladder of refining quadrature grids, each level warm-started from the parameters of the previous one. Each step requires only dense matrix-vector products on mini-batches of collocation rows and network passes, operations that map directly onto GPU hardware. The residual contraction is governed by the empirical neural tangent kernel (NTK), the discrete sample of a single continuum kernel. On a fixed grid, training stalls once the residual concentrates in modes the network contracts slowly, the plateau described by the frequency principle; a spectral analysis explains, and experiments confirm, how refining the quadrature resolves more of the continuum kernel's spectrum and returns these modes to the optimizer's reach. Spectral bias, elsewhere an obstruction to neural network solvers, thus serves as the smoother of a multigrid-type iteration, with quadrature refinement in place of coarse-grid correction. Under a uniform regularity bound on the network, the total work is a constant multiple of the work on the finest grid, and the uniform conditioning of the discrete second-kind operator leaves the NTK as the sole rate-determining spectrum while converting the training residual into an a posteriori error bound. Experiments on interior Dirichlet Laplace/Poisson problems and exterior Neumann Helmholtz problems, using both parametric and signed-distance surface representations, demonstrate the effectiveness and efficiency of the proposed method compared with GMRES at comparable tolerances.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Referee Report

2 major / 2 minor

Summary. The paper develops a multilevel stochastic-gradient neural solver for second-kind boundary integral equations. The unknown density is represented by a multilayer perceptron trained by minimizing the Nyström-discretized residual on successively refined quadrature grids, with each level warm-started from the previous parameters. Convergence is analyzed via the empirical neural tangent kernel (NTK); spectral bias is repurposed as a multigrid smoother via quadrature refinement. Under an assumed uniform regularity bound on the network that is preserved across levels, the method is claimed to achieve total work that is a constant multiple of the work on the finest grid, with the training residual serving as an a posteriori error bound due to uniform conditioning of the discrete operator. Experiments on interior Dirichlet Laplace/Poisson and exterior Neumann Helmholtz problems are presented, comparing favorably to GMRES.

Significance. If the central claims hold, the work offers a novel GPU-mapped neural approach to BIEs that converts spectral bias into a convergence accelerator and supplies a residual-based a posteriori bound, potentially competitive with established iterative solvers for moderate-accuracy regimes. The explicit mapping of dense matrix-vector products and network evaluations to hardware is a practical strength.

major comments (2)
  1. [§4] §4 (complexity analysis) and the paragraph containing the uniform-regularity claim: the linear-work statement and the conversion of training residual to a posteriori error bound both rest on the assumption that a uniform regularity bound on the network is preserved across levels. No derivation is supplied showing that the warm-start from the previous level plus the stochastic-gradient updates on the refined Nyström system maintain this bound; without it the claimed O(work on finest grid) complexity and the residual-to-error equivalence do not follow.
  2. [§3.2] §3.2 (NTK spectral analysis): the argument that quadrature refinement returns stalled modes to the optimizer's reach relies on the empirical NTK being a faithful discrete sample of the continuum kernel. The manuscript does not quantify the approximation error between the finite-sample NTK and its continuum limit, nor does it show that this error remains controlled under the multilevel refinement schedule.
minor comments (2)
  1. [§3.1] Notation for the empirical NTK is introduced without an explicit definition of the mini-batch sampling measure; a one-line equation would remove ambiguity.
  2. [Figure 4] Figure 4 caption states “error versus iterations” but the plotted quantity is the training residual; relabeling or an added panel for true error would improve clarity.

Simulated Author's Rebuttal

2 responses · 0 unresolved

We thank the referee for the careful reading and constructive comments on our manuscript. We respond point-by-point to the major comments below. Where the comments identify gaps in justification or quantification, we commit to revisions that add discussion and clarification without altering the core claims.

read point-by-point responses

  1. Referee: [§4] §4 (complexity analysis) and the paragraph containing the uniform-regularity claim: the linear-work statement and the conversion of training residual to a posteriori error bound both rest on the assumption that a uniform regularity bound on the network is preserved across levels. No derivation is supplied showing that the warm-start from the previous level plus the stochastic-gradient updates on the refined Nyström system maintain this bound; without it the claimed O(work on finest grid) complexity and the residual-to-error equivalence do not follow.

    Authors: We acknowledge that the uniform regularity bound is stated as an assumption rather than derived from the warm-start and update procedure. The manuscript explicitly flags this assumption before stating the complexity result. In the revision we will expand §4 with a heuristic justification: because each level is warm-started from the converged parameters of the previous level and the NTK spectrum of the second-kind operator changes only mildly under quadrature refinement, the network parameters remain in a regime where the Lipschitz constant (hence regularity) is controlled by the same constants as the prior level. We will also add a remark that a complete proof of invariance would require tracking the evolution of the network's Sobolev norm under stochastic gradient flow on successively refined systems, which lies beyond the present scope but is consistent with the observed behavior in the experiments. The complexity claim therefore remains conditional on the assumption, which we will make more prominent. revision: partial

  2. Referee: [§3.2] §3.2 (NTK spectral analysis): the argument that quadrature refinement returns stalled modes to the optimizer's reach relies on the empirical NTK being a faithful discrete sample of the continuum kernel. The manuscript does not quantify the approximation error between the finite-sample NTK and its continuum limit, nor does it show that this error remains controlled under the multilevel refinement schedule.

    Authors: The spectral analysis in §3.2 is performed directly on the empirical NTK matrix constructed from the collocation points at each level; the continuum kernel is invoked only for interpretive context with the frequency principle. In the revised manuscript we will insert a short paragraph after the definition of the empirical NTK that recalls standard concentration bounds (e.g., from the NTK literature) showing that the empirical kernel converges to its continuum counterpart at rate O(1/√N) in operator norm when the network width is large and N is the number of collocation points. Because the multilevel schedule increases N at each level, the approximation error decreases, ensuring that the spectral properties used to explain mode recovery remain valid at every stage. We will also note that all convergence statements in the section are already stated for the empirical kernel, so the refinement argument does not depend on the continuum limit. revision: yes

Circularity Check

0 steps flagged

No significant circularity; claims rest on explicit assumptions and standard external facts

full rationale

The paper states its strongest claims conditionally ('Under a uniform regularity bound on the network...'). The NTK governing residual contraction is the standard empirical NTK arising from the network linearization, not a fitted quantity renamed as a prediction. Conversion of training residual to a posteriori error bound follows from the known uniform conditioning property of second-kind integral operators, an external mathematical fact. No derivation step reduces by construction to its own inputs, no self-citation is load-bearing, and no ansatz or uniqueness result is smuggled via prior work. The derivation chain is self-contained against the enumerated circularity patterns.

Axiom & Free-Parameter Ledger

0 free parameters · 0 axioms · 0 invented entities

Abstract-only review prevents identification of specific free parameters, axioms, or invented entities; the uniform regularity bound on the network is the only explicit modeling assumption stated.

pith-pipeline@v0.9.1-grok · 5798 in / 1219 out tokens · 21508 ms · 2026-07-02T08:09:33.441011+00:00 · methodology

discussion (0)

Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.

Reference graph

Works this paper leans on

19 extracted references

  1. [1]

    P. M. Anselone.Collectively Compact Operator Approximation Theory and Applications to Integral Equations. Prentice-Hall, Englewood Cliffs, NJ, 1971

    1971

  2. [2]

    Atanasov, B

    A. Atanasov, B. Bordelon, and C. Pehlevan. Neural networks as kernel learners: The silent alignment effect. InInternational Conference on Learning Representations (ICLR), 2022

    2022

  3. [3]

    K. E. Atkinson.The Numerical Solution of Integral Equations of the Second Kind. Cambridge Uni- versity Press, Cambridge, 1997

    1997

  4. [4]

    F. A. Bornemann and P. Deuflhard. The cascadic multigrid method for elliptic problems.Numerische Mathematik, 75(2):135–152, 1996

    1996

  5. [5]

    A. J. Burton and G. F. Miller. The application of integral equation methods to the numerical solution of some exterior boundary-value problems.Proceedings of the Royal Society of London. Series A, 323(1553):201–210, 1971

    1971

  6. [6]

    Chen and R

    C. Chen and R. Tsai. Implicit boundary integral methods for the Helmholtz equation in exterior do- mains.Research in the Mathematical Sciences, 4(19), 2017

    2017

  7. [7]

    Cheng, L

    H. Cheng, L. Greengard, and V . Rokhlin. A fast adaptive multipole algorithm in three dimensions. Journal of Computational Physics, 155(2):468–498, 1999

    1999

  8. [8]

    Colton and R

    D. Colton and R. Kress.Integral Equation Methods in Scattering Theory. SIAM, Philadelphia, 2013

    2013

  9. [9]

    Greengard and V

    L. Greengard and V . Rokhlin. A fast algorithm for particle simulations.Journal of Computational Physics, 73(2):325–348, 1987

    1987

  10. [10]

    Hackbusch.Multi-Grid Methods and Applications

    W. Hackbusch.Multi-Grid Methods and Applications. Springer, Berlin, 1985

    1985

  11. [11]

    Hackbusch

    W. Hackbusch. A sparse matrix arithmetic based onH-matrices. Part I: Introduction toH-matrices. Computing, 62(2):89–108, 1999

    1999

  12. [12]

    F. Izzo, O. Runborg, and R. Tsai. Corrected trapezoidal rules for singular implicit boundary integrals. Journal of Computational Physics, 461:111193, 2022

    2022

  13. [13]

    D. P. Kingma and J. Ba. Adam: A method for stochastic optimization. InInternational Conference on Learning Representations (ICLR), 2015

    2015

  14. [14]

    Kress.Linear Integral Equations

    R. Kress.Linear Integral Equations. Springer, New York, 3rd edition, 2014

    2014

  15. [15]

    Kublik, N

    C. Kublik, N. M. Tanushev, and R. Tsai. An implicit interface boundary integral method for Poisson’s equation on arbitrary domains.Journal of Computational Physics, 247:279–311, 2013

    2013

  16. [16]

    Rahaman, A

    N. Rahaman, A. Baratin, D. Arpit, F. Draxler, M. Lin, F. Hamprecht, Y . Bengio, and A. Courville. On the spectral bias of neural networks. InProceedings of the 36th International Conference on Machine Learning (ICML), pages 5301–5310, 2019

    2019

  17. [17]

    Saad and M

    Y . Saad and M. H. Schultz. GMRES: A generalized minimal residual algorithm for solving nonsym- metric linear systems.SIAM Journal on Scientific and Statistical Computing, 7(3):856–869, 1986

    1986

  18. [18]

    Sitzmann, J

    V . Sitzmann, J. N. P. Martel, A. W. Bergman, D. B. Lindell, and G. Wetzstein. Implicit neural repre- sentations with periodic activation functions. InProceedings of the 34th International Conference on Neural Information Processing Systems, NIPS ’20, Red Hook, NY , USA, 2020. Curran Associates Inc. 32 Neural BIE solver

    2020

  19. [19]

    Z.-Q. J. Xu, Y . Zhang, T. Luo, Y . Xiao, and Z. Ma. Frequency principle: Fourier analysis sheds light on implicit regularization of deep neural networks.Communications in Computational Physics, 28(5):1746–1767, 2020. 33

    2020