My name is Hanzhang Yin, and I have developed this website as a resource to facilitate the review of key concepts in abstract algebra, and it is not fully completed. Feel free to email me at hanyin@ku.edu if there are any errors or suggestions for improvement.

Binary

$$x=\pm (0.a_1\cdots a_n\cdots {)}_2\times 2^m$$ where $a_1=1$

$\mathbf{\text{Note. }}$For Marc-32, $n=24$, $-126\le m\le 127$. For Marc-64, $n=53$, $-1021\le m\le 1024$.

Unit Roundoff

$$fl(1+\epsilon )>1$$ $$\epsilon =2^{-24}{10}^{-8}$$ $$\epsilon =2^{-53}{10}^{-16}$$

Roundoff Error

Let $\delta =\frac{x-fl(x)}{x}$. define $$\left|\delta \right|\le \epsilon$$ Suppose that $x=\pm (0.a_1\cdots a_n\cdots {)}_2\times 2^m$, then we pick two points such that $x_{-}$ is less than $x$ and $x_{+}$ is greater than $x$. Then we have $$x_{-}=(0,a_1\cdots a_n{)}_2\times 2^m$$ $$x_{+}=(0,a_1\cdots a_n+0.000\dots 1_n{)}_2\times 2^m$$ $$fl(x)=x_{-}\text{ if }{a}_{n+1}=0\text{(x is closer to x\_-than x\_+) and }fl(x)=x_{+}\text{ if }{a}_{n+1}=1\text{ (x is closer tox\_+ than x\_-)}$$ $$|fl(x)-x|\le \frac{1}{2}\times (x_{+}-x_{-})=\frac{1}{2}\cdot 2^{-n}\cdot 2^m$$ $$\left|\frac{fl(x)-x}{x}\right|\le \frac{\frac{1}{2}{2}^{-n}{2}^m}{(0,a_1\cdots a_n{)}_2\times 2^m}=2^{-n}$$

Number System with base $\beta$

The unit roundoff is $$\epsilon =\frac{1}{2}{\beta }^{1-t}$$ where $t$ is the number of digits in the mantissa.

Floating point error analysis

Let $x,y$ be two machine numbers. Let $\circ$ is a binary operation. $$fl(x\circ y)=(x\circ y)(1+\delta )$$ where $|\delta |\le \epsilon$ is the roundoff error. $$fl(x)=x(1+\delta )$$

$\mathbf{\text{Example.}}$

Orders of Convergence

Let $\{x_n\}$ be a sequence of real numbers tending to a limit $x^{\ast }$. We say that the rate of convergence is at least linear if there are a constant $c<1$ and an integer $N$ such that $$|x_{n+1}-x^{\ast }|\le c|x_n-x^{\ast }|(n\ge N)$$ We say that the rate of convergence is at least superlinear if there exist a sequence ${ϵ}_n$ tending to 0 and an integer $N$ such that $$|x_{n+1}-x^{\ast }|\le {ϵ}_n|x_n-x^{\ast }|(n\ge N)$$ The convergence is at least quadratic if there are a constant $C$ (not necessarily less than 1) and an integer $N$ such that $$|x_{n+1}-x^{\ast }|\le C|x_n-x^{\ast }{|}^2(n\ge N)$$ In general, if there are positive constants $C$ and $\alpha$ and an integer $N$ such that $$|x_{n+1}-x^{\ast }|\le C|x_n-x^{\ast }{|}^{\alpha }(n\ge N)$$ we say that the rate of convergence is of order $\alpha$ at least.

"Numerical Analysis: Mathematics of Scientific Computing" page. 17.

Newton's Method

Definition. Newton's method, also known as the Newton-Raphson method, is an iterative root-finding algorithm used to find successively better approximations to the roots (or zeros) of a real-valued function. Given a function $f(x)$ and its derivative $f^{\prime }(x)$, Newton's method starts with an initial guess $x_0$ and iteratively refines this guess using the following formula: $$x_{n+1}=x_n-\frac{f(x_n)}{f^{\prime }(x_n)},$$ where:

• $x_n$ is the current approximation.
• $x_{n+1}$ is the next approximation.
• $f(x_n)$ is the value of the function at $x_n$.
• $f^{\prime }(x_n)$ is the value of the derivative of the function at $x_n$.

Convergence Criteria. The method typically converges if the initial guess $x_0$ is sufficiently close to the actual root, and $f(x)$ is continuously differentiable in the neighborhood of the root.

Example. Consider finding the root of the function $f(x)=x^2-2$. The derivative of the function is $f^{\prime }(x)=2x$. Firstly, we start with an initial guess $x_0$. Secondly, we apply the Newton's iteration: $$x_{n+1}=x_n-\frac{x_n^2-2}{2x_n}=\frac{2x_n-\frac{x_n^2-2}{x_n}}{2}=\frac{x_n+\frac{2}{x_n}}{2}$$ So, the iterative formula becomes: $$x_{n+1}=\frac{1}{2}(x_n+\frac{2}{x_n})$$ By repeatedly applying this formula, the sequence $\{x_n\}$ will converge to the square root of 2.

