A 2D Advancing-Front Delaunay Mesh Refinement Algorithm

Abstract

I present a generalization of Chew's first algorithm for Delaunay mesh refinement. I split the line segments of an input planar straight line graph (PSLG) such that the lengths of split segments are asymptotically proportional to the local feature size at their endpoints. By employing prior algorithms, I then refine the truly or constrained Delaunay triangulation of the PSLG by inserting off-center Steiner vertices of "skinny" triangles while prioritizing triangles with shortest edges first. This technique inserts Steiner vertices in an advancing front manner such that we obtain a size-optimal, truly or constrained Delaunay mesh if the desired minimum angle is less than 30° (in the absence of small input angles). This is an improvement over prior algorithms that produce size-optimal meshes with minimum angles of about 26.4°and 28.6°for truly and constrained Delaunay meshes, respectively. Even in the presence of small input angles, the upper bound on the maximum angle is an angle strictly greater than 120° (an improvement from about 137°). The lower bound on the minimum angle in the presence of small angles is identical to prior bounds.

Figure A.19:

How the value of R affects the lengths of segments from a vertex on a PSLG segment to a vertex on its adjacent segment. (a) The length of a₁b₁ must be greater than the length of oa₁, so ob₁ > 2|a₁b₁| cos ϕ, which happens when cos ϕ < 1/(2R). (b) The length of a₂m should be greater than a₁a₂. If the ratio of lengths of oa₁ and a₁a₂ is small enough, this condition is easily satisfied.

Figure 1:

The difference between truly Delaunay and constrained Delaunay meshes. Only parts of the meshes are shown for clarity. The thick horizontal lines are input segments to which the mesh should conform. A vertex a is not visible from a point b if the line segment joining a and b intersects an input line segment. (a) All triangles are Delaunay because no vertex is inside their circumcircles. (b) A vertex from the other triangle is inside the circumcircle of both triangles, but the vertex is not visible from the interior of the triangles. (c) Similar to (b), but one of the triangles is not on the input segment. (d) This case is not allowed because an end point of the input segment is inside the circumcircle of a triangle, and the end point is visible from the interior of the triangle.

Figure 2:

The relationship between the minimum angle in a triangle and the ratio of the radius of its circumcircle and the length of the shortest edge. Let ab be the shortest edge of $\Delta abc$. Then, the angle at c is its shortest angle. The triangle’s circumcircle is shown. Note that $\angle c^{\prime }=\angle c$. Also, $\angle aob=2\angle c^{\prime }$, and $\angle mob=\angle c^{\prime }$. Clearly, $\sin \angle c=\sin \angle c^{\prime }=\sin \angle mob=\frac{\left|ab\right|}{2\left|ob\right|}=\frac{l}{2r}$, where $l=\left|ab\right|$ and r is the length of the radius of the circumcircle.

Figure 3:

The off-center vertex. Instead of the circumcenter of a skinny triangle pqo, Üngör and Erten [10] add an off-center vertex into the Delaunay triangulation. The off-center vertex r lies on the perpendicular bisector of the shortest edge pq of the skinny triangle such that the edge subtends the minimum desired angle at the point. Note that r should be on the same side of pq as o is. As in their algorithm, I use the same point for Delaunay refinement if it is closer to the shortest edge than to the circumcenter.

Figure 4:

Encroachment of a Steiner vertex. A possible circumcircle of a skinny triangle is shown. Since the triangulation obeys the Delaunay property, the circumcircle cannot enclose any vertex of a PSLG subsegment pq. The circumcenter o and the skinny triangle have to be on opposite sides of pq. Thus, the skinny triangle has to be in the shaded gray region, which is fully inside the diametral circle of the PSLG subsegment.

Figure 5:

The various distance functions associated with a line segment pq in the PSLG. The blue, thick line segment pq is horizontal and can be considered as part of the x axis with vertex p being at the origin. The distance to the red feature(s) is function of x, and the function is plotted as a thin black curve. (a) The distance to a PSLG vertex is plotted as a function of x. The domain of the function is from p to q. (b) The distance to some other PSLG line segment (red) is plotted. The distance varies linearly, and the domain of the linear distance function is limited. The dashed lines are perpendicular to the PSLG line segment (red), and they limit the domain of the distance function. Beyond the domain, the distance to c or d (whichever is closer) defines the distance function on pq. Those parts of the distance functions look like the distance function in (a). (c) As p and q are not adjacent, the distance to p or q, whichever is larger, also limits the feature size at any point on the line segment pq. A disk centered at a point on pq with the radius equal to the greater of xp or xq contains both p and q. Thus, this piecewise linear function is also considered to compute the local feature size.

Figure 6:

The worst case when computing the local feature size. When there are O(n) adjacent features as shown above and O(n) nonadjacent features interleaved among the adjacent features, there are O(n²) possible pairs of features for which the distance functions have to be computed, whose lower envelope computes the local feature size.

