Size-optimal Steiner points for Delaunay-refinement on curved surfaces

7.9 MB·English

by  Darren Engwirda

(cid:73)
Off-centre Steiner points for Delaunay-refinement on curved surfaces
Darren Engwirdaa,∗∗, David Iversa
aSchool of Mathematics and Statistics F07, University of Sydney, NSW 2006, Australia
Abstract
An extension of the restricted Delaunay-refinement algorithm for surface mesh generation is described, where a new
point-placement scheme is introduced to improve element quality in the presence of mesh size constraints. Specifically,
it is shown that the use of off-centre Steiner points, positioned on the faces of the associated Voronoi diagram, typically
leadstosignificantimprovementsintheshape-andsize-qualityoftheresultingsurfacetessellations. Thenewalgorithm
6
can be viewed as a Frontal-Delaunay approach – a hybridisation of conventional Delaunay-refinement and advancing-
1
fronttechniquesinwhichnewverticesarepositionedtosatisfybothelementsizeandshapeconstraints. Theperformance
0
2of the new scheme is investigated experimentally via a series of comparative studies that contrast its performance with
that of a typical Delaunay-refinement technique. It is shown that the new method inherits many of the best features
n
of classical Delaunay-refinement and advancing-front type methods, leading to the construction of smooth, high quality
u
Jsurface triangulations with bounded radius-edge ratios and convergence guarantees. Experiments are conducted using a
7 range of complex benchmarks, verifying the robustness and practical performance of the proposed scheme.
2
Keywords: Surface mesh generation, Delaunay-refinement, Advancing-front, Frontal-Delaunay, Off-centre Steiner

]points, Element quality
G
C
.1. Introduction 1.1. Related Work
s
c Mesh generation is a broad and evolving area of re-
Surface mesh generation is a key component in a va-
[ search. A number of algorithms have been developed to
riety of computational modelling and simulation tasks,

3 tacklevariousmeshingtasks,includinggrid-overlay meth-
including many forms of computational engineering and
v ods[2,3]inwhichameshisgeneratedbywarpingandre-
numerical modelling, a range of problems in computer
2 fining a semi-structured background grid, advancing-front
0graphics and animation, and various data visualisation
techniques [4, 5, 6, 7] which grow a mesh incrementally in
0applications. Given a geometric domain described by a
a layer-wise fashion, and Delaunay-based strategies [8, 9,
4bounding surface Σ embedded in R3, the surface trian-
0 10,11,12,13,14,15,16,17]whichincrementallyrefinean
gulation problem consists of tessellating Σ into a mesh
. initially coarse Delaunay triangulation.
1of non-overlapping triangular elements, such that all geo-
TheDelaunaytriangulation[18]isimbuedwithanum-
0metrical, topological and user-defined constraints are sat-
5 ber of attractive theoretical properties. In the context of
isfied. While each use-case contributes its own set of spe-
1 mesh generation, for triangulations in R2 and surface tes-
cific considerations, it is typical to require that surface
v:tessellations: (i) consist of elements of high shape-quality, sellationsembeddedinR3,itisknownthatsuchstructures
i are topologically optimal [19], maximising the worst-case
X(ii) provide good geometrical and topological approxima-
element quality in the resulting mesh. The reader is re-
tions to the underlying surface Σ, and (iii) satisfy a set
r ferredto[19]foradditionaldetails. Asthenamesuggests,
aof user-specified element sizing constraints. While various
Delaunay-refinementschemesarebasedontheincremental
strategies have been developed to solve the surface mesh-
refinement of an initially coarse bounding Delaunay tes-
ing problem in the past, a new algorithm is developed in
sellation. At each step of the algorithm, elements that vi-
this study with the aim of improving the quality of the
olate a set of geometrical, topological or user-defined con-
resulting triangulations.
straints are identified and the worst offending element is
eliminated. This is achieved through the insertion of an
(cid:73)A short version of this paper appears in the proceedings of the additional Steiner vertex at the refinement point of the
23rdInternationalMeshingRoundtable[1].
element – either the circumcentre of the element itself, or
∗Correspondingauthor. Tel.: [phone].
∗∗Currentaddress: DepartmentofEarth,AtmosphericandPlane- a point on an adjacent segment of the bounding geome-
tarySciences,Room54-1517,MassachusettsInstituteofTechnology, try. The original planar Delaunay-refinement methods of
77MassachusettsAvenue,Cambridge,MA02139-4307. Chew [8] and Ruppert [9, 10] have since been extended to
Email addresses: [email protected];
support volumetric tessellations [12] and surface triangu-
[email protected](DarrenEngwirda),
[email protected](DavidIvers) lations [15, 16, 14].
Preprint submitted to Computer-Aided Design June 28, 2016
Delaunay-basedmethodshavepreviouslybeenapplied and numerical simulation.
tothesurfacemeshingproblem,viatheso-calledrestricted AbriefreviewofthetraditionalsurfaceDelaunay-refine-
Delaunay paradigm originally introduced by Edelsbrun- ment scheme is presented in Section 2, along with a dis-
ner and Shah [20]. Given a surface Σ embedded in R3 cussion of the underlying theoretical constructs, including
and a set of points X ∈ Σ, the restricted Delaunay sur- the restricted Delaunay tessellation. The new Frontal-
face triangulation Del| (X) is a sub-complex of the full- Delaunay algorithm is presented in Section 3, focusing
Σ
dimensional Delaunay tessellation Del(X), containing the in-detail on the off-centre point-placement scheme pro-
set of 2-faces f ∈ Del(X) that approximate the underly- posed in this work. In Section 5, a detailed comparison
ing surface. Boissonnat and Oudot [15, 16] have shown between the conventional Delaunay-refinement and pro-
that when the point-wise sampling X is sufficiently dense posed Frontal-Delaunay schemes is presented, contrasting
andwell-distributed,therestrictedDelaunaytriangulation output quality and computational performance.
Del| (X)isanaccuratepiecewiseapproximationofΣ,in-
Σ
corporatingseveraltheoreticalguaranteesoffidelity. Given
2. Restricted Delaunay-refinement Techniques
such behaviour, a Delaunay-refinement algorithm for sur-
facemeshgenerationproceedsaspertheplanarcase,with Delaunay-refinement algorithms operate by incremen-
an initially coarse restricted Delaunay triangulation, in- tally adding new Steiner vertices to an initially coarse re-
duced by a sparse sampling X ∈ Σ, being incrementally stricted Delaunay triangulation of the underlying surface.
refined via the introduction of new Steiner vertices posi- Contrarytotypicalplanaralgorithms[8,9,10],refinement
tionedonthesurfaceΣ. Thisprocessisdiscussedindetail schemes for surface meshing are designed not only to en-
in Section 2. sure that the resulting mesh satisfies element shape and
Advancing-fronttechniquesareanalternativeapproach size constraints, but that the geometry and topology of
tomeshgeneration, inwhichelementsarecreatedone-by- the mesh itself is an accurate piecewise approximation of
one, starting from a set of frontal facets initialised on the the underlying surface.
bounding geometry. The mesh is marched inwards layer- BeforediscussingthedetailsofDelaunay-refinemental-
by-layer, with new elements created by carefully position- gorithms for surface meshing problems, a number of im-
ingverticesadjacenttofacetsinthefrontalset. Following portant theoretical concepts are introduced:
the placement of new elements, the set of frontal facets
is updated to incorporate the new mesh topology. Ele- Definition 1. Let Σ be a smooth surface embedded in
mentsareaddedincrementallyinthisfashionuntilacom- R3, enclosing a volume Ω ⊂ R3. Let Del(X) be a full-
pletetessellationofthedomainisobtained. Suchmethods dimensional Delaunay tessellation of a point-wise sample
are applicable to surface meshing problems [5, 6, 7], al- X ⊆ Σ and Vor(X) be the Voronoi diagram associated
though several non-trivial issues of robustness are known with Del(X). The restricted Delaunay surface triangula-
to exist [7]. Fundamentally, advancing-front techniques tion Del| (X) is a sub-complex of Del(X) including any
Σ
are heuristic in nature, and do not typically offer theo- 2-face f ∈ Del(X) associated with an edge v ∈ Vor(X)
f
retical guarantees of convergence or bounds on element suchthatv ∩Σ(cid:54)=∅. Therestricted Delaunay volume tri-
f
quality. Nonetheless, when such methods are successful, angulation Del| (X)isasub-complexofDel(X)including
Ω
they are known to produce very high-quality tessellations. any 3-simplex τ ∈Del(X) associated with an internal cir-
Due to their theoretical robustness, Delaunay-refine- cumcentre c∈Ω.
ment schemes are a popular choice for practical meshing
algorithms. On the other hand, it is well known that It has been shown by a number of authors, including
optimised advancing-front type methods can often pro- Amenta and Bern [21], Cheng, Dey, Levine [22], Cheng,
duce meshes with superior element quality characteris- Dey and Ramos [14], Cheng, Dey, Ramos and Ray [23],
tics, though their heuristics can sometimes break-down BoissonnatandOudot[15,16]andCheng,DeyandShewchuk
in practice. In this study, a new Frontal-Delaunay al- [19], that given a sufficiently dense point-wise sampling
gorithm is presented which seeks to combine the bene- of the surface, X ∈ Σ, the restricted Delaunay tessella-
fits of classical Delaunay-refinement and advancing-front tionDel|Σ(X)isguaranteedtobebothgeometricallyand
typeapproaches. Whilesuchstrategieshavebeenpursued topologically representative of the underlying surface Σ.
previously in the case of planar and/or volumetric refine- Specifically, it is known that, when the restricted tessel-
ment algorithms, it is believed that the present study is lation is a so-called loose (cid:15)-sample [15, 16], Del|Σ(X) is
the first application of such ideas to the restricted surface homeomorphic to the surface Σ, the Hausdorff distance
meshing problem directly. The new algorithm is designed H(Del|Σ(X),Σ) is small and that Del|Σ(X) provides a
to produce smooth, high-quality triangulations consistent good piecewise approximation of the geometrical prop-
with an advancing-front scheme, whilst also inheriting the erties of Σ, including its normals, area and curvature.
theoretical robustness of Delaunay-based techniques. It Similartheoreticalguaranteesextendtotheassociatedre-
is expected that this new algorithm may be of interest to stricted volume tessellation Del|Ω(X). The properties of
userswhoplaceahighpremiumonmeshquality,including restricted Delaunay tessellations are well documented in
those operating in the areas of computational engineering the literature, and a full account is not given here. The
2
Figure 1: Restricted tessellations for a curved domain in R2, showing (i) the curve Σ and enclosed area Ω, (ii) the Delaunay tessellation
Del(X)andVoronoidiagramVor(X),and(iii)therestrictedcurveandareatessellationsDel|Σ(X)andDel|Ω(X).
(i) (ii) (iii)
reader is referred to the detailed expositions presented in where R is the radius of the diametric ball of τ and (cid:107)e (cid:107)
0
[15, 16] and [19] for further details and proofs. is the length of its shortest edge.
Definition 2. LetDel| (X)bearestrictedDelaunaytri-
Σ The radius-edge ratio is a measure of element shape-
angulation of a smooth surface Σ. Given a 2-simplex √
quality. It achieves a minimum, ρ(τ) = 1/ 3, for equilat-
f ∈ Del| (X), any circumball of f centred on the sur-
Σ eral triangles and increases toward +∞ as elements tend
face Σ is known as a surface Delaunay ball of f, denoted
toward degeneracy. For 2-simplexes, the radius-edge ratio
SDB(f).
is robust and can be related to the minimum plane angle
Surface Delaunay balls are centred at the intersection θminbetweenadjacentedges,suchthatρ(τ)= 12(sin(θmin))−1.
oftheVoronoidiagramVor(X)withthesurfaceΣ. Specif- Duetothesummationofanglesinatriangle,givenamin-
ically, each 2-face f ∈ Del|Σ(X) is associated with an or- imumangleθmin thelargestangleθmax isclearlybounded,
thogonaledgeintheVoronoidiagramvf ∈Vor(X),which such that θmax ≤180◦−2θmin.
isguaranteed,bydefinition,tointersectwiththesurfaceΣ
at least once. Insomecases, especiallywhenthesampling 2.1. An Existing Delaunay-refinement Algorithm
X ∈ Σ is relatively coarse, there may be multiple surface Thedevelopmentofprovably-good Delaunay-refinement
intersectionsassociatedwithagiven2-facef ∈Del|Σ(X). schemes for surface-based mesh generation is an ongoing
In such cases, it is typical to consider only the largest as- area of research. An algorithm for the meshing of closed
sociated surface ball. See Figure 2 for an illustration of 2-manifolds embedded in R3, adapted largely from the
the surface Delaunay ball for a general facet f. methods presented by Boissonnat and Oudot in [15, 16]
is presented here. This method is largely equivalent to
Definition 3. LetDel| (X)bearestrictedDelaunaysur-
Σ the cgalmesh algorithm, available as part of the cgal
face triangulation of a smooth surface Σ. Given a 2-
package, and summarised by Jamin, Alliez, Yvinec and
simplex f ∈ Del| (X), the surface discretisation error
Σ
Boissonnat in [17]. A similar algorithm is also described
(cid:15)(f) is the Euclidean distance between the centres of the
by Cheng, Dey and Shewchuk in [19] and Dey and Levine
largest surface Delaunay ball of f and its diametric ball.
in [24], while a variant developed for the reconstruction
The surface discretisation error is a measure of how ofmedicalimagesispresentedbyFoteinos,Chernikovand
well the restricted triangulation Del|Σ(X) approximates Chrisochoidesin[25]. Thealgorithmpresentedinthissec-
the underlying surface Σ geometrically. Considering a 2- tionisreferredtoasthe‘conventional’Delaunay-refinement
face f ∈Del|Σ(X), the surface discretisation error (cid:15)(f) is approachthroughout,duetoitsdirectuseofcircumcentre-
a measure of the distance between the face f and the fur- based Steiner vertices.
thest adjacent point on the surface Σ. This measure can AsperJaminetal.[17],thedevelopmentoftheconven-
be thought of as a one-sided discrete Hausdorff distance, tionalDelaunay-refinementalgorithmisgeometry-agnostic,
(cid:15)(f) = H(Del|Σ(X),Σ), defined from a set of representa- being independent of the specific description used for the
tivepointsonthetriangulationDel|Σ(X)tothesurfaceΣ. underlying surface Σ. It is required only that the geome-
Clearly, if the triangulation Del|Σ(X) is a good piecewise trysupportaso-calledoracle predicatethatcanbeusedto
representation of the surface Σ, the surface discretisation compute the intersection of a given line segment with the
error (cid:15)(f) should be small. See Figure 2 for a description surface Σ. The Frontal-Delaunay algorithm presented in
ofthesurfacediscretisationerrorforafacetf ∈Del|Σ(X). Section3,additionallyrequiresthattheoraclecomputein-
tersectionsbetweenadiskofagivenradiusandthesurface
Definition 4. Given a d-simplex τ, its radius-edge ratio,
Σ. Suchconstructionswillbediscussedinfurtherdetailin
ρ(τ), is given by:
Section3. Whileaverygeneralclassofsurfacedescription
ρ(τ)=R/(cid:107)e0(cid:107), (1)
3
Figure2: ThesurfaceDelaunayballforarestricted2-facef ∈Del|Σ(X),showing(i)placementofthesurfaceballattheintersectionofthe
dualedgev ∈Vor(X)andthesurfaceΣ,and(ii)associatedradiusr(f)andsurfacediscretisationerror(cid:15)(f).
f
(i) (ii)
is supported at the theoretical level, in this study, atten- hedral representation of the surface Σ. Care is taken to
tionisrestrictedtothedevelopmentofso-calledremeshing ensure that each connected component in P is sampled.
operations, in which the underlying surfaces Σ are speci- Alternative initialisation methods making use of robust
fied as existing 2-manifold triangular complexes P. This persistent-triangles techniques are described in [16, 19].
restriction is made in the present study for convenience In the next step, the initial triangulation objects are
only, and future work is intended to focus on more gen- formed. In the present work, the full-dimensional De-
eral surface descriptions, including implicit and analytic launay tessellation, Del(X), is built using an incremental
functions, in addition to those that contain sharp creases. Delaunay triangulation algorithm, based on the Bowyer-
Following Jamin et al. [17], the Delaunay-refinement Watson technique [26]. The restricted surface and vol-
algorithm takes as input a surface domain, described by umetric triangulations, Del| (X) and Del| (X), are de-
Σ Ω
a smooth 2-manifold Σ, an upper bound on the allowable rived from Del(X) by explicitly testing for intersections
element radius-edge ratio ρ¯, a mesh size function h¯(x) de- between edges in the associated Voronoi diagram Vor(X)
finedatallpointsonthesurfaceΣandanupperboundon and the surface Σ. These queries are computed efficiently
the allowable surface discretisation error (cid:15)¯(x). The input by storing the surface definition P in a spatial-tree.
surface Σ encloses a bounded volume Ω. The algorithm The main loop of the algorithm proceeds to incremen-
returns a triangulation T| of the surface Σ, where T| tally refine any 2-faces f ∈ Del| (X) that violate either
Σ Σ Σ
is a restricted Delaunay surface triangulation of a point- the radius-edge, element size or surface discretisation re-
wise sampling X ∈ Σ, such that T| = Del| (X). As a quirements. The refinement process is priority scheduled,
Σ Σ
by-product, the algorithm also returns a coarse triangu- with triangles f ∈ Del| (X) ordered according to their
Σ
lation T| of the enclosed volume Ω, where T| is a re- radius-edge ratios ρ(f), ensuring that the element with
Ω Ω
stricted Delaunay volume triangulation T| = Del| (X). the worst ratio is refined at each iteration. Individual el-
Ω Ω
Both Del| (X) and Del| (X) are sub-complexes of the ements are refined based on their surface Delaunay balls,
Σ Ω
full-dimensional Delaunay tessellation Del(X). Note that with a triangle f ∈ Del(X) eliminated by inserting the
Del| (X) is a triangular complex, while Del| (X) and centreofthelargestballB(c,r) =SDB(f)intothetes-
Σ Ω max
Del(X) are tetrahedral complexes. The method is sum- sellationDel(X). Thisprocessis adirectgeneralisation of
marised in Algorithm 2.1. the circumcentre-based insertion method associated with
TheDelaunay-refinementalgorithmguarantees,firstly, Ruppert’s algorithm for planar domains [9, 10].
that all elements in the output surface triangulation f ∈ As a consequence of changes to the full-dimensional
T| satisfy element shape constraints, ρ(f) ≤ ρ¯, element tessellation following the insertion of a new Steiner ver-
Σ
size constraints h(f) ≤ h¯(x ) and surface discretisation tex,correspondingupdatestotherestrictedtriangulations
f
bounds (cid:15)(f) ≤ (cid:15)¯(x ). Furthermore, for sufficiently small Del| (X) and Del| (X) are instigated, ensuring that all
f Σ Ω
mesh size functions h¯(x), the surface triangulation T| tessellation objects remain valid throughout the refine-
Σ
is guaranteed to be a good piecewise approximation of ment process. The Delaunay-refinement algorithm termi-
the underlying surface Σ, exhibiting both geometrical and nates when all 2-faces f ∈ Del| (X) satisfy all radius-
Σ
topological convergence as h¯(x) → 0, consistent with the edge,sizeandsurfacediscretisationthresholds,suchthat1
characteristicsofrestrictedDelaunaytessellationsoutlined ρ(f) ≤ ρ¯, h(f) ≤ αh¯(x ) and (cid:15)(f) ≤ (cid:15)¯(x ), respectively,
f f
in Definition 1. where the element size h(f) is proportional to the radius
TheDelaunay-refinementalgorithmbeginsbycreating of the associated surface ball B(c,r) = SDB(f), such
max
an initial point-wise sampling of the surface X ∈ Σ. Ex-
0
ploiting the discrete representation available for Σ, in this
study the initial sampling is obtained as a well-distributed 1Thecoefficientα=4/3,ensuringthatthemeanelementsizedoes
not,onaverage,undershootthedesiredtargetsizeh¯(x ).
subsetoftheexistingverticesY ∈P,whereP isthepoly- f
4
√
that2 h(f)= 3r,andthetargetsizeh¯(x )issampledat in [28], developed a planar Frontal-Delaunay algorithm in
f
the centre of the surface ball x =c. which new vertices are positioned along edges of the as-
f
sociated Voronoi diagram. Rebay showed that new ver-
ticescanbepositionedonVor(X)accordingtoana priori
3. Restricted Frontal-Delaunay Methods
meshsizefunctionh¯(x),astrategybroadlyconsistentwith
Frontal-Delaunayalgorithmsareahybridisationofadv- conventional advancing-front techniques. While his algo-
ancing-frontandDelaunay-refinementtechniques,inwhich rithm maintains a Delaunay triangulation T = Del(X)
a Delaunay triangulation is used to define the topology of the current vertex set X, it is still fundamentally an
of a mesh while new Steiner vertices are inserted in a advancing-front scheme – bereft of guarantees on the ele-
mannerconsistentwithadvancing-frontmethodologies. In ment shape quality ρ(τ). Rebay reported that his scheme
practice, such techniques have been observed to produce produced very high-quality output, typically outperform-
very high-quality meshes, inheriting the smooth, semi- ing conventional Delaunay-refinement in practice.
structuredvertexplacementofpureadvancing-frontmeth- U¨ng¨or and Erten have approached the problem from
odsandtheoptimalmeshtopologyofDelaunay-basedap- the flip-side, introducing the notion of generalised Steiner
proaches. While Frontal-Delaunay methods have previ- vertices for planar Delaunay-refinement methods. Their
ouslybeenusedbyarangeofauthorsinthecontextofpla- ‘off-centre’ vertices are points lying along edges in the as-
nar, volumetric and parametric surface meshing, includ- sociated Voronoi diagram, as per Rebay. In [34], U¨ngor
ing, for example studies by U¨ng¨or and Erten [27], Rebay showedthatwhenrefiningatriangleτ,itsoff-centreshould
[28],Mavriplis[29],Frey,BorouchakiandGeorge[30],and bepositioned,ifpossible,suchthatthenewtriangleσ ad-
Remacle, Henrotte, Carrier-Baudouin, B´echet, Marchan- jacenttotheshortestedgeinτ satisfiesthedesiredshape-
dise,GeuzaineandMouton[31],theauthorsarenotaware quality target, such that ρ(σ) ≤ ρ¯. Importantly, U¨ng¨or
ofanypreviousinvestigationsdescribingtheapplicationof demonstrated that such a strategy typically leads to an
such techniques to the surface meshing problem directly. improvement in the performance of Delaunay-refinement
The conventional advancing-front method, on the other in practice, reducing the size of the output |T|. U¨ng¨or
hand, has been generalised to support surface meshing, extended the guarantees derived for Ruppert’s Delaunay-
as, for example, outlined in studies by Rypl [32, 6] and refinementalgorithmtohisoff-centretechniqueandshowed
Schreiner, Scheidegger, Fleishman and Silva [7, 33]. that such bounds are typically improved upon in practice.
Inthesepreviousstudies,theconventionalplanaradva- Off-centrerefinementisbasedonRuppert’sDelaunay-refine-
ncing-front methodology is extended directly to facilitate ment framework directly, involving modifications to the
surface operations. Meshing proceeds via the incremen- procedure used to refine low-quality triangles only.
tal introduction of a well-distributed set of vertices X The use of generalised Delaunay-refinement strategies
positioned on the surface X ∈ Σ. It should be noted has also been explored by Chernikov, Chrisochoides and
that, contrary to the restricted Delaunay techniques in- Foteinosin [35,36], inwhich Steinerpoints arepositioned
troduced in previous sections, advancing-front methods within a set of so-called selection-balls adjacent to ele-
typically maintain a partial, and possibly non-Delaunay mentcircumcentres. In[35], ChernikovandChrisochoides
manifold triangular complex T| only – they do not con- show that a family of ‘provably-good’ generalised two-
Σ
struct a full-dimensional tessellation T| . It should also and three-dimensional Delaunay-refinement schemes ex-
Ω
be noted that, in addition to the usual limitations associ- ist, and can be realised via the specification of appropri-
ated with advancing-front strategies, serious issues of ro- ate selection-ball radii parameters. It is unclear whether
bustness often afflict these techniques in practice due to this selection-ball formalism incorporates the full set of
the difficulties associated with the reliable evaluation of Voronoi-type off-centres currently in use, including those
therequisitegeometricintersectionandoverlappredicates introducedbyRebay[28],andErtenandU¨ng¨or[27],inad-
for sets of non-planar elements. Additionally, for highly dition to the strategies defined in the present study. Such
curved and/or poorly separated surface definitions it is a determination may form the basis for future work.
known to be difficult to ensure the topological correctness
of the output mesh T| . Schreiner et al. [7, 33] introduce 3.2. Point-placement Preliminaries
Σ
a number of heuristic techniques in an effort to overcome In this study, a generalisation of the ideas introduced
these difficulties. by Rebay and U¨ng¨or is formulated for the surface mesh-
ing problem – using off-centre Steiner vertices to simulate
3.1. Off-centres thevertexplacementstrategyofaconventionaladvancing-
The new Frontal-Delaunay algorithm is based primar- front approach, while also preserving the framework of
ily on a generalisation of ideas introduced by Rebay, who, a restricted Delaunay-refinement technique. The aim of
such a strategy is to recover the high element qualities
√ and smooth, semi-structured point-placement generated
2The coefficient 3 represents the mapping between the edge by frontal methods, while inheriting the theoretical guar-
length and diametric ball radius for an equilateral element. Such
antees of Delaunay-refinement methods. Advancing-front
scalingensuresthatsizeconstraintsareappliedwithrespecttomean
edgelength. algorithmstypicallyincorporateameshsizefunctionh¯(x),
5
Algorithm 2.1 Surface mesh generation by restricted Delaunay-refinement
1: functionDelaunaySurface(Σ,Ω,ρ¯,(cid:15)¯(x),h¯(x),T|Σ,T|Ω) 1: functionRefineSimplex(f) (cid:46){surface refinement}
2: FormaninitialsamplingX∈ΣsuchthatX 2: CallSurfaceDelaunayBall(f,B(c,r)max).
iswell-distributed onΣ. 3: FormnewSteinervertexpaboutB(c,r)max.
3: FormDelaunaytessellationDel(X). 4: InsertSteinervertexX←p,updateDel(X)←X.
4: FormrestrictedobjectsDel|Σ(X)andDel|Ω(X). 5: endfunction
5: Enqueue all restricted 2-simplexes Q|Σ ← 1: functionSurfaceDelaunayBall(f) (cid:46){surface ball}
f ∈ Del|Σ(X). A 2-simplex f is enqueued 2: FormVoronoiedgevf orthogonalto2-simplexf.
if BadSimplex(f)returnstrue. 3: FormthesetofassociatedsurfaceDelaunay
6: while(Q|Σ(cid:54)=∅)do (cid:46){main refinement sweeps} ballsB(c,r)ifortherestricted2-simplexf ∈
7: CallRefineSimplex(f ←Q|Σ) Del|Σ(X). Ballsarecentredaboutthesetof
8: Update the restricted Delaunay tessel- surfaceintersectionsvf ∩Σ(cid:54)=∅.
lationsDel|Σ(X)andDel|Ω(X). 4: returnB(c,r)max,wherermax ismaximal.
9: UpdateQ|Σ toreflectchangestoDel|Σ(X). 5: endfunction
10: endwhile 1: functionBadSimplex(f) (cid:46){termination criteria}
11: return T|Σ=Del|Σ(X)andT|Ω=Del|Ω(X) 2: return ρ(f)>ρ¯or(cid:15)(f)>(cid:15)¯(xf)orh(f)>h¯(xf)
12: endfunction 3: endfunction
a function f : R3 → R+ defined over the domain to be generalised off-centreframeworkintroducedbyU¨ng¨or,the
meshed,whereh¯(x)representsthedesirededgelength(cid:107)e(cid:107) ‘ideal’ location of the size-optimal off-centre c(2) is based
at any point x ∈ Σ. This mesh size function typically in- onaconsiderationoftheisoscelestriangleσ formedabout
corporates size constraints dictated by both the user and theshort‘frontal’edgee ∈f. Thepointc(2)ispositioned
0
the geometry of the domain to be meshed. The construc- to ensure that σ satisfies local size constraints.
tion of appropriate mesh size functions will be discussed Given a refinable 2-simplex f ∈ Del| (X), the size-
Σ
insubsequentsections,butfornow,itsufficestonotethat optimalTypeIIvertexc(2) ispositionedatanintersection
h¯(x) is a g-Lipschitz function defined at all points on Σ of the surface Σ and a plane V, where V is aligned with
with 0<g <1. the local face of the Voronoi diagram Vor(X) associated
The proposed Frontal-Delaunay algorithm is an ex- with the frontal edge e ∈ f. The plane is positioned
0
tension of the restricted Delaunay-refinement algorithm such that it passes through three local points on Vor(X):
presented in Section 2, modified to use off-centre rather the midpoint of the frontal edge e ∈ f, the centre of the
0
than circumcentre-based refinement strategies. The basic diametric ball of f and the centre of the surface Delaunay
framework is consistent with the Delaunay-refinement al- ball B(c,r) = SDB(f). The vertex c(2) is positioned
max
gorithm described previously, in which an initially coarse such that it forms an isosceles triangle candidate σ about
restricted Delaunay triangulation of a surface Σ is refined the frontal edge e , such that its size h(σ) satisfies local
0
throughtheintroductionofadditionalSteinerverticesX ∈ constraints. Specifically, the altitude of σ is computed
Σ until all constraints are satisfied. A restricted surface from local mesh-size information, such that
triangulation T| = Del| (X) is constructed as a sub-
complex of the fuΣll-dimensΣional tessellation Del(X) and a aσ =min(cid:0)(cid:0)h¯2σ−(cid:107)12e0(cid:107)2(cid:1)21,√3/2h¯σ(cid:1) (2)
coarse restricted volumetric tessellation T|Ω = Del|Ω(X) h¯ = 1(cid:0)h¯(m )+h¯(m )(cid:1) (3)
is also available as a by-product. The constraints satis- σ 2 1 2
fitheodsebydetshceribFerdonptarle-vDioeulasulynafyoraltghoeritrhesmtriacrteedidDenetlaicuanlatyo- where the mi’s are the edge midpoints and √3/2h¯σ is the
altitude for an ‘ideal’ element.
refinement scheme, with upper bounds on the radius-edge
The position of thepoint c(2) is calculatedbycomput-
ratio ρ¯, surface discretisation error (cid:15)¯(x ) and element size
f ingtheintersectionofthesurfaceΣwithacircleofradius
h¯(x ) all required to be satisfied for convergence. See Al-
f a , centred at the midpoint of the frontal edge e ∈ f
gorithm 2.1 for a detailed outline of the method. σ 0
and inscribed on the plane V. In the case of multiple in-
tersections, the candidate point c (2) of closest alignment
3.3. Point-placement Strategy i
to the frontal direction vector d is selected. Specifi-
f
Given a surface facet f ∈ Del|Σ(X) marked for re- cally, the point ci(2) that maximises the scalar product
finement, the new Steiner vertex introduced to eliminate (cid:0)c (2) −m (cid:1)·d is chosen, where m is the midpoint of
i 0 f 0
f is an off-centre, constructed to adhere to local size and the short edge e ∈ f and the frontal direction vector d
0 f
shapeconstraints. Theoff-centresintroducedinthisstudy istakenfromthemidpointm tothecentreofthesurface
0
involve the placement of two distinct kinds of Steiner ver- ball B(c,r) =SDB(f).
tices. Type I vertices, c(1), are equivalent to conventional For nonm-uanxiform h¯(x), expressions for the position of
element circumcentres, and are used to satisfy constraints the point c(2) are weakly non-linear, with the altitude a
σ
on the element radius-edge ratios. Type II vertices, c(2), dependingonanevaluationofthemeshsizefunctionatthe
are size-optimal points, designed to satisfy element siz- edge midpoints h¯(m ) and visa-versa. In practice, since
i
ing constraints in a locally optimal fashion. Adopting the h¯(x) is guaranteed to be Lipschitz smooth, a simple it-
6
Figure 3: Placement of off-centre vertices, showing (i) the surface ball associated with a facet f ∈Del|Σ(X), (ii) the plane V aligned with
thelocalfacetofVor(X)associatedwiththefrontaledgee0∈f,(iii)placementofthesize-optimalvertexc(2).
(i) (ii) (iii)
erative predictor-corrector procedure is sufficient to solve un-converged elements in the mesh is a common feature
these expressions approximately. The positioning of size- of Frontal-Delaunay algorithms, with similar approaches
optimal Type II Steiner vertices is illustrated in Figure 3. used in [28, 29, 30, 31].
Usingthesize-optimalTypeIIpointc(2)andtheTypeI
point c(1), the final position of the refinement point c for 3.5. Theoretical Guarantees
the facet f is calculated. The point c is selected to satisfy Boissonnat and Oudot [15, 16] have previously shown
the limiting local constraints, setting that the conventional restricted Delaunay-refinement al-
(cid:40) c(2), if (cid:0)d(2) ≤d(1)(cid:1) and (cid:0)d(2) ≥ 1(cid:107)e (cid:107)(cid:1), gorithim presented in Section 2 is guaranteed to: (i) ter-
c= 2 0 (4) minate in a finite number of steps, (ii) produce elements
c(1), otherwise of bounded radius-edge ratios ρ(f) ≤ ρ¯, (iii) satisfy non-
uniform user-defined mesh-sizing and surface error con-
where d(1) = (cid:107)c(1) − m (cid:107) and d(2) = (cid:107)c(2) − m (cid:107) are straints h(f) ≤ h¯(x ), (cid:15)(f) ≤ (cid:15)¯(x ), and (iv) provide
0 0 f f
distances from the midpoint of the frontal edge e to the a good geometrical and topological approximation to the
0
Type I and Type II vertices, respectively. underlying surface Σ when the mesh size function h¯(x) is
The selection criteria is designed to ensure that the sufficiently small.
refinement scheme smoothly degenerates to that of a con- The behaviour of the new Frontal-Delaunay algorithm
ventional circumcentre-based Delaunay-refinement strat- can be assessed using a similar framework.
egy in limiting cases, while using a locally shape-optimal
approach where possible. 3.6. Termination & Convergence
Specifically, the condition d(2) ≤ d(1) guarantees that
Firstly, the termination and convergence of the algo-
c lies no further from the frontal edge e than the cen-
0 rithm is analysed:
tre of the surface Delaunay ball SDB(f). Such behaviour
ensures that a conventional circumcentre-based scheme is Proposition 1. Given a smooth surface Σ, a radius-edge
selected when the element size becomes sufficiently small. threshold ρ¯ ≥ 1 and positive mesh-size and surface-error
Additionally, the condition d(2) ≥ 1(cid:107)e (cid:107) ensures that the functions h¯(x) ≥ h¯ , (cid:15)¯(x) ≥ (cid:15)¯ , h¯ ,(cid:15)¯ ∈ R+ the Frontal-
2 0 0 0 0 0
size-optimalschemeisonlychosenwhentheedgee issuf- Delaunay algorithm is guaranteed to terminate in a finite
0
ficiently small compared to the local mesh size function. number of steps.
This behaviour ensures that e is a good ‘frontal’ edge
0
Before seeking to prove Proposition 1 itself, a number
candidate.
of intermediate results are first established.
3.4. Refinement Order Lemma 1. Let D ⊂ Rd be a bounded domain. Given a
In addition to the use of ‘off-centre’ point-placement subsetX ⊂D,satisfying(cid:107)u−v(cid:107)≥γ forallpairsu,v∈X
schemes, the Frontal-Delaunay algorithm also introduces and some scalar γ ∈ R+, the size of X is bounded, such
changes to the order in which elements are refined. To that |X|≤µ for some constant µ∈Z+.
better mimic the behaviour of an advancing-front type
This so-called packing-lemma states that any point-
method, elements are refined only if they are adjacent to
set X in a bounded domain Rd for which there exists a
anexisting‘frontal’entity. Inthecaseofsurfacefacetsf ∈
positive minimum separation distance between the points
Del| (X),theshortfrontaledgee ∈f mustbesharedby
Σ 0
in X must, consequently, be finite. As such, if it can be
at least one adjacent facet f ∈ Del| (X) that has ‘con-
j Σ
shownthatameshingalgorithmpreservessuchaminimum
verged’ – satisfying its associated radius-edge, mesh-size
separation length between its vertices, Lemma 1 can be
andsurface-errorconstraints. Theideaofdefining‘frontal’
invoked to prove that such a vertex-set is finite, and, as
entities as a dynamic boundary between converged and
7
a result, that termination of the algorithm is guaranteed. where H is the desired edge-length, computed from lo-
Lemma 1 is stated here without proof – interested readers cal mesh-size constraints (2)–(3). This expressions can be
are referred to [19] for full details. further simplified, leading to a simple inequality for the
BeforeexaminingthebehaviouroftheFrontal-Delaunay selection of Type II vertices:
algorithm as a whole, the impact of a single refinement
operation is analysed. Firstly, the behaviour of the vari- H ≤rk (9)
ous off-centre point-placement schemes introduced in Sec-
It can be seen that the Type II point-placement scheme
tion 3.3 is examined:
is declined when the surface ball associated with a given
Proposition 2. Let the threshold on radius-edge ratios be facet f is sufficiently small compared to the local target
ρ¯ ≥ 1 and Del|Σ(Xk) be the restricted surface triangula- edge-length H.
tion induced after k refinement operations. Given a low-
quality surface facet f ∈ Del| (X ), for which ρ(f) ≥ ρ¯, Using Propositions 2 and 3, the off-centre point-place-
Σ k
the minimum edge length (cid:107)e (cid:107) ∈ Del| (X ) is not de- mentstrategiesutilisedintheFrontal-Delaunayalgorithm
0 Σ k
creased following the insertion of a new Type I Steiner are shown to constitute a ‘safe’ set of refinement opera-
vertex at the centre of SDB(f). tions,withtheinsertionofTypeIverticesleadingtoanon-
decreasing minimum edge-length and use of the Type II
Proof. Recalling that the surface facet f is locally Delau- scheme declined once the mesh is sufficiently well refined.
nay,theinsertionofanewTypeIvertexck atthecentreof Using these results, termination of the full Frontal-Del-
the associated surface Delaunay ball B(ck,rk) = SDB(f) aunay algorithm can be re-visited:
isguaranteedtocreateedgesnoshorterthantheradiusof
the ball of f: Proof of Proposition 1. The Frontal-Delaunay algo-
rithm refines any surface facet f ∈Del| (X) if: (i) it is of
Σ
(cid:107)e (cid:107) ≥r . (5)
0 k+1 k poor shape quality, such that ρ(f)≥ρ¯, (ii) it is too large,
such that h(f) ≥ αh¯(x ), or (iii) it violates the local sur-
Here (cid:107)e (cid:107) is the length of the shortest edge created f
0 k+1 facediscretisationerrorthreshold, suchthat(cid:15)(f)≥(cid:15)¯(x ).
by the insertion of the new vertex c . Dividing by the f
k
current minimum edge length, and using the definition of (i) Using Proposition 3, it is known that the introduc-
the radius-edge ratio: tion of Type II off-centres is declined once the mesh
is sufficiently well refined. Type I points alone can
(cid:107)e (cid:107) r
0 k+1 ≥ k ≥ρ¯. (6) therefore be relied upon to satisfy element shape-
(cid:107)e (cid:107) (cid:107)e (cid:107)
0 k 0 k constraints.
Clearly, when ρ¯≥1 the existing minimum length (cid:107)e (cid:107) is
0 k Given a radius-edge threshold ρ¯ ≥ 1, Proposition 2
preserved by the insertion of the new vertex c .
k states that refinement does not lead to a reduction
Proposition 3. Let h¯(x) ≥ h¯ , h¯ ∈ R+ be a positive in minimum edge length. Shape-based refinement
0 0
therefore preserves the existing minimal edge length
mesh-size function and Del| (X ) be the restricted sur-
Σ k
(cid:107)e (cid:107) ,associatedwitheithertheinitialsamplingX ∈
face triangulation induced after k refinement operations. 0 k 0
Σ, or the insertion of some previous vertex x ∈ X
Given a refinable surface facet f ∈ Del| (X ) use of the k
Σ k
due to size- or surface-error-driven refinement.
TypeIIpoint-placementschemeis‘declined’whenf issuf-
ficiently small. Specifically, selection of the Type II point- (ii) Notingthath¯(x)ispositive,withh¯(x)≥h¯ forsome
0
placement scheme requires that rk ≥ H, where rk is the h¯0 ∈ R+, size-driven refinement is clearly declined
radius of the surface ball associated with the facet f and once the local element size is sufficiently small. Re-
√
H is the length of the new frontal edge candidates defined calling that h(f) = 3r, where r is the radius of
via expressions (2)–(3). the surface ball associated with a given facet f, the
mesh-sizing constraints h(f) ≤ αh¯(x ) are satisfied
Proof. Givenasurfacefacetf, theoff-centreselectioncri- f
teria expressed in (4) requires that the Type II off-centre once all r ≤(α/√3)h¯(xf).
c(2) lies along a ‘safe’ region of the Voronoi diagram – (iii) Noting that (cid:15)¯(x) is positive, with (cid:15)¯(x)≥(cid:15)¯ for some
0
bounded by the centre of the ball B(ck,rk) = SDB(f) (cid:15)0 ∈ R+, it is clear that surface-error driven refine-
such that: ment is declined once the local element size is suffi-
ciently small. Recalling that (cid:15)(f) is the orthogonal
(cid:107)c(2)−m (cid:107)≤(cid:107)c −m (cid:107) (7)
0 k 0 distance between a given facet f and the centre of
its associated surface ball B(c,r) = SDB(f), it is
where m is the midpoint of the base edge (cid:107)e (cid:107). Express-
0 0
clear that (cid:15)(f) ≤ r, as SDB(f) must circumscribe
ing(7)intermsofthetargetedgelengthsandsurfaceball
the facet. In the worst-case, surface-error driven re-
radii:
finement is declined once all r ≤(cid:15)(x ).
(cid:113) (cid:113) f
H2−(cid:107)1e (cid:107)2 ≤ r2−(cid:107)1e (cid:107)2 (8)
2 0 k 2 0
8
Using(i)–(iii),theprogressionoftheFrontal-Delaunay SDB(f)isguaranteedtocreateedgesnoshorterthan
algorithm can be examined. Starting with an initial ver- the radius of the ball of f:
tex distribution X ∈ Σ, both Type I and Type II point-
0
(cid:107)e (cid:107) ≥r . (10)
placementschemesareusedtosatisfymesh-sizeandsurface- 0 k+1 k
error constraints – with refinement continuing until all
Here (cid:107)e (cid:107) is the length of the shortest edge cre-
0 k+1
surface balls are sufficiently small, as per (ii)–(iii). Any
ated by the insertion of the new vertex c . Dividing
k
remaining low-quality elements in violation of the radius-
by the current minimum edge length, and using the
edge constraints are subsequently refined through the in-
definition of the radius-edge ratio:
sertion of Type I vertices alone. This final phase is guar-
anteed to preserve minimum edge-length. (cid:107)e0(cid:107)k+1 ≥ rk =ρ(f). (11)
Using Lemma 1, it is clear that the Frontal-Delaunay (cid:107)e0(cid:107)k (cid:107)e0(cid:107)k
algorithm is guaranteed to produce a point-wise sampling
Clearly, the minimum edge length is decreased when
X ∈Σoffinitesizeandmustterminateinafinitenumber
any element with ρ(f) < 1 is refined. Noting that
of steps as a consequence. √
ρ(f) achieves a minimum of 1/ 3 when f is equi-
lateral, it can be seen that insertion of a Type I
Finite termination leads directly to a number of useful
Steiner vertex leads to a reduction in the minimum
auxiliaryguaranteesonboththenatureandqualityofthe
√
output tessellation. Adopting the conventional terminol- edge length by a factor of 1/ 3 in the worst case.
ogy, surface meshes generated using the Frontal-Delaunay (ii) ATypeIIvertexc positionedonafacetV ⊂Vor(X)
k
algorithm can be considered to be ‘provably-good’: associated with the short edge e , is equidistant to
0
theverticesx ,x ∈e . SinceV isalocalfacetofthe
i j 0
Corollary 1. Given a smooth surface Σ, a radius-edge
Voronoi complex, no other vertex x ∈X lies closer
m
threshold ρ¯ ≥ 1 and positive mesh-size and surface-error
to the point c .
functions h¯(x)>0, (cid:15)¯(x)>0, all surface facets f ∈T| in k
Σ
Recalling that the Type II refinement strategy re-
the restricted surface tessellation T| = Del| (X) gener-
Σ Σ
quiresthediametricballoftheedgee remainempty,
ated by the Frontal-Delaunay algorithm satisfy: (i) ρ(f)< 0
ρ¯, (ii) h(f)<αh¯(x ), and (iii) (cid:15)(f)<(cid:15)¯(x ). a worst case reduction in (cid:107)e0(cid:107) occurs when ck is po-
f f
sitioned at the intersection of V and the diametric
Proof. The Frontal-Delaunay maintains a queue of ‘bad’ ballofe ,formingarighttriangle. NotingthatV in-
0
elements,containinganysurfacefacetsf ∈Del| (X)that tersectse atitsmidpoint,theminimumedgelength
Σ 0
are in violation of one or more local constraints. Specif- is decreased by a factor of 1/√2 in such cases.
ically, a given 2-face f ∈ Del| (X) is enqueued if: (i)
Σ
ρ(f) ≥ ρ¯, (ii) h(f) ≥ αh¯(x ), or (iii) (cid:15)(f) ≥ (cid:15)¯(x ). Given Overall, Type I point-placement is a limiting case, re-
f f √
that termination is assured, it is clear that the refinement ducing minimum edge length by an O(1) factor of 1/ 3 in
the worst case.
queue must become empty after a finite number of steps
andthatallfacetsintheresultingtriangulationDel| (X)
Σ 3.7. Quality of Approximation
satisfy the requisite radius-edge ratio, element size and
While the preceding results address the behaviour and
surface error constraints as a consequence.
performance of the underlying Delaunay-refinement algo-
In addition to proofs of termination and convergence, rithm, the approximation quality of the resulting mesh is
it is also necessary to ensure the point-set X ∈Σ remains also a key consideration. Using the framework introduced
well-distributed throughout the refinement process. Such byBoissonnatandOudot[15,16],avarietyofgeometrical
behaviour requires that no spurious short edges be intro- and topological guarantees are sought. Firstly, a number
ducedwithintheinitialrefinementphaseasthealgorithm of results from [15, 16] are summarised:
seeks to satisfy mesh-size and surface-error constraints.
Definition 5. Given a bounded domain Ω enclosed by a
Proposition 4. Let Del| (X ) be the current restricted surface Σ, the medial-axis of Ω is the topological closure
Σ k
surface triangulation. Given a large refinable surface facet of the set of centres of empty balls that touch the surface
f ∈ Del| (X ), the minimum edge length (cid:107)e (cid:107) is de- Σ at more than one point.
Σ k 0 k
√
creased by an O(1) factor of 1/ 3 in the worst-case fol-
Noting that the medial-axis lies at the ‘centre’ of any
lowing the insertion of a new Type I or Type II Steiner
geometrical features in the domain induced by Σ, it is
vertex.
clear that the distance to the medial axis, denoted d (x)
M
Proof. ConsideringtheinsertionofnewTypeIandTypeII throughout, isameasureofthelocal‘thickness’ofthedo-
vertices separately: main. Using such considerations, Boissonnat and Oudot
introduce the notion of loose (cid:15)-samples to describe point-
(i) Recallingthatthesurfacefacetf islocallyDelaunay, wisesamplingsX ∈Σthatinducerestrictedsurfacetrian-
the insertion of a new Type I vertex ck at the centre gulations Del|Σ(X) exhibiting favourable geometrical and
of the associated surface Delaunay ball B(ck,rk) = topological guarantees.
9
Definition 6. Let X ∈ Σ be a point-wise sampling of a 4. Mesh Size Functions
smooth surface Σ. Let E (X) ⊆ Vor(X) be the set of
V
edgesintheassociatedVoronoicomplex. ThesamplingX The construction of high-quality mesh-size functions
iscalledaiscalledaloose(cid:15)-sample ofΣif: (i)forallinter- is an important aspect of both the restricted Delaunay-
sections x ∈ {E (X)∩Σ} the associated (cid:15)-balls, B(c,r), refinement and Frontal-Delaunay algorithms presented in
are non-empty, sVuch that X ∩B(x,(cid:15)d (x)) (cid:54)= ∅ for some Sections2and3. Agoodmeshsizefunctionh¯(x)incorpo-
M
(cid:15)∈R+, and (ii) the restricted triangulation Del| (X) in- rates sizing constraints imposed by both the user and the
Σ
cludes at least one vertex on all connected components of geometryofthedomaintobemeshed. Thesecontributions
the surface Σ. can be considered via two separate size functions, where
h¯ (x)representsuser-definedsizinginformationandh¯ (x)
u g
Boissonnat and Oudot additionally establish the fol- encapsulatesgeometrydrivenconstraints. Inthisstudy,it
lowinggeometricalandtopologicalguaranteesforrestricted is required that both h¯ (x) and h¯ (x) be piecewise linear
u g
Delaunay tessellations Del| (X) induced by such loose (cid:15)- functions defined on a supporting tetrahedral complex S.
Σ
samples: Construction of appropriate user-defined functions h¯ (x)
u
ishighlyproblemdependentanddetailedformulationsare
Proposition 5. Let X ∈ Σ be a loose (cid:15)-sampling of a
not considered here.
smooth surface Σ, with (cid:15)≤0.08. The associated restricted
As per the discussions presented in Section 3.5, it is
surfacetriangulationDel|Σ(X)exhibitsthefollowingprop- known that the geometric mesh size function h¯ (x) must
erties: (i) the triangulation Del| (X) is homeomorphic to g
Σ besufficientlysmallrelativetothelocal-feature-size ofthe
the surface Σ, (ii) for any 2-face f ∈ Del| (X) the angle
Σ domain, denoted lfs(x), to ensure that the resulting sur-
between the facet normal nˆ and the surface normal nˆ ,
f Σ face triangulation Del| (X) is both topologically and ge-
sampled at the vertices x ∈ f, is O((cid:15)), and (iii) the one- Σ
ometrically correct. Specifically, to guarantee that the re-
sided Hausdorff distance H(Del| (X),Σ), sampled at the
Σ sulting tessellation is a loose (cid:15)-sample, it is required that
centre of the surface Delaunay balls B(c,r) = SDB(f)
is O(cid:0)(cid:15)2(cid:1). max hg(x) ≤ 0.08lfs(x), where lfs(x) is the distance to the
medial-axis of the bounded domain Ω. In practice, such
The reader is referred to [15, 16] or [19] for additional constraintsaretypicallyfoundtobeoverlyrestrictive[19],
details and proofs. Boissonnat and Oudot go on to show allowing for a relaxation of the coefficient (cid:15). In this study,
that restricted Delaunay refinement, consistent with the (cid:15) = 1/2 is used throughout. The local-feature-size can be
algorithmpresentedinSection2,isa‘provably-good’tech- expressed as the distance from any point xi ∈ Σ to the
nique, being guaranteed to produce a loose (cid:15)-sample of a closestpointonthemedial-axis ofthedomain. Giventhat
surface Σ in a finite number of steps, and to generate a Σ is represented as an underlying triangulation P in this
high-quality surface tessellation Del| (X) as a result. work, a discrete approximation to lfs(x) can be calculated
Σ
TheFrontal-Delaunayalgorithmpresentedinthisstudy directly via the method of Amenta and Bern [21, 37, 38].
derivesfromthesameDelaunay-refinementframeworkde- For the sake of brevity, we do not describe these methods
veloped in [15, 16], and inherits much of its theoretical in full here, but simply note that they provide an efficient
robustness as a result. Using the properties of loose (cid:15)- means to estimate lfs(x) at the vertices of the supporting
samples, the surface triangulations generated using the complex S.
Frontal-Delaunay algorithm are guaranteed to be a good Given an estimate of lfs(x) on the boundary of S, it
geometricalandtopologicalapproximationtotheunderly- is possible to exact a degree of user-defined control on the
ingsurface–providedthemeshsizefunctionissufficiently resulting size function via a Lipschitz smoothing process.
small: FollowingtheworkofPersson[39],agradient-limited func-
tion h˜(x) can be constructed by limiting the variation in
Proposition 6. Given a smooth surface Σ and an (cid:15)-con- h¯(x)overtheelementsofS. Inthisstudy,ascalarsmooth-
formingsizefunction(α/√3)h¯(x)≤(cid:15)dM(x), thepoint-wise ingparameterg ∈R+ isusedtoglobally,andisotropically
sampling X ∈ Σ generated by the Frontal-Delaunay algo-
limit variation, such that
rithm is a loose (cid:15)-sample. When (cid:15) ≤ 0.08 the associated
triangulation Del|Σ(X) is a tight topological and geomet- h˜(xi)≤h˜(xj)+g(cid:107)xi−xj(cid:107) (12)
rical approximation of the surface Σ.
for all vertex pairs x ,x ∈ S. The gradient-limited size
i j
Proof. The Frontal-Delaunay algorithm refines any large
function h˜(x) becomes more uniform as g →0.
surface facets, ensuring that h(f) ≤ αh¯(x ) for all f ∈
f
Del| (X). Rearranging the previous expression, and sub-
Σ
stitutingforthegivenmeshsizeconstraintsleadstoh(f)= 5. Results & Discussions
√ (cid:0)√ (cid:1)
3r ≤α 3/α (cid:15)dM(x),whichsimplifiestor ≤0.08dM(x)
when (cid:15)≤0.08. Recalling Definition 6, it is clear that X is TheperformanceoftheDelaunay-refinementandFrontal-
a loose (cid:15)-sample of Σ, and, via Proposition 5, Del| (X) is Delaunay surface meshing algorithms presented in Sec-
Σ
a tight geometrical and topological approximation of Σ as tions 2 and 3 was investigated experimentally, with both
a result. techniques used to mesh a series of benchmark problems.
10