Contents

Stability of polynomial interpolation

9.3. Stability of polynomial interpolation

With barycentric interpolation available in the form of Function 9.2.3, we can explore polynomial interpolation using a numerically stable algorithm. Any remaining sensitivity to error is due to the conditioning of the interpolation process itself.

Demo 9.3.1

We choose a function over the interval \([0,1]\).

f = x -> sin(exp(2*x));

Here is a graph of \(f\) and its polynomial interpolant using seven equally spaced nodes.

plot(f,0,1,label="function",legend=:bottomleft)
t = (0:6)/6
y = f.(t)
scatter!(t,y,label="nodes")

p = FNC.polyinterp(t,y)
plot!(p,0,1,label="interpolant",title="Equispaced interpolant, n=6")
Equispaced interpolant, n=6

This looks pretty good. We want to track the behavior of the error as \(n\) increases. We will estimate the error in the continuous interpolant by sampling it at a large number of points and taking the max-norm.

n = 5:5:60;  err = zeros(size(n))
x = range(0,1,length=2001)      # for measuring error
for (i,n) in enumerate(n) 
  t = (0:n)/n                   # equally spaced nodes
  y = f.(t)                     # interpolation data
  p = FNC.polyinterp(t,y)
  err[i] = norm( (@. f(x)-p(x)), Inf )
end
plot(n,err,m=:o,title="Interpolation error for equispaced nodes",
    xaxis=(L"n"),yaxis=(:log10,"max error"),)
max error Interpolation error for equispaced nodes

The error initially decreases as one would expect but then begins to grow. Both phases occur at rates that are exponential in \(n\), i.e., \(O(K^n\)) for a constant \(K\), appearing linear on a semi-log plot.

Runge phenomenon

The disappointing loss of convergence in Demo 9.3.1 is a sign of ill conditioning due to the use of equally spaced nodes. We will examine this effect using the error formula (9.1.5) as a guide:

\[ f(x) - p(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!} \Phi(x), \qquad \Phi(x) = \prod_{i=0}^n (x-t_i). \]

While the dependence on \(f\) is messy here, the error indicator \(\Phi(x)\) can be studied as a function of the nodes only.

Demo 9.3.2

We plot \(|\Phi(x)|\) over the interval \([-1,1]\) with equispaced nodes for different values of \(n\).

plot(xaxis=(L"x"),yaxis=(:log10,L"|\Phi(x)|",[1e-25,1]),legend=:bottomleft)

x = range(-1,1,length=2001)
for n in 10:10:50
    t = range(-1,1,length=n+1)
    Φ = [ prod(xₖ.-t) for xₖ in x ]
    scatter!(x,abs.(Φ),m=(1,stroke(0)),label="n=$n")
end
title!("Error indicator for equispaced nodes")
Error indicator for equispaced nodes