Fast measure modification of orthogonal polynomials via matrices with displacement structure

TL;DR

This work addresses efficiently modifying and factorizing Gram matrices arising in orthogonal polynomial measure modification by exploiting displacement structure. The authors develop fast Cholesky methods based on a Sylvester-type displacement equation with rank-2 generators, enabling $O(n^2)$ work for principal blocks and $O(b\,n)$ when the Gram matrix is banded. They present two practical routes to obtain modified moments—recurrences from weight differential equations and simple-function approximations—plus specialized treatment for unbounded domains, and they validate the approach with numerical experiments across single/multiple algebraic factors and Laguerre-type weights. The paper further introduces hierarchical Cholesky factorization for $\mathcal{H}_\ell(r)$ matrices, with randomized low-rank techniques yielding $O(n\log^4 n)$-type precomputation and $O(r n\log^2 n)$ matvecs, and shows that the Chebyshev Gram matrix has a Toeplitz-plus-Hankel structure enabling FFT-based $O(N\log N)$ subblock operations. Overall, the framework offers substantial speedups for computing connection coefficients between original and modified OP families and opens avenues for multivariate extensions, Sobolev variants, and further exploitation of structured polynomial transforms.

Abstract

It is well known that matrices with low Hessenberg-structured displacement rank enjoy fast algorithms for certain matrix factorizations. We show how $n\times n$ principal finite sections of the Gram matrix for the orthogonal polynomial measure modification problem has such a displacement structure, unlocking a collection of fast algorithms for computing connection coefficients (as the upper-triangular Cholesky factor) between a known orthogonal polynomial family and the modified family. In general, the ${\cal O}(n^3)$ complexity is reduced to ${\cal O}(n^2)$, and if the symmetric Gram matrix has upper and lower bandwidth b, then the ${\cal O}(b^2n)$ complexity for a banded Cholesky factorization is reduced to ${\cal O}(b n)$. In the case of modified Chebyshev polynomials, we show that the Gram matrix is a symmetric Toeplitz-plus-Hankel matrix, and if the modified Chebyshev moments decay algebraically, then a hierarchical off-diagonal low-rank structure is observed in the Gram matrix, enabling a further reduction in the complexity of an approximate Cholesky factorization powered by randomized numerical linear algebra.

Figures (6)

• Figure 1: Comparison of performance between the ${\@fontswitch{}{\mathcal{}} O}(b^2n)$ algorithms of Gutleb-Olver-Slevinsky-1-24 (solid lines) and the ${\@fontswitch{}{\mathcal{}} O}(bn)$ algorithms in the present work (dashed lines). In both plots, the dashed lines correspond to the same $\delta$ value as the solid lines of the same colour. Left: the relative error measured in the Frobenius norm $\left\|{P_nW_PP_n^\top - (P_nRP_n^\top)^\top P_nRP_n^\top}\right\|_F/\left\|{P_nW_PP_n^\top}\right\|_F$. Right: calculation times of the construction of the Legendre--Gram matrix and its Cholesky factorization for four different values of $\delta$. The difference in the bandwidth scalings results in gaps between the solid lines being twice as wide as those between the dashed lines.
• Figure 1: Left: Numerical ranks of the hierarchical approximation of the principal finite section of the Chebyshev--Gram matrix for the log-Chebyshev weight $w(x) = \log\left(\dfrac{2}{1-x}\right)\dfrac{1}{\sqrt{1-x^2}}$. Right: numerical ranks of the corresponding Cholesky factor. In both panels, green blocks indicate low-rank approximations and red indicates that dense blocks are used. The opacity $\alpha\in(0,1)$ of a green block of size $m\times n$ and numerical rank $r$ is proportional to its data-sparsity compared to a dense filling: $\alpha mn = (m+n)r$.
• Figure 2: A weight with multiple algebraic factors on $(-1,1)$. Left: the relative error in the modified Chebyshev moments. This is calculated by comparison with an extended precision calculation using MPFR Fousse-et-al-33-13-1-07 following the same steps. Right: synthesis of $q_{500}(x)$ and the Szegő evelope Szego-75.
• Figure 2: Conversion of a degree-$n-1$ expansion in log-Chebyshev OPs with standard normally distributed pseudorandom coefficients to Chebyshev polynomials. Left: $2$-norm and $\infty$-norm relative error in the forward and backward transformation. Right: calculation times of the hierarchical construction of the Chebyshev--Gram matrix, its hierarchical Cholesky factorization, and the matrix-vector product with the hierarchical Cholesky factor.
• Figure 3: Modified Laguerre polynomials. Left: the eigenvalues of the connection matrix $P_{512}RP_{512}^\top$. Right: three of the modified orthogonal polynomials multiplied by $\sqrt{w(x)}$. At low degree, this function is highly localized to the interval $[0,4]$, whereas at high degree, the function is not so different from $e^{-x/2}L_n(x)$, with a phase shift in the tail.

Theorems & Definitions (36)

• Definition 1.1
• Proposition 1.2
• Proof 1
• Definition 2.1
• Definition 2.2
• Definition 2.3
• Theorem 2.4: Gautschi Gautschi-24-245-70
• Proof 2
• Remark 2.5
• Lemma 2.6