Multiple Hermite polynomials and simultaneous Gaussian quadrature

TL;DR

This paper studies multiple Hermite polynomials, their properties, zero distributions, and uses their zeros for simultaneous Gaussian quadrature to approximate multiple integrals with Gaussian weights.

Contribution

It introduces properties and asymptotic zero distributions of multiple Hermite polynomials and applies their zeros to develop a new simultaneous quadrature method for Gaussian-weighted integrals.

Findings

Zeros may accumulate on multiple disjoint intervals depending on parameters.

Asymptotic zero distribution is characterized for different parameter regimes.

A new simultaneous Gaussian quadrature method is proposed for multiple integrals.

Abstract

Multiple Hermite polynomials are an extension of the classical Hermite polynomials for which orthogonality conditions are imposed with respect to $r>1$ normal (Gaussian) weights $w_j(x)=e^{-x^2+c_jx}$ with different means $c_j/2$, $1 \leq j \leq r$. These polynomials have a number of properties, such as a Rodrigues formula, recurrence relations (connecting polynomials with nearest neighbor multi-indices), a differential equation, etc. The asymptotic distribution of the (scaled) zeros is investigated and an interesting new feature happens: depending on the distance between the $c_j$, $1 \leq j \leq r$, the zeros may accumulate on $s$ disjoint intervals, where $1 \leq s \leq r$. We will use the zeros of these multiple Hermite polynomials to approximate integrals of the form $\displaystyle \int_{-\infty}^{\infty} f(x) \exp(-x^2 + c_jx)\, dx$ simultaneously for $1 \leq j \leq r$ for the case $r=3$ and the situation when the zeros accumulate on three disjoint intervals. We also give some properties of the corresponding quadrature weights.

Tables

Table 1. Table 1: The quadrature weights λ k , 30 ( 1 ) superscript subscript 𝜆 𝑘 30 1 \lambda_{k,30}^{(1)} for the first integral ( c ^ = 4.7434 ^ 𝑐 4.7434 \hat{c}=4.7434 ).

$k$	${\lambda }_{k,30}^{\left(1\right)}$
1	${\mathrm{6.887653865\; 10}}^{-9}$
2	${\mathrm{4.384111578\; 10}}^{-6}$
3	${\mathrm{0.3591983034\; 10}}^{-3}$
4	${\mathrm{0.8149617619\; 10}}^{-2}$
5	${\mathrm{0.6836500666\; 10}}^{-1}$
6	$0.2410330694$
7	$0.3725933960$
8	$0.2452710131$
9	${\mathrm{0.6041135610\; 10}}^{-1}$
10	${\mathrm{0.3809098858\; 10}}^{-2}$
11	${\mathrm{6.755525278\; 10}}^{-6}$
12	$-{\mathrm{5.189883715\; 10}}^{-6}$
13	${\mathrm{3.848392520\; 10}}^{-6}$
14	$-{\mathrm{2.434636570\; 10}}^{-6}$
15	${\mathrm{1.261797315\; 10}}^{-6}$

Equations

$${\int }_{\mathbb{R}}f(x)w_j(x)dx,1\le j\le r$$

$$\sum _{k=1}^N{\lambda }_{k,N}^{(j)}f(x_{k,N})$$

$${\int }_{-\infty }^{\infty }f(x)w_R(x)dx,{\int }_{-\infty }^{\infty }f(x)w_G(x)dx,{\int }_{-\infty }^{\infty }f(x)w_B(x)dx.$$

$$\int x^k{P}_{\vec{n}}(x)d{\mu }_j(x)=0,0\le k\le n_j-1,$$

$$\begin{gathered} w_1,x{w}_1,x^2{w}_1,\dots ,x^{n_1-1}{w}_1,w_2,x{w}_2,x^2{w}^2,\dots ,x^{n_2-1}{w}_2,\dots , \\ {w}_r,x{w}_r,x^2{w}_r,\dots ,x^{n_r-1}{w}_r \end{gathered}$$

$$w_j(x)=e^{-x^2+c_{j}x},x\in \mathbb{R},$$

$$e^{-x^2}{H}_{\vec{n}}(x)=(-1{)}^{|\vec{n}|}{2}^{-|\vec{n}|}\left(\prod _{j=1}^r{e}^{-c_{j}x}\frac{d^{n_j}}{d{x}^{n_j}}{e}^{c_{j}x}\right)e^{-x^2},$$

$$H_{\vec{n}}(x)=(-1{)}^{|\vec{n}|}{2}^{-|\vec{n}|}\sum _{k_1=0}^{n_1}\cdots \sum _{k_r=0}^{n_r}(\genfrac{}{}{0ex}{}{n_1}{k_1})\cdots (\genfrac{}{}{0ex}{}{n_r}{k_r})c_1^{n_1-k_1}\cdots c_r^{n_r-k_r}(-1{)}^{|\vec{k}|}{H}_{|\vec{k}|}(x).$$

