# Developing a meshfree contouring algorithm

I would like to find the contours of a scalar function $$u(x,y)$$ available as a discrete set of function values $$u_i = u(x_i,y_i)$$ over a scattered set of points $$(x_i,y_i), i=1,...,N$$. In my case, the scattered data comes from a numerical solution to a PDE using the RBF finite difference method, however interpolation of any type of scattered data is just as relevant.

In a related thread, the recommendation was to form a triangulation of the scattered points first, and then using a common contour algorithm such as marching triangles in 2D. While this is a tested approach, building the triangulation defeats one of the primary strengths of the RBF-FD method, which is exactly the fact it can avoid polygonal meshes. One of the answers in that thread suggested a RBF approach could work better.

Apart from the point coordinates and function values, suppose I also have the adjacency graph for the RBF-FD stencils of size $$s$$, and the RBF-FD differentation weights $$w_i$$:

$$L(u(\boldsymbol{x}))|_{\boldsymbol{x}_c=\boldsymbol{x}_1} \approx \sum_{j=1}^s w_j u_j$$

which are constructed in a way, that they can provide approximations to various (linear) spatial differential operators $$\partial_x, \partial_y, ...$$.

Returning to the initial problem, finding a contour level $$h$$ can be recast to finding the level set (zero contour) of the implicit function,

$$\phi(x,y) = u(x,y) - h = 0,$$ which can be seen as a signed distance-like function (except for the fact it doesn't really return the distances).

Given some initial guess $$\boldsymbol{p} = (p_x, p_y)$$ of a point near the contour, we can try and find a projection $$\Delta\boldsymbol{p}$$ so that

$$\phi(\boldsymbol{p} + \Delta\boldsymbol{p}) = 0, \\ \Delta \boldsymbol{p} \,|| \, \nabla\phi(\boldsymbol{p}+\Delta \boldsymbol{p}),$$

where the projection is to be orthogonal to the gradient at the projected point. In the thesis of Per-Olof Persson (pp. 47-48) a first order approximation of the projection is found using Lagrange multipliers and a truncated Taylor expansion. The approximate projection formula is

$$\Delta \boldsymbol{p} = \frac{\phi}{|\nabla \phi|^2} \nabla \phi.$$

An example of the approximate projection result is shown in the figure below. First, I identified the scattered points near the contour by looking for stencils, where the sign of $$\phi$$ changes. Using the RBF-FD spatial operators to construct the gradient $$\nabla = (\partial_x, \partial_y)$$, I could project the nearby points to the target contour:

What I am missing now is a robust way to connect the projected points into a single contour. Is their some graph approach that could be used (construct the shortest path)?

In view of the marching squares algorithm, I am also anticipating the path-finding will have disambiguation difficulties near saddle points. Is their any way to detect such saddle regions and if yes, how to make the path-finding robust?

• Your introduction is quite long. Do you want to find an algorithm that connects the blue/green points with the shortest path? If yes, the problem is a NP hard problem and can be solved with the traveling salesman algorithm. Nov 28 '20 at 7:47
• Precisely. I was hoping that by using the adjacency graph of the points near the contour (e.g. to zigzag between +/x symbols) it might be possible to solve this more elegantly. Nov 28 '20 at 11:14
• I think for this you have define your problem more precisly. Most of the stuff you mentioned is not needed (my point of view) to describe your actual problem. Are the adjacency points ordered in some kind? For what kind of application do you need the exact contour line? Nov 28 '20 at 15:00
• I’m not sure I understand the problem. Couldn’t you just fit a curve to the points? May 1 at 13:06

Maybe this is not a full answer to your question. However, I am currently developing my own unstructured mesh generator and found this toolbox quite helpful. Perhaps there are some algorithms which will help you.

Short (not complete) overview of the genetic algorithms for the TPS problem:

• Solves the classic Traveling Salesman Problem (TSP)
• Solves an open variation of the TSP with no start/end constraints
• Solves an open variation of the TSP with fixed start & end constraints
• Solves an open variation of the TSP with only a fixed start constraint
• Solves a variation of the TSP that maximizes the total distance
• Nearest Neighbor (NN) solution to the TSP
• Nearest Neighbor (NN) solution to the open variation of the TSP (no start/end constraints)
• Farthest Neighbor (FN) solution to the TSP
• Random Search (RS) solution to the TSP (provides interesting performance comparisons)
• Solves the TSP with a variation of the GA that uses a cross-over operator
• Solves the TSP with a variation of the GA that uses a hybrid set of operators
• Variation of the TSP to find the shortest tour through clusters of points
• Solves the TSP with a variation of the GA that increases the mutation rate on the best route
• Solves the Open TSP with a variation of the GA that increases the mutation rate on the best route
• Script to compare tsp_ga performance with an approach that seeds the algorithm with the NN solution
• Script to compare tsp_ga_turbo performance with an approach that seeds the algorithm with the NN solution
• Solves the classic Multiple Traveling Salesmen Problem (M-TSP)
• Random Search (RS) solution to the M-TSP (provides interesting performance comparisons)
• Solves a variation of the M-TSP where all salesmen start and end at the first city
• ...