Generalized Pseudospectral Method and Zeros of Orthogonal Polynomials

TL;DR

This paper extends the pseudospectral method to derive new algebraic identities for the zeros of a broad class of orthogonal polynomials, including nonclassical types, by modifying matrix representations of differential operators.

Contribution

It introduces a generalized pseudospectral approach that depends on two sets of nodes, enabling the derivation of identities for zeros of diverse orthogonal polynomial families.

Findings

Derived nonlinear algebraic identities for polynomial zeros.

Generalized identities apply to classical and nonclassical orthogonal polynomials.

Presented new relations for zeros of Sonin-Markov (generalized Hermite) polynomials.

Abstract

Via a generalization of the pseudospectral method for numerical solution of differential equations, a family of nonlinear algebraic identities satisfied by the zeros of a wide class of orthogonal polynomials is derived. The generalization is based on a modification of pseudospectral matrix representations of linear differential operators proposed in the paper, which allows these representations to depend on two, rather than one, sets of interpolation nodes. The identities hold for every polynomial family $\{p_\nu(x)\}_{\nu=0}^\infty$ orthogonal with respect to a measure supported on the real line that satisfies some standard assumptions, as long as the polynomials in the family satisfy differential equations $\mathcal{A} p_\nu(x) =q_\nu(x) p_\nu(x)$, where $\mathcal{A}$ is a linear differential operator and each $q_\nu(x)$ is a polynomial of degree at most $n_0 \in \mathbb{N}$; $n_0$ does not depend on $\nu$. The proposed identities generalize known identities for classical and Krall orthogonal polynomials, to the case of the nonclassical orthogonal polynomials that belong to the class described above. The generalized pseudospectral representations of the differential operator $\mathcal{A}$ for the case of the Sonin-Markov orthogonal polynomials, also known as generalized Hermite polynomials, are presented. The general result is illustrated by new algebraic relations satisfied by the zeros of the Sonin-Markov polynomials.

Equations

$$\mathcal{A}{p}_{\nu }(x)=q_{\nu }(x)p_{\nu }(x),$$

$$\mathcal{A}{\mathbb{P}}^{\nu }\subseteq {\mathbb{P}}^{\nu +n_0}$$

$$\mathcal{A}u=f,$$

$${\pi }_{N}g=\left(g\left(x_1^{(N)}\right),\dots ,g\left(x_N^{(N)}\right)\right).$$

$${\pi }_N^{-1}{\vec{g}}^{(N)}=\sum _{j=1}^N{g}_j^{(N)}{\ell }_{N-1,j}(x).$$

$${\pi }_N\mathcal{A}g={\tilde{A}}^c{\pi }_{N}g$$

$$\tilde{A}^{c}_{kj}=\mathcal{A}\ell_{N-1,j}(x)\big{|}_{x=x_{k}^{(N)}},$$

$${\tilde{A}}^c{\vec{u}}^{(N)}={\pi }_{N}f$$

$${\ell }_{\rho -1,j}(x)=\frac{{\psi }_{\rho }(x)}{{\psi }_{\rho }^{\prime }\left(x_j^{(\rho )}\right)\left(x-x_j^{(\rho )}\right)},$$

$$A_{kj}^c=(\mathcal{A}{\ell }_{N-1,j})\left(x_k^{(N+n0)}\right),1\le k\le N+n_0,1\le j\le N,$$

$${\pi }_{N+n_0}\mathcal{A}g=A^c{\pi }_{N}g$$

$$A^c{\vec{u}}^{(N)}={\pi }_{N+n_0}f$$

$$A_{kj}^{\tau }=\frac{⟨\mathcal{A}{p}_{j-1},p_{k-1}⟩}{\Vert p_{k-1}{\Vert }^2},1\le k\le N+n_0,1\le j\le N.$$

$$L^{(N+n_0-1)}{A}^c=A^{\tau }{L}^{(N-1)},$$

$$L_{mj}^{(\rho -1)}=\frac{p_{m-1}(x_j^{(\rho )})}{\Vert p_{m-1}{\Vert }^2}{\lambda }_j^{(\rho -1)},1\le m,j\le \rho ,$$

$${\lambda }_j^{(\rho -1)}=\int {\ell }_{\rho -1,j}(x)d\omega ,1\le j\le \rho ,$$

$${\left\{{\left[L^{(\rho -1)}\right]}^{-1}\right\}}_{mj}=p_{j-1}\left(x_m^{(\rho )}\right),1\le m,j\le \rho .$$

$$\sum _{k=1}^N{A}_{mk}^c\left({\vec{x}}^{(N)},{\vec{x}}^{(N+n_0)}\right)p_{n-1}\left(x_k^{(N)}\right)$$

$$=\sum _{k=\max \{n-n_0,1\}}^{\min \{n+n_0,N+n_0\}}{A}_{kn}^{\tau }{p}_{k-1}\left(x_m^{(N+n_0)}\right),$$

$${\left[{\vec{e}}_m^{(N+n_0)}\right]}^T{A}^c{\vec{v}}_n^{(N)}={\vec{w}}_m^{(N+n_0)}{A}^{\tau }{\vec{e}}_n^{(N)}$$

$$\sum _{k=1}^N{A}_{mk}^c\left({\vec{x}}^{(N)}\right)p_{n-1}\left(x_k^{(N)}\right)=q_{n-1}{p}_{n-1}(x_m^{(N)})$$

$$p_{2n}^{\beta }(x)=c_n{L}_n^{\alpha }(x^2),$$

$$p_{2n+1}^{\beta }(x)=d_{n}x{L}_n^{\alpha +1}(x^2),$$

$$\alpha =(\beta -1)/2$$

$$c_n=(-1{)}^n\sqrt{\frac{n!}{\Gamma (n+1+\alpha )}},d_n=(-1{)}^n\sqrt{\frac{n!}{\Gamma (n+2+\alpha )}}$$

$$K_{2n}={\left[n!\Gamma (n+1+\alpha )\right]}^{-1/2},$$

$$K_{2n+1}={\left[n!\Gamma (n+2+\alpha )\right]}^{-1/2}.$$