$$x{H}_{\vec{n}}(x)=H_{\vec{n}+{\vec{e}}_k}(x)+\frac{c_k}{2}{H}_{\vec{n}}(x)+\frac{1}{2}\sum _{j=1}^r{n}_j{H}_{\vec{n}-{\vec{e}}_j}(x),1\le k\le r,$$

$${\left(e^{-x^2+c_{j}x}{H}_{\vec{n}-{\vec{e}}_j}(x)\right)}^{\prime }=-2e^{-x^2+c_{j}x}{H}_{\vec{n}}(x),1\le j\le r.$$

$$H_{\vec{n}}^{\prime }(x)=\sum _{j=1}^r{n}_j{H}_{\vec{n}-{\vec{e}}_j}(x).$$

$$\left(\prod _{j=1}^r{D}_j\right)D{H}_{\vec{n}}(x)=-2\left(\sum _{j=1}^r{n}_j\prod _{i\ne j}{D}_j\right)H_{\vec{n}}(x),$$

$$D=\frac{d}{dx},D_j=e^{x^2-c_{j}x}D{e}^{-x^2+c_{j}x}=D+(-2x+c_j)I.$$

$$w_1(x)=e^{-x^2-cx},w_2(x)=e^{-x^2},w_3(x)=e^{-x^2+cx}.$$

$$I_1=[-\frac{c}{2}-\sqrt{4n+1},-\frac{c}{2}+\sqrt{4n+1}],I_2=[-\sqrt{4n+1},\sqrt{4n+1}],$$

$$I_3=[\frac{c}{2}-\sqrt{4n+1},\frac{c}{2}+\sqrt{4n+1}],$$

$${\int }_{-\infty }^{\infty }{H}_{n,n,n}(x){\pi }_m(x)e^{-x^2+cx}dx=0.$$

$${\int }_{\mathbb{R}\setminus I_3}|H_{n,n,n}(x){\pi }_m(x)|e^{-x^2+cx}dx<{\int }_{I_3}{H}_{n,n,n}(x){\pi }_m(x)e^{-x^2+cx}dx,$$

$${\int }_{\mathbb{R}\setminus I_3}{H}_{n,n,n}(x){\pi }_m(x)e^{-x^2+cx}dx>-{\int }_{I_3}{H}_{n,n,n}(x){\pi }_m(x)e^{-x^2+cx}dx$$

$$\begin{gathered} {\int }_{-\infty }^{\infty }{H}_{n,n,n}(x){\pi }_m(x)e^{-x^2+cx}dx \\ ={\int }_{I_3}{H}_{n,n,n}(x){\pi }_m(x)e^{-x^2+cx}dx+{\int }_{\mathbb{R}\setminus I_3}{H}_{n,n,n}(x){\pi }_m(x)e^{-x^2+cx}dx>0, \end{gathered}$$

$$\underset{n\to \infty }{\lim }\frac{1}{3n}\sum _{j=1}^{3n}f\left(\frac{x_{j,3n}}{\sqrt{n}}\right)=\int f(x)v(x)dx$$

$$y^{\prime \prime \prime \prime }-6x{y}^{\prime \prime \prime }+(12x^2-c^2-6)y^{\prime \prime }+[-8x^3+(2c^2+12)x]y^{\prime }=-2n[3y^{\prime \prime }-12x{y}^{\prime }+(12x^2-c^2-6)y].$$

$$w_1(x)=e^{-n(x^2+\hat{c}x)},w_2(x)=e^{-n{x}^2},w_3(x)=e^{-n(x^2-\hat{c}x)}.$$

$$S_n(z)=\frac{1}{\sqrt{n}}\frac{H_{n,n,n}^{\prime }(\sqrt{n}z)}{H_{n,n,n}(\sqrt{n}z)}=\frac{1}{n}\sum _{j=1}^{3n}\frac{1}{z-\frac{x_{j,3n}}{\sqrt{n}}}=\int \frac{d{\mu }_n(x)}{z-x},$$

$${\mu }_n=\frac{1}{n}\sum _{j=1}^{3n}{\delta }_{x_{j,3n}/\sqrt{n}}.$$

$$S_{n_k}\to S,S_{n_k}^{\prime }\to S^{\prime },S_{n_k}^{\prime \prime }\to S^{\prime \prime },S_{n_k}^{\prime \prime \prime }\to S^{\prime \prime \prime }.$$

$$S(z)=3\int \frac{d\mu (x)}{z-x}dx,$$

$$H_{n,n,n}^{\prime }(\sqrt{n}z)=\sqrt{n}{S}_n{H}_{n,n,n}(\sqrt{n}z),$$

$$H_{n,n,n}^{\prime \prime }(\sqrt{n}z)$$

$$H_{n,n,n}^{\prime \prime \prime }(\sqrt{n}z)$$

$$H_{n,n,n}^{\prime \prime \prime \prime }(\sqrt{n}z)$$

Full text

1 Simultaneous Gaussian quadrature

