# How does Interpolation really work?

I'm looking for some explanation or advice, not help in solving something. Recently I finished my program and my supervisor said "Ok, now it's time for your first paper: write a scientific text about how your program works". If I understand him correctly it means I have to describe the algorithm I used instead of writing in a manner like "for this purpose I use the built-in Interpolation function, and for this purpose I use NDSolve" etc. I know how FindRoot and NDSolve work because there is an explanation in the doc pages about the methods they use, but I did not find detailed information about Interpolation. The only thing I know it fits with polynomial curves.

So my question is: What exactly does the Interpolation function do? How does it work? How does it determinate (partial) derivatives? And why do 3D data points have to be situated in rectangle order to interpolate a surface?

If there is some literature I may read and then reference it would be great too.

• My paper is about numerically implementation of some theoretical smooth model which deals with differential equations on a given surface. So the way surface is given is not highly important and polynomials (built-in interpolation) work fine. I don't think I have to re-invnet the wheel and create my own version of Interpolation function but I think I have to know how it works exactly. – ddd Apr 13 '12 at 4:49
• I just don't know: if " fitting polynomial curves between successive data points" is speciefed well-known method like Runge-Kutta for NDSolve or Newton for FindRoot, then my question is useless. But I guess (I don't know) there are many of them, so I want to know which one is implemented in Mathematica – ddd Apr 13 '12 at 4:53
• @lcanix, Interpolation uses several different methods. You could try different methods and see how they affect your solution. If no affect is seen, then the method may not be that important. If you see a influence it might be helpful if you could add a simple example in your text to see where the Interpolation algorithm might be going. – user21 Apr 13 '12 at 6:05
• Were you able to find enough research material for your paper? Not to detract from mathematica in any way, but since your goal is to understand splines in general, start off with beta splines because they are so common and you can see them applied in a variety of contexts. For example, they are popular in animation programs for drawing character paths. If playing with parameters in mathematica is not inspiring, look at the animation program Maya on the autodesk website and look up path generation using splines. Hope that helps, Iceberg – user994 Apr 16 '12 at 6:36

# Interpolation function methods

Interpolation supports two methods:

• Hermite interpolation (default, or Method->"Hermite")
• B-spline interpolation (Method->"Spline")

## Hermite method

I really can't find any good reference to Hermite method within Mathematica's documentation. Instead, I recommend you to take a look at this Wikipedia article.

The benefits of Hermite interpolation are:

1. You can compute them locally at the time of evaluation. No global system solving required. So construction time is shorter, and the resulting InterpolatingFunction is smaller.
2. Multi-level derivatives can be specified at each point.

One problem is that the resulting function is not continuously differentiable ($C^1$ or higher), even if InterpolationOrder->2 or higher is used. See the following example:

## Spline method

To be specific, we are using B-spline interpolation with certain knot configuration--depending on the distribution of sample points. I could not find a good web source to describe the method (the Wikipedia article is not great). Although, you can find a step-by-step description of the method in 1D case within Mathematica's documentation (BSplineCurve documentation, Applications -> Interpolation section). Multi-dimension is simply tensor product version.

The benefits:

1. InterpolationOrder->d always guarantees a smooth function of $C^{d-1}$ class.

2. Evaluation/derivative computation is very fast.

3. You can take BSplineFunction out of the resulting InterpolatingFunction (it's the 4th part), which is compatible with BSplineCurve and BSplineSurface for fast rendering.

The problems (of current implementation in V8):

1. It is machine precision only--although, it is not hard to implement it manually for arbitrary precision using BSplineBasis.
2. It does not support derivative specification.
3. Initially it solves global linear system and store the result. So the resulting function is much larger than Hermite method (this is not implementation problem).

# Other functions

Some plot functions such as ListPlot3D have their own methods. Sometimes they call the B-spline method, sometimes they use a method based on distance field (for unorganized points), etc. But probably it is not useful here since they are only supported as a visual measure.

• I understand cubic spline interpolation has a few different varieties depending what the derivative does at the end points. I thought InterpolationOrder->3 used a cubic spline interpolation. Is a cubic spline interpolation a specific case of a B-Spline? – Ted Ersek Apr 13 '12 at 21:57
• Yes, so called cubic spline interpolation is a special case of B-spline interpolation. Now the problem is that the current Mathematica implementation uses something called "clamped" knot configuration, where as the cubic spline interpolation uses "unclamped" or "natural" configuration. The difference lies in how the end points (derivatives) are being treated. This link has a brief explanation. – Yu-Sung Chang Apr 14 '12 at 1:27
• @Yu-Sung What is the relation between InterpolationOrder->n (Method->"Spline") and commonly used term "n-point spline interpolation"? – Alexey Popkov Oct 5 '12 at 1:43
• @Yu-Sung Chang: Do we know that Interpolation uses only those two methods? All that ref/Interpolation says is, "Interpolation supports a Method option. Possible settings include 'Spline' for spline interpolation and 'Hermite' for Hermite interpolation." [emph added] – murray Apr 21 '13 at 15:06

## Method -> "Spline"

Let us see how interpolation using Method -> "Spline" works in Mathematica.

A spline of degree m is a piecewise polynomial with (m-1) continuous derivatives (although the term "spline" sometimes is used in more general sence). The m'th derivative will (generally) have discontinuities at the points of splicing of individual polynomials of degree m from which the spline is formed. Let us see where these points are located (the green line is (m-1)'th derivative and the red line is m'th derivative):

SeedRandom[60]
n = 6;
data = Transpose[{Range[0, n - 1]/(n - 1), RandomReal[{-1, 1}, n]}];
Do[sp[ord] =
Interpolation[data, Method -> "Spline",
InterpolationOrder -> ord], {ord, n - 2}];
Grid[Partition[
Table[Show[
Plot[{D[sp[ord][x], {x, ord - 1}], D[sp[ord][x], {x, ord}]}, {x,
0, 1}, PlotRange -> All, PlotStyle -> {{Green, Thick}, Red},
Exclusions ->
If[OddQ[ord], data[[(ord + 3)/2 ;; -(ord + 3)/2, 1]],
MovingAverage[data[[ord/2 + 1 ;; -ord/2 - 1, 1]], 2]],
PlotLabel ->
Row[{Method -> "\"Spline\", ", InterpolationOrder -> ord}],
Filling -> {2 -> Axis}, Evaluated -> True],
ListPlot[data, PlotStyle -> Directive[Blue, PointSize[Medium]]],
ImageSize -> Medium], {ord, n - 2}], 2, 2, 1, ""]]


It is easy to see that for odd-degree splines the discontinuities are located at the data points while for even-degree splines they are located in the middle of them. And the spline consists of (n - m) individual polynomials where n is number of data points and m is degree of spline.

P.S. I'm looking for readable reference where this kind of splines is described...