Pith. sign in

REVIEW 2 major objections 5 minor 22 references

Newtonian potentials of Legendre polynomials on rectangles have displacement structure

T0 review · 2 major / 5 minor · reviewed 2026-08-07 · deepseek-v4-flash

Pith's one-line read The paper proves that the Newtonian potential and its gradient, applied to tensor Legendre polynomials on a square, are entries of infinite matrices with displacement structure: a rank-3 Sylvester equation for the potential and a rank-1…

desk verdict Solid exact-recurrence result for Newtonian potentials on rectangles, with a minor boundary-case gap in the theorem statement and an honest but unproven stability caveat. read the letter →

arxiv 2506.03003 v1 pith:522VI56Q submitted 2025-06-03 math.NA cs.NA

classification math.NAcs.NA MSC 65D3065R10
keywords NewtonianpotentialPoissonequationLegendrepolynomialsdisplacementstructureSylvestersingularintegralscomplexlogarithmicintegralStieltjes
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

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

The reading

The paper aims to show that the Newtonian potential and its gradient, applied to tensor-product Legendre polynomials on a rectangle, are the real and imaginary parts of complex integrals that satisfy simple recurrence relations with low-rank inhomogeneous terms. Because the inhomogeneous terms have rank at most three for the potential and rank one for its gradient, the infinite arrays of these integrals have displacement structure. This means every entry can be reconstructed exactly from the first row and column by forward substitution, bypassing the logarithmic singularity entirely and avoiding specialized singular quadrature. If true, building particular solutions of Poisson's equation becomes a matter of solving a cheap recurrence, accurate for low-degree polynomial approximations on or near the domain and usable to moderate degree with high-precision arithmetic.

What carries the argument

The central object is the pair of infinite matrices $L(z)$ and $S(z)$, whose entries are the complex logarithmic and Stieltjes integrals of tensor Legendre polynomials over the square. The carrying mechanism is displacement structure: each matrix obeys a Sylvester equation of the form $AX - iXB^\top = G$ in which the right-hand side $G$ is low rank, rank three for $L$ and rank one for $S$. The tridiagonal matrices $A$ and $B$ are the multiplication-by-$x$ matrices for Legendre and ultraspherical polynomials, so the recurrences are just three-term polynomial recurrences applied in each index. The low-rank right-hand side is what allows all entries to be filled from boundary data: subtracting the two Sylvester equations for $L$ isolates the inhomogeneity into a rank-3 matrix $F(z)=F^{(2)}(z)-F^{(1)}(z)$, which is constructed explicitly from ultraspherical polynomials.

What would settle it

For a fixed interior point $z$, assemble a truncated $F(z)$ from the paper's formulas, fill $L_p(z)$ by direct high-precision quadrature of the integrals, and compute the singular values of $C L_p(z) - i L_p(z) C^\top$: the rank-three displacement claim stands only if exactly three singular values exceed round-off level.

Watch

Extended reading notes

Core claim

The central claim is that the infinite matrix $L(z)$ of complex logarithmic integrals $L_{kj}(z)=\iint_\Omega P_k(s)P_j(t)\log(z-(s+it))\,ds\,dt$ satisfies the Sylvester equation $C L(z) - i L(z) C^\top = F(z)$, where $C$ is a fixed tridiagonal matrix and $F(z)$ has rank at most three. A companion statement holds for the Stieltjes matrix $S(z)$: its entries satisfy $(zI - B)S(z) - i S(z)B^\top = 4e_0e_0^\top$, a rank-one equation. Since the real part of $L$ is the Newtonian potential and the real and imaginary parts of $S$ form its gradient, these identities reduce singular two-dimensional integration to algebraic recurrences: the first row and column are computed from one-dimensional integrals, and each subsequent row or column follows by a three-term recurrence. The recurrences compute the integrals exactly in exact arithmetic, so the only numerical error is round-off, which grows with polynomial degree.

Load-bearing premise

