MTH 451/551
Lecture 1

0 What is numerical analysis?

This is a course on numerical analysis (specifically numerical calculus, which focuses specifically on calculus topics such as integration and differentiation). Numerical analysis is concerned with finding (or approximating) solutions to problems numerically, rather than symbolically. To give a simple example that illustrates this distinction, we could ask:

What are the roots of this polynomial?

$$x^2+x-1=0$$

We could solve this problem “symbolically” or “analytically” using standard techniques such as the quadratic formula or completing the square. We would find that the solution is

$$x=\frac{-1\pm \sqrt{5}}{2}.$$ (1)

We could then use, for example, a calculator to see that the roots are approximately

$$x_1\approx -1.61803\mathrm{\dots },x_2\approx 0.618034\mathrm{\dots }.$$ (2)

Numerical analysis takes a different approach to questions like these. There are many questions in mathematics, physics, and engineering that arise quite naturally, and that do not admit a closed-form solution. This means that we could not write a formula for the exact roots as in (1). In fact, it is a famous theorem that a general formula (in terms of basic arithmetic operations and roots) for the roots of polynomials of degree five or higher does not exist. So, instead of searching for the exact solution, we try to find the approximate solution (2) directly. It turns out that there are some very good computer algorithms for finding these approximate solutions for very many mathematics problems, and that these algorithms give us solutions to many problems for which finding a symbolic solution would be impossible. Numerical analysis transforms large classes of problems from “come up with a clever substitution or apply complicated transformations” to “run a deterministic computer algorithm”. The utility of such methods is enormous in many fields of science and engineering.

Numerical analysis is about developing and understanding computer algorithms to numerically solve mathematics problems of relevance to science and engineering.

1 Polynomial Approximation and Taylor’s Theorem

[Textbook reference: §1.1]

As described above, the goal of numerical calculus is to numerically approximate the solution to calculus problems, such as computing derivatives or integrals of given functions. One of the most useful tools to this end will be polynomial approximation. We know how to differentiate and integrate polynomials (just apply the power rule). So, if, given a function $f(x)$, we can find a polynomial $p(x)$ that is very close to $f(x)$ (in some sense of “close”), then differentiating or integrating $p(x)$ could give us a good approximation of the derivatives or integrals of $f(x)$. Here, differentiation and integration are essentially stand-ins for many other operations that we might want to perform on $f(x)$, and which we know how to do for polynomials like $p(x)$. Recall that a polynomial is a function $p:ℝ\to ℝ$ given by

$$p(x)=a_0+a_{1}x+a_2{x}^2+\mathrm{\cdots }+a_n{x}^n=\sum _{k=0}^n{a}_k{x}^k.$$

The largest number $n$ such that $a_n\ne 0$ is called the degree of the polynomial.

It turns out that (under some conditions), functions (and their derivatives and integrals) can be approximated very well by polynomials. The key result is called Taylor’s Theorem, which we will state shortly. We consider a real-valued function $f:[a,b]\to ℝ$. Here, $[a,b]$ is the domain of $f$, which just means that it’s the interval on which the function is defined. We allow $a=-\mathrm{\infty }$ and $b=\mathrm{\infty }$, so depending on the function, $f$ may be defined on all of $ℝ$. We are interested in continuous functions. A function is continuous at $x_0$ if

$$\underset{x\to x_0}{lim}f(x)=f(x_0).$$

A function is continuous on $[a,b]$ if it is continuous at $x_0$ for all $x_0\in [a,b]$. This means that $f(x)$ doesn’t have any “jumps” or “gaps”. Its graph can be drawn without lifting the pencil off of the page. We also assume that $f$ is differentiable, meaning that its derivative exists. We use the notation

$$f^{\prime }(x)=\frac{df}{dx}(x).$$

If $f^{\prime }(x)$ is also differentiable, then we say $f$ is twice differentiable, and write

$$f^{\prime \prime }(x)=f^{(2)}(x)=\frac{d^{2}f}{d{x}^2}(x).$$

