Make barycentric interpolation retain the bounding polygon

[Huite]

A current problem of the BarycentricInterpolator when called with face associated data is that it shrinks the resulting mesh, since it cannot easily interpolate beyond centroids. To illustrate:

The black lines are the original (triangular) grid. The orange lines form the centroidal voronoi tesselation. The red lines are the voronoi rays, which creates new vertices through projecting the centroids on the exterior edges. During interpolation, we compute the Wachspress barycentric coordinates, giving weights for each centroid. Newly created vertices (4, 5, 6, 7) are not centroids, but can be directly associated with them, providing not too much issue during interpolation.

The problems start with the maintained exterior vertices (1, 2, 3). 1 and 3 are not so much a problem, like the projected vertices they can be associated with a single centroid. 2 forms the problem, since it creates a convex face in which we cannot use Wachspress barycentric coordinates; alternative barycentric coordinates won't do much good either since they may be negative which I'm pretty sure would directly result in extrapolation.

I've though about simply padding in #149, but it's not very satisfying.

Perhaps a more robust scheme would use skip the exterior vertices for vertices like 2, but keep them for 1 and 3. We associate al the new vertices (projections, exterior vertices) with the centroids. Any "concave" vertex would be inside. However, this doesn't work if the vertex is exterior and convex:

We could use barycentric coordinates just fine since it's convex, but it doesn't have a clear value associated with it.

[Huite]

Some options here:

1. Use the nearest centroid
2. IDW weighted value of centroids
3. Ghost cells / centroids

The ideal option would be one that generates a grid with exclusively convex faces such that the barycentric interpolation is straightforward. The problem with e.g. ghost centroids from mirroring of the centroid across the edge is that they may end up in other cells. Then, we cannot simply search the grid using the celltree because it will conflicting answers.

If everything is strictly convex and larger than the original grid, we could compute weights, then discard the points that fall outside of the original grid.

[Huite]

Alternative scheme: "walk" along the edges, from points 5 and 6 to point 2, measure the distance. Use the distances to assign a weighted mean to point 2. Make the result always larger than the original: reject concave original vertices, keep convex original vertices. Then run a barycentric interpolation.

[Huite]

This last scheme seems to perform well. In case of quads, it does exactly the right thing.

Instead of:

We get:

To explain, the centroidal voronoi tesselation generates the projections of the centroid q and r on the exterior edges. The original vertex p is always located "between" them. We measure the distance p -> q and the distance p -> r.

p isn't associated with any centroid, but we can associate it with nodes q and r, which themselves are directly associated with a centroid. r with c1, and q with c2. Then during barycentric interpolation, we compute barycentric weights as usual, including p. But afterwards, we redistribute the weight for p to c1 and c2.

The effect is the same as first linearly interpolating along the exterior to find a value for p, and then to find all values inside of (c1, r, p, q, c2).