Let $w_{1},\ldots,w_{r}$ be $r\geq 1$ weight functions on $\mathbb{R}$ and $f:\mathbb{R}\to\mathbb{R}$. Simultaneous quadrature is a numerical method to approximate $r$ integrals

[TABLE]

by sums

[TABLE]

at the same $N$ points $\{x_{k,N},1\leq k\leq N\}$ but with weights $\{\lambda_{k,N}^{(j)},1\leq k\leq N\}$ which depend on $j$. This was described by Borges [3] in 1994 but was originally suggested by Aurel Angelescu [1] in 1918, whose work seems to have gone unnoticed. The past few decades it became clear that this is closely related to multiple orthogonal polynomials, in a similar way as Gauss quadrature is related to orthogonal polynomials. The motivation in [3] was involving a color signal $f$, which can be transmitted using three colors: red-green-blue (RGB). For this we need the amount of R-G-B in $f$, given by

[TABLE]

A natural question is whether this can be done with $N$ function evaluations and maximum degree of accuracy. If we choose $n$ points for each integral, and then use Gaussian quadrature, then this would require $3n$ function evaluations for a degree of accuracy $2n-1$. A better choice is to choose the zeros of the multiple orthogonal polynomials $P_{n,n,n}$ for the weights $(w_{R},w_{G},w_{B})$ and then use interpolatory quadrature. This again requires $3n$ function evaluations, but the degree of accuracy is increased to $4n-1$. We will call this method based on zeros of multiple orthogonal polynomials the simultaneous Gaussian quadrature method. Some interesting research problems for simultaneous Gaussian quadrature are

• •

To find the multiple orthogonal polynomials when the weights $w_{1},\ldots,w_{r}$ are given.

• •

To study the location and computation of the zeros of the multiple orthogonal polynomials.

• •

To study the behavior and computation of the weights $\lambda_{k,N}^{(j)}$.

• •

To investigate the convergence of the quadrature rules.

Part of this research has already been started in [4], [5], [8], [9], [13], but there is still a lot to be done in this field.

2 Multiple orthogonal polynomials

Definition 2.1**.**

Let $\mu_{1},\ldots,\mu_{r}$ be $r$ positive measures on $\mathbb{R}$ and let $\vec{n}=(n_{1},\ldots,n_{r})$ be a multi-index in $\mathbb{N}^{r}$. The (type II) multiple orthogonal polynomial $P_{\vec{n}}$ is the monic polynomial of degree $|\vec{n}|=n_{1}+n_{2}+\cdots+n_{r}$ that satisfies the orthogonality conditions

[TABLE]

for $1\leq j\leq r$.

Such a monic polynomial may not exist, or may not be unique. One needs conditions on the (moments of) the measures $(\mu_{1},\ldots,\mu_{r})$. Two important cases have been introduced for which all the multiple orthogonal polynomials exist and are unique. The measures $(\mu_{1},\ldots,\mu_{r})$ are an Angelesco system if $\textup{supp}(\mu_{j})\subset\Delta_{j}$, where the $\Delta_{j}$ are intervals which are pairwise disjoint: $\Delta_{i}\cap\Delta_{j}=\emptyset$ whenever $i\neq j$.

Theorem 2.1** (Angelesco).**

For an Angelesco system the multiple orthogonal polynomials exist for every multi-index $\vec{n}$. Furthermore $P_{\vec{n}}$ has $n_{j}$ simple zeros in each interval $\Delta_{j}$.

For a proof, see [10, Ch. 4, Prop. 3.3] or [6, Thm. 23.1.3]. The behavior of the quadrature weights for simultaneous Gaussian quadrature is known for this case (see [10, Ch. 4 4, Prop. 3.5], [8, Thm. 1.1]).

Theorem 2.2**.**

The quadrature weights $\lambda_{k,n}^{(j)}$ are positive for the $n_{j}$ zeros on $\Delta_{j}$. The remaining quadrature weights have alternating sign, with those for the zeros closest to the interval $\Delta_{j}$ positive.

Another important case is when the measures form an AT-system. The weight functions $(w_{1},\ldots,w_{r})$ are an algebraic Chebyshev system (AT-system) on $[a,b]$ if

[TABLE]

are a Chebyshev system on $[a,b]$ for every $(n_{1},\ldots,n_{r})\in\mathbb{N}^{r}$.

Theorem 2.3**.**

For an AT-system the multiple orthogonal polynomials exist for every multi-index $(n_{1},\ldots,n_{r})$. Furthermore $P_{\vec{n}}$ has $|\vec{n}|$ simple zeros on the interval $[a,b]$.

For a proof, see [10, Ch. 4, Corr. of Thm. 4.3] or [6, Thm. 23.1.4].

3 Multiple Hermite polynomials

We will consider the weight functions

[TABLE]