Error Analysis for Newton's Method

Let $r$ be a solution of $f(x)=0$ (i.e. $f(r)=0$). Suppose that we have already computed $x_n$ and the error in $x_n$ is $e_n=|x_n-r|$. We now derive a formula that relates the error after the next step, $|x_{n+1}-r|$, to $|x_n-r|$. According to Taylor's theorem, there is a $c$ between $x_n$ and $r$ such that $$0=f(r)=f(x_n)+f^{\prime }(x_n)(r-x_n)+\frac{1}{2}{f}^{″}(c)(r-x_n{)}^2(1)$$ By the definition of Newton's method, we have $$0=f(x_n)+f^{\prime }(x_n)(x_{n+1}-x_n)(2)$$ Subtracting (2) from (1). $$\begin{array}{rl}0& =f^{\prime }(x_n)(r-x_{n+1})+\frac{1}{2}{f}^{″}(c)(r-x_n{)}^2\\ x_{n+1}-r& =\frac{f^{″}(c)}{2f^{\prime }(x_n)}(x_n-r{)}^2\\ e_{n+1}& =|x_{n+1}-r|=\left|\frac{f^{″}(c)}{2f^{\prime }(x_n)}\right||x_n-r{|}^2\end{array}$$ If $x_n$ is close to $r$, then $c$, which must be between $x_n$ and $r$, is also close to $r$ and $$e_{n+1}=|x_{n+1}-r|\approx \left|\frac{f^{″}(r)}{2f^{\prime }(r)}\right||x_n-r{|}^2=\left|\frac{f^{″}(r)}{2f^{\prime }(r)}\right|e_n^2=C{e}_n^2$$

Thus, we can see that Newton's Method is quadratically convergent.

However, the disadvantage of Newton's Method is that if the initial value is closed enough to the root.

Theorems on Newton's Method

Theorem. Let $f^{″}$ be continuous and let $r$ be a simple zero of $f$. Then there is a neighborhood of $r$ and a constant $C$ such that if Newton's method is started in that neighborhood, the successive points become steadily closer to $r$ and satisfy $$|x_{n+1}-r|\le C(x_n-r{)}^2(n\ge 0)$$

