Gaussian Quadrature

In Newton-Cotes formula, we choose $n+1$ evenly spaced points on the interval to do Lagrange interpolation and deduce a formula from the integral of that polynomial. To improve the precision, we will choose $n$ special points in the interval; with the same interpolation–integration approach, we can achieve a degree of exactness of $2n-1$. (Recall: When an algorithm works exactly right for all polynomials of degree $\leq m$ but not for degree $m+1$, we say that it has a degree of exactness of $m$. )

One of the typical solutions is to choose the $n$ roots of the Legendre polynomial $P_n(x)$. The Legendre polynomials form a complete orthogonal system on $[-1,1]$ with weight function $w(x)=1$ (where $P_n(x)$ is a polynomial of degree $n$), i.e.

$$\int_{-1}^{1}{P_m(x)P_n(x)\,\mathrm{d}x}=0\quad \text{if}\, m\neq n.$$

With the standardization condition $P_n(1)=1$ for all $n$, all the polynomials can be uniquely determined. In fact, by Rodrigues' Formula,

$$P_n(x)=\frac{1}{2^{n}n!}\frac{\mathrm{d}^{n}}{\mathrm{d}x^n}{(x^2-1)^n}.$$

Let $x_1,x_2,\cdots,x_n$ be the roots of $P_n(x)$. We squeeze our target interval $[a,b]$ to $[-1,1]$ and do Lagrange interpolation with points $(x_i,g(x_i))$, where $g(x)$ is the transformed $f(x)$. Gaussian Quadrature is then

$$\int_{a}^{b}{f(x)\,\mathrm{d}x}\approx \frac{b-a}{2}\sum_{i=1}^{n}{w_{i}f(\frac{b-a}{2}x_i+\frac{a+b}{2})}$$

where the weights $w_i$ can be found similarly to the Newton–Cotes $\alpha_i$ values, or via the closed form

$$w_i=\frac{2}{(1-x_i^2)[P_n'(x_i)]^2}.$$

We have known the procedure of Gaussian quadrature. But why are Legendre polynomials so useful that our exactness is doubled? It all comes from the orthogonality. Say we have a polynomial $T(x)$ whose degree is no more than $2n-1$. Divide $T(x)$ by $P_n(x)$and we have $T(x)=Q(x)P_n(x)+R(x)$, where $deg\,Q\leq n-1,\, deg\,R\leq n-1$. We write $Q(x)$ as a linear combination of lower-order Legendre polynomials: $Q(x)=\sum_{i=0}^{n-1}{k_{i}P_i(x)}$. Then by the definition of Legendre polynomials, $\int_{-1}^{1}{P_i(x)P_n(x)=0}$, thus $\int_{-1}^{1}{Q(x)P_n(x)=0}$. Now we are ready to prove that

$$\int_{-1}^{1}{T(x)\,\mathrm{d}x}=\sum_{i=1}^{n}{w_{i}T(x_i)}.$$

For the left side,

$$\int_{-1}^{1}{T(x)\,\mathrm{d}x}=\int_{-1}^{1}{Q(x)P_n(x)\,\mathrm{d}x}+\int_{-1}^{1}{R(x)\,\mathrm{d}x}=\int_{-1}^{1}{R(x)\,\mathrm{d}x}.$$

For the right side,

$$\forall i=1,2,\cdots,n\quad P_n(x_i)=0$$ $$\begin{aligned} \sum_{i=1}^{n}{w_{i}T(x_i)} &=\sum_{i=1}^{n}{w_{i}[Q(x_i)P_n(x_i)+R(x_i)]} \\ &=\sum_{i=1}^{n}{w_{i}R(x_i)}. \end{aligned}$$

By its construction, an n-point Gaussian quadrature rule is exact for any polynomial of degree up to $n−1$. Since $deg\,R\leq n-1$, the quadrature gives the exact integral:

$$\int_{-1}^{1}{R(x)\,\mathrm{d}x}=\sum_{i=1}^{n}{w_{i}R(x_i)}.$$

Thus our quadrature works exactly right for $T(x)$. On the other hand, the error

$$E_n(f)=\frac{2^{2n+1}(n!)^4}{(2n+1)[(2n)!]^3}f^{(2n)}(\xi)$$

for some $\xi\in(-1,1)$. Therefore, the degree of exactness is $2n-1$. Above is all of the magic of Gaussian quadrature.