Generalizing, $f$ may be $n$ times differentiable, with $n$th derivative $f^{(n)}(x)=\frac{d^{n}f}{d{x}^n}(x)$. We usually want $f$ to be $n$ times continuously differentiable, meaning that, not only does $f^{(n)}$ exist, but that $f^{(n)}$ is also continuous. For example, we would say that $f$ is continuously differentiable if $f^{\prime }(x)$ exists and is continuous. To express that $f$ is $n$ times continuously differentiable, we write

$$f\in C^n[a,b].$$

We are now ready to state Taylor’s Theorem.

Theorem 1 (Taylor’s Theorem).

Suppose that $f$ is $n$ times differentiable at $c$. There is a unique polynomial $p$ of degree $\le n$ for which $p^{(k)}(c)=f^{(k)}(c)$ for $0\le k\le n$, and it can be expressed as

$$p(x)=\sum _{k=0}^n\frac{f^{(k)}(c)}{k!}{(x-c)}^k.$$ (3)

Now suppose that $f\in C^{n+1}[a,b]$ and $c\in [a,b]$. Then for each $x\in [a,b]$ there is a $\xi =\xi (x)$ between $x$ and $c$ such that

$$f(x)-p(x)=\frac{f^{(n+1)}(\xi (x))}{(n+1)!}{(x-c)}^{n+1}\text{for all }x\in [a,b].$$ (4)

Furthermore, $\xi (x)$ may be chosen so that $f^{(n+1)}(\xi (x))$ is continuous on $[a,b]$.

The polynomial $p$ is called the Taylor polynomial of order $n$ for $f$ at $c$. The degree of $p$ is at most $n$ (but it may sometimes be less than $n$, so “degree” and “order” may be different).

Example 1 (Approximation by Taylor polynomials).

Let $f(x)=\ln (1+x)$, which is infinitely differentiable on $(-1,\mathrm{\infty })$. The derivatives of $f$ can be computed:

$$f^{(n)}(x)=\frac{{(-1)}^{n-1}(n-1)!}{{(1+x)}^n}\text{for }n\ge 1.$$

The Taylor polynomials of order $n$ for $f$ at $c=1$, for $n=0$, $n=1$, and $n=2$, are

$$p_0(x)=\ln (2),p_1(x)=\ln (2)+\frac{x-1}{2},p_2(x)=\ln (2)+\frac{x-1}{2}-\frac{{(x-1)}^2}{8}.$$

Using (4) with $n=2$, we see that, for example,

$$f(x)-p_2(x)=\frac{2/{(1+\xi )}^3}{6}{(x-1)}^3$$

for some $\xi =\xi (x)$ between $x$ and $c=1$. This gives us an expression for the error. It tells us how far the polynomial approximation $p_2(x)$ can be from $f(x)$. Graphs of these polynomials over $[0,2]$ are pictured below, together with $f$, as well as the errors $f-p_n$ and the absolute errors $|f-p_n|$.

The plots above show Taylor polynomials of order 0, 1, and 2, centered at $c=1$, for the function $f(x)=\ln (1+x)$ (left), as well as errors (center) and absolute errors (right). $\mathrm{◆}$

In the next example, we will use the following theorem.

Theorem 2 (Integral Mean Value Theorem).

Let $f,g\in C(a,b)$, with $f$ bounded on $(a,b)$ and . If $g$ does not change sign on $(a,b)$, then there is an $\eta \in (a,b)$ such that

$${\int }_a^{b}f(x)g(x)𝑑x=f(\eta ){\int }_a^{b}g(x)𝑑x.$$

Example 2 (The midpoint rule for approximating ${\int }_a^{b}f(x)𝑑x$).

Suppose that $f\in C^2[a,b]$ and $c=(a+b)/2$. Let $p(x)=f(c)+f^{\prime }(c)(x-c)$. This is the first-order Taylor polynomial of $f$ at the midpoint of $[a,b]$. We call this function the linearization of $f$ at $c$, since $p(x)$ is a linear approximation to $f(x)$. Using (4), we obtain