$$\mathcal{D}{p}_{\nu }(x)=({\mu }_{\nu }{x}^2+{\eta }_{\nu })p_{\nu }(x),$$

$$\mathcal{D}=x^2\frac{d^2}{d{x}^2}+x(-2x^2+2\alpha +1)\frac{d}{dx}$$

$${\mu }_{\nu }=-2\nu ,$$

$$\displaystyle\eta_{\nu}=\left\{\begin{array}[]{l}0\mbox{ if }\nu=2n,\\ 2\alpha+1\mbox{ if }\nu=2n+1.\end{array}\right.$$

$$Z^{(k)}_{mn}=\left[\frac{d^{k}}{dx^{k}}\ell_{N-1,n}(x)\right]\Bigg{|}_{x=x_{m}^{(N+2)}}.$$

$$D_{mn}^c$$

Full text

Keywords: Pseudospectral methods; Spectral methods; Zeros of polynomials; Nonclassical orthogonal polynomials; Sonin-Markov polynomials; Generalized Hermite polynomials.

MSC2010 33C45; 65L60; 33C47; 26C10.

1 Introduction and Main Results

1.1 Summary of results

In this paper we identify a class of algebraic relations satisfied by the zeros of a wide class of orthogonal polynomials. To prove the identities, we generalize the notion of pseudospectral matrix representations of linear differential operators, by allowing these representations to depend on two, rather than one, sets of interpolation nodes. The identities hold for all polynomials $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$ orthogonal with respect to a measure satisfying some standard assumptions, as long as they satisfy the differential equations

[TABLE]

where $\mathcal{A}$ is a linear differential operator and each $q_{\nu}(x)$ is a polynomial of degree at most $n_{0}\in\mathbb{N}$; $n_{0}$ does not depend on $\nu$. This includes classical orthogonal polynomials [1, 2], polynomials in the Askey scheme [3] and Krall polynomials [4], as well as additional classes of nonclassical orthogonal polynomials. If applied to classical or Krall orthogonal polynomials, the proposed result reduces to known identities [5, 6, 7], see Subsection 1.4.

The motivation of this study stems from understanding that zeros of orthogonal polynomials play an important role in mathematical physics, numerical analysis and related areas. For example, zeros of some orthogonal polynomials are equilibria of important N-body problems [8, 9, 10, 11]. They transpire as building blocks of remarkable isospectral matrices [12, 13, 14, 15] and play an important role in construction of highly accurate approximation schemes for numerical integration [16, 17, 18, 19].

To prove algebraic identities satisfied by the zeros of the polynomials $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$, we generalize and relate the notions of spectral and pseudospectral matrix representations of linear differential operators used in the corresponding numerical methods for solving differential equations. The standard pseudospectral matrix representations of linear differential operators are based on Lagrange collocation on the real line. These representations were proposed by Calogero in the context of numerical solving of eigenvalue and boundary value problems for linear ODEs [20, 21], see also [22, 23, 24, 8, 25]. Mitropolsky, Prykarpatsky and Samoylenko set up a general algebraic-projection framework based on Calogero’s method and furthered its applications to solution of evolution equations in Mathematical Physics [26]. The convergence analysis of Calogero’s method was studied in [27].

The standard pseudospectral method was utilized in [7] to prove new properties of the zeros Krall polynomials. While Krall polynomials are eigenfunctions of linear differential operators, the polynomial families considered in this paper satisfy differential equations (1) with $q_{n}(x)$ being polynomials (as opposed to eigenvalues) of degree $n_{0}$ that does not depend on $n$. To prove new properties of the zeros of the latter polynomial families, we propose a generalization of the standard pseudospectral method.

Because, in general, the differential operator $\mathcal{A}$ in (1) raises the degree of polynomials by a summand of $n_{0}$, the standard pseudospectral method does not allow exact discretization of differential equations (1). The main idea of the proposed generalization is to construct Lagrange collocation type matrices that exactly represent linear differential operators acting between spaces of polynomials of different degrees. In this paper, the last goal is achieved by allowing the matrix representations to depend on two rather than one set of interpolation nodes. We thus find exact discretizations of the differential equations (1) satisfied by the polynomials $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$; the discretizations are constructed using the zeros of these polynomials as the nodes. By comparing the generalized pseudospectral and the generalized spectral matrix representations of the differential operator $\mathcal{A}$, we derive a family of algebraic identities satisfied by the zeros of polynomials $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$.

The proposed generalization of the pseudospectral matrix representations of linear differential operators has applications beyond those outlined in this paper; one such application allows to simplify the process of incorporation of initial or boundary conditions into linear systems that discretize certain ODEs, see the discussion in Section 4.

We illustrate the general result of the main Theorem 1.1 by applying it to the case of the Sonin-Markov polynomials, also known as generalized Hermite polynomials, see  [28, 29, 16] and references therein. These zeros play an important role in the computation of integrals of singular or oscillatory functions [17, 18] as well as extended Lagrange interpolation on the real line [19].

To prove the identities satisfied by the zeros of the Sonin-Markov polynomials stated in Theorem 2.1, we compute the generalized pseudospectral matrix representations of the differential operators $d/dx$ and $d^{2}/dx^{2}$ as well as the differential operator $\mathcal{D}$ associated with the Sonin-Markov family, see (20b) of (20b) and Subsection 2.2, assuming that the interpolation nodes are zeros of the Sonin-Markov polynomials. The formulas for these matrix representations as well as the identities of Theorem 2.1 have been verified using programming environment MATLAB, for several particular values of the relevant paramenets (the degree $N$ and the parameter $\beta$, see Theorem 2.1 and definitions (16), (17)).

1.2 The orthogonal polynomial family $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$

Let $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$ be a sequence of polynomials orthogonal with respect to a measure $\omega$ and the corresponding inner product $\langle f,g\rangle=\int fg\;d\omega$. We denote the norm associated with this inner product by $\|\cdot\|$, that is, $\|f\|^{2}=\int f^{2}\;d\omega$. Assume that $\omega$ is a Borel measure with support on the real line satisfying the following three conditions:

• (a)

$\omega$ is positive;

• (b)

all its moments $\int x^{\nu}\,d\omega$ exist and are finite;

• (c)

$\omega$ has infinitely many points in its support $I=\operatorname{supp}\omega$.

Under the above assumptions on the measure $\omega$, the zeros of each polynomial $p_{\nu},\;\nu\geq 1,$ are real, simple and belong to the convex hull of the support of $\omega$, see for example [16].

Notation 1*.*

Here and throughout the rest of the paper $N$ denotes a fixed integer strictly larger than 1, while $n_{0}$ is a fixed nonnegative integer. The small Greek letter $\rho$ denotes an index that may take values $N$ or $N+n_{0}$. The small Greek letter $\nu$ denotes an integer index that usually takes values $0,1,2\ldots$, unless otherwise indicated. The small Latin letters $n,m,j,k$ etc. denote integer indices that usually run from $1$ to $N$ or from $1$ to $N+n_{0}$, see (2) and (1), thus we indicate the range of the indices each time they are used. We reserve the letter $\ell$ to denote polynomials in Lagrange interpolation bases.

Let $\mathbb{P}^{\nu}$ denote the space of all algebraic polynomials with real coefficients of degree at most $\nu$. Assume that for each $\nu$, the polynomials $\{p_{j}(x)\}_{j=0}^{\nu}$ form a basis of $\mathbb{P}^{\nu}$. Let $\mathcal{A}$ be a linear differential operator acting on functions of one variable. Assume that $\mathcal{A}$ has the property

[TABLE]

for all $\nu$. Recall that $n_{0}$ is a fixed nonnegative integer; it does not depend on $\nu$. For example, the differential operator $\mathcal{D}=a_{0}+a_{1}(x)\frac{d}{dx}+\ldots+a_{q}(x)\frac{d^{q}}{dx^{q}}$ with $q\in\mathbb{N}$ and $a_{j}(x)\in\mathbb{P}^{j}$ for all $j=1,2,\ldots,q$ has property (2) with $n_{0}=0$, while the operator $x^{2}\mathcal{D}$ has property (2) with $n_{0}=2$.

Suppose that the orthogonal polynomials $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$ satisfy differential equations (1).

1.3 Generalized pseudospectral and spectral matrix representations of linear differential operators

In this Subsection we introduce the notions of a generalized pseudospectral and a generalized spectral matrix representations of the linear differential operator operator $\mathcal{A}$ introduced in the previous Subsection. Note that, in general, the definitions of these generalized matrix representations hold for any linear differential operator $\mathcal{D}$, which may or may not satisfy property (2), and for $n_{0}$ being a (positive or negative) integer such that $N+n_{0}\geq 1$.

The definition of the standard pseudospectral $N\times N$ matrix representation $\tilde{A}^{c}$ of $\mathcal{A}$, see [30, 31, 7, 25], is motivated by a search for an exact discretization of a differential equation

[TABLE]

under the assumption that it possesses a polynomial solution $u\in\mathbb{P}^{N-1}$. More precisely, choose a vector $\vec{x}^{(N)}=\left(x^{(N)}_{1},\ldots,x^{(N)}_{N}\right)$ of distinct real nodes and define the isomorphism $\pi_{N}:\mathbb{P}^{N-1}\to\mathbb{R}^{N}$ by

[TABLE]

The inverse of $\pi_{N}$ is of course given in terms of the standard Lagrange interpolation basis $\{\ell_{N-1,j}(x)\}_{j=1}^{N}$ of $\mathbb{P}^{N-1}$ constructed using the nodes $\vec{x}^{(N)}$: for every vector $\vec{g}^{(N)}=\left(g^{(N)}_{1},\ldots,g^{(N)}_{N}\right)\in\mathbb{R}^{N}$,

[TABLE]

The standard pseudospectral $N\times N$ matrix representation $\tilde{A}^{c}$ of the differential operator $\mathcal{A}$ is defined as the unique $N\times N$ matrix that satisfies the condition

[TABLE]

for all $g\in\mathbb{P}^{N-1}$. Note that the superscript ‘$c$’ in the notation of the matrix $\tilde{A}^{c}$ stands for “collocation” in the spectral collocation method for numerical solving of differential equations, also known as the pseudospectral method [25]. It is not difficult to conclude that the matrix $\tilde{A}^{c}$ is given componentwise by

[TABLE]

where $1\leq k,j\leq N$, see [20, 21, 26, 27, 31, 30, 7]. The last definition implies that a vector $\vec{u}^{(N)}\in\mathbb{R}^{N}$ solves the linear system

[TABLE]

if and only if the polynomial $u(x)=(\pi_{N})^{-1}\vec{u}^{(N)}=\sum_{j=1}^{N}{u}^{(N)}_{j}\ell_{N-1,j}(x)$ solves ODE (3). Of course, if ODE (3) does not possess a polynomial solution $u\in\mathbb{P}^{N-1}$, linear system (5) approximates ODE (3) and its solution $\vec{u}^{(N)}$, if it exists, allows to construct an approximate solution of ODE (3) given by $u(x)=(\pi_{N})^{-1}\vec{u}^{(N)}$.

Let us now generalize this notion of pseudospectral matrix representation of the differential operator $\mathcal{A}$ to take advantage of its property (2).

In addition to the interpolation nodes $\vec{x}^{(N)}$, consider another vector of distinct real nodes $\vec{x}^{(N+n_{0})}=\left(x^{(N+n_{0})}_{1},\ldots,x^{(N+n_{0})}_{N+n_{0}}\right)$. In short, we will work with two vectors of nodes $\vec{x}^{(\rho)}$, where $\rho=N$ or $\rho=N+n_{0}$, see Notation 1, and the respective Lagrange interpolation bases $\{\ell_{\rho-1,j}(x)\}_{j=1}^{\rho}$. Recall that for each $j\in\{1,\ldots,\rho\}$,

[TABLE]

where $\psi_{\rho}(x)=\left(x-x_{1}^{(\rho)}\right)\left(x-x_{2}^{(\rho)}\right)\cdots\left(x-x_{\rho}^{(\rho)}\right)$ is the node polynomial.

We define the generalized pseudospectral $(N+n_{0})\times N$ matrix representation $A^{c}\left(\vec{x}^{(N)},\vec{x}^{(N+n_{0})}\right)\equiv A^{c}$ of the linear differential operator $\mathcal{A}$ componentwise by

[TABLE]

where, as before, the superscript ‘$c$’ stands for “collocation”. This definition is motivated by the following relation:

[TABLE]

for all $g\in\mathbb{P}^{N-1}$, where the isomorphisms $\pi_{\rho}:\mathbb{P}^{\rho-1}\to\mathbb{R}^{\rho}$ are defined by (4) with $N$ replaced by $\rho$, $\rho\in\{N,N+n_{0}\}$. In other words, if ODE (3) has a polynomial solution $u\in\mathbb{P}^{N-1}$, then the vector $\vec{u}^{(N)}\in\mathbb{R}^{N}$ solves the system of linear equations

[TABLE]

if and only if the polynomial $u(x)=(\pi_{N})^{-1}\vec{u}^{(N)}=\sum_{j=1}^{N}{u}^{(N)}_{j}\ell_{N-1,j}(x)$ solves ODE (3).

Using analogous motivation, we define the generalized spectral $(N+n_{0})\times N$ matrix representation $A^{\tau}$ of the linear differential operator $\mathcal{A}$ componentwise by

[TABLE]

Here, the superscript ‘$\tau$’ indicates that the $\tau$-variant of the spectral method is used [25].

We prove that the $(N+n_{0})\times N$ matrices $A^{\tau}$ and $A^{c}$ satisfy the following property:

[TABLE]

where each of the two $\rho\times\rho$ matrices $L^{(\rho-1)}$ with $\rho=N$ or $\rho=N+n_{0}$ is the transition matrix from the orthogonal polynomial basis $\{p_{j}(x)\}_{j=0}^{\rho-1}$ to the Lagrange interpolation basis $\{\ell_{\rho-1,j}\}_{j=1}^{\rho}$, see Theorem 3.1 in Section 3.

1.4 Main result: algebraic identities satisfied by the zeros of the polynomials $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$

Let us now assume that $x_{1}^{(N)},\ldots,x_{N}^{(N)}$ are the zeros of the polynomial $p_{N}(x)$ from the orthogonal family $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$, while $x_{1}^{(N+n_{0})},\ldots,x_{N+n_{0}}^{(N+n_{0})}$ are the zeros of the polynomial $p_{N+n_{0}}(x)$. We therefore use these two sets of zeros as the nodes in the definition of the pseudospectral matrix representation $A^{c}$ of the differential operator $\mathcal{A}$, see (7) and (6). In this case, each of the two matrices $L^{(\rho-1)}$ in the relation (9), where $\rho=N$ or $\rho=N+n_{0}$, can be expressed in terms of the values $p_{m-1}(x_{j}^{(\rho)})$ and the Christoffel numbers $\lambda_{j}^{(\rho-1)}$:

[TABLE]

where the Christoffel numbers are defined by

[TABLE]

see Theorem 3.1.

Recall that Christoffel numbers arise in the Gaussian quadrature numerical integration formulas; they are always positive [16]. Christoffel numbers play an important role in the proof of the main identity (13) presented in this paper, although they are eliminated from that identity in the process of inversion of the matrices $L^{(\rho-1)}$: the inverses of $L^{(\rho-1)}$ are given componentwise by

[TABLE]

see (45) in the proof of Theorem 3.1.

Using property (9) of the matrix representations $A^{\tau}$ and $A^{c}$, together with the neat formulas (12) for the matrices $\left[L^{(\rho-1)}\right]^{-1}$ in the case where $\vec{x}^{(\rho)}$ are the zeros of $p_{\rho}(x)$, $\rho\in\{N,N+n_{0}\}$, we prove the following algebraic identities satisfied by the zeros of the polynomials in the family $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$.

Theorem 1.1**.**

The zeros $\vec{x}^{(N)}=\left(x_{1}^{(N)},\ldots,x_{N}^{(N)}\right)$ of the polynomial $p_{N}(x)$ and the zeros $\vec{x}^{(N+n_{0})}=\left(x_{1}^{(N+n_{0})},\ldots,x_{N+n_{0}}^{(N+n_{0})}\right)$ of the polynomial $p_{N+n_{0}}(x)$ in the orthogonal polynomial family $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$ of generalized eigenfunctions of the linear differential operator $\mathcal{A}$, see (1), satisfy the following algebraic relations for all integer $m,n$ such that $1\leq m\leq N+n_{0}$, $1\leq n\leq N$:

[TABLE]

where the $(N+n_{0})\times N$ pseudospectral and spectral matrix representations $A^{c}=A^{c}\left(\vec{x}^{(N)},\vec{x}^{(N+n_{0})}\right)\equiv A^{c}$ and $A^{\tau}$, respectively, are defined by (7) and (8).

Remark 1*.*

For every pair of integers $n,m$ such that $0\leq n\leq N-1$ and $1\leq m\leq N+n_{0}$, identity (13) relates the zeros of the polynomials $p_{N}(x)$, $p_{N+n_{0}}(x)$ and $p_{n}(x)$ with the zeros of all the polynomials $p_{k}(x)$ such that the index $k\in\{0,1,\ldots,N+n_{0}-1\}$ satisfies $n-n_{0}\leq k\leq n+n_{0}$.

Remark 2*.*

The main identity of Theorem 1.1 may be recast as follows:

[TABLE]

for all $1\leq m\leq N+n_{0}$ and $1\leq n\leq N$, where each $\vec{v}_{n}^{(N)}\in\mathbb{R}^{N}$ is a column-vector with the components $\left[\vec{v}_{n}^{(N)}\right]_{k}=p_{n-1}\left(x_{k}^{(N)}\right)$, each $\vec{w}_{m}^{(N+n_{0})}\in\mathbb{R}^{N+n_{0}}$ is a row-vector with the components $\left[\vec{w}_{m}^{(N+n_{0})}\right]_{k}=p_{k-1}\left(x_{m}^{(N+n_{0})}\right)$ and the column-vectors $\{\vec{e}^{(\rho)}_{k}\}_{k=1}^{\rho}$ form the standard basis of $\mathbb{R}^{\rho}$.

Theorem 1.1 is proved in Section 3.

Identities (13) or, equivalently, (14), are remarkable in the sense that they reveal a deeper structure and relation between the linear operators $A^{c}$ and $A^{\tau}$, the generalized spectral and pseudospectral representations, respectively, of the differential operator $\mathcal{A}$, for the case where these operators are constructed using two sets of zeros of orthogonal polynomials as the nodes.

Let us compare the main result of Theorem 1.1 with other results of this kind. By setting $n_{0}=0$ in identities (13), (14), we obtain that

[TABLE]

for all $m,n$ such that $1\leq m,n\leq N$. In this case, of course, $q_{\nu}$ are constants, see (1), and $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$ are either classical or Krall orthogonal polynomials. In this special case where $n_{0}=0$, identities of Theorem 1.1 reduce to those reported in [7]. It was shown in [7] that, if applied to the classical Jacobi, Hermite or Laguerre polynomials, identities (15) reduce to known identities for the zeros of these polynomials reported in [5, 6].

Therefore, identities (13), (14) generalize similar results for classical orthogonal polynomials proved in [5, 6] and for Krall polynomials proved in [7]. These identities may be considered as analogues of the properties of the zeros of the Askey scheme and generalized hypergeometric polynomials proved in [12, 13, 14, 15], for the case of the polynomial families considered in this paper. An application of the identities proved in this paper is related to the study of the asymptotic behavior of algebraic expressions involving the zeros of orthogonal polynomials of degree $N$ as $N\to\infty$, see [32, 33].

In the next Section 2 we apply Theorem 1.1 to prove new identities satisfied by the zeros of the nonclassical Sonin-Markov orthogonal polynomials. In Section 3, “Proofs”, we elaborate on the proofs of most of the theorems of this paper, except for those that are straightforward consequences of another theorem. In Section 4 titled “Discussion and Outlook” we summarize the results proposed in this paper and discuss their importance, possible applications and further developments.

2 Application: Properties of the Zeros of the Sonin-Markov Orthogonal Polynomials

In this Section 2 we illustrate Theorem 1.1 by applying it to the case of the Sonin-Markov orthogonal polynomials, which are generalized eigenfunctions of a certain linear differential operator $\mathcal{D}$, see (20b). This application requires a computation of the generalized spectral and the generalized pseudospectral matrix representations of the differential operator $\mathcal{D}$ associated with the Sonin-Markov polynomials. We thus proceed as follows. In Subsection 2.1 we define the Sonin-Markov polynomials and state their basic properties. In Subsections 2.2 and 2.3, respectively, we compute the generalized pseudospectral and the generalized spectral representations, respectively, of the differential operator $\mathcal{D}$. Finally, in Subsection 2.4 we provide a family of new algebraic properties satisfied by the zeros of the Sonin-Markov polynomials.

2.1 Definition and Basic Properties of the Sonin-Markov Polynomials

The Sonin-Markov polynomials $\{p^{\beta}_{\nu}(x)\}_{\nu=0}^{\infty}$, which are also known as generalized Hermite polynomials, are orthogonal on $I=(-\infty,\infty)$ with respect to the weight $w(x)=|x|^{\beta}e^{-x^{2}}$ (see [28, 16] and references therein). These polynomials can be expressed in terms of the generalized Laguerre polynomials as follows:

[TABLE]

where

[TABLE]

and the coefficients

[TABLE]

are chosen to make the polynomials orthonormal and the leading coefficients positive. Note that the leading coefficients $K_{\nu}$ of $p_{\nu}^{\beta}(x)$ are given by

[TABLE]

The zeros of the Sonin-Markov polynomials are distinct, real and symmetric with respect to the origin.

The Sonin-Markov polynomials satisfy the following differential equations:

[TABLE]

and

[TABLE]

These differential equations can be derived from the corresponding differential equations satisfied by the generalized Laguerre polynomials stated in [3], see also [1] and formula (3.5) in [34]. Clearly, $\mathcal{D}\mathbb{P}^{\nu}\subseteq\mathbb{P}^{\nu+2}$, so we can apply Theorem 1.1 with $n_{0}=2$ to these polynomials.

2.2 The Generalized Pseudospectral Matrix Representation of the Sonin-Markov Differential Operator (20b)

Let $x_{1}^{(N)},\ldots,x_{N}^{(N)}$ be the zeros of the Sonin-Markov polynomial $p_{N}(x)\equiv p_{N}^{\beta}(x)$ defined by (16) and let $x_{1}^{(N+2)},\ldots,x_{N+2}^{(N+2)}$ be the zeros of the Sonin-Markov polynomial $p_{N+2}(x)\equiv p_{N+2}^{\beta}(x)$. In this subsection we find the generalized pseudospectral matrix representation $D^{c}\left(\vec{x}^{(N)},\vec{x}^{(N+2)}\right)\equiv D^{c}$ of the Sonin-Markov differential operator (20b) with respect to the two vectors of nodes $\vec{x}^{(N)}=\left(x_{1}^{(N)},\ldots,x_{N}^{(N)}\right)$ and $\vec{x}^{(N+2)}=\left(x_{1}^{(N+2)},\ldots,x_{N+2}^{(N+2)}\right)$. To compute this representation $D^{c}$, we use definition (7) with $n_{0}=2$. In the following, we use the notation of Subsection 2.1, note the appropriate definitions of $\alpha$, $c_{n},d_{n}$ and $\mu_{\nu},\eta_{\nu}$.

Let $\{\ell_{N-1,j}(x)\}_{j=1}^{N}$ be the Lagrange interpolation basis with respect to the nodes $\{x_{j}^{(N)}\}_{j=1}^{N}$ defined by (6) with $\rho=N$ and $\psi_{N}(x)$ replaced by $p_{N}(x)$. Let $Z^{(k)}\left(\vec{x}^{(N)},\vec{x}^{(N+2)}\right)\equiv Z^{(k)}$ be the $(N+2)\times N$ generalized pseudospectral matrix representation of the differential operator $\frac{d^{k}}{dx^{k}}$ with respect to the two vectors of nodes $\vec{x}^{(N)}$ and $\vec{x}^{(N+2)}$. By definition (7), its components $Z^{(k)}_{mn}$ are given by

[TABLE]

The generalized pseudospectral $(N+2)\times N$ matrix representation $D^{c}$ of the differential operator $\mathcal{D}$ is given componentwise by

[TABLE]

see definition (7).

By using the fact that the Sonin-Markov polynomial $p_{N}(x)$ satisfies differential equation (20b) and the differential equation obtained by differentiation of (20b) with respect to $x$ (with $\nu=N$) , we simplify the formulas for $Z^{(k)}_{mn}$, where $k=1,2$, to read

[TABLE]

and

[TABLE]

where $B_{mn}$ are defined by

[TABLE]

By using the last two expressions for the components of the matrices $Z^{(1)}$ and $Z^{(2)}$ in (23), we obtain the following formulas for the components of the generalized pseudospectral matrix representation $D^{c}$:

[TABLE]

where, again, $B_{mn}$ is given by definition (26).

Remark 3*.*

If $\rho$ is even, the roots $\{x_{j}^{(\rho)}\}_{j=1}^{\rho}$ of the Sonin-Markov polynomial $p_{\rho}(x)\equiv p_{\rho}^{\beta}(x)$ are all distinct from zero because the roots of every generalized Laguerre polynomial are all distinct from zero, see definition (16a). Therefore, case (27b) does not apply to even values of $N$.

Remark 4*.*

It is known that two generalized Laguerre polynomials $L_{n}^{\alpha}(x)$ and $L_{n+1}^{\alpha}(x)$ do not have common roots  [16, 2]. But then, by definition (16), two Sonin-Markov polynomials $p_{N}^{\beta}(x)$ and $p_{N+2}^{\beta}(x)$ have a common root $\hat{x}$ if and only if $N$ is odd and $\hat{x}=0$. This is why formulas (24) and (25) include the case where $x_{m}^{(N+2)}=x_{n}^{(N)}=0$ (wich is possible for odd $N$ but not for even $N$), but do not include the case where $x_{m}^{(N+2)}=x_{n}^{(N)}\neq 0$ (the latter case is not possible).

2.3 The Generalized Spectral Matrix Representation of the Sonin-Markov Differential Operator (20b)

The generalized spectral $(N+2)\times N$ matrix representation $D^{\tau}$ of the operator $\mathcal{D}$ is defined by formula (8). To find an explicit formula for the components of $D^{\tau}$, we will employ recurrence relations satisfied by the Sonin-Markov polynomials.

Recall that the generalized Laguerre polynomials satisfy certain three-term recurrence relations [16, 2], which imply

[TABLE]

with

[TABLE]

where $p_{k}^{\beta}(x)$ are assumed to equal zero if $k<0$. Using recurrence relations (28), we obtain

[TABLE]

where $n=1,2,\ldots,N+2$, $j=1,2,\ldots,N$, and the coefficients $\mu_{\nu},\eta_{\nu}$ and $\alpha_{\nu},\beta_{\nu},\gamma_{\nu}$ are given by (21) and (29).

Remark 5*.*

The coefficients $\alpha_{\nu}$ and $\beta_{\nu}$ defined in (29c) and (29d), respectively, must not be confused with the parameters $\alpha$ and $\beta$ in the definition of the Sonin-Markov polynomials, see (16) and (17).

For completeness, let us mention that as a consequence of their orthogonality, Sonin-Markov polynomials satisfy the standard three-term recurrence relation

[TABLE]

where

[TABLE]

and $p_{j}(x)=0$ if $j<0$. Of course, the three-term recurrence relation (28) is a consequence of the recurrence relation (31).

2.4 Algebraic Properties of the Zeros of the Sonin-Markov Polynomials

Having found the generalized pseudospectral and the generalized spectral matrix representations $D^{c}$ and $D^{\tau}$, respectively, of the differential operator $\mathcal{D}$, we apply Theorem 1.1 to the case of the Sonin-Markov polynomials. We thus obtain the following algebraic identities satisfied by their zeros.

Theorem 2.1**.**

For every pair of integers $n,m$ such that $1\leq n\leq N$ and $1\leq m\leq N+2$, the zeros $x_{1}^{(N+2)},x_{2}^{(N+2)},\ldots,x_{N+2}^{(N+2)}$ of the Sonin-Markov polynomial $p_{N+2}(x)\equiv p_{N+2}^{\beta}(x)$ and the zeros $x_{1}^{(N)},\ldots,x_{N}^{(N)}$ of the Sonin-Markov polynomial $p_{N}(x)\equiv p_{N}^{\beta}(x)$ satisfy the following relations.

If $x_{m}^{(N+2)}=0$, then

[TABLE]

where $p_{j}(x)=0$ if $j<0$.

If $x_{m}^{(N+2)}\neq 0$, then

[TABLE]

where $B_{mn}$ are defined by (26) and, as before, $p_{j}(x)=0$ if $j<0$.

3 Proofs

The proof of Theorem 1.1 is based on the following result.

Theorem 3.1**.**

Let $\mathcal{A}$ be a linear differential operator that satisfies condition (2). Let $\{\ell_{N-1,j}(x)\}_{j=1}^{N}$ be the Lagrange interpolation basis of $\mathbb{P}^{N-1}$ constructed using the $N$ distinct real nodes $x_{1}^{(N)},\ldots,x_{N}^{(N)}$ and let $\{\ell_{N+n_{0}-1,j}(x)\}_{j=1}^{N+n_{0}}$ be the Lagrange interpolation basis of $\mathbb{P}^{N+n_{0}-1}$ constructed using the $(N+n_{0})$ distinct real nodes $x_{1}^{(N+n_{0})},\ldots,x_{N+n_{0}}^{(N+n_{0})}$, see (6). If the $(N+n_{0})\times N$ matrices $A^{c},A^{\tau}$ are defined by (7) and (8), respectively, while the two $\rho\times\rho$ matrices $L^{(\rho-1)}$ with $\rho=N$ or $\rho=N+n_{0}$ are defined componentwise by

[TABLE]

then

[TABLE]

Moreover, if the interpolation nodes $x_{1}^{(\rho)},\ldots,x_{\rho}^{(\rho)}$ are the distinct real zeros of the polynomial $p_{\rho}(x)$ from the orthogonal family $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$ introduced in Section 1, then the transition matrix $L^{(\rho-1)}$ and its inverse $\left[L^{(\rho-1)}\right]^{-1}$ are given by (10) and (12), respectively, where $\rho\in\{N,N+n_{0}\}$.

Remark 6*.*

Let us note that $L^{(N-1)}$ is the transition matrix from the polynomial basis $\{p_{m}(x)\}_{m=0}^{N-1}$ to the basis $\{\ell_{N-1,m}(x)\}_{m=1}^{N}$ of $\mathbb{P}^{N-1}$, while $L^{(N+n_{0}-1)}$ is the transition matrix from the polynomial basis $\{p_{m}(x)\}_{m=0}^{N+n_{0}-1}$ to the basis $\{\ell_{N+n_{0}-1,m}(x)\}_{m=1}^{N+n_{0}}$ of $\mathbb{P}^{N+n_{0}-1}$.

The proof of this theorem is similar to the proof of Theorem 1.1 in [7]. It is provided below for the convenience of the reader.

Proof of Theorem 3.1.

First, let us prove property (36). Let $u$ be a polynomial of degree $N-1$. Then

[TABLE]

where the coefficients

[TABLE]

respectively, are the components of the column-vectors $\vec{u}^{c}$ and $\vec{u}^{\tau}$, respectively, and $1\leq j\leq N$. To prove relation (36), we will show that $\vec{u}^{\tau}=L^{(N-1)}\vec{u}^{c}$ and $L^{(N+n_{0}-1)}A^{c}\vec{u}^{c}=A^{\tau}\vec{u}^{\tau}$.

Let us expand

[TABLE]

where the coefficients $L_{mj}^{(N-1)}$ are given by (35). Upon a substitution of (41) into (37), we obtain

[TABLE]

To obtain the equality $L^{(N+n_{0}-1)}A^{c}\vec{u}^{c}=A^{\tau}\vec{u}^{\tau}$, we first notice that because $u\in\mathbb{P}^{N-1}$ and the operator $\mathcal{A}$ satisfies $\mathcal{A}\mathbb{P}^{N-1}\subseteq\mathbb{P}^{N+n_{0}-1}$, we have

[TABLE]

On the other hand,

[TABLE]

By comparing expansions (43) and (44), we obtain $L^{(N+n_{0}-1)}A^{c}\vec{u}^{c}=A^{\tau}\vec{u}^{\tau}$. Because $\vec{u}^{\tau}=L^{(N-1)}\vec{u}^{c}$, we conclude that $L^{(N+n_{0}-1)}A^{c}=A^{\tau}L^{(N-1)}$ and so finish the proof of relation (36).

Second, let us assume that $x_{1}^{(\rho)},\ldots,x_{\rho}^{(\rho)}$ are the zeros of the polynomial $p_{\rho}(x)$, where $\rho=N$ or $\rho=N+n_{0}$. Let us prove that the transition matrix $L^{(\rho-1)}$ is given componentwise by (10). The Gaussian rule for approximate integration with respect to the measure $\omega$ based on these nodes $x_{1}^{(\rho)},\ldots,x_{\rho}^{(\rho)}$ has degree of exactness $2\rho-1$, see, for example, Theorem 5.1.2 of [16]. Therefore, for the polynomial $\ell_{\rho-1,j}(x)p_{m-1}(x)$ of degree $\rho-1+m-1\leq 2\rho-2$ we have

[TABLE]

which implies (10). By applying the Gaussian rule to the polynomial $p_{m-1}(x)p_{n-1}(x)$, where $0\leq m,n\leq N$, we obtain

[TABLE]

which implies (12). ∎

Proof of Theorem 1.1.

Theorem 1.1 is a straightforward consequence of identity (36) rewritten as $\left\{A^{c}\left[L^{(N-1)}\right]^{-1}\right\}_{mn}=\left\{\left[L^{(N+n_{0}-1)}\right]^{-1}A^{\tau}\right\}_{mn}$, where $1\leq m\leq N+n_{0}$ and $1\leq n\leq N$, and formula (12) in Theorem 3.1.

We note that the index $k\in\{1,\ldots,N+n_{0}\}$ in the sum from the right hand side of identity (13) ranges from $n-n_{0}$ to $n+n_{0}$ rather than from $1$ to $N+n_{0}$. This is due to the following consideration. Because the polynomial family $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$ is orthogonal and has the property that $\{p_{j}(x)\}_{j=0}^{\nu}$ is a basis of $\mathbb{P}^{\nu}$ for each $\nu$, this family must satisfy a three-term recurrence relation

[TABLE]

where $a_{\nu,k}$ are constants that equal zero if $k<0$ or $k$ is outside of the range $\{\nu-1,\nu,\nu+1\}$ [16]. From the last recurrence relation we derive

[TABLE]

where $b_{\nu,k}$ are constants that equal zero if $k<0$. Thus, for all integer $n,k$ such that $1\leq k\leq N+n_{0}$ and $1\leq n\leq N$,

[TABLE]

if the index $k$ is outside of the range $\{n-n_{0},\ldots,n+n_{0}\}$. ∎

4 Discussion and Outlook

The identities of the main Theorem 1.1 are derived using a very special exact discretization of the differential equation

[TABLE]

compare with (1). This discretization is constructed using the fact that ODE (49) has a polynomial solution $p_{N}(x)$, where $p_{N}(x)$ is a member of the nonclassical orthogonal polynomial family $\{p_{\nu}(x)\}_{\nu=0}^{\infty}$. Because the differential operator $\mathcal{A}$ maps $\mathbb{P}^{N}$ to $\mathbb{P}^{N+n_{0}}$, see property (2), we choose to represent the operator by a $(N+n_{0})\times N$ matrix $A^{c}\left(\vec{x}^{(N)},\vec{x}^{(N+n_{0})}\right)\equiv A^{c}$ defined by (7) in terms of two vectors of interpolation nodes, $\vec{x}^{(N)}$ and $\vec{x}^{(N+n_{0})}$. We thus generalize the notion of pseudospectral matrix representations of a linear differential operator.

This generalization can be utilized to resolve the issues related to the incorporation of the initial or boundary conditions into pseudpsectral methods for solving a linear ODE $\mathcal{L}u=f$, where $\mathcal{L}$ is a linear differential operator, see [30, 31]. For example, if the differential operator $\mathcal{L}$ has the property $\mathcal{L}\mathbb{P}^{\nu}\subseteq\mathbb{P}^{\nu-k_{0}}$, where $k_{0}$ is a positive integer, it is beneficial to use the generalized pseudospectral $N\times(N-k_{0})$ matrix representation $L^{c}$ of the differential operator $\mathcal{L}$ with respect to two vectors of nodes $\vec{x}^{(N)}\in\mathbb{R}^{N}$ and $\vec{x}^{(N-k_{0})}\in\mathbb{R}^{N-k_{0}}$. Of course, $L^{c}$ is defined componentwise by $L^{c}_{kj}=(\mathcal{L}\ell_{N-1,j})\left(x_{k}^{(N-k0)}\right),\;\;1\leq k\leq N-k_{0},1\leq j\leq N$, compare with (7), where $\{\ell_{N-1,j}(x)\}_{j=1}^{N}$ is the standard Lagrange interpolation basis with respect to the nodes $\vec{x}^{N}$.

We can approximate the ODE $\mathcal{L}u=f$ by the $N\times(N-k_{0})$ system of linear equations $L^{c}\vec{u}^{(N)}=\vec{f}^{(N-k_{0})}$, where $\vec{f}^{(N-k_{0})}=\left(f\left(x_{1}^{(N-k_{0})}\right),\ldots,f\left(x_{N-k_{0}}^{(N-k_{0})}\right)\right)$, compare with (4). The linear system $L^{c}\vec{u}^{(N)}=\vec{f}^{(N-k_{0})}$ for the $N$ components of $\vec{u}^{(N)}$ consists of $N-k_{0}$ linear equations, which allows to easily incorporate $k_{0}$ boundary or initial conditions into the system. Once the solution $\vec{u}^{(N)}$ of this augmented linear system is found, an approximate solution of the ODE $\mathcal{L}u=f$ is given by $u(x)\approx(\pi_{N})^{-1}\vec{u}^{(N)}=\sum_{j=1}^{N}{u}^{(N)}_{j}\ell_{N-1,j}(x)$, where $\{\ell_{N-1,j}(x)\}_{j=1}^{N}$ is the standard Lagrange interpolation basis with respect to the nodes $\vec{x}^{N}$.

The ODE (49) may also be discretized using only one vector of nodes $\vec{x}^{(N+n_{0})}$ and the standard pseudospectral $(N+n_{0})\times(N+n_{0})$ matrix representation $\hat{A}^{c}$ of the differential operator $\mathcal{A}$ defined componentwise by

[TABLE]

see [7] and compare with (7). Indeed, in the proof of Theorem 3.1 we may choose to represent the polynomial $u\in\mathbb{P}^{N-1}$ in terms of the Lagrange interpolation basis $\left\{\ell_{N+n_{0}-1,j}(x)\right\}_{j=1}^{N+n_{0}}$ or the basis $\{p_{j-1}(x)\}_{j=1}^{N+n_{0}}$ , see (37) and (38). This modification allows to prove the relation

[TABLE]

where the standard spectral $(N+n_{0})\times(N+n_{0})$ matrix representation $\hat{A}^{\tau}$ is defined by

[TABLE]

Let us assume again that $\vec{x}^{(N+n_{0})}$ are the zeros of the polynomial $p_{N+n_{0}}$. Then relation (51) rewritten as $\left\{\hat{A}^{c}\left[L^{(N+n_{0}-1)}\right]^{-1}\right\}_{mn}=\left\{\left[L^{(N+n_{0}-1)}\right]^{-1}\hat{A}^{\tau}\right\}_{mn}$ together with expression (12) for the inverse of the transition matrices $L^{(N+n_{0}-1)}$ yield the identities

[TABLE]

which hold for all integer $m,n$ such that $1\leq m,n\leq N+n_{0}$. Given a fixed integer $n$ such that $0\leq n\leq N+n_{0}-1$, the last identity relates the zeros of the polynomials $p_{N+n_{0}}(x)$ and $p_{n}(x)$ with the zeros of all the polynomials $p_{j}(x)$ such that the integer index $j\in\{0,1,\ldots,N+n_{0}-1\}$ satisfies $n-n_{0}\leq j\leq n+n_{0}$.

The set of identities (53) is a straighforward generalization of Theorem 1.1 in [7], thus it is not the main focus of this paper. We invite the reader to apply identities (53) to the Sonin-Markov polynomials, replacing $\hat{A}^{c}$ and $\hat{A}^{\tau}$ with $\hat{D}^{c}$ and $\hat{D}^{\tau}$, respectively, the latter being matrix representations of the differential operator $\mathcal{D}$ associated with the Sonin-Markov family, see (20b). The pseudospectral matrix representation $\hat{D}^{c}$ can be found using definition (50) and the techniques outlined in Appendix A of [7]. The spectral matrix representation $\hat{D}^{\tau}$ can be found using definition (52) and recurrence relations (28).

In summary, the method of generalized pseudospectral matrix represenetations proposed in this paper can be utilized to prove new identities satisfied by nonclassical orthogonal polynomials that are generalized, rather than standard, eigenfunctions of linear differential operators. The proposed generalized matrix representations of linear differential operators are quite useful in the process of incorporation of initial or boundary conditions into a discretization of a linear ODE via a pseudospectral method. The last application deserves further exploration involving numerical tests, which the authors may pursue in the near future. Other interesting developments of the method proposed in this paper may involve its generalization to exceptional orthogonal polynomials [35, 36, 37].

Acknowledgements

This work was supported by a grant awarded to O. Bihun by the Committee on Research and Creative Works (CRCW) of the University of Colorado, Colorado Springs.

The authors declare that there is no conflict of interest regarding the publication of this paper.