The practical claim rests on the forward recurrence being numerically reliable for the chosen polynomial degree and evaluation point, but the paper gives no rigorous rounding-error bound and its figures show errors growing roughly exponentially with the degree, with the usable degree dropping sharply for points far from the square.

Editorial extensions

If this is right

  • Newtonian potentials of polynomial densities on a square can be evaluated on or inside the domain without any quadrature rule for the logarithmic singularity, at a cost of $O(p^2)$ operations for degree $p$.
  • The same recurrences supply the gradient of the potential, because its components are the real and imaginary parts of the Stieltjes matrix $S(z)$.
  • For low-degree polynomial approximations, results are accurate to roughly machine precision on and near the square; double-word arithmetic extends machine-precision use to about degree 90--100 inside the square.
  • The recurrence approach can be combined with traditional quadrature away from the square, where the kernel is smooth and the recurrence loses accuracy.
  • Displacement structure opens the way to fast matrix-vector products and low-rank approximation of the potential matrices, and to solvers that exploit the Sylvester form directly.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

  • A natural next test is whether the same displacement-structure argument extends to triangles and parallelograms, and whether a different derivation could lift it to cubes and higher-dimensional power-law kernels.
  • Because the recurrence error grows exponentially with degree and with distance from the square, an a priori error estimate for the forward substitution would be a direct next step; without one, the practical switch-over point between recurrence and quadrature remains empirical.
  • The rank-three inhomogeneity suggests that the potential matrix is numerically low rank away from the domain boundary; a quantitative singular-value decay bound would let users truncate the matrix without computing all entries.
  • If the Sylvester structure survives for power-law kernels, then equilibrium measures and fractional integrals on squares could inherit the same recurrence machinery.
Share X Bluesky LinkedIn Reddit HN

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

2 major / 5 minor

Summary. The paper develops recurrences and Sylvester equations for the complex logarithmic integral L_kj(z) and the Stieltjes integral S_kj(z) of tensor products of Legendre polynomials on the unit square. The real part of L and the real/imaginary parts of S give the Newtonian potential and its gradient, so the recurrences offer a way to evaluate these singular integrals without singular quadrature. The main algebraic results are Lemma 4.1, giving a rank-1 Sylvester equation for S(z), and Theorem 4.4, giving a displacement equation C L(z) - i L(z) C^T = F(z) whose right-hand side has rank at most three. The paper also proposes Algorithms 5.1 and 5.2, which fill L_p(z) and S_p(z) from first-row and first-column data by forward substitution, and reports numerical experiments in double and double-word arithmetic showing exponential round-off growth in p.

Significance. If the stated algebraic structure is correct, this is a valuable contribution: it replaces singular quadrature for a useful class of volume potentials by explicit recurrences with low-rank displacement structure. The derivations are explicit and self-contained, the rank statements are parameter-free, and the paper ships an experimental Julia implementation. The main theorem is elegant and should interest researchers in numerical PDEs, orthogonal polynomials, and structured matrix computations. However, the practical usefulness depends on numerical stability of the forward recurrences, which is only demonstrated empirically, and there is a branch-cut omission that makes one of the stated exact recurrences incorrect on a whole line segment.

