Why do equi-spaced points behave badly?

Question

Experiment description:

In Lagrange interpolation, the exact equation is sampled at $N$ points (polynomial order $N - 1$) and it is interpolated at 101 points. Here $N$ is varied from 2 to 64. Each time $L_1$, $L_2$ and $L_\infty$ error plots are prepared. It is seen that, when the function is sampled at equi-spaced points, the error drops initially (it happens till $N$ is less than about 15 or so) and then the error goes up with further increase in $N$.

Whereas, if the initial sampling is done at Legendre-Gauss (LG) points (roots of Legendre polynomials), or Legendre-Gauss-Lobatto (LGL) points (roots of Lobatto polynomials), the error drops to machine level and doesn't increase when $N$ is further increased.

My questions are,

What exactly happens in the case of equi-spaced points?

Why does increase in polynomial order cause the error to rise after a certain point?

Does this also mean that if I use equi-spaced points for WENO / ENO reconstruction (using Lagrange polynomials), then in the smooth region, I would get errors? (well, these are only hypothetical questions (for my understanding), it is really not reasonable to reconstruct polynomial of the order of 15 or higher for WENO scheme)

Additional details:

Function approximated:

$f(x) = \cos(\frac{\pi}{2}~x)$, $x \in [-1, 1]$

$x$ divided into $N$ equispaced (and later LG) points. The function is interpolated at 101 points each time.

Results:

1. a) Equi-spaced points (interpolation for $N = 65$):

1. b) Equi-spaced points (error plot, log scale):

2. a) LG points (Interpolation for $N = 65$):

3. b) LG points (error plot, log scale):

Answer 1

The problem with equispaced points is that the interpolation error polynomial, i.e.

$$ f(x) - P_n(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!} \prod_{i=0}^n (x - x_i),\quad \xi\in[x_0,x_n] $$

behaves differently for different sets of nodes $x_i$. In the case of equispaced points, this polynomial blows up at the edges.

If you use Gauss-Legendre points, the error polynomial is significantly better behaved, i.e. it doesn't blow up at the edges. If you use Chebyshev nodes, this polynomial equioscillates and the interpolation error is minimal.

Answer 2

This is a really interesting question, and there are a lot of possible explanations. If we are attempting to use a polynomial interpolation, then note that polynomial satisfy the following annoying inequality

Given a polynomial $P$ of degree not exceeding $N$ we have

$$ |P^{\prime}(x)| \leq \frac{N}{\sqrt{1-x^2}}\max _x |P(x) | $$

for every $x \in (-1,1)$. This is known as Bernstein's inequality, note the singularity in this inequality. This can be bounded by the Markov inequality

$$\max _x |P^{\prime} (x) | \leq N^2 \max _x |P(x) | $$

and note that this is sharp in the sense that Chebysehv polynomials make this an equation. So in other words we have the following combined bound.

$$ |P^{\prime}(x)| \leq \min \left(\frac{N}{\sqrt{1-x^2}},N^2\right)\max _x |P(x) |$$

What this means: Gradients of polynomials grow linearly in their order everywhere except in small neighborhoods of the interval boundaries. At the boundaries they grow more like $N^2$. It is no accident that stable interpolation nodes all have a $1/{N^2}$ clustering near boundaries. The clustering is necessary to control the gradients of the basis, whereas near the midpoint one can be a bit more relaxed.

It turns out however that this is not necessarily a polynomial phenomena, I suggest the following paper:

http://math.la.asu.edu/~platte/pub/prevised.pdf

It says loosely: If you have the same approximation power of the polynomial basis, then you can't use equally spaced points in a stable way.

Answer 3

It is not the equally spaced points that are the problem. It's the global support of the basis functions along with equally spaced points that is the problem. A perfectly well conditioned interpolant using equally spaced points is described in Kress's Numerical Analysis, using cubic-b spline basis functions of compact support.

Answer 4

What exactly happens in the case of equi-spaced points?

Why does increase in polynomial order cause the error to rise after a certain point?

This is similar to the Runge's phenomenon where, with equi-spaced nodes, the interpolation error goes to infinity with the increase of the polynomial degree, i.e. the number of points.

One of the roots of this problem can be found in the Lebesgue's constant as noted by @Subodh's comment to @Pedro answer. This constant relates the interpolation with the best approximation.

---

Some notations

We have a function $f \in C([a,b])$ to interpolate over the nodes $x_k$. In the Lagrange interpolation are defined the Lagrange polynomials:

$$ L_k(x) = \prod_{i=0, i\neq j}^{n}\frac{x-x_i}{x_k-x_i} $$

with this is defined the interpolation polynomial $p_n \in P_n$ over the couples $(x_k, f(x_k))$ for light notation $(x_k, f_k)$

$$ p_n(x) = \sum_{k=0}^n f_kL_k(x) $$

---

Now consider a perturbation over the data, this can be for example for rounding, so we have got $\tilde{f}_k$. With this the new polynomial $\tilde{p}_n$ is:

$$ \tilde{p}_n(x) = \sum_{k=0}^n \tilde{f}_k L_k(x) $$

The error estimates are:

$$ p_n(x) - \tilde{p}_n(x) = \sum_{k=0}^n (f_k - \tilde{f}_k) L_k(x) $$

