Gaussian quadrature

Each monomial $x^j$ with $j\le 2k$ is integrable by [F1]. Expanding the finite coefficient sums in [F2] expresses $P_c{P}_d$ as a finite real linear combination of such monomials, so [F5] makes the product integrable. Distributivity and [F5] then give bilinearity in each coefficient vector.

Real multiplication is commutative by [F7], so $P_c{P}_d=P_d{P}_c$ pointwise and the pairing is symmetric.

If $c\mathrm{\ne }0$, its coefficient list defines a nonzero formal polynomial $P_c$ of degree at most $k$. By [F3], its zero set has at most $k$ points. Since the support in [F1] has at least $m+1\ge k+1$ distinct points, choose a support point $s$ with $P_c(s)\mathrm{\ne }0$.

Put $\alpha =P_c(s{)}^2>0$. Continuity in [F6] gives a relative open neighborhood $U$ of $s$ on which $P_c(x{)}^2>\alpha \mathrm{/}2$. This set is Borel and has $\mu (U)>0$ by the support definition. Monotonicity, homogeneity and the indicator identity in [F4] imply $⟨c,c{⟩}_{\mu }=\int P_c^{2}d\mu \ge \frac{\alpha }{2}\mu (U)>0.$

The pairing is well-defined by step 1.1, real-bilinear and symmetric by steps 1.1–1.2, nonnegative on the diagonal, and zero on the diagonal only for $c=0$ by step 2.1. These are the real inner product axioms in [F8].

By [F2], the pairing on the coefficient space ${\mathbb{R}}^m$ for degrees at most $m-1$ is an inner product. The list $e_0,\dots ,e_{m-1}$ is independent by [F4], so [F3] gives orthonormal coefficient vectors $u_0,\dots ,u_{m-1}$ whose first $j+1$ vectors span the same space as $e_0,\dots ,e_j$. The polynomial $P_{u_j}$ therefore has degree exactly $j$: it lies in the span of $1,x,\dots ,x^j$ but not in the span through degree $j-1$. Dividing by its nonzero leading coefficient, which is invertible by [F8], gives a monic $p_j$, orthogonal to every polynomial of degree less than $j$.

Put $G_{ik}=\int x^{i+k}d\mu$ for $0\le i,k<m$, and $b_i=\int x^{i+m}d\mu$. All entries are finite by [F1]. The matrix $G$ represents the inner product of [F2] in the standard basis. If $Gc=0$, then $c^{\mathsf{T}}Gc=⟨c,c{⟩}_{\mu }=0$, so [F2] gives $c=0$. Thus the linear map $G:{\mathbb{R}}^m\to {\mathbb{R}}^m$ has trivial kernel; [F5] makes it surjective. There is therefore a unique coefficient vector $c$ solving $Gc=-b$.

Define $p_m(X)=X^m+{\sum }_{k=0}^{m-1}{c}_k{X}^k$. Its leading coefficient is $1$, and the $i$th equation of $Gc=-b$ gives $\int p_m(x)x^{i}d\mu (x)=0$ for each $i<m$. Linearity [F6] extends this equality to every polynomial of degree less than $m$. If another monic polynomial had the same orthogonality, its difference from $p_m$ would have degree at most $m-1$ and pair to zero with itself; positive definiteness [F2] makes that difference zero.

Step 1.1 gives the unique monic polynomials below degree $m$, and steps 1.2–2.1 give the unique degree-$m$ polynomial using only moments through $2m-1$. This proves the assertion for every $j=0,\dots ,m$.

For $j<m$, $p_j^2$ and $x{p}_j^2$ have degrees at most $2j$ and $2j+1$, both at most $2m-1$. Thus the moments in $h_j$ and ${\alpha }_j$ are finite by [F1]. Since $p_j$ is nonzero, [F2] gives $h_j>0$, so ${\alpha }_j$ is defined; for $j\ge 1$, the same positivity for $p_{j-1}$ gives ${\beta }_j>0$.

The polynomials $p_0,\dots ,p_j$ form a basis for the polynomials of degree at most $j$: their distinct degrees and unit leading coefficients permit successive subtraction of the top term. Because $x{p}_j-p_{j+1}$ has degree at most $j$ by cancellation of the leading term, write it in this basis.

For $k\le j-2$, symmetry of the real pairing gives $⟨x{p}_j,p_k{⟩}_{\mu }=⟨p_j,x{p}_k{⟩}_{\mu }=0$, since $x{p}_k$ has degree at most $j-1$. As the basis in step 1.2 is orthogonal and each squared norm is positive by [F2], the expansion has no $p_k$ term for $k\le j-2$. Its $p_j$ coefficient is $⟨x{p}_j,p_j{⟩}_{\mu }\mathrm{/}{h}_j={\alpha }_j$.