major comments (2)
  1. [Proposition 2.3 / Theorem 2.4 / Corollary 2.5] The case split omits the line x = 0 with -1 < y < 1. Proposition 2.3 covers x > 0 or y ≥ 1, x ≤ 0 and y ≤ -1, and x < 0 with -1 < y < 1; the pure-imaginary segment is in no case, and Theorem 2.4 and Corollary 2.5 inherit the omission. This is a genuine correctness gap, not a round-off issue. For example, at z = 0 one has M_1(0) = ∫_{-1}^{1} t log(-i t) dt = -iπ/2, whereas the uncorrected base relation would give L_1(0) = iπ/2, an error of -iπ. The omission is not harmless for the square: Algorithms 5.1 and 5.2 evaluate M_j(z-1) and M_j(z+1) to seed the first row and column, and for z on the right edge z = 1 + iy these arguments are exactly on the omitted line, so the boundary data are wrong and the recurrences in Corollaries 5.1 and 5.2 propagate the error through all entries. Including x = 0 in the third case (x ≤ 0) or defining M_k at x = 0 by the one-sided limit from x < 0 restores the formulas; as stated, the exact recurrences are not valid for all z ∈ Ω.
  2. [Section 5 and Figures 2 and 4] The practical claim that Algorithms 5.1 and 5.2 provide an accurate way to compute Newtonian potentials and their gradients rests on the forward recurrence from first-row and first-column data, but no rigorous round-off bound is supplied. The paper's own experiments show errors growing roughly exponentially with p, with the usable degree dropping sharply as z moves away from the square. This does not affect the algebraic displacement structure, but the abstract and Section 5 currently promise accurate moderate-degree evaluation; please either provide a stability analysis or explicitly restrict the claim to the tested regime and describe the method as heuristically stable.
minor comments (5)
  1. [Theorem 3.3 / Theorem 4.4] The polynomial C_k^{(-3/2)} appears in the definition of β_kj and in the rank argument, but the paper defines ultraspherical polynomials only for λ > -1/2 with an extension to λ = -1/2. Please give the definition used for λ = -3/2, for example via the DLMF hypergeometric expression, so the displayed formulas are unambiguous.
  2. [Remark 2.6] Remark 2.6 addresses signed zero in floating-point implementations of the logarithm, but it does not address the mathematical omission of x = 0, -1 < y < 1. Please clarify that the one-sided limit is the correct mathematical definition, independent of any floating-point artifact.
  3. [Theorem 4.4] The rank analysis for z outside the square is only indicated by 'similarly'. Please spell out the four exterior regions and verify that the rank bound remains at most three in each of them.
  4. [Algorithm 5.1] Algorithm 5.1 says that L_0(z ± i), ..., L_{p+1}(z ± i) can be computed via the recurrence in Corollary 2.5, but Corollary 2.5 is the recurrence for M_k, not L_k. The one-dimensional L_k values should be computed with Theorem 2.2; please correct this citation.
  5. [Various] There are several small typographical issues: 'leeds' should be 'leads' in Section 6, 'Stietljes' should be 'Stieltjes' in Section 1, and the expression in Proposition 2.3 written as '∫_y^{-1} -3∫_y^1 f(t)dt' is hard to parse and should be typeset as ∫_{-1}^{y} f(t) dt - 3 ∫_{y}^{1} f(t) dt.

Circularity Check

0 steps flagged · score 0.0 of 10

No significant circularity: recurrences and displacement structure are derived from standard orthogonal-polynomial identities and integration by parts, with no fitted inputs or load-bearing self-citations.

full rationale

The paper's central derivation is self-contained. Theorem 4.4 follows by subtracting two Sylvester equations obtained in Lemma 4.2, whose inhomogeneous terms F^(1) and F^(2) are constructed explicitly in Theorem 3.3 from one-dimensional recurrences, integration by parts, and explicit ultraspherical polynomial identities; Lemma 4.1 and Corollaries 3.5 and 5.1 similarly derive from the three-term Legendre recurrence and the definitions of the integrals. No parameter is fitted to data, no prediction is made from a fitted subset, and the cited prior work (Gautschi, DLMF, and even the author's own orthogonal-polynomial survey) supplies standard identities or algorithmic context rather than importing the claimed Sylvester equation. The omitted x=0 branch-cut case noted in the skeptical summary is a correctness gap in a case split, not a circular reduction: the formulas are asserted, not derived from the target result. The Julia package and timing experiments are illustrative implementation, not input to the proof. Hence the circularity score is 0.

Assumptions & free parameters 0 free parameters · 4 assumptions · 0 invented entities

The derivation uses only standard orthogonal polynomial identities, an explicit branch convention for the complex logarithm, and the assumption that low-degree Legendre expansions are adequate. No fitted parameters, no new physical entities, and no circular predictions are introduced.