$$ | p_n(x) - \tilde{p}_n(x) | \leq \sum_{k=0}^n |f_k - \tilde{f}_k| |L_k(x)| \leq \left ( \max_k |f_k - \tilde{f}_k| \right) \sum_{k=0}^n |L_k(x)| $$

Now it is possible define the Lebesgue's constant $\Lambda_n$ as:

$$ \Lambda_n = \max_{x \in [a,b]} \sum_{k=0}^n |L_k(x)| $$

With this the final estimates is:

$$ || p_n - \tilde{p}_n ||_{\infty} \leq \left ( \max_k |f_k - \tilde{f}_k| \right) \Lambda_n $$

(marginal note, we look only $\infty$ norm also because we are over a space of finite measure so $L^{\infty} \subseteq \dots \subseteq L^1 $)

From the above calculation we have got that $\Lambda_n$ is:

• independent from the date:
• depends only from the nodes distribution;
• an indicator of stability (the smaller it is, the better it is).

It is also the norm of the interpolation operator respect the $|| \cdot||_\infty$ norm.

Withe the follow theorem we con have got an estimate of the interpolation error with the Lebesgue's constant:

Let $f$ and $p_n$ as above we have $$ || f - p_n ||_{\infty} \leq (1 + \Lambda_n) d_n(f) $$ where $$ d_n(f) = \inf_{q_n \in P_n} || f - q_n ||_{\infty} $$ is the error by the best uniform approximation polynomial

I.e. if $\Lambda_n$ is small the error of the interpolation is not far from the error of the best uniform approximation and the theorem compares the interpolation error with the smallest possible error which is the error of best uniform approximation.

For this the behavior of the interpolation depends by the nodes distribution. There is a lower bounds about $\Lambda_n$ that given a node distribution exist a constant $c$ such that: $$ \Lambda_n \geq \frac{2}{\pi} \log(n) - c $$ so the constant grows, but how it grow is importan.

For equi-spaced nodes $$\Lambda_n \approx \frac{2^{n+1}}{en \log(n)} $$ I omitted some details, but we see that the grow is exponential.

For Chebyshev nodes $$\Lambda_n \leq \frac{2}{\pi} \log(n) + 4 $$ also here I omitted some details, there are more accurate and complicate estimate. See [1] for more details. Note that nodes of Chebyshev family have got logarithmic grow and from the previous estimates is near the best you can obtain.

For other nodes distributions see for example table 1 of this article.

---

There are a lot of reference on book about interpolation. On-line these slides are nice as resume.

Also this open article ([1])

A Numerical Seven Grids Interpolation Comparison of for polynomial on the Interval for various comparisons.

Answer 5

It's good to be aware of Floater-Hormann interpolants when you have to (or want to) work with equidistant points $\{x_i\}_{i=1\ldots n}$.

Given the integer $d$ with $0 \le d \le n$, let $p_i$ be the polynomial interpolant of $\{ x_i, \ldots x_{i+d} \}$. Then the FH interpolant of a function $f$ at $\{x_i\}_{i=1\ldots n}$ has the form

$$ r_n(x) := \frac{\sum_{i=0}^{n-d} \lambda_i(x) \, p_i(x)}{\sum_{i=0}^{n-d} \lambda_i(x)} $$

with the "blending functions"

$$ \lambda_i(x) = \frac{(-1)^i}{(x-x_i) \ldots (x - x_{i+d})} $$

Some properties of these interpolants:

• they are barycentric rational interpolants with no real poles;
• achieve arbitrary approximation orders ${\cal O}(h^{d+1})$ for $f \in C^{d+2}[a,b]$, regardless of the distribution of points;
• are somewhat similar to splines, in that they blend (local) polynomial interpolants $p_0, \ldots p_{n-d}$ with the $\lambda$'s acting as the blending functions;
• they reproduce polynomials of degree at most $d$ (or $d+1$ if $n-d$ is odd);
• can be written in barycentric form (see section 4 of Floater and Hormann's paper).

Caveat emptor: As expected (see the paper referenced by @Reid.Atcheson), increasing $d$ quickly degrades the conditioning of the approximation process.

There is some fairly recent work done by Klein to alleviate this problem. He modified the original Floater-Hormann approach by adding $2 d$ new data values corresponding to points outside the original interpolation interval $[a,b]$ constructed from a smooth extension of $f$ outside $[a,b]$ using only the given data $f_0, \ldots f_n$. This "global" data set is then interpolated by a new FH rational function $r_{n+2d}$ and evaluated only inside $[a,b]$.

The details are nicely laid out in Klein's paper (linked below), where it is shown that these extended rational interpolants have Lebesgue constants that grow logarithmically with $n$ and $d$ (whereas for the original FH scheme, said growth is exponential in $d$, see Bos et al.).

The Chebfun library uses FH interpolants when building chebfuns out of equispaced data, as explained here.

References:

M. S. Floater and K. Hormann, Barycentric rational interpolation with no poles and high rates of approximation, Numerische Mathematik 107 (2007).

G. Klein, An Extension of the Floater–Hormann Family of Barycentric Rational Interpolants, Mathematics of Computation, 82 (2011) - preprint

L. Bos, S. De Marchi, K. Hormann, and G. Klein, On the Lebesgue constant of barycentric rational interpolation at equidistant nodes, Numer. Math. 121 (2012)