11
$\begingroup$

What specific method does Interpolation use for unstructured multi-dimensional data when we set InterpolationOrder -> All? Documentation links are welcome.

Example 2D data:

data = RandomReal[1, {20, 3}];

When the data points are not on a grid, the only allowed settings for InterpolationOrder are 1 and All, according to the error message issued when trying something else.

With 1, it is clear how it works: a Delaunay triangulation is computed and linear interpolation is done over each triangle.

But how does All work, and what determines the actual order that is chosen?

if = Interpolation[data, InterpolationOrder -> All];

if["InterpolationOrder"]
(* 5 *)

Show[
 Plot3D[if[x, y], {x, 0, 1}, {y, 0, 1}],
 Graphics3D[{PointSize[Large], Point[data]}]
]

enter image description here

|improve this question
$\endgroup$
  • $\begingroup$ Dunno, but the return value of if["InterpolationOrder"] that I get is {9223372036854775806, 9223372036854775806}. Oo $\endgroup$ – Henrik Schumacher Apr 16 '19 at 8:49
  • 1
    $\begingroup$ @HenrikSchumacher Oops ... It seems I tried this with M12.0 (it's available in the cloud). $\endgroup$ – Szabolcs Apr 16 '19 at 8:54
  • 1
    $\begingroup$ Anyways, very good questions. I am also curious what works there in the background. $\endgroup$ – Henrik Schumacher Apr 16 '19 at 9:00
  • 1
    $\begingroup$ That sounds as if they were using straight-forward global interpolation by a polynomial of degree up to n. Then you have Binomial[n, 2] basis functions. In that case, this should become nasty for higher point counts due to Runge's phenomenon and ill-conditioned linear systems (for solving for the coefficients). So I presume, that they will switch to another method when the point count becomes larger... $\endgroup$ – Henrik Schumacher Apr 16 '19 at 9:06
  • 4
    $\begingroup$ In this context, I just want to point to the "Obtuse" package for doing multidimensional interpolation on scattered data points. You can download it (for free) from familydahl.se/mathematica. The package implements five different methods: the Delaunay interpolation method, the Voronoi interpolation method, the Shepard interpolation method, the RBF interpolation method and the ObtuseAngle interpolation method. For the ObtuseAngle method the InterpolationOrder can be set to 0, 1, 2, and 3, but not to All. The number of dimensions is not limited. $\endgroup$ – Ingolf Dahl Jun 12 '19 at 19:54
13
$\begingroup$

This is code that has been written many moons ago... here is an example:

d = {{0.4138352728412389, 0.02365673668161028}, {0.5509946389658635, 
    0.7254061374370833}, {0.14521595926324116, 
    0.6528630823305817}, {0.48768962246740544, 
    0.22066264105073286}, {0.8309710560928056, 
    0.3496966364384875}, {0.4553589220242207, 
    0.9383446951847001}, {0.2126873262146789, 
    0.017512080396716145}, {0.967248982535015, 
    0.6211273372083488}, {0.3548669163916416, 
    0.737108322193581}, {0.6919974835480842, 0.9322403408098401}};
f = {{0.9953617542392983}, {0.14070666511222818}, \
{0.285662339441511}, {0.7988192898854105}, {0.3592646208757597}, \
{0.565455746009103}, {0.22110814761432618}, {0.2735048548887764}, \
{0.08792348530403005}, {0.4202942851818514}};
data = Join[d, f, 2];
if = Interpolation[data, InterpolationOrder -> All];
if[0.5, 0.5]

0.268157

And here is roughly what it does:

dt = Transpose[d];
temp = Join[{ConstantArray[1., {Length[d]}]}, dt, dt[[{1}]]^2, 
   dt[[{1}]]*dt[[{2}]], dt[[{2}]]^2, dt[[{1}]]^3, 
   dt[[{1}]]^2*dt[[{2}]], dt[[{1}]]*dt[[{2}]]^2, dt[[{2}]]^3];
p = Transpose[temp];
ls = LinearSolve[p];
vals = ls[Flatten[f]];
System`Private`EvaluateListPolynomial[vals, {0.5, 0.5}]

0.268157
|improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.