Exercise 1

$4.2-1$

Use a divided difference table to show that the following data can be represented by a polynomial of degree 3:

x -2 -1 0 1 2 3
y 1 4 11 16 13 -4
In [1]:
#computes the set of divided differences
def DividedDifferences(points):
    #we compute n from the size of the array
    n = len(points)-1
    
    #create an array n+1 by n+1
    divDiff = [[0 for i in range(n+1)] for i in range(n+1)]
    
    #first column will contain the y values
    #this is column 0!
    for row in [0,1,2,..,n]:
        divDiff[row][0] = points[row][1]

    #compute, recursively the div differences
    #one column at a time
    #in second column will have f[.]
    #in third column will have f[.,.], etc
    for column in [1,2,..,n]:
        #for a column, compute the value in each row
        for row in [0,1, 2, .., n-column]:
            divDiff[row][column] = ((divDiff[row+1][column-1] - divDiff[row][column-1])
                / (points[row+column][0] - points[row][0]))
    return divDiff

def printDividedDifferences(points):
    #we compute n from the size of the array
    n = len(points)-1
    divDiff = DividedDifferences(points)
    for row in [0,1, 2, .., n]:
        print points[row][0],"\t|", [divDiff[row][column] for column in [0,1,2,..,n-row]]


#compute the Newton interpolation polynomial
def Newton_Interpolation(x, points):
    #we compute n from the size of the array
    n = len(points) - 1   
    #the very first step is computed here
    #p0(x)=y0
    p(x) = points[0][1]
    #new term will be x-x0
    #it will eventually store (xk-x0)(xk-x1)...(xk-xk-1)
    newTerm(x) = x - points[0][0]
    for k in range(1,n+1):
        #we repeatedly add the new data point
        #first compute the new c
        c = ((points[k][1] - p(points[k][0])) / newTerm(points[k][0]))
        #and then add it to the interpolating function
        p(x) = p(x) + c*newTerm(x)
        #prepare for next iteration
        newTerm(x) = newTerm(x) * (x-points[k][0])
    return p
In [2]:
points = [[-2, 1],[-1, 4],[0, 11],[1,16],[2,13],[3,-4]]
printDividedDifferences(points)
print
print "Create Newton interpolation polynomial and check results"
print
#call the function and find the final answer
f(x) = Newton_Interpolation(x, points)
print "Interpolating polynomial f(x) = ", f(x)
print "f([points]) = ", [f(pt[0]) for pt in points]
print
-2 	| [1, 3, 2, -1, 0, 0]
-1 	| [4, 7, -1, -1, 0]
0 	| [11, 5, -4, -1]
1 	| [16, -3, -7]
2 	| [13, -17]
3 	| [-4]

Create Newton interpolation polynomial and check results

Interpolating polynomial f(x) =  -(x + 2)*(x + 1)*x + 2*(x + 2)*(x + 1) + 3*x + 7
f([points]) =  [1, 4, 11, 16, 13, -4]

Analysis

We get a cubic interpolating polynomial from the difference table.
This divided difference table when augmented with two additional data points, shows 0 coefficients for higher power terms
and so the interpolating polynomial remains unchanged as a cubic.

Exercise 2

$4.2-6$

Linear interpolation in a table of function values means the following: if $y_0 = f(x_0)$ and $y_1 = f(x_1)$ are tabulated values, and if $x_0 < x < x_1$, then an interpolated value of $f(x)$ is \begin{equation*} y_0 + \frac{y_1 − y_0}{ x_1 − x_0} (x − x_0). \end{equation*}

How accurately can we determine $\sin x$ by linear interpolation, given a table of $\sin x$ to 10 decimal places for $x \in [0, 2]$, with $h = 0.01$?

In [3]:
n = 1
M = 1
h = 0.01

1/(4*(n+1))*M*h^(n+1)
Out[3]:
0.0000125000000000000

Analysis