with real parameters $c_{1},\ldots,c_{j}$ such that $c_{i}\neq c_{j}$ whenever $i\neq j$. These weights are proportional to normal weights with means at $c_{j}/2$ and variance $\sigma^{2}=\frac{1}{2}$. They form an AT-system and the corresponding multiple orthogonal polynomials are known as multiple Hermite polynomials $H_{\vec{n}}$. They can be obtained by using the Rodrigues formula

[TABLE]

from which one can find the explicit expression

[TABLE]

See [6, §23.5]. Multiple orthogonal polynomials satisfy a system of recurrence relations connecting the nearest neighbors. For multiple Hermite polynomials one has

[TABLE]

where $\vec{e}_{1}=(1,0,0,\ldots,0),\vec{e}_{2}=(0,1,0,\ldots,0),\ldots,\vec{e}_{r}=(0,0,\ldots,0,1)$. They also have interesting differential relations, such as $r$ raising operators

[TABLE]

and a lowering operator

[TABLE]

See [6, Eqs. (23.8.5)–(23.8.6)]. Combining these raising operators and the lowering operator gives a differential equation of order $r+1$

[TABLE]

where

[TABLE]

From now on we deal with the case $r=3$ and weights $c_{1}=-c,c_{2}=0,c_{3}=c$:

[TABLE]

3.1 Zeros

Let $x_{1,3n}<\ldots<x_{3n,3n}$ be the zeros of $H_{n,n,n}$. First we will show that, for $c$ large enough, the zeros of $H_{n,n,n}$ lie on three disjoint intervals around $-c/2$, [math] and $c/2$.

Proposition 3.1**.**

For $c$ sufficiently large (e.g., $c>4\sqrt{4n+1}$) the zeros of $H_{n,n,n}$ are on three disjoint intervals $I_{1}\cup I_{2}\cup I_{3}$, where

[TABLE]

[TABLE]

and each interval contains $n$ simple zeros.

Proof.

Suppose $x_{1},x_{2},\ldots,x_{m}$ are the sign changes of $H_{n,n,n}$ on $I_{3}$ and that $m<n$. Let $\pi_{m}(x)=(x-x_{1})(x-x_{2})\cdots(x-x_{m})$, then $H_{n,n,n}(x)\pi_{m}(x)$ does not change sign on $I_{3}$. By the multiple orthogonality one has

[TABLE]

Suppose that $H_{n,n,n}(x)\pi_{m}(x)$ is positive on $I_{3}$, then by the infinite-finite range inequalities (see, e.g., [7, Ch. 4, Thm. 4.1], where we take $Q(x)=x^{2}-cx$, $p=1$, $t=4n+1$, so that $\Delta_{t}=I_{3}$) one has

[TABLE]

so that

[TABLE]

and

[TABLE]

which is in contradiction with (3.1). This mean that our assumption that $m<n$ is false and hence $m\geq n$. We can repeat the reasoning for $I_{2}$ and $I_{1}$, and since $H_{n,n,n,}$ is a polynomial of degree $3n$, we must conclude that each interval contains $n$ zeros of $H_{n,n,n}$, which are all simple. Clearly the three intervals are disjoint when $c>4\sqrt{4n+1}$. ∎

This result shows that, for large $c$, the multiple Hermite polynomials behave very much like an Angelesco system, i.e., multiple orthogonal polynomials for which the orthogonality conditions are on disjoint intervals. Some results for simultaneous Gauss quadrature for Angelesco systems were proved earlier in [10, Ch. 4, Prop. 3.4 and 3.5] and [8]. In this paper we will show that similar results are true for multiple Hermite polynomials when $c$ is large.

The intervals $I_{1},I_{2},I_{3}$ are in fact a bit too large, because they were obtained by using the infinite-finite range inequalities for one weight only, and not for the three weights simultaneously. In order to study the zeros in more detail, we will take the parameter $c$ proportional to $\sqrt{n}$ and scale the zeros by a factor $\sqrt{n}$ as well. This amounts to investigating the polynomials $H_{n,n,n}(\sqrt{n}x)$ with $c=\sqrt{n}\hat{c}$. In order to find for which values of $c$ the zeros are accumulating on three disjoint intervals as $n\to\infty$, we will investigate the asymptotic distribution of the zeros. Our main theorem in this section is

Theorem 3.2**.**

There exists a $c^{*}>0$ such that for the zeros of $H_{n,n,n}$ with $c=\sqrt{n}\hat{c}$ one has

[TABLE]

where $v$ is a probability density supported on three intervals $[-b,-a]\cup[-d,d]\cup[a,b]$ $(0<d<a<b)$ when $\hat{c}>c^{*}$ and $v$ is supported on one interval $[-b,b]$ when $\hat{c}<c^{*}$. The numerical value is $c^{*}=4.10938818$.

Such a phase transition when the zeros cluster on one interval when the parameters are close together or on two intervals when the parameters are far apart, was first observed and proved for $r=2$ by Bleher and Kuijlaars [2].

Proof.

The differential equation for $y=H_{n,n,n}(x)$ becomes

[TABLE]

