Polynomial Interpolation

It is often needed to estimate the value of a function at certan point based on the known values of the function at a set of node points in the interval . This process is called interpolation if or extrapolation if either or . One way to carry out these operations is to approximate the function by an nth degree polynomial:

 (1)

where the coefficients can be obtained based on the given points. Once is available, any operation applied to , such as differentiation, intergration, and root finding, can be carried out approximately based on . This is particulary useful if is non-elementary and therefore difficult to manipulate, or it is only available as a set of discrete samples without a closed-form expression.

Specifically, to find the coefficients of , we require it to pass through all node points , i.e., the following linear equations hold:

 (2)

Now the coefficients can be found by solving these linear equations, which can be expressed in matrix form as:

 (3)

where , , and

 (4)

is known as the Vandermonde matrix. Solving this linear equation system, we get the coefficients . Here polynomials can be considered as a set of polynomial basis functions that span the space of all nth degree polynomials (which can also be spanned by any other possible bases). If the node points are distinct, i.e., has a full rank and its inverse exists, then the solution of the system is unique and so is .

In practice, however, this method is not widely used for two reasons: (1) the high computational complexity for calculating , and (2) matrix becomes more ill-conditioned as increases. Instead, other methods to be discussed in the following section are practically used.

The error of the polynomial interpolation is defined as

 (5)

which is non-zero in general, except if at which . In other words, the error function has zeros at , and can therefore be written as

 (6)

where is an unknown function and is a polynomial of degree defined as

 (7)

To find , we construct another function of variable with any as a parameter:

 (8)

which is zero when :

 (9)

We therefore see that has zeros at and . According to Rolle's theorem, which states that the derivative function of any differentiable function satisfying must have at least a zero at some point at which , has at least zeros each between two consecutive zeros of , has at least zeros, and has at least one zero at some :
 (10)

The last equation is due to the fact that and are respectively an nth and (n+1)th degree polynomials of . Solving the above we get

 (11)

Now the error function can be written as

 (12)

where is a point located anywhere between and dependending on . The error can be quantitatively measured by the 2-normal of :

 (13)

In particular, if is a polynomial of degree , then , and it can be exactly interpolated by . But if , the interpolation has a non-zero error term . In particular, if is a polynomial of degree , such as , then and the error term becomes:

 (14)

Due to the uniqueness of the polynomial interpolation, the error analysis above also applies to all other methods to be considered in the following sections, such as the Lagrange and Newton interpolations.

Example:

Approximate function by a polynomial of degree , based on the following node points:

We first find the Vandermonde matrix:

and get the coefficients:

and then the interpolating polynomial can be found as a weighted sum of the first power functions used as the basis functions to span the polynomial space:

This interpolation polynomial is plotted in the figure below, in comparison to the orginal function , together with the basis polynomials, the power functions . The error can be approximated by a set of discrete samples of the function and the interpolating polynomial :

The Matlab code that implements this method is listed below.

function [v P]=PI(u,x,y)
% vectors x and y contain n+1 points and the corresponding function values
% vector u contains all discrete samples of the continuous argument of f(x)
n=length(x);     % number of interpolating points
k=length(u);     % number of discrete sample points
v=zeros(1,k);    % polynomial interpolation
P=zeros(n,k);    % power function basis polynomials
V=zeros(n);      % Vandermonde matrix
for i=1:n
for j=1:n
V(i,j)=x(i)^(j-1);
end
end
c=inv(V)*y;      % coefficients
for i=1:n
P(i,:)=u.^(i-1);
v=v+c(i)*P(i,:);
end
end