$\textbf{Theorem on Newton's Method for a Convex Function. }$ If $f$ belongs to $C^2(\mathbb{R})$, is increasing, is convex, and has a zero, then the zero is unique, and the Newton iteration will converge to it from any starting point.

$\textbf{Convergence of Newton's method.}$ Assume that $f(x)$, $f^{\prime }(x)$, and $f^{″}(x)$ are continuous for all $x$ in some neighborhood of $\alpha$ and $f(\alpha )=0$ and $f^{\prime }(\alpha )\ne 0$. If $x_0$ is sufficiently close to $\alpha$, then $x_n\to \alpha$ as $n\to \mathrm{\infty }$. Moreover, $$\underset{n\to \mathrm{\infty }}{lim}\frac{|x_{n+1}-\alpha |}{|x_n-\alpha {|}^2}=\frac{1}{2}\left|\frac{f^{″}(\alpha )}{f^{\prime }(\alpha )}\right|.$$

Horner's Algorithm

This method is also known as Horner's method, nested multiplication, or synthetic division.

$\mathbf{\text{Definition. }}$ Horner's algorithm is a method for evaluating a polynomial at a point. Given a polynomial $p(x)=a_0+a_{1}x+a_2{x}^2+\cdots +a_n{x}^n$, Horner's algorithm evaluates $p(x)$ at a point $x=c$ as follows: $$p(c)=a_0+c(a_1+c(a_2+\cdots +c(a_{n-1}+c{a}_n)\cdots )).$$

Given a polynomial $p(x)=a_n{x}^n+a_{n-1}{x}^{n-1}+\cdots +a_{1}x+a_0$. Then $$\begin{array}{rl}{b}_{n-1}& =a_n\\ b_{n-2}& =a_{n-1}+z_0{b}_{n-1}\\ & ⋮\\ b_0& =a_1+z_0{b}_1\\ p(z_0)& =a_0+z_0{b}_0\end{array}$$ Horner's algorithm evaluates $p(x)$ at a point $x=z_0$ as follows: $$\begin{array}{cccccc}& a_n& a_{n-1}& a_{n-2}& \cdots & a_0\\ z_0& & z_0{b}_{n-1}& z_0{b}_{n-2}& \cdots & z_0{b}_0\\ & b_{n-1}& b_{n-2}& b_{n-3}& \cdots & \boxed{{\displaystyle b_{-1}}}\end{array}$$ The boxed number satisfies $p(z_0)=b_{-1}$.

$\mathbf{\text{Example. }}$ Given the polynomial $p(x)=4x^3-3x^2+2x-1$ and $z_0=2$, evaluate $p(z_0)$. $$\begin{array}{ccccc}& 4& -3& 2& -1\\ 2& & 8& 10& 24\\ & 4& 5& 12& \boxed{{\displaystyle 23}}\end{array}$$ Hence, we have $p(2)=23$.

$\mathbf{\text{Theorem on Horner's Method. }}$ Let $p(x)=a_n{x}^n+\cdots +a_{1}x+a_0$. Define pairs $({\alpha }_j,{\beta }_j)$ for $j=n,n-1,\dots ,0$ by the algorithm $$\{\begin{array}{l}({\alpha }_n,{\beta }_n)=(a_n,0)\\ ({\alpha }_j,{\beta }_j)=(a_j+x{\alpha }_{j+1},{\alpha }_{j+1}+x{\beta }_{j+1})& (n-1\ge j\ge 0)\end{array}$$ Then ${\alpha }_0=p(x)$ and ${\beta }_0=p^{\prime }(x)$.

Contractive Mapping

Newton's method and Steffensen's method are examples of procedures whereby a sequence of points is computed by a formula of the form $$\begin{array}{}\text{(1)}& x_{n+1}=F(x_n)(n\ge 0)\end{array}$$ The algorithm defined by such an equation is called functional iteration.

"Numerical Analysis: Mathematics of Scientific Computing" page. 100.

$\mathbf{\text{Definition. }}$ A mapping $F$ of a metric space $(X,d)$ into itself is called a contraction if there exists a constant $k$ with $0\le k<1$ such that $$d(F(x),F(y))\le kd(x,y)$$ for all $x,y$ in $X$.

"Numerical Analysis: Mathematics of Scientific Computing" page. 101.

Contractive Mapping Theorem

Theorem (Contractive Mapping Theorem). Let $C$ be a closed subset of the real line. If $F$ is a contractive mapping of $C$ into $C$, then $F$ has a unique fixed point. Moreover, this fixed point is the limit of every sequence obtained from Equation (1) with a starting point $x_0\in C$.

"Numerical Analysis: Mathematics of Scientific Computing" page. 102.

Polynomial Interpretation

Theorem on Polynomial Interpolation Error

$\mathbf{\text{Theorem. }}$Let $f$ be a function in $C^{n+1}[a,b]$, and let $p$ be the polynomial of degree at most $n$ that interpolates the function $f$ at $n+1$ distinct points $x_0,x_1,\dots ,x_n$ in the interval $[a,b]$. To each $x$ in $[a,b]$ there corresponds a point ${\xi }_x$ in $(a,b)$ such that $$f(x)-p(x)=\frac{1}{(n+1)!}{f}^{(n+1)}({\xi }_x)\prod _{i=0}^n(x-x_i).$$

Lagrange's Form

The Lagrange form of the interpolation polynomial is a method for constructing a polynomial that passes through a given set of points. Specifically, given a set of $n+1$ distinct data points $(x_0,y_0),(x_1,y_1),\dots ,(x_n,y_n)$, the Lagrange interpolation polynomial $P(x)$ is defined as: $$P(x)=\sum _{i=0}^n{y}_i{\ell }_i(x)$$ where ${\ell }_i(x)$ are the Lagrange basis polynomials given by: $${\ell }_i(x)=\prod _{\begin{array}{c}0\le j\le n\\ j\ne i\end{array}}\frac{x-x_j}{x_i-x_j}$$ Each ${\ell }_i(x)$ is a polynomial of degree $n$ that is 1 at $x=x_i$ and 0 at $x=x_j$ for $j\ne i$.

$\mathbf{\text{Example. }}$ Given points: $(1,1),(2,4),(3,9)$. Here $x_0=1,x_1=2,x_2=3$, and $y_0=1,y_1=4,y_2=9$.

Compute the Lagrange basis polynomials ${\ell }_i(x)$: $${\ell }_0(x)=\frac{(x-2)(x-3)}{(1-2)(1-3)}=\frac{(x-2)(x-3)}{2}$$ $${\ell }_1(x)=\frac{(x-1)(x-3)}{(2-1)(2-3)}=-(x-1)(x-3)$$ $${\ell }_2(x)=\frac{(x-1)(x-2)}{(3-1)(3-2)}=\frac{(x-1)(x-2)}{2}$$

Construct the interpolation polynomial $P(x)$: $$P(x)=1\cdot {\ell }_0(x)+4\cdot {\ell }_1(x)+9\cdot {\ell }_2(x).$$ $$P(x)=1\cdot \frac{(x-2)(x-3)}{2}-4\cdot (x-1)(x-3)+9\cdot \frac{(x-1)(x-2)}{2}.$$ Simplifying, we get: $$P(x)=\frac{(x-2)(x-3)}{2}-4(x-1)(x-3)+\frac{9(x-1)(x-2)}{2}.$$

Divided Differences

The divided differences of a set of data points are a sequence of numbers that are used to construct the Newton form of the interpolation polynomial.

$$f[x_i,x_{i+1},\dots ,x_{i+j}]=\frac{f[x_{i+1},x_{i+2},\dots ,x_{i+j}]-f[x_i,x_{i+1},\dots ,x_{i+j-1}]}{x_{i+j}-x_i}$$ The Following table shows the divided differences for a set of data points: $$\begin{array}{ccccc}{x}_0& f[x_0]& f[x_0,x_1]& f[x_0,x_1,x_2]& f[x_0,x_1,x_2,x_3]\\ x_1& f[x_1]& f[x_1,x_2]& f[x_1,x_2,x_3]& \\ x_2& f[x_2]& f[x_2,x_3]& & \\ x_3& f[x_3]& & & \end{array}$$

$\mathbf{\text{Example. }}$ Given the following table of data points: $$\begin{array}{ccccc}x& 3& 1& 5& 6\\ f(x)& 1& -3& 2& 4\end{array}$$ Solution: We arrange the given table vertically and compute divided differences by the use of Formula above, arriving at $$\begin{array}{ccccc}3& 1& 2& -3/8& 7/40\\ 1& -3& 5/4& 3/20& \\ 5& 2& 2& & \\ 6& 4& & & \end{array}$$

Theorem on Permutations in Divided Differences. The divided difference is a symmetric function of its arguments. Thus, if $(z_0,z_1,\dots ,z_n)$ is a permutation of $(x_0,x_1,\dots ,x_n)$, then $$f[z_0,z_1,\dots ,z_n]=f[x_0,x_1,\dots ,x_n].$$

Theorem on Derivatives and Divided Differences. If $f$ is $n$ times continuously differentiable on $[a,b]$ and if $x_0,x_1,\dots ,x_n$ are distinct points in $[a,b]$, then there exists a point $\xi$ in $(a,b)$ such that $$f[x_0,x_1,\dots ,x_n]=\frac{1}{n!}{f}^{(n)}(\xi ).$$

Newton's Form

The Newton form of the interpolation polynomial is a method for constructing a polynomial that passes through a given set of points. Specifically, given a set of $n+1$ distinct data points $(x_0,y_0),(x_1,y_1),\dots ,(x_n,y_n)$, the Newton interpolation polynomial $P(x)$ is defined as: $$P(x)=a_0+a_1(x-x_0)+a_2(x-x_0)(x-x_1)+\cdots +a_n(x-x_0)(x-x_1)\cdots (x-x_{n-1})$$ where $a_0,a_1,\dots ,a_n$ are the divided differences of the data points.

$\mathbf{\text{Remark. }}$Here, $a_n=f[x_0,x_1,\dots ,x_n]$.

$\mathbf{\text{Example. }}$ Given the following table of data points: $$\begin{array}{ccccc}x& 3& 1& 5& 6\\ f(x)& 1& -3& 2& 4\end{array}$$ We already found the divided differences in the previous example. The Newton form of the interpolation polynomial is $$p(x)=1+2(x-3)-\frac{3}{8}(x-3)(x-1)+\frac{7}{40}(x-3)(x-1)(x-5).$$

Theorem on Polynomial Interpolation Error

$\mathbf{\text{Theorem. }}$Let $f$ be a function in $C^{n+1}[a,b]$, and let $p$ be the polynomial of degree at most $n$ that interpolates the function $f$ at $n+1$ distinct points $x_0,x_1,\dots ,x_n$ in the interval $[a,b]$. To each $x$ in $[a,b]$ there corresponds a point ${\xi }_x$ in $(a,b)$ such that $$f(x)-p(x)=\frac{1}{(n+1)!}{f}^{(n+1)}({\xi }_x)\prod _{i=0}^n(x-x_i).$$

Hermite Interpolation

Hermite Interpolation Theorem. Let $n\ge 0$, and suppose that $x_i$, $i=0,\dots ,n$, are distinct real numbers. Then, given two sets of real numbers $y_i$, $i=0,\dots ,n$, and $z_i$, $i=0,\dots ,n$, there is a unique polynomial $p_{2n+1}$ in ${\mathcal{P}}_{2n+1}$ such that $$p_{2n+1}(x_i)=y_i,p_{2n+1}^{\prime }(x_i)=z_i,i=0,\dots ,n.(6.13)$$

Remark. In general, if there is a set of $m$ distinct points, there is a unique polynomial of degree $2m-1$ that interpolates the function at those points for hermite interpolation.

"Introduction to Numerical Analysis", page. 188.

Theorem. Suppose that $n\ge 0$ and let $f$ be a real-valued function, defined, continuous and $2n+2$ times differentiable on the interval $[a,b]$, such that $f^{(2n+2)}$ is continuous on $[a,b]$. Further, let $p_{2n+1}$ denote the Hermite interpolation polynomial of $f$ defined by (6.16). Then, for each $x\in [a,b]$ there exists $\xi =\xi (x)$ in $(a,b)$ such that $$f(x)-p_{2n+1}(x)=\frac{f^{(2n+2)}(\xi )}{(2n+2)!}[{\pi }_{n+1}(x){]}^2,(6.17)$$ where ${\pi }_{n+1}$ is as defined in (6.9). Moreover, $$|f(x)-p_{2n+1}(x)|\le \frac{M_{2n+2}}{(2n+2)!}[{\pi }_{n+1}(x){]}^2,(6.18)$$ where $M_{2n+2}=\underset{\zeta \in [a,b]}{max}|f^{(2n+2)}(\zeta )|$.

"Introduction to Numerical Analysis", page. 190.

Theorem on Hermite Interpolation Error Estimate Let $x_0,x_1,\dots ,x_n$ be distinct nodes in $[a,b]$ and let $f\in C^{2n+2}[a,b]$. If $p$ is the polynomial of degree at most $2n+1$ such that $$p(x_i)=f(x_i)\text{and}{p}^{\prime }(x_i)=f^{\prime }(x_i)(0\le i\le n)$$ then to each $x$ in $[a,b]$ there corresponds a point $\xi$ in $(a,b)$ such that $$f(x)-p(x)=\frac{f^{(2n+2)}(\xi )}{(2n+2)!}\prod _{i=0}^n(x-x_i{)}^2$$

"Numerical Analysis: Mathematics of Scientific Computing" page. 344.

Change of Intervals

Suppose that a numerical integration formula is given: $${\int }_c^{d}f(t)dt\approx \sum _{i=0}^n{A}_{i}f(t_i).$$ Let us assume that it is exact for all polynomials of degree $m$ or less. If we want to integrate the function $f$ over the interval $[a,b]$, we can use the change of interval formula. We first define a linear function $\lambda$ of $t$ such that if $t$ traverses $[c,d]$, $\lambda (t)$ will traverse $[a,b]$. The function $\lambda$ is given explicitly by $$\lambda (t)=\frac{b-a}{d-c}t+\frac{ad-bc}{d-c}.$$ Now in the integral $${\int }_a^{b}f(x)dx$$ we make the change of variable $x=\lambda (t)$. Then $dx={\lambda }^{\prime }(t)dt=(b-a)(d-c{)}^{-1}dt$, and so we have $${\int }_a^{b}f(x)dx=\frac{b-a}{d-c}{\int }_c^{d}f(\lambda (t))dt$$ $$\approx \frac{b-a}{d-c}\sum _{i=0}^n{A}_{i}f(\lambda (t_i))$$ Hence, we have $${\int }_a^{b}f(x)dx\approx \frac{b-a}{d-c}\sum _{i=0}^n{A}_{i}f(\frac{b-a}{d-c}{t}_i+\frac{ad-bc}{d-c})$$ Observe that, because $\lambda$ is linear, $f(\lambda (t))$ is a polynomial in $t$ if $f$ is a polynomial, and the degrees are the same. Hence, the new formula is exact for polynomials of degree $m$. For Simpson’s rule, this procedure yields Equation (6) as a consequence of Equation (5) using $\lambda (t)=(b-a)t+a$.

"Numerical Analysis: Mathematics of Scientific Computing," page. 486.

Gaussian Quadrature

The theory can be formulated for quadrature rules of a slightly more general form; namely, $${\int }_a^{b}f(x)w(x)dx\approx \sum _{i=0}^n{A}_{i}f(x_i),$$ where $w$ is a fixed positive weight function.

The Gauss quadrature rule: $$\begin{array}{}\text{(10.6)}& {\int }_a^{b}w(x)f(x)dx\approx G_n(f)=\sum _{k=0}^n{W}_{k}f(x_k),\end{array}$$ where the quadrature weights are $$\begin{array}{}\text{(10.7)}& W_k={\int }_a^{b}w(x)[L_k(x){]}^{2}dx,\end{array}$$ and the quadrature points $x_k,k=0,\dots ,n$, are chosen as the zeros of the polynomial of degree $n+1$ from a system of orthogonal polynomials over the interval $(a,b)$ with respect to the weight function $w$. Since this quadrature rule was obtained by exact integration of the Hermite interpolation polynomial of degree $2n+1$ for $f$, it gives the exact result whenever $f$ is a polynomial of degree $2n+1$ or less.

"Introduction to Numerical Analysis", page. 279.

Error Analysis for Gaussian Quadrature

Given the Gaussian quadrature formula $${\int }_a^{b}f(x)dx=\sum _{i=1}^n{w}_{i}f(x_i)+R_n(f),$$ where $$R_n(f)=\frac{(b-a{)}^{2n+1}(n!{)}^4}{(2n+1)[(2n)!{]}^3}{f}^{(2n)}(c).$$

Theorem. Suppose that $w$ is a weight function, defined, integrable, continuous and positive on $(a,b)$, and that $f$ is defined and continuous on $[a,b]$; suppose further that $f$ has a continuous derivative of order $2n+2$ on $[a,b]$, $n\ge 0$. Then, there exists a number $\eta$ in $(a,b)$ such that $$\begin{array}{}\text{(10.18)}& {\int }_a^{b}w(x)f(x)dx-\sum _{k=0}^n{W}_{k}f(x_k)=K_n{f}^{(2n+2)}(\eta ),\end{array}$$ and $$K_n=\frac{1}{(2n+2)!}{\int }_a^{b}w(x)[{\pi }_{n+1}(x){]}^{2}dx.$$ Consequently, the integration formula (10.6), (10.7) will give the exact result for every polynomial of degree $2n+1$.

"Introduction to Numerical Analysis," page. 282.

Theorem on Gaussian Formula with Error Term Consider a Gaussian formula with error term: $${\int }_a^{b}f(x)w(x)dx=\sum _{i=0}^{n-1}{A}_{i}f(x_i)+E$$ For an $f\in C^{2n}[a,b]$, we have $$E=\frac{f^{(2n)}(\xi )}{(2n)!}{\int }_a^b{q}^2(x)w(x)dx$$ where $a<\xi <b$ and $q(x)=\prod _{i=0}^{n-1}(x-x_i)$.

Proof. By Hermite interpolation (Section 6.3, p. 338), a polynomial $p$ of degree at most $2n-1$ exists such that $$p(x_i)=f(x_i),p^{\prime }(x_i)=f^{\prime }(x_i)(0\le i\le n-1)$$ The error formula for this interpolation is $$f(x)-p(x)=\frac{f^{(2n)}(\zeta (x))q^2(x)}{(2n)!}(11)$$ as given in Theorem 2 of Section 6.3 (p. 344). Here $q(x)$ is the polynomial of degree $n$ which is $w$-orthogonal to the polynomials of degree at most $n-1$. It follows that $${\int }_a^{b}f(x)w(x)dx-{\int }_a^{b}p(x)w(x)dx=\frac{1}{(2n)!}{\int }_a^b{f}^{(2n)}(\zeta (x))q^2(x)w(x)dx$$ Now use the fact that $p$ is of degree at most $2n-1$ to see that $${\int }_a^{b}p(x)w(x)dx=\sum _{i=0}^{n-1}{A}_{i}p(x_i)=\sum _{i=0}^{n-1}{A}_{i}f(x_i)$$ Furthermore, the Mean-Value Theorem for Integrals can be used to write $${\int }_a^b{f}^{(2n)}(\zeta (x))q^2(x)w(x)dx=f^{(2n)}(\xi ){\int }_a^b{q}^2(x)w(x)dx$$ This requires continuity of $f^{(2n)}(\zeta (x))$, which can be inferred from Equation (11). The desired error formula results from making the substitution described. $◼$

"Numerical Analysis: Mathematics of Scientific Computing" page. 498.

Theorem 8. Let $w(x)$ be a weight function and $p_n(x)$ be the monic polynomial of degree $n$ orthogonal to all polynomials of smaller degrees. Let $x_j$, $j=1,2,\dots ,n$, be zeros of $p_n(x)$. Let $Q(f)$ be the quadrature rule defined in Eqs. (36) and (37). Suppose $f\in C^{2n}[a,b]$. Then one can find $\lambda \in (a,b)$ such that $${\int }_a^{b}f(x)w(x)dx=Q(f)+{\gamma }_n\frac{f^{(2n)}(\lambda )}{(2n)!},\text{where}{\gamma }_n={\int }_a^b{p}_n^2(x)w(x)dx.(43)$$

Proof. We use the Hermite interpolation to prove the error estimate. There exists a unique polynomial $h_{2n-1}(x)$ of degree $2n-1$ such that $$h_{2n-1}(x_j)=f(x_j)\text{and}{h}_{2n-1}^{\prime }(x_j)=f^{\prime }(x_j).$$ In addition, there exists $\zeta \in (a,b)$ depending on $x$ such that $$f(x)=h_{2n-1}(x)+\frac{f^{(2n)}(\zeta )}{(2n)!}(x-x_1{)}^2\cdots (x-x_n{)}^2.$$ Note that $(x-x_1)\cdots (x-x_n)=p_n(x)$ as $p_n$ is monic and $x_j$, $j=1,\dots ,n$, are its roots. Therefore, $${\int }_a^{b}f(x)w(x)dx={\int }_a^b{h}_{2n-1}(x)w(x)dx+{\int }_a^b\frac{f^{(2n)}(\zeta (x))}{(2n)!}{p}_n^2(x)w(x)dx.$$ Since the quadrature is exact for polynomials of degree $2n-1$, it is exact for $h_{2n-1}(x)$. Hence $${\int }_a^b{h}_{2n-1}(x)w(x)dx=Q(h_{2n-1})=\sum _{j=1}^n{h}_{2n-1}(x_j)w_j=\sum _{j=1}^{n}f(x_j)w_j=Q(f).$$ On the other hand, because $p_n^2(x)w(x)\ge 0$ on $[a,b]$, we can apply the mean value theorem and get $${\int }_a^b\frac{f^{(2n)}(\zeta (x))}{(2n)!}{p}_n^2(x)w(x)dx=\frac{f^{(2n)}(\lambda )}{(2n)!}{\int }_a^b{p}_n^2(x)w(x)dx.$$ for some $\lambda \in (a,b)$. This completes the proof. $◻$

Composite Formulae

Definition(Composite trapezium rule). $${\int }_a^{b}f(x)dx\approx h[\frac{1}{2}f(x_0)+f(x_1)+\cdots +f(x_{m-1})+\frac{1}{2}f(x_m)].$$

Introduction to Numerical Analysis page. 209.

Definition (Composite Simpson rule). $${\int }_a^{b}f(x)dx\approx \frac{h}{3}[f(x_0)+4f(x_1)+2f(x_2)+4f(x_3)+\cdots +2f(x_{2m-2})+4f(x_{2m-1})+f(x_{2m})].$$

Introduction to Numerical Analysis page. 210.

Definition (Composite Midpoint Rule). The composite midpoint rule is the composite Gauss formula with $w(x)\equiv 1$ and $n=0$ defined by $${\int }_a^{b}f(x)dx\approx h\sum _{j=1}^{m}f(a+(j-\frac{1}{2})h).$$ This follows from the fact that when $n=0$ there is one quadrature point ${\xi }_0=0$ in $(-1,1)$, which is at the midpoint of the interval, and the corresponding quadrature weight $W_0$ is equal to the length of the interval $(-1,1)$, i.e., $W_0=2$. It follows from (10.24) with $n=0$ and $$C_0={\int }_{-1}^1{t}^{2}dt=\frac{2}{3}$$ that the error in the composite midpoint rule is $$E_{0,m}=\frac{(b-a{)}^3}{24m^2}{f}^{″}(\eta ),$$ where $\eta \in (a,b)$, provided that the function $f$ has a continuous second derivative on $[a,b]$.

"Introduction to Numerical Analysis", page. 286.

$\theta$-Method

Theta methods. We consider methods of the form $$y_{n+1}=y_n+h[\theta f(t_n,y_n)+(1-\theta )f(t_{n+1},y_{n+1})],n=0,1,\dots ,$$ where $\theta \in [0,1]$ is a parameter: - If $\theta =1$, we recover Euler's method. - If $\theta \in (0,1)$, then the theta method (3.3) is implicit: Each time step requires the solution of $N$ (in general, nonlinear) algebraic equations for the unknown vector $y_{n+1}$. - The choices $\theta =0$ and $\theta =\frac{1}{2}$ are known as Backward Euler: $$y_{n+1}=y_n+hf(t_{n+1},y_{n+1}),$$ Trapezoidal rule: $$y_{n+1}=y_n+\frac{1}{2}h[f(t_n,y_n)+f(t_{n+1},y_{n+1})].$$ Solution of nonlinear algebraic equations can be done by iteration. For example, for backward Euler, letting $y_{n+1}^{[0]}=y_n$, we may use $\textit{Direct iteration}$: $$y_{n+1}^{[j+1]}=y_n+hf(t_{n+1},y_{n+1}^{[j]});$$ $\textit{Newton–Raphson}$: $$y_{n+1}^{[j+1]}=y_{n+1}^{[j]}-{[I-h\frac{\partial f(t_n,y_n)}{\partial y}]}^{-1}[y_{n+1}^{[j]}-y_n-hf(t_{n+1},y_{n+1}^{[j]})];$$ $\textit{Modified Newton–Raphson}$: $$y_{n+1}^{[j+1]}=y_{n+1}^{[j]}-{[I-h\frac{\partial f(t_n,y_n)}{\partial y}]}^{-1}[y_{n+1}^{[j]}-y_n-hf(t_{n+1},y_{n+1}^{[j]})].$$

Error

In numerically solving a differential equation, several types of errors arise. These are conveniently classified as follows:

• Local truncation error
• Local roundoff error
• Global truncation error
• Global roundoff error
• Total error

Local Truncation Error. The local truncation error is the error made in a single step of a numerical method.

{Local Roundoff Error. The local roundoff error is the error made in a single step of a numerical method due to roundoff.

Global Truncation Error. The global truncation error is the accumulation of the local truncation errors in the previous steps.

Global Roundoff Error. The global roundoff error is the accumulation of the local roundoff errors in the previous steps.

Total Error. The total error is the sum of the global truncation error and the global roundoff error.

"Numerical Analysis: Mathematics of Scientific Computing" page. 533.

Richardson Extrapolation

$\mathbf{\text{Definition. }}$ Richardson extrapolation is a method for improving the accuracy of a numerical method. Given a numerical method that approximates a quantity with error $E(h)$ for a step size $h$, Richardson extrapolation constructs a new approximation with error $E(h^2)$ by combining two approximations with step sizes $h$ and $h/2$.

Initial value problems for ODEs

Our model is an initial-value problem written in the form $$\begin{array}{}\text{(1)}& \{\begin{array}{l}{x}^{\prime }=f(t,x)\\ x(t_0)=x_0\end{array}\end{array}$$ Here $x$ is an unknown function of $t$ that we hope to construct from the information given in Equations (1), where $x^{\prime }=\frac{dx(t)}{dt}$. The second of the two equations in (1) specifies one particular value of the function $x(t)$. The first equation gives the slope.

"Numerical Analysis: Mathematics of Scientific Computing" page.524.

Existence

$\textbf{First Existence Theorem, Initial-Value Problem}$ If $f$ is continuous in a rectangle $R$ centered at $(t_0,x_0)$, say $$\begin{array}{}\text{(4)}& R=\{(t,x):|t-t_0|\le \alpha ,|x-x_0|\le \beta \}\end{array}$$ then the initial-value problem above has a solution $x(t)$ for $|t-t_0|\le min(\alpha ,\beta /M)$, where $M$ is the maximum of $|f(t,x)|$ in the rectangle $R$.

"Numerical Analysis: Mathematics of Scientific Computing" page 525.

$\textbf{Second Existence Theorem, Initial-Value Problem.}$ If $f$ is continuous in the strip $a\le t\le b,-\mathrm{\infty }<x<\mathrm{\infty }$ and satisfies there an inequality $$\begin{array}{}\text{(5)}& |f(t,x_1)-f(t,x_2)|\le L|x_1-x_2|\end{array}$$ then the initial-value problem (1) has a unique solution in the interval $[a,b]$.

"Numerical Analysis: Mathematics of Scientific Computing" page.526.

Lipschitz Condition

Definition: A function $f:{\mathbb{R}}^n\to \mathbb{R}$ satisfies the Lipschitz condition on a domain $D\subset {\mathbb{R}}^n$ if there exists a constant $L\ge 0$ such that for all $x,y\in D$: $$\Vert f(x)-f(y)\Vert \le L\Vert x-y\Vert$$ where $\Vert \cdot \Vert$ denotes a norm (commonly the Euclidean norm). And $L$ is called the Lipschitz constant.

Note. Hence, inequality (5) is the Lipschitz condition.

Uniqueness

Uniqueness Theorem, Initial-Value Problem. If $f$ and $\partial f/\partial x$ are continuous in the rectangle $R$ defined by $\textbf{(4)}$, then the initial-value problem (1) has a unique solution in the interval $|t-t_0|<min(\alpha ,\beta /M)$.

$\textbf{Theorem (Picard's Theorem).}$ Suppose that the real-valued function $(x,y)\mapsto f(x,y)$ is continuous in the rectangular region $D$ defined by $x_0\le x\le X_M$, $y_0-C\le y\le y_0+C$; that $|f(x,y_0)|\le K$ when $x_0\le x\le X_M$; and that $f$ satisfies the Lipschitz condition: there exists $L>0$ such that $$|f(x,u)-f(x,v)|\le L|u-v|\text{for all}(x,u)\in D,(x,v)\in D.$$ Assume further that $$\begin{array}{}\text{(12.3)}& C\ge \frac{K}{L}(e^{L(X_M-x_0)}-1).\end{array}$$ Then, there exists a unique function $y\in C^1[x_0,X_M]$ such that $y(x_0)=y_0$ and $y^{\prime }=f(x,y)$ for $x\in [x_0,X_M]$; moreover, $$|y(x)-y_0|\le C,x_0\le x\le X_M.$$

"Introduction to Numerical Analysis", page. 311.

One-step methods

Definition. Generally, a one-step method may be written in the form $$\begin{array}{}\text{(12.13)}& y_{n+1}=y_n+h\mathrm{\Phi }(x_n,y_n;h),n=0,1,\dots ,N-1,y(x_0)=y_0,\end{array}$$ where $\mathrm{\Phi }(\cdot ,\cdot ;\cdot )$ is a continuous function of its variables.