Figure 7:

An example of reference-to-PSLG mappings M_i(t) from a reference segment T_i to a PSLG segment L_i. The reference segment is uniformly split into n subsegments, and the corresponding splits are made in the PSLG segment. The mapping function is defined such that uniform splits in the reference segment correspond to asymptotically proportional (to the local feature size) splits in the PSLG segment. Note that the reference segment and the PSLG segment may be of different lengths.

Figure 8:

Images depicting the problem with small angles in the PSLG. The thick lines are part of the PSLG. The black dots are the vertices added to the PSLG (not all are shown). (a) When the angle is obtuse, the line segment pq is longer than subsegments at p and q. (b) When the angle is acute, but greater than 60°, pq might be longer than subsegments at p and q, but it is not guaranteed unless the segments are adequately refined (see the appendix for a detailed explanation). (c) When the angle is very small, there is a good chance that pq might be shorter than the threshold for size optimality.

Figure 9:

The diametral semicircle of a PSLG subsegment pq is shown. If the condition in the Lemma 4.12 is satisfied, no vertex is placed within the area bounded by the arcs centered at p and q and the diametral semicircle. Vertices are allowed only in the shaded region. If ${\theta }^{*}>0$, the maximum possible distance between two points in the shaded region is half the length of subsegment pq.

Figure 10:

The shortest edge of a skinny triangle is pq. The arc is the locus of points at which pq subtends an angle ${\theta }^{*}$, which is the threshold for a triangle to be considered skinny. The third vertex should be outside the dashed arc pq. Since pq is the shortest side, the third vertex should be outside the circles centered at p and q with the radii equal to the length of pq. The length of the longest edge is, therefore, at least the distance between p (or q) and the point of intersection of the arc and one of the circles.

Figure 11:

The diametral circles on PSLG subsegments at a small angle ϕ. The length of the subsegment oa₁, which is closest to o on the horizontal PSLG segment oa, has been normalized and set to 1. (a) In the hypothetical limiting case when R = 1, all subsegments have an equal length. In this case, the diametral circles do not contain vertices from the adjacent segment. (b) When R > 1, there is a window for each vertex within which the vertices have to lie. The windows are shown as filled gray circles, but the vertices have to lie on the line segment. Consequently, the diametral circles have windows, too. The windows are shown with a pattern of diagonal lines. When R is sufficiently small, the windows of the vertices on a segment and the windows of diametral circles of subsegments on an adjacent segment do not overlap.

Figure 12:

A part of a PSLG is shown as thick lines. The thin, dashed lines are for reference. We assume $\left|xp\right|\le \left|xq\right|$ and that all segments in the PSLG have been split using at least one vertex.

Figure 13:

The bound on the maximum angle in the triangle across a small angle. The thick lines segment are part of the input PSLG. (a) This is a generic case that will be seen when the PSLG is not highly refined. Notice that $\left|pq\right|$ is less than any segment adjacent to p or q. (b) When the PSLG is highly refined, the end points of mesh edges tend to vertex positions as shown, i.e., they will move toward positions such that $\left|ps\right|=\left|qr\right|$. The bound on the angle is $\pi /2+\varphi /2$. In contrast to the other diagrams in this paper, no PSLG segment is horizontal in this diagram because there is symmetry along the horizontal line, which is easy to observe when the PSLG segments are rotated.

Figure 14:

A comparison of meshes produced by Triangle [23] with the existing algorithms and with my algorithm for a trapezoidal domain, where one of the parallel sides is much shorter than the other. The minimum angle in the meshes is 35°. (a) Results with no modification to PSLG (poorly graded). (b) Result after PSLG segments were split into a few subsegments (well graded). (c) Results after PSLG segments were split into more subsegments (well graded).

Figure 15:

A zoomed-in view of the meshes produced in Fig. 14 above. (a) The short segment is not split. (b) The short segment is split into four subsegments using my algorithm. (c) The short segment is split into five subsegments using my algorithm.

Figure 16:

The input PSLG. Note that there are tiny segments at the two corners. The dimension of the domain is 1000 units by 500 units. The tiny segments’ length is 1 unit.

Figure 17:

The vertex cloud produced by my algorithm for the domain in Fig. 16 above. (a) The short segment is split into one subsegment using my algorithm. (b) The short segment is split into two subsegments using my algorithm. (c) The short segment is split into three subsegments using my algorithm. (d) The short segment is split into four subsegments using my algorithm.

Figure 18:

The meshes produced by my algorithm for the domain in Fig. 16 above. (a) The short segment is split into one subsegment using my algorithm. (b) The short segment is split into two subsegments using my algorithm. (c) The short segment is split into three subsegments using my algorithm. (d) The short segment is split into four subsegments using my algorithm.