MATH 417 Lecture 10

« previous | Thursday, February 13, 2014 | next »

Section 4.3: Numerical Integration

If we have a function $f(x)$ that is impossible (or really hard) to integrate exactly on an interval $\left[a,b\right]$, we can estimate the integral based on the values of $f$ at $n+1$ points $x_{0},x_{1},\ldots ,x_{n}$.

1. Choose a partition of the interval $\left[a,b\right]$ → $\left\{I_{j}\right\}_{j=1}^{N}$
2. Pick $\xi _{j}\in I_{j}$
3. Form the riemann sum $\sum _{j=1}^{N}f(\xi _{j})\,\left|I_{j}\right|$.

We usually define partition intervals by a diameter $h=\max _{j}\left|I_{j}\right|$. The Riemann integral, if it exists, is the limit of the sum as $h\to 0$.

In practice, we approximate $f(x)$ with a polynomial $p_{n}(x)$, and then integrate the polynomial:

$\int _{a}^{b}f(x)\,\mathrm {d} x=\int _{a}^{b}p_{n}(x)\,\mathrm {d} x+\int _{a}^{b}{\frac {f^{(n+1)}(\xi (x))}{(n+1)!}}\,\left(x-x_{0}\right)\dots \left(x-x_{n}\right)\,\mathrm {d} x$

Our error is going to be 0 if $f$ is a polynomial of degree $n$.

We can get an upper bound on this error by approximating $x-x_{0}\leq b-a$. Hence

$\left|\int _{a}^{b}{\frac {f^{(n+1)}(\xi (x))}{(n+1)!}}\,\left(x-x_{0}\right)\dots \left(x-x_{n}\right)\,\mathrm {d} x\right|\leq {\frac {C}{(n+1)!}}\,\left|b-a\right|^{n+1}$

For equidistant points, this bound can be better estimated by

$\left|\int _{a}^{b}{\frac {f^{(n+1)}(\xi (x))}{(n+1)!}}\,\left(x-x_{0}\right)\dots \left(x-x_{n}\right)\,\mathrm {d} x\right|\leq {\frac {C}{(n+1)!}}\,h^{n}$

1. Construct a numerical rule which is exact for degree at most $n$
2. Find the degree of accuracy of the rule

Quadrature Rule

Given $x_{0},x_{1},\ldots ,x_{n}$,

$\int _{a}^{b}f(x)\,\mathrm {d} x\approx \int _{a}^{b}p_{n}(x)\,\mathrm {d} x=\int _{a}^{b}\sum _{j=0}^{n}f(x_{j})\,\prod _{i=0;i\neq j}^{n}{\frac {x-x_{i}}{x_{j}-x_{i}}}\,\mathrm {d} x=\sum _{j=0}^{n}f(x_{j})\,\int _{a}^{b}p_{n,j}(x)\,\mathrm {d} x$

Thus

$\int _{a}^{b}f(x)\,\mathrm {d} x\approx Q_{n}(f)=\sum _{j=0}^{n}A_{j}\,f(x_{j})$, where $A_{j}=\int _{a}^{b}p_{n,j}(x)\,\mathrm {d} x=\int _{a}^{b}{\frac {(x-x_{0})\dots (x-x_{n})}{(x_{j}-x_{0})\dots (x_{j}-x_{n})}}\,\mathrm {d} x$

In a simple case where $x_{0}=a$ and $x_{1}=b$, we get the trapezoid rule.

The degree of accuracy is between $n$ and $2n+1$.

To find $A_{i}$'s, we can set up a linear system instead of integrating the Lagrange polynomials

${\begin{aligned}b_{0}=\int _{a}^{b}1\,\mathrm {d} x&=A_{0}\cdot 1+A_{1}\cdot 1+\dots +A_{n}\cdot 1\\b_{1}=\int _{a}^{b}x\,\mathrm {d} x&=A_{0}\cdot x_{0}+A_{1}\cdot x_{1}+\dots +A_{n}\cdot x_{n}\\b_{2}=\int _{a}^{b}x^{2}\,\mathrm {d} x&=A_{0}\cdot x_{0}^{2}+A_{1}\cdot x_{1}^{2}+\dots +A_{n}\cdot x_{n}^{2}\\\vdots &=\vdots \\b_{n}=\int _{a}^{b}x^{n}\,\mathrm {d} x&=A_{0}\cdot x_{0}^{n}+A_{1}\cdot x_{1}^{n}+\dots A_{n}\cdot x_{n}^{n}\end{aligned}}$

Change of Variables

Suppose we have a rule for $f(x)$ when $x\in \left[-1,1\right]$. We can scale the rule to fit an arbitrary interval $y\in \left[a,b\right]$ as follows:

$-1\leq x\leq 1$

$a\leq y(x)\leq b$

$y(x)=\alpha \,x+\beta ={\frac {b-a}{2}}x+{\frac {a+b}{2}}$

$f(x)=g(y(x))$, where $g$ is a polynomial of the same degree:

Performing a substitution in the integral $\int _{-1}^{1}g(x)\,\mathrm {d} x$ gives Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \int_{-1}^1 g(y(x)) \, y'(x) \,\mathrm{d}x}

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle y'(x) = \frac{b-a}{2}}

This integral can now be approximated by Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle B_0 \, g(y(-1)) + B_1 \, g(y(0)) + B_2 \, g(y(1))} , where

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle B_j = \frac{b-a}{2} \, A_j}

Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \begin{align} f(x) &= 1 + x + x^2 \\ g(y) &= 1 + \frac{y+1}{2} + \left( \frac{y+1}{2} \right)^2 \end{align}}

If Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \int_{a}^{b} f(x) \,\mathrm{d}x = \sum_{j=0}^n A_j \, f(x_j)} for Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle f(x) = 1, x, x^2, \ldots, x^4} , then it is also exact for Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle p(x)} of degree Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \le n} .

Take Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle \int_{a}^{b} f(x) \,\mathrm{d}x = \int_{a}^{b} \sum_{k=0}^n c_k \, x^k = \sum_{k=0}^n c_k \, \int_{a}^{b} x^k \,\mathrm{d}x = \sum_{k=0}^n c_k \, \sum_{j=0}^n A_j \, x_j^k = \sum_{k=0}^n \sum_{j=0}^n A_j \, c_k \, x_j^k = \sum_{j=0}^n A_j \sum_{k=0}^n c_k \, x_j^k = \sum_{j=0}^n A_j \, f(x_j)}