For $j\ge 1$, the $p_{j-1}$ coefficient is $\frac{⟨x{p}_j,p_{j-1}{⟩}_{\mu }}{h_{j-1}}=\frac{⟨p_j,x{p}_{j-1}{⟩}_{\mu }}{h_{j-1}}=\frac{h_j}{h_{j-1}}={\beta }_j\mathrm{.}$ Indeed, $x{p}_{j-1}-p_j$ has degree at most $j-1$, so it is orthogonal to $p_j$ by [F3].

The expansion in step 1.2 therefore has only the displayed $p_j$ and, when $j\ge 1$, $p_{j-1}$ terms. Rearranging it yields the stated recurrence for every $0\le j<m$.

Call an interior root sign-changing if $p_m$ has opposite nonzero signs on the two adjacent intervals containing no root. Let $R$ be the set of such roots; it is finite by [F3]. Suppose, for contradiction, that $R$ has $k<m$ elements, and put $q(X)={\prod }_{s\in R}(X-s)$, with the empty product equal to $1$. Between consecutive roots of $p_m$, the product $h=p_{m}q$ is continuous and nonzero, so its image is an interval avoiding $0$ and its sign is constant by [F7]–[F8]. Crossing a root in $R$ flips both factors; crossing any other root of $p_m$ flips neither product sign. Thus $h$ has one weak sign on $[a,b]$.

Since $\mathrm{deg}q=k<m$, orthogonality [F2] gives $\int hd\mu =0$.

The polynomial $p_m$ has at most $m$ distinct roots by [F3]. The support condition [F1] therefore supplies a support point $y$ with $p_m(y)\mathrm{\ne }0$; as every root of $q$ is a root of $p_m$, also $h(y)\mathrm{\ne }0$. Continuity of $h$ gives a relative open neighborhood $U$ of $y$ on which $\mathrm{\mid }h\mathrm{\mid }\ge c$ for some $c>0$ and the sign of $h$ is unchanged. This neighborhood is Borel and $\mu (U)>0$ by the definition of support.

If $h\ge 0$ on $[a,b]$, [F4] gives $\int hd\mu \ge c\mu (U)>0$. If $h\le 0$, apply [F4] to $-h$ and use [F5] to obtain $\int hd\mu <0$. Both contradict step 2.1, so $k\ge m$.

By [F3], $p_m$ has at most $m$ distinct roots, and step 3.1 gives at least $m$, so $k=m$. Write them as $x_1,\dots ,x_m\in (a,b)$ and let $s(X)={\prod }_{i=1}^m(X-x_i)$. Both $s$ and $p_m$ are monic of degree $m$ by [F2, F6]. Their degree-$m$ coefficients cancel, so their difference is zero or has degree at most $m-1$ by [F9]. It vanishes at all $m$ distinct $x_i$; if nonzero, [F3] is contradicted. Thus $p_m=s$, so each root occurs as exactly one linear factor and is simple.

For each $i$, set $h_i(t)=(1-2{\mathrm{\ell }}_i^{\mathrm{\prime }}(x_i)(t-x_i)){\mathrm{\ell }}_i(t{)}^2,k_i(t)=(t-x_i){\mathrm{\ell }}_i(t{)}^2\mathrm{.}$ At $x_i$, ${\mathrm{\ell }}_i(x_i)=1$, so $h_i(x_i)=1$, $h_i^{\mathrm{\prime }}(x_i)=0$, $k_i(x_i)=0$, and $k_i^{\mathrm{\prime }}(x_i)=1$ by the product rule [F3]. At $x_j$ with $j\mathrm{\ne }i$, ${\mathrm{\ell }}_i(x_j)=0$, so both ${\mathrm{\ell }}_i^2$ and its derivative vanish there; hence $h_i(x_j)=h_i^{\mathrm{\prime }}(x_j)=k_i(x_j)=k_i^{\mathrm{\prime }}(x_j)=0$. Therefore the displayed sum defining $H$ satisfies every value and first-derivative condition.

Each ${\mathrm{\ell }}_i$ has degree at most $m-1$ by [F1]. Thus ${\mathrm{\ell }}_i^2$ has degree at most $2m-2$, and multiplication by an affine polynomial gives degree at most $2m-1$ by [F8]. The finite sum $H$ also has degree at most $2m-1$.

Suppose $x$ is not a node. By [F1], $\omega (x{)}^2\mathrm{\ne }0$, so put $c=(f(x)-H(x))\mathrm{/}\omega (x{)}^2$ and define $g(t)=f(t)-H(t)-c\omega (t{)}^2$. At every node, step 1.1 and $\omega (x_i)=0$ give $g(x_i)=0$ and $g^{\mathrm{\prime }}(x_i)=0$, since $({\omega }^2{)}^{\mathrm{\prime }}=2\omega {\omega }^{\mathrm{\prime }}$ by [F3]; by the choice of $c$, also $g(x)=0$.

Sort the $m$ nodes and $x$ as $z_0<\cdots <z_m$. Rolle's theorem on each consecutive interval gives a zero of $g^{\mathrm{\prime }}$ strictly between its endpoints; these $m$ zeros are distinct and are not nodes. Along with $g^{\mathrm{\prime }}(x_i)=0$ at the $m$ nodes from step 2.1, this gives $2m$ distinct ordered zeros of $g^{\mathrm{\prime }}$.

