How to interpolate multipoint data to the cell centres of an unstructured mesh?

I have sets of multipoint field data, each point data set relates to a single cell of an unstructured mesh. The goal is to interpolate the data to the cell centre, directly or indirectly, in the most accurate way.

If I use Inverse Distance Weighted interpolation, in the case when the distance between the source and the target (cell centre) is very small, I may end up with a floating point exception.

For this kind of interpolation on a structured mesh, a volume weighted interpolation is used. This does not translate directly to an arbitrarily shaped mesh cell.

Introducing a tolerance for a IDW interpolation to circumvent the SIGFPE makes sense only if I do not introduce any tests that could render the interpolation inefficient. Is adding a sufficiently small $\delta$ to the denominator for every weight a possible option with the IDW interpolation? What interpolation methods suitable for this problem do you know?

For the interpolation from the mesh to the points, I am using an interpolation based on the barcycentric coordinates. Each polyhedral cell of the mesh is decomposed into tetrahedra. Cell centred field is interpolated to the cell points using IDW interpolation. A search is conducted for each point to find the tetrahedron within which it lies, and the values are interpolaed using the barycentric interpolation.

For the interpolation from the points to the mesh, this is not possible. The cell centred values are unknown. There is no way to assemble a tetrahedral composition that would enforce $\sum_p W_{PC} = 1$, where $W_{PC}$ is the weight related to a point P and a cell centre C. This comes from the fact that the point configuration is arbitrary. So, I am currently using IDW for this, making sure that I don't get a floating point exeption. Are there any better suited interpolation methods for this problem?

• Can you be a bit more specific regarding the weighting function you are using? There are several interpolation approaches which use polynomial weighting functions which do not have singularities at the endpoints. – Pedro May 21 '12 at 10:23
• If the function to interpolate is smooth, have you thought of using radial basis functions to globally interpolate from the mesh? You could built a quick prototype using python and scipy, see here: docs.scipy.org/doc/scipy/reference/generated/… – fcruz May 21 '12 at 10:40
• The library I'm working in actually supports RBF (both compact and global support), but I have very little experience with this. Thanks for the tip.. :) – tmaric May 21 '12 at 11:48
• If you do try the SciPy RBF library, be aware of this outstanding bug – n00b May 21 '12 at 16:50
• Do you have values of gradient of dependent variable at cell vertices? – Johntra Volta May 25 '12 at 19:21

Links to diverse software packages for scattered data interpolation are on my web page http://www.mat.univie.ac.at/~neum/stat.html#fit

The book
G.E. Fasshauer, Meshfree Approximation Methods using MATLAB, World Scientiﬁc 2007.
gives a comprehensive state of the art (as of 2006).

A few more recent papers on scattered data interpolation:
http://www.stanford.edu/group/uq/pdfs/journals/jcp_scattered_2010.pdf
http://www.math.auckland.ac.nz/~waldron/Preprints/Box-splines/box-splines.pdf

Which method to use depends a lot on the use made of the resulting interpolant. Kriging methods are based on a stochastic model, hence are good if the data to be interpolated are somewhat noisy. Radial basis functions are to be preferred if (implemented stably) and a visually pleasing result is wanted (low curvature interpolation).

• Professor Neumaier, of the methods you've gathered in the links on your web site and in your answer, do you have any specific recommendations for the application described in the question? – Geoff Oxberry May 22 '12 at 2:45
• I guess I have some learning to do... as always, the amount of things that need learning grow exponentially with time, and the amount of things that I manage to learn, grow linearily at best. :) Thanks. I'll use IDW interpolation for the initial implementation, and read about scattered data interpolation (RBF seems to be the trendy choice). :) – tmaric May 22 '12 at 10:49
• @GeoffOxberry: I collected the links over time, without having tried the packages. Thus I cannot recommend a particular one. But I added a comment on qualitative recommendations. – Arnold Neumaier May 22 '12 at 14:11

Below I will give an example how do I interpolate from one set of points to another one, on finite volume mesh.

I have collocated arrangement of variables - the data I store in memory represent values at cell-centers. I store field variables and their gradients. Gradients are found from surrounding values solving a least-squares problem (with QR via Householder reflections).

Your arrangement may differ but the principle is the same.

Then if I'm looking for $\phi_f$ - a value at cell face center I may get it from:

$\phi_{nb1} + \nabla\phi_{nb1}\mathbb{r}_{nb1,f} = \phi_f$

$\phi_{nb2} + \nabla\phi_{nb2}\mathbb{r}_{nb2,f} = \phi_f$

...

$\phi_{nbn} + \nabla\phi_{nbn}\mathbb{r}_{nbn,f} = \phi_f$

where, $nb$ designates a neighboring cell center,they go from 1 to n (very often just 1 suffices, I use 2, i.e. I use neighboring cells that share the face). $\mathbb{r}_{nbn,f}$ is a distance vector from n-th neighbors cell's center to face center $f$.

Then I write

$\phi_f = \frac{1}{n}( \sum_{i=1}^{n} \phi_{nbi} + \sum_{i=1}^n (\nabla \phi_{nbi}\mathbb{r}_{nbi,f}) )$

So you need one set of field values and gradients at those points. You need to decide which surrounding points will contribute to your interpolated point, as well as distance vectors from these points to point to which we interpolate.

For example: If one stores data representative of values at cell vertices you use this equation to find cell-center values, etc., all depending of what situation you have.

So this is based on Taylor series around the point. One can use also second derivatives to derive a more accurate expression.

• Hi, didn't have time to reply until now... the ewuation for $\phi_f$ with a $\frac{1}{n}$ multiplier works for cell-point communication on a uniform mesh. I'm dealing with general cell shapes and general distribution of points, so I need weighted interpolation, not averaging.... at least I think so.. the 1/n prefactor will make it behave as if the values are averaged, right? What I could do is a LSQ : make a square of your bracket and minimize... – tmaric Jun 3 '12 at 17:48
• I wrote this in my code on the fly. Here it is only assumed that data varies linearly in space, that's why only first derivatives of Taylor series are included. There are no assumptions about the mesh, it works for any mesh. – Johntra Volta Jun 3 '12 at 21:32
• However, you may add weights. Like in least-square gradient reconstruction, where we have both weighted and unweighted version. – Johntra Volta Jun 3 '12 at 21:38
• Well, thanks, but I need this for general fields, and especially for fields involving sharp jumps (two phase flow)... if I just do an average, this will smooth out the fields. This seems like a reverse LSQ problem, instead of finding the gradient, I can find the value, but I will have to use minimization, not $\frac{1}{n}$ at the end. I'll try it, if it works, I'll accept your answer! Thanks! – tmaric Jun 4 '12 at 8:55
• You may contact me by email (you can find it on my user profile) for further discussion. I'm interested to see how will you solve the problem. – Johntra Volta Jun 4 '12 at 14:36