The scaling amounts to studying zeros of $H_{n,n,n}(\sqrt{n}x)$ and these are multiple orthogonal polynomials for the weight functions

[TABLE]

Consider the rational function

[TABLE]

where $\mu_{n}$ is the discrete measure with mass $1/n$ at each scaled zero $x_{j,3n}/\sqrt{n}$:

[TABLE]

The sequence $(S_{n})_{n\in\mathbb{N}}$ is a family of analytic functions which is uniformly bounded on every compact subset of $\mathbb{C}\setminus\mathbb{R}$, hence by Montel’s theorem there exists a subsequence $(S_{n_{k}})_{k}$ that converges uniformly on compact subsets of $\mathbb{C}\setminus\mathbb{R}$ to an analytic function $S$, and also its derivatives converge uniformly on these compact subsets:

[TABLE]

Since each $S_{n}$ is a Stieltjes transform of a positive measure (with total mass 3), the limit is of the form

[TABLE]

with $\mu$ a probability measure on $\mathbb{R}$, which describes the asymptotic distribution of the scaled zeros, and $\mu_{n}$ converges weakly to the measure $3\mu$ alongs the chosen subsequence. This function $S$ may depend on the selected subsequence $(n_{k})_{k}$, but we will show that every convergent subsequence has the same limit $S$. Observe that

[TABLE]

from which we can find

[TABLE]

Put this in the differential equation (with $x=\sqrt{n}z$ and $c=\sqrt{n}\hat{c}$), then as $n=n_{k}\to\infty$ one finds

[TABLE]

This is an algebraic equation of order 4 and hence it has four solutions $S_{(1)},S_{(2)},S_{(3)},S_{(4)}$. A careful analysis of these solutions and equation (3.2) near infinity shows that for $z\to\infty$

[TABLE]

[TABLE]

We are therefore interested in $S_{(1)}(z)$ since it gives the required Stieltjes transform

[TABLE]

The algebraic equation is independent of the selected subsequence, which implies that every subsequence $(S_{n_{k}})_{k}$ has the same limit, which in turn implies that the full sequence $(S_{n})_{n}$ converges to this limit $S$. The measure $\mu$ can be retrieved by using the Stieltjes-Perron inversion theorem

[TABLE]

If $\mu$ has no mass points, then the density $v$ of $\mu$ is given by

[TABLE]

hence the support of the density $v$ is given by the set on $\mathbb{R}$ where $S$ has a jump discontinuity. This can be analyzed by investigating the discriminant of the algebraic equation:

[TABLE]

This is a polynomial of degree $6$ in the variable $z$. The support of $v$ is where this polynomial is negative. There is a phase transition from one interval to three intervals when the $z$-polynomial (3.3) has two double roots.

This happens when the discriminant of the $z$-polynomial (3.3) is zero

[TABLE]

The only positive real zero is the positive real root of

[TABLE]

and this is $c^{*}=4.10938818$. ∎

Observe that the phase transition $c^{*}$ is at a smaller value than the one suggested by Proposition 3.1, which would give the value $8$. As mentioned before, this is because in Proposition 3.1 we used the infinite-finite range inequalities for one single weight and not of the three weights simultaneously.

4 Some potential theory

From now one we assume $\hat{c}>c^{*}=4.10938818$. The Stieltjes transform of the asymptotic zero distribution is

[TABLE]

The measures $\nu_{1},\nu_{2},\nu_{3}$ are unit measures that are minimizing the expression

[TABLE]

over all unit measures $\mu_{1},\mu_{2},\mu_{3}$ supported on $\mathbb{R}$, with

[TABLE]

and

[TABLE]

This is the vector equilibrium problem for an Angelesco system [10, Ch. 5, §6]. Define the logarithmic potential

[TABLE]

The variational conditions for this vector equilibrium problem are

[TABLE]

[TABLE]

[TABLE]

where $\ell_{1},\ell_{2},\ell_{3}$ are constants (Lagrange multipliers). As an example, we have plotted these functions in Figure 3 for $\hat{c}=6$.

The measures $\nu_{1},\nu_{2},\nu_{3}$ give the asymptotic distribution of the (scaled) zeros of $H_{n,n,n}$ on the intervals $[-b,-a]$, $[-d,d]$ and $[a,b]$ respectively. They are absolutely continuous and their densities can be found as the jumps of an algebraic function $\xi$ on the real line. The function $\xi$ satisfies the algebraic equation

[TABLE]

which has four solutions $\xi_{1},\xi_{2},\xi_{3},\xi_{4}$ which behave near infinity as

[TABLE]

[TABLE]

The densities $\nu_{1}^{\prime},\nu_{2}^{\prime},\nu_{3}^{\prime}$ are given by

[TABLE]

[TABLE]

The relation between the algebraic function $S$ from (3.2) is given by

[TABLE]

The Stieltjes transforms of $\nu_{1},\nu_{2},\nu_{3}$ are related to the solutions of (3.2) by