The function $g^{\mathrm{\prime }}$ is continuous on the closed span of those zeros and is $(2m-1)$-times differentiable on its interior by [F2]–[F4]. Apply [F5] to those $2m$ zeros. There is ${\xi }_x\in (a,b)$ with $g^{(2m)}({\xi }_x)=0$.

By step 1.2, $H^{(2m)}=0$. By [F6], ${\omega }^2$ is monic of degree $2m$, so [F7] gives $({\omega }^2{)}^{(2m)}=(2m){!}_{\mathbb{R}}$. Differentiating the definition of $g$ in step 2.1 at ${\xi }_x$ therefore gives $0=f^{(2m)}({\xi }_x)-c(2m){!}_{\mathbb{R}}$.

By [F9], $(2m){!}_{\mathbb{R}}\mathrm{\ne }0$, so step 5.1 gives $c=f^{(2m)}({\xi }_x)\mathrm{/}(2m){!}_{\mathbb{R}}$. Substituting this into the definition of $c$ in step 2.1 proves the remainder formula. At a node, step 1.1 gives zero interpolation error and $\omega =0$.

By [F1] and [F5], the normalized polynomials $q_j:=p_j\mathrm{/}\sqrt{h_j}$ satisfy $x{q}_j=\sqrt{{\beta }_{j+1}}{q}_{j+1}+{\alpha }_j{q}_j+\sqrt{{\beta }_j}{q}_{j-1}$ for $0\le j<m-1$, with $q_{-1}=0$. For the final index, the monic recurrence gives $x{p}_{m-1}=p_m+{\alpha }_{m-1}{p}_{m-1}+{\beta }_{m-1}{p}_{m-2}$; at a root $x_i$ the $p_m$ term vanishes. After dividing by $\sqrt{h_{m-1}}$, this is the final row of $J_m{z}_i=x_i{z}_i$ for $z_i=(q_0(x_i),\dots ,q_{m-1}(x_i){)}^{\mathsf{T}}$. When $m=1$, the identity $x{p}_0=p_1+{\alpha }_0{p}_0$ gives the same conclusion at the sole root.

The vector $z_i$ is nonzero because $q_0=1\mathrm{/}\sqrt{h_0}$ and $h_0=\mu ([a,b])>0$. Thus each of the $m$ distinct roots in [F2] is an eigenvalue of $J_m$.

The matrix $J_m$ is real symmetric by its definition, so [F4] makes its endomorphism self-adjoint on ${\mathbb{R}}^m$. By [F3], it has an orthonormal eigenbasis of $m$ vectors; hence it has exactly $m$ eigenvalues counted with multiplicity. The $m$ distinct eigenvalues already exhibited in step 2.1 exhaust this list, so each is simple and there are no others. By [F6], this eigenvalue set is $\mathrm{spec}(J_m)$.

The reduced-order uses in [F1] are valid: the original support has at least $m+1\ge k+1$ points and the absolute moments through $2k-1$ are among those assumed. By [F2], every $p_k$ has the stated simple roots and monic factorization. If $m=1$, there is no pair to prove. For the remaining case $m\ge 2$, take $j=1$; if $x_1$ is the root of $p_1$, [F3] gives $p_2(x_1)=-{\beta }_1<0$. At both $a$ and $b$, the degree-two monic factorization [F2] gives $p_2(a)>0$ and $p_2(b)>0$. By [F4] there is one root on each side of $x_1$ inside $(a,b)$, proving the first pair when $m\ge 2$.

For $2\le j<m$, assume the roots of $p_{j-1}$ strictly interlace those $x_1<\cdots <x_j$ of $p_j$. At $x_i$, exactly $j-i$ roots of the monic polynomial $p_{j-1}$ lie to its right and none equals it. Its factorization [F2] therefore makes its sign $(-1{)}^{j-i}$. By [F3], the sign of $p_{j+1}(x_i)=-{\beta }_j{p}_{j-1}(x_i)$ is $(-1{)}^{j-i+1}$. These signs alternate, so [F4] supplies a root in each $(x_i,x_{i+1})$. Also [F2] gives sign $(-1{)}^{j+1}$ for $p_{j+1}(a)$ and positive sign for $p_{j+1}(b)$; these are opposite to the signs $(-1{)}^j$ at $x_1$ and $-1$ at $x_j$. Thus [F4] supplies an additional root in each of $(a,x_1)$ and $(x_j,b)$.

These $j+1$ disjoint intervals contain $j+1$ distinct roots, and [F2] says there are exactly that many. They therefore exhaust the roots of $p_{j+1}$ in the required order. Starting with step 1.1, finite induction through step 1.2 proves every pair through $p_m$. By [F6], the same ordered roots are the eigenvalues of the leading principal truncations, since their recurrence coefficients are identical prefixes. This proves the spectral assertion as well; when $m=1$ the quantified pair set is empty.