A Brute-Force Constrained Delaunay Triangulation?

Question

I need to create a triangle mesh from a set of points. The set has very few points so it doesn't need to be fast or optimised (I will deal with 100 points maximum). The mesh needs to be a constrained "delaunay triangulation". In the image below I showed (on the left) the set of points I start from (blue and red dots). I also know the connections between these points (the outline in black). The mesh needs to look like the example on the right (including the edges in grey that form outside and inner triangles).

I can't use libraries.

I looked at many different algorithms. They are many and it's easy to be confused. I would like to know if there is a naive and thus hopefully simpler algorithm I can use in order to produce the mesh on the right? Brute force approach is fine (ps: I can do a Delaunay triangulation).

Sourdine

Answer 1

Thanks for all the answers.

I went through the process of developing a solution to this problem so thought I would share my own experience, hoping people facing the problem will find the insight useful.

So from my own experience implementing an algorithm I came to the conclusion:

1. There is not really quick way to this problem. It is not reasonable to think it can be achieved in just 50 lines of code. In fact the routine that I wrote (C++) is about 400 to 500 lines (hard to tell with comments). So reasonably compact yet challenging and it took me a 2 to 3 days to get the logic right (it can be tricky).

2. I found the algorithm propose by Sloan in "A FAST ALGORITHM FOR GENERATING CONSTRAINED DELAUNAY TRIANGULATIONS" to be perfectly well suited for the problem at hand. The reality when it comes to Delaunay triangulation which was a new subject for me, is that there seems to be a lot of different algorithms approach and this research is pretty old. So for a new comer it's really hard to know where to start.

   2.1. It's hard to know which algorithm is recent, simple in its comprehension and fast and simple to implement.

   2.2 Generally once you have understood the principle it's mostly a matter of coding the logic in the most efficient way (and it seems that this what most of the algorithms/papers are fighting above).

   2.3 I found the paper from Sloan to be understandable, very well explained. If you follow the logic and the instructions, then anyone can really implement a constrained Delaunay triangulation.

So in conclusion:

1. I recommend the Sloan paper because it includes an explanation on how to create a Delaunay triangulation followed by a constrained triangulation if necessary.

2. To answer my own question there is not really brute force to this problem. Implementing this technique just requires to implement the full logic and most implementation must require more or less the same amount of work

3. With the nuance, that I wasn't looking after much optimisations because my point sets are really small. So I am sure many algorithms are better than the one described by Sloan; they probably propose optimized data structures and algorithms optimized to minimize steps such point insertion in triangulation etc.

So anyway Sloan worked. A small image to illustrate the answer and make it more attractive;-)

EDIT

This is production code so helas I can't share that... I could lead me to be fired. The process is very simple though. You look for the intersection between a segment (your constraint) and all edges in the model. Then for each intersected edge, you swap the diagonal between the 2 triangles that this edges belongs to. If the new diagonal intersects the segment still, then add the new diagonal back onto the stack of intersected edges for this segment. If the new diagonal doesn't intersect the segment then add it to the stack of newly created edges. Keep processing the stack of intersected edges until it's empty.

Then once this is finished you need to process the list of new added edges. For each one of them, check that the Delaunay triangulation criterion is respected. If not swap the diagonal of the triangle this edge belongs to. Simple ...

This is just the paper ...

Point Set

26.9375 10.6875
32.75 9.96875
31.375 4.875
27.6562 2.0625
23.9375 -0.75
18.1562 -0.75
10.875 -0.75
6.60938 3.73438
2.34375 8.21875
2.34375 16.3125
2.34375 24.6875
6.65627 29.3125
10.9688 33.9375
17.8438 33.9375
24.5 33.9375
28.7188 29.4062
32.9375 24.875
32.9375 16.6562
32.9375 16.1562
32.9062 15.1562
8.15625 15.1562
8.46875 9.6875
11.25 6.78125
14.0312 3.875
18.1875 3.875
21.2812 3.875
23.4687 5.5
25.6562 7.125
8.46875 19.7812
27 19.7812
26.625 23.9688
24.875 26.0625
22.1875 29.3125
17.9062 29.3125
14.0312 29.3125
11.3906 26.7188
8.75 24.125

These are x/y/z coordinates (z=0)

Segments:

0 1
1 3
3 5
5 7
7 9
9 11
11 13
13 15
15 17
17 19
19 20
20 22
22 24
24 26
26 0
28 29
29 31
31 33
33 35
35 28

Indices start at 0 (0 -> first vertex in vertex list)