$${\int }_a^{b}f(x)𝑑x={\int }_a^{b}p(x)𝑑x+\frac{1}{2}{\int }_a^b{f}^{\prime \prime }(\xi (x)){(x-c)}^2𝑑x.$$

We can compute

$${\int }_a^{b}p(x)𝑑x={\int }_a^{b}f(c)𝑑x+{\int }_a^b{f}^{\prime }(c)(x-c)𝑑x.$$

The first term is given by

$${\int }_a^{b}f(c)𝑑x=f(c)(b-a).$$

The second term is

${\displaystyle {\int }_a^b}{f}^{\prime }(c)(x-c)𝑑x$ $={f^{\prime }(c){\displaystyle \frac{x^2}{2}}|}_a^b-f^{\prime }(c)c(b-a)$

$=f^{\prime }(c){\displaystyle \frac{b^2-a^2}{2}}-f^{\prime }(c){\displaystyle \frac{b^2-a^2}{2}}=0.$

We could also see that this integral is zero since the function $x-c$ is antisymmetric (odd) about the midpoint of the interval $[a,b]$.

Using the Integral Mean Value Theorem (Theorem 2), we can compute

${\displaystyle \frac{1}{2}}{\displaystyle {\int }_a^b}{f}^{\prime \prime }(\xi (x)){(x-c)}^2𝑑x$ $={\displaystyle \frac{1}{2}}{f}^{\prime \prime }(z){\displaystyle {\int }_a^b}{(x-c)}^2𝑑x$

$={{\displaystyle \frac{1}{2}}{f}^{\prime \prime }(z){\displaystyle \frac{1}{3}}{(x-c)}^3|}_a^b$

$={\displaystyle \frac{1}{24}}{f}^{\prime \prime }(z){(b-a)}^3,$

for some $z=\xi (\eta )\in (a,b)$. This gives the approximation

$${\int }_a^{b}f(x)𝑑x=f(c)(b-a)+\frac{1}{24}{f}^{\prime \prime }(z){(b-a)}^3.$$

The quantity $f(c)(b-a)$ is the basic midpoint rule approximation of ${\int }_a^{b}f(x)𝑑x$. If the interval we’re integrating over is very small, then the error will also be small (as long as the second derivative of $f$ is not too large). $\mathrm{◆}$

The midpoint rule for approximating integrals is the first numerical analysis algorithm we have seen in this class. If we want to approximate the integral of $f(x)$ over the interval $[a,b]$, then the midpoint rule gives us the approximation ${\int }_a^{b}f(x)𝑑x\approx f(c)(b-a)$, which can be computed easily as long as we have some way to evaluate $f$. This is an example of numerical quadrature, which is the term used for methods to approximate integrals. We will discuss much more about quadrature later in this course.

In the above example, the error scaled like ${(b-a)}^3$. The width of the interval is $b-a$, and we considered what happens as the intervals become smaller. We often use the symbol $h$ to denote a small parameter, such as the interval width. That makes the following form of Taylor’s Theorem often convenient.

Theorem 3 (Taylor’s Theorem with $h$).

Suppose that $f\in C^{n+1}[c-\delta ,c+\delta ]$ for some $\delta >0$. For all $h$, $|h|\le \delta$, there is $z$ between $c$ and $c+h$ such that

$$f(c+h)=\sum _{k=0}^n\frac{f^{(k)}(c)}{k!}{h}^k+\frac{f^{(n+1)}(z)}{(n+1)!}{h}^{n+1}.$$

Example 3 (Difference approximations of $f^{\prime }(c)$).

Suppose that $f\in C^2[c-\delta ,c+\delta ]$ for some $\delta >0$, and let . Applying Taylor’s Theorem with $h$ and $n=1$, we find that there are $y_{+}\in (c,c+h)$ and $y_{-}\in (c-h,c)$ such that

