Math 441

Interpolation: Function Approximation

Evaluating Newton’s Polynomials

"""
    Horner(a::Array, z::Array)

Evaluates the polynomial p(z₁,z₂,…,zₙ) = a₁ + a₂z₁ + a₃z₁z₂ + ⋯ + aₙ₊₁z₁z₂⋯zₙ.

The length of `a` must be at least one more than the length of `z`!
"""
function Horner(a :: Array, z :: Array)
    n = length(z)
    p = a[n+1]
    for j in n:-1:1
        p = muladd(p, z[j], a[j])
    end
    p
end
struct NewtonPoly <: Function
    coeffs :: Array{Float64}
    nodes :: Array{Float64}
end
function (p::NewtonPoly)(x)
    Horner(p.coeffs, x .- p.nodes)
end

Generating Newton’s Polynomials

function divided_diffs(x::Array, y::Array)
    m = length(x)
    if length(y) != m
        throw(ArgumentError("x and y must have the same length."))
    end
    a = copy(y)
    for i in 2:m
        for k in m:-1:i
            a[k] = (a[k] - a[k-1])/(x[k] - x[k - i + 1])
        end
    end
    a
end
function Newtons_poly(x::Array, y::Array)
    a = divided_diffs(x, y)
    NewtonsPoly(a, x[1:end-1])
end

Function Interpolation

Given a function \(f\colon [a,b] \to \mathbb{R}\), we want to find a polynomial of degree at most \(n\) that will approximate \(f\).

Theorem (Weierstrass): Given \(f\in \mathcal{C}([a,b])\), then for every \(\varepsilon > 0\) there exists a polynomial \(p\) such that \[\left\lvert f(x) - p(x)\right\rvert < \varepsilon\] for every \(x \in [a,b]\).

Problem: we do not know what the degree of \(p\) is.

Idea: Pick \(n+1\) points \(x_0, x_1, \dots, x_n\) in the interval \([a,b]\), define \(y_i = f(x_i)\) for \(i = 0, 1, \dots, n\), and interpolate this data using a Newton’s polynomial \(p\).

Then the degree of \(p\) is at most \(n\), and \(p(x_i) = f(x_i)\) for each \(i = 0, 1, 2, \dots, n\).

Question: How close is \(p(x)\) to \(f(x)\) if \(x \in [a,b]\) is not one of the \(x_i\)’s?

Interpolation Error

Theorem: Let \(x_0, x_1, \dots, x_n\) be distinct numbers in the interval \([a,b]\), and let \(f\in \mathcal{C}^{n+1}([a,b])\). Let \(p\) be the unique polynomial of degree at most \(n\) that interpolates the data \(\{(x_i, f(x_i))\mid i = 0, 1, \dots, n\}\). Then for each \(x \in [a,b]\), there is a \(\xi \in (a,b)\) such that \[f(x) - p(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!}\prod_{i=0}^n (x - x_i).\]

Proof:

  • Define \(R(x) = f(x) - p(x)\).
  • Fix an \(x \in [a,b]\).
  • Define \(M(t) = \displaystyle\prod_{i=0}^n(t - x_i)\)
  • Define \(Y(t) = R(t) - \frac{R(x)}{M(x)}M(t)\)
  • \(M^{(n+1)}(t) = (n+1)!\)
  • \(R^{(n+1)}(t) = f^{(n+1)}(t)\)
  • \(Y^{(n+1)}(t) = f^{(n+1)}(t) - \frac{R(x)}{M(x)}(n+1)!\)
  • Suppose \(Y^{(n+1)}(\xi) = 0\) for some \(\xi \in (a,b)\):
  • \(f^{(n+1)}(\xi) - \frac{R(x)}{M(x)}(n+1)! = 0\), and \(R(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!}M(x)\).