Let $f$ be a function such that $f^{(n+1)}$ is continuous on $[a, b]$ and satisfies $|f^{(n+1)}(x)| \leq M$. Let $p$ be the polynomial of degree $\leq n$ that interpolates $f$ at $n + 1$ equally spaced nodes in $[a, b]$, including the endpoints. Then on $[a, b]$, \begin{align} |f(x) − p(x)| \leq \frac{1}{4(n + 1)} Mh^{n+1} \end{align} where h = (b − a)/n is the spacing between nodes.

Linear interpolation is used for each segment therefore we have $n+1=2$ points and with $h=0.01$ interval width the accuracy is determined by \begin{align*} \frac{1}{4(n+1)} M h^{(n+1)} = \frac{1}{8}1\cdot0.01^2 = 1.25\times10^{-5}\ \end{align*}

Setting $M = 1$ is justified by trigonometric analysis since all derivatives of $\sin$ or $\cos$ are themselves $\sin$ or $\cos$ trig functions and hence bounded by 1.
The maximum accuracy of $\sin x$ given by this linear interpolation is $1.25\times10^{-5}$.

Exercise 3

$4.2-8$

A table of values of $\cos x$ is required so that the linear interpolation will yield five decimal places accuracy for any value of $x$ in $[0, \pi]$. Assume that the tabular values are equally spaced, and determine the minimum number of entries needed in this table.

In [4]:
n = 1
M = 1
h = sqrt(8 * .5 * 10^-5.)
N = ceil(pi.n(100) / h); N
Out[4]:
497

Analysis

This analysis is based on the same error estimation formula used in Exercise 2.
Five decimal places of accuracy means the error must be $\leq \frac{1}{2} \times 10^{-5}$.
Linear interpolation is used for each segment so we have $n+1=2$ points.
Therefore,

\begin{align*} \frac{1}{4(n+1)} M h^{(n+1)} &= 0.5 * 10^{-5} \\ \frac{1}{8}\cdot1\cdot h^2 &= 0.5 * 10^{-5} \\ h &= \sqrt{8 * 0.5 * 10^{-5.}} \\ h &= 0.00632455532033676 \\ N &= \lceil\pi / h\rceil = 497 \end{align*}

The minimum number of entries needed is 497.

Exercise 4

$4.3-2$

Using Taylor series, establish the error term for the formula \begin{equation*} f'(0) \approx \frac{1}{2h}[f(2h) − f(0)]. \end{equation*}

