Numerical methods. Note 12
Home ] Interpolation and Extrapolation. Numerical Methods. Note  « 12 »

Polynomial interpolation. Lagrange interpolating polynomial Pn-1(x) is the polynomial of the order n-1 that passes through the n given points (xi,yi), i=1..n

P(x) = ∑i yik≠i (x-xk)/(xi-xk) .
Might show erratic behaviour for higher degree polynomials. Not useful for n larger than, say, 5.

Rational function interpolation. The rational function R(x) is a ratio of two polynomials Pm(x) and Qk(x),

R(x) = Pm(x)/Qk(x) .
It has potentially larger range of conversion than polynomial. Still dangerous for higher orders.

Spline interpolation. At the interval between the tabulated points i and i+1 the function is approximated by a polynomial:

y[i](x) = yi + bi (x-xi) + ci (x-xi)2 + di (x-xi)3 + ... .
The coefficients bi, ci, di,... are calculated from the continuity equations for the function and its derivatives at the boundarys
y[i-1](xi)=y[i](xi),  y[i-1]'(xi)=y[i]'(xi),  y[i-1]''(xi)=y[i]''(xi), ...
For the quadratic spline the equations become
yi+bih+cih2 = yi+1 (continuity of the function at the (i+1)st point),
bi + 2cih = bi+1 (continuity of the derivative)
The most popular spline is the cubic spline. Here are the fortran routines from netlib.org: spline.f - for constructing the spline (calculating coefficients b,c,d) and seval.f - for the very interpolation.

Integration and differentiation of tabulated functions. First interpolate the tabulated data and then integrate or differentiate the intepolating function.

Problems

  1. Make a function locate(x,xt[]) which, given an ordered array xt[] and a number x, finds index j such that xt[j] < x < xt[j+1].
  2. Make a function polinterp(x,xt[],yt[],k) which does a polynomial interpolation of a table (xt[],yt[]) at a point x using k nearest tabulated points.
  3. (Non-oblig.) Make a function for the quadratic spline interpolation.
  4. Make a Lagrange polynomial interpolation (and a spline interpolation) of the following table:
    {sinh(1.5 i/5), 1/cosh2(1.5 i/5)}, i=-5..5
    Make a plot of the result.

"Copyleft" © 2004 D.V.Fedorov (fedorov at phys dot au dot dk)