[TABLE]

[TABLE]

5 The quadrature weights

Recall that for polynomials $f$ of degree $\leq 4n-1$

[TABLE]

Here $x_{k,3n}$ are the zeros of $H_{n,n,n}(x)=p_{n}(x)q_{n}(x)r_{n}(x)$, where $p_{n}$ has its zeros on $[-b,-a]$, $q_{n}$ on $[-d,d]$, and $r_{n}$ on $[a,b]$. Take $f(x)=\pi_{2n-1}(x)q_{n}(x)r_{n}(x)$, with $\pi_{2n-1}$ of degree $\leq 2n-1$, then (5.1) gives

[TABLE]

This is the Gauss quadrature formula for the weight function $q_{n}(x)r_{n}(x)e^{-n(x^{2}+\hat{c}x)}$ with quadrature nodes at the zeros of $p_{n}$. So we have

Lemma 5.1**.**

The first $n$ quadrature weights for the first integral (5.1) are

[TABLE]

where $\lambda_{k,n}(q_{n}r_{n}\,d\mu_{1})$ are the usual Christoffel numbers of Gaussian quadrature for the weight $q_{n}(x)p_{n}(x)e^{-n(x^{2}+\hat{c}x)}$ on $\mathbb{R}$.

For the middle $n$ quadrature weights and the last $n$ quadrature weights we have a weaker statement. By taking $f(x)=\pi_{n-1}(x)p_{n}^{2}(x)r_{n}(x)$, with $\pi_{n-1}$ of degree $\leq n-1$, the quadrature formula (5.1) gives

[TABLE]

This is not a Gauss quadrature rule but the Lagrange interpolatory rule for the weight function $p_{n}^{2}(x)r_{n}(x)e^{-n(x^{2}+\hat{c}x)}$, with quadrature nodes at the zeros of $q_{n}$. So now we have

Lemma 5.2**.**

The middle $n$ quadrature weights for the first integral are

[TABLE]

where $w_{k,n}(q_{n})$ are the quadrature weights for the Lagrange interpolatory quadrature at the zeros of $q_{n}$ and weight function $p_{n}^{2}(x)r_{n}(x)e^{-n(x^{2}+\hat{c}x)}$.

In a similar way, we take $f(x)=\pi_{n-1}(x)p_{n}^{2}(x)q_{n}(x)$, with $\pi_{n-1}$ of degree $\leq n-1$, so that (5.1) becomes

[TABLE]

We then have

Lemma 5.3**.**

The last $n$ quadrature weights for the first integral are

[TABLE]

where $w_{k,n}(r_{n})$ are the quadrature weights for the Lagrange interpolatory quadrature at the zeros of $r_{n}$ and weight function $p_{n}^{2}(x)q_{n}(x)e^{-n(x^{2}+\hat{c}x)}$.

Of course similar results are true for the quadrature weights $\lambda_{k,3n}^{(2)}$ for the second integral (5.2) and the quadrature weights $\lambda_{k,3n}^{(3)}$ for the third integral (5.3).

The weight function $q_{n}(x)r_{n}(x)e^{-n(x^{2}+\hat{c}x)}$ is not a positive weight on the whole real line, but it is positive on $[-b,-a]$ since the zeros of $q_{n}$ and $r_{n}$ are on $[-d,d]$ and $[a,b]$ respectively, at least when $n$ is large. We can prove the following result.

Theorem 5.4**.**

Let $\hat{c}$ be sufficiently large**$\hat{c}>8$ certainly works, but we conjecture that $\hat{c}>c^{*}$ is sufficient.. For the quadrature weights of the first integral (5.1) one has*

[TABLE]

and

[TABLE]

Proof.

For the first $n$ weights we use $f(x)=p_{n}^{2}(x)q_{n}(x)r_{n}(x)/(x-x_{k,3n})^{2}$ in (5.1) to find (we write $x_{k}=x_{k,3n}$)

[TABLE]

Clearly $[p_{n}^{\prime}(x_{k})]^{2}q_{n}(x_{k})r_{n}(x_{k})>0$ since $x_{k}\in[-b,-a]$ and the zeros of $q_{n}$ and $r_{n}$ are on $[-d,d]$ and $[a,b]$ respectively. So we need to prove that the integral is positive. Let $I_{1}=[-\frac{\hat{c}}{2}-\sqrt{4+1/n},-\frac{\hat{c}}{2}+\sqrt{4+1/n}]$, then by Proposition 3.1 all the zeros of $p_{n}$ are in $I_{1}$ and hence $[-b,-a]\subset I_{1}$. If $\hat{c}$ is large enough, then $q_{n}r_{n}$ is positive on $I_{1}$ and by the infinite-finite range inequality (see Proposition 3.1)

[TABLE]

so that

[TABLE]