assumptions (4)
  • standard math Legendre polynomials satisfy the three-term recurrence x P_k = k/(2k+1) P_{k-1} + (k+1)/(2k+1) P_{k+1} with P_{-1}=0.
    Used throughout Sections 2 and 3 to derive recurrences for Stieltjes and logarithmic integrals.
  • standard math The ultraspherical polynomial family C^(lambda)_n extends to lambda = -1/2 through the explicit hypergeometric and DLMF formulas cited in Section 2.
    Proposition 2.1 and Theorem 2.2 use this extension to connect logarithmic integrals to Stieltjes integrals.
  • domain assumption The complex logarithm branch convention in Proposition 2.3 and Corollary 2.5, including the modified branch log_minus, is the correct branch for representing the Newtonian potential as Re L.
    Proposition 2.3 and Remark 2.6 depend on this branch split; floating-point negative zero can silently select the wrong branch.
  • domain assumption In practical use, the right-hand side f can be approximated well by a low or moderate total degree tensor Legendre expansion.
    Section 5 shows accuracy degrades as p grows, so the practical value of the recurrence depends on low or moderate degree approximations.

how reviews work

0 comments
Cite this review

Pith. "Pith review of Newtonian potentials of Legendre polynomials on rectangles have displacement structure." pith.science (2026). https://pith.science/paper/522VI56Q

@misc{pith2026250603003,
  author       = {Pith},
  title        = {Pith review of: Newtonian potentials of Legendre polynomials on rectangles have displacement structure},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/522VI56Q}},
  note         = {Machine review of arXiv:2506.03003}
}
read the original abstract

Particular solutions of the Poisson equation can be constructed via Newtonian potentials, integrals involving the corresponding Green's function which in two-dimensions has a logarithmic singularity. The singularity represents a significant challenge for computing the integrals, which is typically overcome via specially designed quadrature methods involving a large number of evaluations of the function and kernel. We present an attractive alternative: we show that Newtonian potentials (and their gradient) applied to (tensor products of) Legendre polynomials can be expressed in terms of complex integrals which satisfy simple and explicit recurrences that can be utilised to exactly compute singular integrals, i.e., singular integral quadrature is completely avoided. The inhomogeneous part of the recurrence has low rank structure (its rank is at most three for the Newtonian potential) and hence these recurrences have displacement structure. Using the recurrence directly is a fast approach for evaluation on or near the integration domain that remains accurate for low degree polynomial approximations, while high-precision arithmetic allows accurate use of the approach for moderate degree polynomials.

Discussion (0). Sign in to comment.

Reference graph

Works this paper leans on

22 extracted references · 21 canonical work pages

  1. [21]

    Shen and K

    Z. Shen and K. Serkh , Rapid evaluation of Newtonian potentials on planar domains , SIAM Journal on Scientific Computing, 46 (2024), pp. A609–A628

  2. [1]

    af Klinteberg and A.-K

    L. af Klinteberg and A.-K. Tornberg , Adaptive quadrature by expansion for layer po- tential evaluation in two dimensions , SIAM Journal on Scientific Computing, 40 (2018), pp. A1225–A1249

  3. [2]

    T. G. Anderson, M. Bonnet, L. M. F aria, and C. P ´erez-Arancibia, Fast, high-order numerical evaluation of volume potentials via polynomial density interpolation , Journal of Computational Physics, 511 (2024), p. 113091

  4. [3]

    Atkinson, The numerical evaluation of particular solutions for Poisson ’s equation , IMA J

    K. Atkinson, The numerical evaluation of particular solutions for Poisson ’s equation , IMA J. Numer. Anal., 5 (1985), pp. 319–338

  5. [4]

    Ballew, T

    C. Ballew, T. Trogdon, and H. Wilber , The akhiezer iteration and an inverse-free solver for sylvester matrix equations , arXiv preprint arXiv:2503.17496, (2025)

  6. [5]

    A. H. Barnett, Evaluation of layer potentials close to the boundary for laplace and helmholtz problems on analytic planar domains , SIAM Journal on Scientific Computing, 36 (2014), pp. A427–A451

  7. [6]

    Beckermann and A

    B. Beckermann and A. Townsend , On the singular values of matrices with displacement structure, SIAM Journal on Matrix Analysis and Applications, 38 (2017), pp. 1227–1248

  8. [7]

    Beckermann and A

    B. Beckermann and A. Townsend , Bounds on the singular values of matrices with displace- ment structure, SIAM Review, 61 (2019), pp. 319–344

