5
$\begingroup$

I've seen a few other posts asking about the "unstructured grid" error in Interpolation. Here's a minimal (non)working example of the issue:

z = RandomReal[{-1, 1}, 50];
xy = RandomPoint[Rectangle[{0, 0}, {1, 1}], 50];
Interpolation[Thread[{xy, z}], InterpolationOrder -> 2]

enter image description here

One way around this is to use InterpolatingPolynomial like this:

data = Thread[{xy, z}];
ip = InterpolatingPolynomial[data, {x, y}];
Show[Plot3D[ip, {x, 0, 1}, {y, 0, 1}], 
 Graphics3D[Point[Flatten /@ data]]]

enter image description here

However, it doesn't work well in practice as the polynomial orders are uncontrollable (too high):

enter image description here

Are there any viable alternatives for successfully getting quadratic InterpolationOrder in Interpolation?

$\endgroup$
1
  • $\begingroup$ A workaround: Perhaps you can create a structured mesh using the FEMfunctionality and ElementMeshInterpolation[] for quadratuic elements? $\endgroup$ – Ulrich Neumann Dec 19 '18 at 11:01
6
$\begingroup$

You can use interpolation based on radial basis functions. To demonstrate the method method see the following example.

Lets start with some data:

data = GeoElevationData[GeoDisk[Entity["Mountain", "MountEverest"], Quantity[3, "Miles"]]];
ListPlot3D[data, MeshFunctions -> {#3 &}, Mesh -> 40] 

enter image description here

To make the computation faster I downsample the data

smallData = ArrayResample[QuantityMagnitude[data], {100, 100}];

Now we select 30 random positions as sample points that we are going to use for the interpolation:

randomx = RandomInteger[{1, 100}, 30];
randomy = RandomInteger[{1, 100}, 30];
stuetzstellen = Table[{randomx[[i]], randomy[[i]], QuantityMagnitude[smallData[[randomy[[i]], randomx[[i]]]]]}, {i,1, Length[randomx]}];
Show[{
  ListPlot3D[smallData, MeshFunctions -> {#3 &}, Mesh -> 40], 
  ListPointPlot3D[stuetzstellen, PlotStyle -> PointSize -> 0.03]}]

Mathematica graphics

Now we can define the RBF. There are many possible choices. You can experiment with other functions.

rbf[x_] /; (x > 0) := N[x^2 Log[x]]
rbf[x_] /; (x == 0) := 0
distances[{x_, y_}, data_] := EuclideanDistance[{x, y}, #] & /@ data[[All, {1, 2}]]

The basis for the interpolation is the distance matrix for our sample points:

distMat = DistanceMatrix[stuetzstellen[[All, {1, 2}]]];
coeffs = LinearSolve[Map[rbf, distMat, {2}], stuetzstellen[[All, 3]]]

Now we can compute the interpolated mountain:

Show[{
  Plot3D[coeffs.rbf /@ distances[{x, y}, stuetzstellen], {x, 0, 100}, {y, 0, 100}], 
  ListPlot3D[smallData, MeshFunctions -> {#3 &}, Mesh -> 40, PlotStyle -> Directive[Red, Opacity[0.5]]], 
  ListPointPlot3D[stuetzstellen, PlotStyle -> PointSize -> Large]}]

Mathematica graphics

The more sample points we have the better will be the approximation. This can also be expanded to higher dimensions.

$\endgroup$
1
  • 3
    $\begingroup$ This is a nice answer. $\endgroup$ – user21 Dec 20 '18 at 6:32
4
$\begingroup$

The reason only first order interpolation is possible on unstructured grids is because it is not at all clear how to construct the mesh from the unstructured data given. Additional information is necessary; if you have that you can construct a higher order interpolating function.

First, we construct a Delaunay triangulation of the coordinates given.

Needs["NDSolve`FEM`"]
m = ToElementMesh[xy]
(*m["Wireframe"]*)

Next, we convert that mesh to a second order mesh.

m2 = MeshOrderAlteration[m, 2]

Now the mesh has coordinates at the mid side nodes and you will need data points at these positions. If you have those you can proceed with:

(* generate some random values for demonstration *)
z2 = RandomReal[{-1, 1}, Length[m2["Coordinates"]]];

Create a second order interpolating function.

if = ElementMeshInterpolation[{m2}, z2]

Then you can visualize the interpolating function and do whatever you would like to do with it.

Plot3D[if[x, y], Element[{x, y}, m2]]

To sum up, the problem with arbitrary coordinates is that no higher order mesh can be generated. If you can generate that mesh and have function values at the higher order nodes, you are good to go.

$\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.