For the middle $n$ quadrature weights we use Lemma 5.2. Clearly $p_{n}^{2}(x_{k})>0$ and $\textup{sign }r_{n}(x_{k})=(-1)^{n}$ since all the zeros of $r_{n}$ are on $[a,b]$ and $x_{k}\in[-d,d]$ for $n+1\leq k\leq 2n$. Furthermore for the Lagrange quadrature nodes one has

[TABLE]

where $\textup{sign }q_{n}^{\prime}(x_{k})=(-1)^{k-2n}$. Observe that for large enough $\hat{c}$ one has $\textup{sign }q_{n}(x)/(x-x_{k})=(-1)^{n-1}$ on $I_{1}$ since all the zeros of $q_{n}$ are on $[-d,d]$, and also $\textup{sign }r_{n}(x)=(-1)^{n}$ on $I_{1}$ since all the zeros of $r_{n}$ are on $[a,b]$. By the infinite-finite range inequality one has

[TABLE]

so that

[TABLE]

This gives $\textup{sign }\lambda_{k,3n}^{(1)}=(-1)^{k-n+1}$ for $n+1\leq k\leq 2n$. In a similar way one finds the sign of $\lambda_{k,3n}^{(1)}$ for $2n+1\leq k\leq 3n$ by using Lemma 5.3. ∎

For the quadrature weights $\lambda_{k,3n}^{(2)}$ one has a similar result, which we state without proof.

Theorem 5.5**.**

Let $\hat{c}$ be sufficiently large. For the quadrature weights of the second integral (5.2) one has

[TABLE]

and

[TABLE]

Observe that the quadrature weights for the nodes outside $[-d,d]$ are alternating, but the weights for the nodes closest to $[-d,d]$ are positive.

For the quadrature nodes $\lambda_{k,3n}^{(3)}$ one has

Theorem 5.6**.**

Let $\hat{c}$ be sufficiently large. For the quadrature weights of the third integral (5.3) one has

[TABLE]

and

[TABLE]

Having positive quadrature weights is a nice property, as is well known for Gauss quadrature. The alternating quadrature weights are not so nice, but we can show that they are exponentially small.

Theorem 5.7**.**

Suppose $\hat{c}$ is sufficiently large (see the footnote in Theorem 5.4). For the positive quadrature weights one has

[TABLE]

whenever $x_{k}\to x\in(-b,-a)$. For the quadrature weights with alternating sign

[TABLE]

whenever $x_{k,3n}\to x\in(-d,d)\cup(a,b)$.

Proof.

Let $x\in(-b,-a)$, then we use Lemma 5.1 to see that $\lambda_{k,3n}^{(1)}q_{n}(x_{k})r_{n}(x_{k})=\lambda_{k,n}$, where $\lambda_{k,n}$ are the Gauss quadrature weights for the weight function $q_{n}(x)r_{n}(x)e^{-nV_{1}(x)}$. We can use the Chebyshev-Markov-Stieltjes inequalities [12, §3.41] for the Gauss quadrature weights to find

[TABLE]

By the mean value theorem, we have

[TABLE]

for some $\xi_{n}\in(x_{k-1},x_{k+1})$. Then, since $x_{k+1}-x_{k-1}\leq b-a$, we find

[TABLE]

whenever $x_{k}\to x\in(-b,-a)$, since

[TABLE]

and

[TABLE]

Let $x\in(-d,d)$, then we use Lemma 5.2 to find

[TABLE]

For the polynomials $p_{n}$ and $r_{n}$ one has

[TABLE]

uniformly in $x\in[-d,d]$, which already gives

[TABLE]

For the integral we use the infinite-finite range inequality (see Proposition 3.1) to find

[TABLE]

For $\hat{c}$ sufficiently large the intervals $I_{1}$, $[-d,d]$ and $[a,b]$ are disjoint, hence for $x\in I_{1}$ and $x_{k}\in[-d,d]$ we have $|x-x_{k}|>\textup{dist}(I_{1},[-d,d])=\delta_{1}$. We thus have

[TABLE]

Observe that the integrand is

[TABLE]

and as $n\to\infty$ the $n$th root thus converges to $\exp(-\ell_{1})$ when $x\in[-b,-a]$ or $\leq\exp(-\ell_{1})$ when $x\notin[-b,-a]$. We thus have (see the third Corollary [10, p. 199] for an Angelesco system)

[TABLE]

The $n$th root behavior of $|q_{n}^{\prime}(x_{k})|$ is more difficult because we evaluate $q_{n}^{\prime}$ at a point in $(-d,d)$, which is on the support of $\nu_{2}$ where the zeros of $q_{n}$ are dense. Clearly $q_{n}^{\prime}$ has $n-1$ zeros between the zeros of $q_{n}$ and the asymptotic distribution of the zeros of $q_{n}^{\prime}$ is the same as that of $q_{n}$, hence $|q_{n}^{\prime}(x)|^{1/n}$ converges to \exp\bigl{(}-U(x;\nu_{2})\bigr{)} whenever $x\notin[-d,d]$. When $x_{k}\to x\in(-d,d)$ one can use the principle of descent [11, Thm. 6.8 in Ch. I] to find

