Ch5.1: Finite Differences & Derivatives

Numerical Differentiation

The derivative of a function is defined as

\[ f'(x) = \lim_{h \rightarrow 0}\frac{f(x+h)-f(x)}{h} \]

title

Numerical Differentiation

This is not a good approximation because of cancellation error in numerator and division by (close to) zero in denominator.
However, it is a good place to start.

\[ f'(x) = \lim_{h \rightarrow 0}\frac{f(x+h)-f(x)}{h} \]

title

Two Jokes



Knock knock.

Who's there?

Dozen.

Dozen who?

Dozen anyone wanna let me in?



It's inappropriate to make a dad joke if you are not a dad.

It's a faux pa.

Finite Difference Formula

The finite difference formula is the natural numerical approximation to the derivative:

\[ f'(x) \cong \frac{f(x+h)-f(x)}{h} \]

title

Forward & Backward Differences

Forward difference
- \( h > 0 \)
Backward difference
- \( h < 0 \)

\[ \small{ \begin{aligned} f'(x) \cong \frac{f(x+h)-f(x)}{h} \\ f'(x) \cong \frac{f(x)-f(x-h)}{h} \end{aligned} } \]

title

Centered (Symmetric) Difference

This formula looks both ahead and behind:

\[ f'(x) \cong \frac{f(x+h)-f(x-h)}{2h} \]

title

Forward/Backward Method using R

R program for finite difference

findiff <- function(f,x,
               h = x*sqrt(.Machine$double.eps))
  {return((f(x + h)-f(x))/h)}

\[ f'(x) \cong \frac{f(x+h)-f(x)}{h} \]

Linear Example: Forward Difference

f<-function(x){4*x-2}

findiff(f,2,h=0.1)

[1] 4

findiff(f,2)

[1] 4

title

Quadratic Example: Forward Difference

f <- function(x) 
  {x^2 + 1}

title

findiff(f,2,h = 0.1)

[1] 4.1

findiff(f,2,h = 0.01)

[1] 4.01

findiff(f,2)

[1] 4

Quadratic Example: Backward Difference

f <- function(x) 
  {x^2 + 1}

title

findiff(f,2,h = -1)

[1] 3

findiff(f,2,h = -0.1)

[1] 3.9

findiff(f,2,h = -0.01)

[1] 3.99

Centered Difference using R

R program for symmetric difference

symdiff <- function(f,x,
             h = x*(.Machine$double.eps)^(1/3))
{return((f(x+h)-f(x-h))/(2*h))}

\[ f'(x) \cong \frac{f(x+h)-f(x-h)}{2h} \]

Linear Example: Forward Difference

f<-function(x){4*x-2}

symdiff(f,2,h=0.1)

[1] 4

symdiff(f,2)

[1] 4

title

Quadratic Example: Symmetric Difference

f <- function(x) 
  {x^2 + 1}

title

symdiff(f,2,h = 1)

[1] 4

symdiff(f,2,h = 0.1)

[1] 4

symdiff(f,2,h = 0.01)

[1] 4

Numerical Analysis & Optimal h

How do we find the right step size h?
- Small h leads to better accuracy
- Too small leads to cancellation error
- Too small leads to division by zero
Trial & error with h values and see what is best?
Use Calculus somehow to minimize error?

\[ \small{ \begin{aligned} f'(x) & \cong \frac{f(x+h)-f(x)}{h} \\ f'(x) &\cong \frac{f(x+h)-f(x-h)}{2h} \end{aligned} } \]

Trial and Error: Forward Difference

f <- function(x){x^2 + 1}; x <- 2

h <- c(3e-4,3e-8,3e-12,3e-16)
findiff(f,2,h)

[1] 4.000300 4.000000 3.999763 5.921189

(h <- x*sqrt(.Machine$double.eps))

[1] 2.980232e-08

findiff(f,x)

[1] 4

Trial and Error: Centered Difference

f <- function(x){x^2 + 1}; x <- 2

h <- c(1e-1,1e-5,1e-12,1e-16)
symdiff(f,x,h)

[1] 4.000000 4.000000 4.000356 0.000000

(h <- x*(.Machine$double.eps)^(1/3))

[1] 1.211091e-05

symdiff(f,x)

[1] 4

Trial and Error: Centered Difference

f <- function(x){sin(x)}; x <- 0.9

h <- c(5e-2,5e-6,5e-10,5e-14)
symdiff(f,x,h)

[1] 0.6213510 0.6216100 0.6216101 0.6217249

(h <- x*(.Machine$double.eps)^(1/3))

[1] 5.449909e-06

symdiff(f,x)

[1] 0.62161

Numerical Analysis & Optimal h

Rather than trial & error with h values to guess what is best, we will use calculus.

Taylor series can provide error bounds.
A separate error formula can be minimized using calculus.

\[ \small{ \begin{aligned} f'(x) & \cong \frac{f(x+h)-f(x)}{h} \\ f'(x) &\cong \frac{f(x+h)-f(x-h)}{2h} \end{aligned} } \]

Forward Difference & Taylor Polynomials

Taylor series expansion of \( f(x) \) about \( x \):

\[ \small{ f(x+h) = f(x) + f'(x)\left((x+h)-x)\right) + \frac{f''(x)}{2}\left((x+h)-x)\right)^2 + \ldots } \]

It follows that

\[ \small{ \begin{aligned} f(x+h) - f(x) & = f'(x)h + \frac{f''(x)}{2}h^2 + \ldots \\ f'(x) & = \frac{f(x+h) - f(x)}{h} - \frac{f''(x)}{2}h + \ldots \\ f'(x) & = \frac{f(x+h) - f(x)}{h} - \frac{f''(c(x))}{2}h, \,\, \,\, c(x) \in [x,x+h] \\ \end{aligned} } \]

Centered Difference & Taylor Polynomials

Taylor series expansion of \( f(x) \) about \( x \): \[ \small{ \begin{aligned} f(x+h) & = f(x) + f'(x)\left((x+h)-x)\right) + \frac{f''(x)}{2}\left((x+h)-x)\right)^2 + \ldots \\ &= f(x) + f'(x)h + \frac{f''(x)}{2}h^2 + \frac{f'''(x)}{6}h^3 + \ldots \\ f(x-h) & = f(x) - f'(x)h + \frac{f''(x)}{2}h^2 - \frac{f'''(x)}{6}h^3 + \ldots \\ \end{aligned} } \]
Then \[ \small{ \begin{aligned} f'(x) & = \frac{f(x+h) - f(x-h)}{2h} - \frac{f'''(x)}{6}h^3 + \ldots \\ f'(x) & = \frac{f(x+h) - f(x-h)}{2h} - \frac{h^2}{6}f'''(c(x)), \,\, \,\, c(x) \in [x-h,x+h] \\ \end{aligned} } \]

Forward & Centered Difference Results

Using Taylor polynomials, we found error terms that go with the forward and centered difference formulas:

\[ \small{ \begin{aligned} f'(x) & = \frac{f(x+h) - f(x)}{h} - \frac{h}{2}f''(c(x)), \,\,\, \, c(x) \in [x,x+h]\\ \\ f'(x) & = \frac{f(x+h) - f(x-h)}{2h} - \frac{h^2}{6}f'''(c(x)), \, \,\,\, c(x) \in [x-h,x+h] \end{aligned} } \]

Forward & Centered Difference Results

Truncation error for each is defined as

\[ \begin{aligned} E_T & = \frac{h}{2}f''(c(x)), \,\,\,\, c(x) \in [x,x+h] \\ E_T & = \frac{h^2}{6}f'''(c(x)), \,\,\,\, c(x) \in [x-h,x+h] \end{aligned} \]

Error bounds are found as follows:

\[ \begin{aligned} |E_T| & \leq \frac{h}{2}M, \,\, M = \max_{u \in [x,x+h]}|f''(u)|, \\ |E_T| & \leq \frac{h^2}{6}M, \,\, M = \max_{u \in [x-h,x+h]}|f'''(u)| \end{aligned} \]

Desmos Illustration

Overview of computation folders.
Desmos link

title

Desmos Illustration: Finite Differences

Forward and centered difference results for \( h=0.1 \).

title

Desmos Illustration: Forward Error

Actual error and error bound estimate.

title

Desmos Illustration: Centered Error

Actual error and error bound estimate.

title

Numerical Analysis & Optimal h

Recall: How do we find the right step size h?
- Small h leads to better accuracy
- Too small leads to cancellation error
- Too small leads to division by zero
Trial & error with h values and see what is best?
Use Calculus somehow to minimize error?

\[ \small{ \begin{aligned} f'(x) & = \frac{f(x+h) - f(x)}{h} - \frac{h}{2}f''(c(x)), \\ f'(x) & = \frac{f(x+h) - f(x-h)}{2h} - \frac{h^2}{6}f'''(c(x)) \end{aligned} } \]

Forward Difference Error

For the forward difference results,

\[ \small{ f'(x) = \frac{f(x+h) - f(x)}{h} - \frac{h}{2}f''(c), \,\, c(x) \in [x,x+h] } \]

Residual (error) term from truncating infinite Taylor polynomial:

\[ \small{ E_T = \frac{h}{2}f''(c(x)) } \]

Next, we look at roundoff error (see next slide).

Forward Difference Error

For roundoff error, denote the computed values by

\[ \begin{aligned} \tilde{f}(x+h) &= f(x+h) + e_1(x) \\ \tilde{f}(x) &= f(x) + e_2(x) \end{aligned} \]

Then the computed forward difference formula is

\[ \small{ \begin{aligned} \tilde{f}'(x) &= \frac{\tilde{f}(x+h) - \tilde{f}(x)}{h} \\ &= \frac{f(x+h) - f(x) }{h} + \frac{e_1(x) - e_2(x)}{h} \\ &= \frac{f(x+h) - f(x) }{h} + E_R \end{aligned} } \]

Total Error: Forward Difference

From previous slide,

\[ \begin{aligned} \tilde{f}'(x) = \frac{f(x+h) - f(x) }{h} + E_R \end{aligned} \]

Here, we have denoted roundoff error by \( E_R \), with

\[ E_R = \frac{e_1(x) - e_2(x)}{h} \]

Thus the total error \( E \) is given by

\[ E = E_R + E_T = \frac{e_1(x) - e_2(x)}{h} + \frac{h}{2}f''(c(x)) \]

Total Error: Forward Difference

Our total error \( E \) is given by

\[ E = \frac{e_1(x) - e_2(x)}{h} + \frac{h}{2}f''(c) \]

Observe that round off increases as \( h \) gets smaller, while truncation error decreases as \( h \) get smaller.
We will use calculus to determine the optimal value of \( h \).

Total Error: Forward Difference

The total error \( E \) is given by

\[ E = \frac{e_1(x) - e_2(x)}{h} + \frac{h}{2}f''(c) \]

First, we observe that

\[ \left|e_1 - e_2\right| \leq |e_1| + |e_2| \leq 2 \max \left\{ |e_1|, |e_2| \right\} = 2e \]

Also, let

\[ M = \max_{x \leq u \leq x+h} |f''(u)| \]

Total Error: Forward Difference

Total error \( E \) is given by

\[ \small{ E = \frac{e_1 - e_2}{h} + \frac{h}{2}f''(c) } \]

From the previous slide, \( E \) is a function of \( h \), with

\[ \small{ E(h) \leq \frac{2e}{h} + \frac{M}{2}h } \]

Worst case scenario is that

\[ \small{ E(h) = \frac{2e}{h} + \frac{M}{2}h} \]

title

Total Error: Forward Difference

From previous slide,

\[ \small{ E(h) = \frac{2e}{h} + \frac{M}{2}h } \]

Set \( E'(h) = 0 \):

\[ \small{ \begin{aligned} E'(h) & = -\frac{2e}{h^2} + \frac{M}{2} \\ 0 & = -\frac{2e}{h^2} + \frac{M}{2} \\ h &= \sqrt{\frac{4e}{M}} \end{aligned} } \]

title

Total Error: Forward Difference

Worst-case scenario:

\[ \small{ E(h) = \frac{2e}{h} + \frac{M}{2}h } \]

Minimum occurs at

\[ \small{ h = \sqrt{\frac{4e}{M}} } \]

Book value:

\[ \small{ h^* = x\sqrt{\epsilon} } \]

title

Total Error: Forward Difference

Optimal h values:

\[ \small{ h = \sqrt{\frac{4e}{M}} , \,\, h^* = x\sqrt{\epsilon} } \]

library(pracma)
c(eps(1), eps(8))

[1] 2.220446e-16 1.776357e-15

\[ \small{ \begin{aligned} M &= \max |f''(u) | \\ e_1 &= \tilde{f}(x+h) - f(x+h) \\ e_2 &= \tilde{f}(x) - f(x) \end{aligned} } \]

title

Optimal Step Sizes in Book

\[ \begin{array}{cllr} f'(x) & = \frac{f(x+h) - f(x)}{h} &\leftrightarrow h^* = x\sqrt{\epsilon} \\ f'(x) & = \frac{f(x+h) - f(x-h)}{2h} & \leftrightarrow h^* = x\sqrt[3]{\epsilon} \end{array} \]

title

Second Derivative & Taylor Polynomials

Recall Taylor series expansion of \( f(x) \) about \( x \):

\[ \small{ \begin{aligned} f(x+h) = f(x) + f'(x)h & + \frac{f''(x)}{2}h^2 + \frac{f'''(x)}{6}h^3 + \frac{f^{(4)}(x)}{24}h^4 + \ldots \\ f(x-h) = f(x) - f'(x)h & + \frac{f''(x)}{2}h^2 - \frac{f'''(x)}{6}h^3 + \frac{f^{(4)}(x)}{24}h^4 + \ldots \end{aligned} } \]

Add together and solve for \( f''(x) \):

\[ \small{ \begin{aligned} f(x+h) + f(x-h) & = 2f(x) + f''(x)h^2 + \frac{f^{(4)}(x)}{12}h^4 + \ldots \\ \\ f''(x) & = \frac{f(x+h) - 2f(x) + f(x-h)}{h^2} - \frac{h^2}{12}f^{(4)}(c(x)) \end{aligned} } \]