Show all 22 references
  1. [8]

    https://dlmf.nist.gov/, Release 1.2.4 of 2025- 03-15, https://dlmf.nist.gov/

    NIST Digital Library of Mathematical Functions. https://dlmf.nist.gov/, Release 1.2.4 of 2025- 03-15, https://dlmf.nist.gov/. F. W. J. Olver, A. B. Olde Daalhuis, D. W. Lozier, B. I. Schneider, R. F. Boisvert, C. W. Clark, B. R. Miller, B. V. Saunders, H. S. Cohl, and M. A. Mc...

  2. [9]

    Fortunato and A

    D. Fortunato and A. Townsend, Fast Poisson solvers for spectral methods , IMA Journal of Numerical Analysis, 40 (2020), pp. 1994–2018

  3. [10]

    Gautschi and J

    W. Gautschi and J. Wimp, Computing the Hilbert transform of a Jacobi weight function , BIT Numerical Mathematics, 27 (1987), pp. 203–215

  4. [11]

    Greengard and J.-Y

    L. Greengard and J.-Y. Lee , A direct adaptive Poisson solver of arbitrary order accuracy , Journal of Computational Physics, 125 (1996), pp. 415–424

  5. [12]

    Greengard, M

    L. Greengard, M. O’Neil, M. Rachh, and F. Vico, Fast multipole methods for the evaluation of layer potentials with locally-corrected quadratures , Journal of Computational Physics: X, 10 (2021), p. 100092

  6. [13]

    Gutleb, J

    T. Gutleb, J. Carrillo, and S. Olver , Computing equilibrium measures with power law kernels, Mathematics of Computation, 91 (2022), pp. 2247–2281

  7. [14]

    T. S. Gutleb, J. A. Carrillo, and S. Olver, Computation of power law equilibrium measures on balls of arbitrary dimension , Constructive Approximation, 58 (2023), pp. 75–120

  8. [15]

    Heinig and K

    G. Heinig and K. Rost , Algebraic methods for Toeplitz-like matrices and operators , vol. 13, Birkh¨ auser Basel, 1984

  9. [16]

    Joldes, J.-M

    M. Joldes, J.-M. Muller, and V. Popescu, Tight and rigorous error bounds for basic building blocks of double-word arithmetic , ACM Transactions on Mathematical Software (TOMS), 44 (2017), pp. 1–27

  10. [17]

    Olver, MultivariateSingularIntegrals.jl

    S. Olver, MultivariateSingularIntegrals.jl. v0.1, 2025

  11. [18]

    Olver, R

    S. Olver, R. M. Slevinsky, and A. Townsend, Fast algorithms using orthogonal polynomials, Acta Numerica, 29 (2020), pp. 573–699

  12. [19]

    Pu and M

    T. Pu and M. F asondini, The numerical solution of fractional integral equations via orthogonal 18 polynomials in fractional powers, Advances in Computational Mathematics, 49 (2023), p. 7

  13. [20]

    Sarnoff et al

    J. Sarnoff et al. , DoubleFloats.jl: math with more good bits . v1.2, 2025

  14. [22]

    R. M. Slevinsky and S. Olver , A fast and well-conditioned spectral method for singular integral equations, Journal of Computational Physics, 332 (2017), pp. 290–315. 19

Pith tools

Reviewed August 7, 2026 · model on record in the stance chip above.