\begin{align*} &\left(1\right)&\quad f(x+h) &= f(x) + hf'(x) + \frac{1}{2}h^2f''(\xi) && \text { use first order Taylor formula } \\ &\left(2\right)&\quad f(x+2h) &= f(x) + 2hf'(x) + \frac{1}{2}(2h)^2f''(\xi) && \text { substitute $2h$ for $h$ } \\ &&2hf'(x) &= f(x+2h) - f(x) - 2h^2f''(\xi) && \text { move $f'(x)$ to LHS } \\ &&f'(x) &= \frac{1}{2h}\left[f(x+2h) - f(x)\right] - hf''(\xi) && \text { cancel 2h from both sides } \\ &&f'(0) &= \frac{1}{2h}\left[f(2h) - f(0)\right] - hf''(\xi) && \text { evaluate $f'(x)$ at $x = 0$ }\\ &&f'(0) &\approx \frac{1}{2h}\left[f(2h) - f(0)\right] \text{ with an error term } = -hf''(\xi) && \text { approximate $f'(0)$ with error term }\\ \end{align*}

Exercise 5

$4.3-3$

Derive the approximation formula

\begin{align*} f'(x) \approx \frac{1}{2h}[4f(x + h) − 3f(x) − f(x + 2h)]. \end{align*}

and show that its error term is of the form \begin{align*} \frac{1}{3}h^2f'''(\xi). \end{align*}

\begin{align*} &\left(1\right)&\quad f(x+h) &= f(x) + hf'(x) + \frac{1}{2}h^2f''(x) + \frac{1}{6}h^3f'''(\xi) && \text { use second order Taylor formula }\\ &\left(2\right)&\quad f(x+2h) &= f(x) + 2hf'(x) + \frac{1}{2}(2h)^2f''(x) + \frac{1}{6}(2h)^3f'''(\xi) &&\text { substitute $2h$ for $h$ }\\ &&4f(x+h) &= 4f(x) + 4hf'(x) + 2h^2f''(x) + \frac{4}{6}h^3f'''(\xi) && \text { multiply (1) by 4 }\\ &&f(x+2h) &= f(x) + 2hf'(x) + 2h^2f''(x) + \frac{8}{6}h^3f'''(\xi) && \text { simplify (2) }\\ &&4f(x+h) - f(x+2h) &= 3f(x) + 2hf'(x) - \frac{4}{6}h^3f'''(\xi) && \text { (1) - (2) }\\ &&2hf'(x) &= 4f(x+h) - f(x+2h) - 3f(x) + \frac{2}{3}h^3f'''(\xi) && \text { rearrange with $f'$ on the LHS }\\ &&f'(x) &\approx \frac{1}{2h}[4f(x + h) − 3f(x) − f(x + 2h)] \text{ with an error term } = \frac{1}{3}h^2f'''(\xi) && \text { approximate $f'(0)$ with error term }\\ \end{align*}

Exercise 6

$4.3-5$

Averaging the forward-difference formula $f'(x) ≈ [f(x + h) − f(x)]/h$ and the backward difference formula $f'(x) ≈ [f(x) − f(x − h)]/h$, each with error term $O(h)$, results in the central-difference formula $f'(x) ≈ [f(x + h) − f(x − h)]/(2h)$, with error term $O(h^2)$. Show why Hint: Determine at least the first term in the error seriesfor each formula.

\begin{align*} &\left(1\right)&\quad f(x+h) &= f(x) + hf'(x) + \frac{1}{2}h^2f''(\xi_1) && \text { use first order Taylor formula }\\ &\left(2\right)&\quad f(x-h) &= f(x) - hf'(x) + \frac{1}{2}h^2f''(\xi_2) &&\text { substitute $2h$ for $h$ }\\ &&-f(x-h) &= -f(x) + hf'(x) - \frac{1}{2}h^2f''(\xi_2) &&\text { multiply (2) by $-1$ }\\ &&f(x+h) - f(x-h) &= 2hf'(x) + \frac{1}{2}h^2f''(\xi_1) - \frac{1}{2}h^2f''(\xi_2)&& \text { (1) - (2) }\\ &&2hf'(x) &= f(x+h) - f(x-h) - \frac{1}{2}h^2f''(\xi_1) + \frac{1}{2}h^2f''(\xi_2)&& \text { rearrange with $f'$ on the LHS }\\ &&f'(x) &\approx [f(x + h) − f(x − h)]/(2h) && \text{ approximate $f'(x)$ with error term $O(h^2)$ }. \end{align*}

Exercise 7

$4.3-6$

Criticize the following analysis. By Taylor’s formula, we have

\begin{align} f(x + h) − f(x) &= hf'(x) + \frac{h^2}{2}f''(x) + \frac{h^3}{6}f'''(\xi)\\ f(x - h) − f(x) &=-hf'(x) + \frac{h^2}{2}f''(x) - \frac{h^3}{6}f'''(\xi)\\ \end{align}

So by adding, it seems that we obtain an $exact(!)$ expression for $f''(x)$: $$ f(x + h) + f(x − h) − 2f(x) = h^2 f''(x). $$ Is this correct?

No. The formulas are incomplete. The correct formulas are: \begin{align} f(x + h) − f(x) &= hf'(x) + \frac{h^2}{2}f''(x) + \frac{h^3}{6}f'''(\xi_1) && \text{ where } \xi_1 \text { is between $x$ and $x+h$ } \ f(x - h) − f(x) &=-hf'(x) + \frac{h^2}{2}f''(x) - \frac{h^3}{6}f'''(\xi_2) && \text{ where } \xi_2 \text { is between $x-h$ and $h$ }\ f''(x)&\approx \left[f(x + h) + f(x − h) − 2f(x)\right] / h^2 && \text { $\mathbf { approximate }$ $f'(x)$ with error term $O(h^3)$. } \ \end{align}