Proof cont.

  • Recap:
    • \(f \in \mathcal{C}^{(n+1)}([a,b])\)
    • \(p\) is a polynomial of degree at most \(n\) such that \(p(x_i) = f(x_i)\) for \(i = 0, 1, \dots, n\).
    • \(x \in [a,b]\), and \(x \neq x_i\) for \(i = 0, 1, \dots, n\).
    • \(R(x) = f(x) - p(x)\) and \(Y(t) = R(t) - \frac{R(x)}{M(x)}M(t)\), where \(M(t) = \displaystyle\prod_{i=0}^n(t - x_i)\).
  • We need to show that \(Y^{(n+1)}(\xi) = 0\) for some \(\xi \in (a,b)\).
  • We know that \(Y(x_i) = 0\) for \(i = 0, 1, \dots, n\), and also \(Y(x) = 0\).
  • So \(Y\) has at least \(n+2\) distinct zeros on \([a,b]\).
  • Rolle: Between any two consecutive zeros of \(Y\), there is a zero of \(Y'\)!
  • So \(Y'\) has at least \(n+1\) distinct zeros on \((a,b)\).
  • So \(Y''\) has at least \(n\) distinct zeros on \((a,b)\).

Proof finish

  • So \(Y^{(k)}\) has at least \(n+2-k\) distinct zeros on \((a,b)\).
  • Finally, \(Y^{(n+1)}\) has at least \(n+2-(n+1) = 1\) zero on \((a,b)\).
  • Therefore, \(Y^{(n+1)}(\xi) = 0\) for some \(\xi \in (a,b)\).

We have this error estimate:

\[f(x) - p(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!}{\color{red}{\prod_{i=0}^n (x - x_i)}}.\]

How large can \(\displaystyle \prod_{i=0}^n (x - x_i)\) be?

Uniformly Spaced Nodes

Let \(h = \frac{b-a}{n}\), and define \(x_i = a + ih\) for \(i = 0, 1, \dots, n\). Then for any \(x \in [a,b]\),

\[\prod_{i=0}^n \left\lvert x - x_i\right\rvert \le \frac{1}{4}h^{n+1} n!\]

Proof: The inequality is true if \(x = x_j\) for some \(j\).

Otherwise, there must be \(j\) such that \(x_j < x < x_{j+1}\).

𝑥 0 𝑥 1 𝑥 2 𝑥 𝑗 𝑥 𝑗+1 𝑥 𝑛−2 𝑥 𝑛−1 𝑥 𝑛 ⋯ ⋯ 𝑥

We will use the following estimates:

  • If \(i < j\), then \(\left\lvert x - x_i\right\rvert < x_{j+1} - x_i = (j + 1 - i)h\).
  • If \(i > j + 1\), then \(\left\lvert x - x_i\right\rvert < x_i - x_{j} = (i-j)h\).
  • We still need to estimate \(\left\lvert x - x_j\right\rvert \left\lvert x - x_{j+1}\right\rvert = (x - x_j)(x_j + h - x)\).

Uniformly Spaced Nodes cont.

What is the maximum of \((x - x_j)(x_{j+1} - x)\) on \([x_j, x_{j+1}]\)?

\[\begin{aligned} \prod_{i=0}^n \left\lvert x - x_i\right\rvert &< \prod_{i=0}^{j-1}(j - i + 1)h \cdot \frac{h^2}{4}\cdot\prod_{i=j+2}^n(i - j)h\\ &= \frac{1}{4}h^{j + 2 + n - j - 1}\prod_{i=0}^{j-1}(j - i + 1)\cdot \prod_{i=j+2}^n(i - j)\\ &= \frac{1}{4}h^{n+1}\prod_{k=2}^{j+1}k\cdot\prod_{k=2}^{n-j}k\\ &= \frac{1}{4}h^{n+1}(j+1)!(n-j)!\\ &= \frac{1}{4}h^{n+1}{\color{green}1\cdot 2\cdot 3 \cdots (j+1)}{\color{blue}\cdot 2\cdot 3 \cdots(n-j)}\\ &\le \frac{1}{4}h^{n+1}{\color{green}1\cdot 2\cdot 3 \cdots (j+1)}{\color{blue}\cdot (j+2)\cdot(j + 3)\cdots(j + n-j)}\\ &= \frac{1}{4}h^{n+1}n! \end{aligned}\]

Uniformly Spaced Nodes: Conclusion

\[\begin{aligned} \left\lvert f(x) - p(x)\right\rvert &= \frac{\left\lvert f^{(n+1)}(\xi)\right\rvert}{(n+1)!}{\color{red}{\prod_{i=0}^n \left\lvert x - x_i\right\rvert}}\\ &\le \frac{\left\lvert f^{(n+1)}(\xi)\right\rvert}{(n+1)!}{\color{red}{\frac{1}{4}h^{n+1}n!}}\\ &= \frac{\left\lvert f^{(n+1)}(\xi)\right\rvert}{4(n+1)}\left(\frac{b-a}{n}\right)^{n+1} \end{aligned}\]