Sourdine

:Hey, can u share ur code? I am interested in sloan in the part where the triangle intersects with a constraint edge. It's stated to change the diagonal of the quadrilateral and repeat it until it fits. Eventually find new vertices Vn and Vm! I didn't understand this part! Thanks in advance! – Concepcion

:Did u tried my algorithm? Can u send me the points from ur example? – Concepcion

: Hey, thank you, I get this:imgur.com/a/DPZUx. Upper is concave dt, lower is dt. IMO the concave algo has a hole because of the alphashape but why is there a red inner border? Well a bit confusing result but it's CONCAVE ALGO. Do you need the constraints edges or just the points for sloan? – Concepcion

Yes this seems like DT only. But good;-) Yes I specified what the constrained edges need to be. Sorry this is not in the data set. Let me add this for you. It will be in the form of indices in the array of vertices I gave you. Will do this later this week. – Sourdine

@bettedev: just added data for segments – Sourdine

:Hey, thank you very much but I have a little problem I am missing a point and edge. Please see my dt:imgur.com/a/OReMN. Care to elaborate? There is also 37 indices but only 36 points!? – Concepcion

Answer 2

I tried it with alpha shapes with good results for a few shapes https://concavehull.codeplex.com/ but it's nowhere near the original constrained delaunay triangulation.

Here is my alpha-shape algorithm:https://alphashape.codeplex.com.

Concepcion

Great thanks for the tip. I am trying the Sloan algorithm for now (1992) and I will post about my results if i get it to work. If I don't I will try/have a look at your solution. – Sourdine

@user18490: I read a bit from Sloan. IMO not worth it but u can try Bowyer-Watson. U can use it with a Hilbert sort. But sloan doesn't solve cdt. – Concepcion

thanks for your advice. I am looking at this paper. researchgate.net/profile/Scott_Sloan/publication/… There is an algorithm for CDT. I will have a look at the other technique you mentioned – Sourdine

(at)Bettedev so the Bowyer-Watson seems to give me a DT but how do I get this CDT? by using the Hilbert sort that you mentioned? Could you be a bit more specific please. I would greatly appreciate. Thx – Sourdine

@user18490: Ear-clipping but it doesn't produce dt. Or try mine :). – Concepcion

@user18490: I read the paper it's useing quadrilaterals:www-cs-students.stanford.edu/~amitp/game-programming/… – Concepcion

Answer 3

A simple approach seems to be to implement a ear clipping algorithm. Without optimisation as in hash grids or quad trees. For ear clipping you just check every three consecutive vertices a,b, and c. If b is convex and no other vertex of the polygon lies inside the triangle abc then you can clip this triangle reducing the boundary of the polygon by one vertex, b.

Additionally you have to store the neighbourhood relations. Thus, reference from each triangle its, at most three, neighbours.

When the triangulation is finished you convert it to the constrained Delaunay triangulation (CDT). This can be done by edge flipping. Therefore you have to check for every triangle the circumcircle. If no vertex of a neighbouring triangle lies inside the triangle is conform to the CDT otherwise flip the edge of the triangle where the violation occurs.

Edit due to @Betterdev in the comments blow: Possible holes in the input polygon can be added to the initial boundary by adding a bridge. As a preprocessing one can connect a vertex of a hole to a vertex of the boundary by a "double" edge. This is always possible and makes each hole part of the main polygon boundary; and works well with ear clipping. Storing the neighbour through these bridges is vital to the flipping however.

Subdeacon

:Is ear clipping delaunay? – Concepcion

@Betterdev No, ear clipping is only a simple triangulation algorithm. That is why I also described the flipping algorithm as a second step. Edge flipping will then convert the triangulation into a CDT. – Subdeacon

:What about holes and crossings? – Concepcion

If the input polygon is simple there should not be a crossing. If you have crossing constraints one can use a sweep line or just brute force to compute the intersections and resolve the crossings. – Subdeacon

:What about the sloan algorithm? Can it triangulated crossings and holes? – Concepcion

@Betterdev yes sloan's algorithm can be used. Holes should be no problem as the algorithm works on PSLG. Crossings I think should also be resoved via some preprocessing. I do not know how complicates sloan's algorithm will be to implement though. – Subdeacon

Answer 4

I previously worked on a vector graphics package, so I can't tell you how many hours I've stared at that exact "e" graphic. I eventually settled on the earcut library for triangulation of point data. It's extremely fast and much simpler compared to libraries such as libtess-2.

Cushion