$f(c+h)$ $=f(c)+f^{\prime }(c)h+{\displaystyle \frac{f^{\prime \prime }(y_{+})}{2}}{h}^2,$

$f(c-h)$ $=f(c)-f^{\prime }(c)h+{\displaystyle \frac{f^{\prime \prime }(y_{-})}{2}}{h}^2.$

Rearranging the first of these and dividing by $h$ gives the forward difference approximation of $f^{\prime }(c)$,

$$\frac{f(c+h)-f(c)}{h}=f^{\prime }(c)+\frac{f^{\prime \prime }(y_{+})}{2}h.$$

Similarly, rearranging the second and dividing by $h$ gives the backward difference approximation,

$$\frac{f(c)-f(c-h)}{h}=f^{\prime }(c)-\frac{f^{\prime \prime }(y_{-})}{2}h.$$

In both cases, the error scales like $\mathrm{h}$. This means that halving $h$ (the distance between the points) gives a result that is about twice as accurate (the error is roughly half).

Instead of picking one direction, we can use the symmetric centered difference approximation,

$$\frac{f(c+h)-f(c-h)}{2h}.$$

To study its error, we need a higher-order Taylor polynomial. Now suppose that $f\in C^3[c-\delta ,c+\delta ]$. Applying Taylor’s Theorem with $h$ and $n=2$, there are $z_{+}\in (c,c+h)$ and $z_{-}\in (c-h,c)$ such that

$f(c+h)$ $=f(c)+f^{\prime }(c)h+{\displaystyle \frac{f^{\prime \prime }(c)}{2}}{h}^2+{\displaystyle \frac{f^{\prime \prime \prime }(z_{+})}{6}}{h}^3,$ (5)

$f(c-h)$ $=f(c)-f^{\prime }(c)h+{\displaystyle \frac{f^{\prime \prime }(c)}{2}}{h}^2-{\displaystyle \frac{f^{\prime \prime \prime }(z_{-})}{6}}{h}^3.$ (6)

Subtracting (6) from (5), the $f(c)$ and $f^{\prime \prime }(c)$ terms cancel, and dividing by $2h$ gives

$$\frac{f(c+h)-f(c-h)}{2h}=f^{\prime }(c)+\frac{h^2}{6}\cdot \frac{f^{\prime \prime \prime }(z_{+})+f^{\prime \prime \prime }(z_{-})}{2}.$$

Now the error scales like ${\mathrm{h}}^2$. Halving $h$ reduces the error by about a factor of four. The centered difference uses the same number of evaluations of $f$ as the forward and backward differences, but it is much more accurate when $h$ is small.

This is the second example of a basic numerical method that we have seen in this course. If we have a function $f$ that we can evaluate, then we can approximate the derivative $f^{\prime }(c)$ using two evaluations of $f$. So far, we have seen a method to approximate integrals (quadrature; specifically, the midpoint method), and a method to approximate derivatives (difference quotients, often called finite differences). $\mathrm{◆}$

Finally, we introduce some notation. Suppose $y$ is a quantity of interest that we wish to approximate by $\widehat{y}$. We write $\widehat{y}\approx y$ to mean that $\widehat{y}$ approximates $y$. The signed error of the approximation is $y-\widehat{y}$ and the absolute error is $|y-\widehat{y}|$. The relative error is $|y-\widehat{y}|/|y|$ (note that this is not defined if $y=0$). The advantage of the relative error is that it is scale invariant, meaning that if both $y$ and $\widehat{y}$ are scaled by the same nonzero number, then the relative error does not change, as can be seen as follows. Let $Y=cy$, and $\widehat{Y}=c\widehat{y}$ for $c\ne 0$. Then,

$$\frac{|Y-\widehat{Y}|}{|Y|}=\frac{|cy-c\widehat{y}|}{|cy|}=\frac{|c||y-\widehat{y}|}{|c||y|}=\frac{|y-\widehat{y}|}{|y|}.$$