[TABLE]

To prove the inequality in the other direction, we look at the quadrature weights $\lambda_{k,3n}^{(2)}$ for the second integral (5.2) corresponding to the nodes on $[-d,d]$ (the zeros of $q_{n}$). These nodes are positive and related to the Gauss quadrature nodes for the orthogonal polynomials with weight function $p_{n}(x)r_{n}(x)e^{-nx^{2}}$, see Theorem 5.5. The result corresponding to (5.4) is

[TABLE]

On the other hand, by taking $f(x)=p_{n}(x)q_{n}^{2}(x)r_{n}(x)/(x-x_{k})^{2}$ in (5.2), we see that the quadrature weight $\lambda_{k,3n}^{(2)}$ satisfies

[TABLE]

Observe that the sign of $r_{n}(x)$ on $[-d,d]$ is $(-1)^{n}$, hence by the infinite-finite range inequalities, one finds

[TABLE]

where $|r_{n}|<1$. On $I_{2}$ we have that $|x-x_{k}|\leq\delta_{2}$ where $\delta_{2}$ is the length of $I_{2}$, hence

[TABLE]

from which we find

[TABLE]

By taking the $n$th root and by using the same reasoning as before, we thus find

[TABLE]

and since $x\in(-d,d)$, it follows from (4.3) that the right hand side is \exp\bigl{(}-2U(x;\nu_{2})\bigr{)}. Combined with (5.6) we then have

[TABLE]

whenever $x_{k}\to x\in(-d,d)$. Combining all these results gives (5.5) for $x_{k,3n}\to x\in(-d,d)$. The proof for $x_{k,3n}\to x\in(a,b)$ is similar, using Lemma 5.3. ∎

The results corresponding the quadrature weights for the second integral (5.2) and the thirs integral (5.3) are:

Theorem 5.8**.**

Suppose $\hat{c}$ is sufficiently large (see the footnote in Theorem 5.4). For the positive quadrature weights one has

[TABLE]

whenever $x_{k}\to x\in(-d,d)$. For the quadrature weights with alternating sign

[TABLE]

whenever $x_{k,3n}\to x\in(-b,-a)\cup(a,b)$.

Theorem 5.9**.**

Suppose $\hat{c}$ is sufficiently large (see the footnote in Theorem 5.4). For the positive quadrature weights one has

[TABLE]

whenever $x_{k}\to x\in(a,b)$. For the quadrature weights with alternating sign

[TABLE]

whenever $x_{k,3n}\to x\in(-b,-a)\cup(-d,d)$.

All the upper bounds in Theorems 5.7–5.9 depend on the logarithmic potentials $U(x;\nu_{1})$, $U(x;\nu_{2}),U(x;\nu_{3})$, and in particular on the linear combination of them that appears in the variational conditions (4.1)–(4.6). Observe that by combining (4.3) with (5.5) we find

[TABLE]

whenever $x_{k}\to x\in(-d,d)$, and by using (4.5) we find

[TABLE]

whenever $x_{k}\to x\in(a,b)$. Hence on $(-d,d)$ the quadrature weights are bounded from above by $e^{-V_{2}(x)}$ times a factor which is small, since $\ell_{2}-U(x;\nu_{2})-\ell_{1}+U(x;\nu_{1})<0$ on $[-d,d]$. On $(a,b)$ the quadrature weights are bounded by $e^{-V_{3}(x)}$ times an even smaller factor, since $\ell_{3}=\ell_{1}$ (by symmetry) and $U(x;\nu_{3})>U(x;\nu_{2})>U(x;\nu_{1})$ for $x\in(a,b)$, see Figure 5. This makes the alternating quadrature weights exponentially small.

6 Numerical example

In Table 1 and Figure 4 we give the quadrature weights $\lambda_{k,3n}^{(1)}$ for the zeros of $H_{10,10,10}$ with $c=15$, which after scaling by $\sqrt{10}$ corresponds to $\hat{c}=4.734$.

This clearly shows that the first 10 zeros are positive and the remaining 20 zeros are alternating in sign and very small in absolute value. The zeros and the quadrature weights behave in a similar way as in an Angelesco system (see [8]) when $\hat{c}$ is sufficiently large. Our scaling and the use of the weight functions

[TABLE]

means that we are using the densities of normal distributions with means $-\hat{c}/2,0,\hat{c}/2$ and variance $\sigma^{2}=1/2n$. In such case we can ignore the alternating weights and only use the positive quadrature weights $\{\lambda_{k,3n}^{(1)}:1\leq k\leq 3n\}$ to approximate the first integral (5.1). In a similar way, when we approximate the second integral (5.2) we can ignore the alternating weights and only use the positive weights $\{\lambda_{k,3n}^{(2)}:n+1\leq k\leq 2n\}$, and for approximating the third integral, one can only use $\{\lambda_{k,3n}^{(3)}:2n+1\leq k\leq 3n\}$.