BIB-VERSION:: CS-TR-v2.0
ID:: CORNELLCS//TR92-1267
ENTRY:: 1993-10-14
ORGANIZATION:: Cornell University, Computer Science Department
LANGUAGE:: English
TITLE:: Quality Mesh Generation in Three Dimensions
AUTHOR:: Mitchell, Scott A. 
AUTHOR:: Vavasis, Stephen A.   
DATE:: February 1992
PAGES:: 62
ABSTRACT::
We show how to triangulate a three dimensional polyhedral region with holes.
Our triangulation is optimal in the following two senses: First, our 
triangulation achieves the best possible aspect ratio up to a constant. 
Second, for any other triangulation of the same region into $m$ triangles with 
bounded aspect ratio, our triangulation has size $n$ = $O$($m$). 
Such a triangulation is desired as an initial mesh for a finite element mesh 
refinement algorithm. Previous three dimensional triangulation schemes either 
worked only on a restricted class of input, or did not guarantee well-shaped 
tetrahedra, or were not able to bound the output size. We build on some of 
the ideas presented in previous work by Bern, Eppstein and Gilbert, who have 
shown how to triangulate a two dimensional polyhedral region with holes, with 
similar quality and optimality bounds.
END:: CORNELLCS//TR92-1267
BODY::
Quality Mesh Generation
in Three Dimensions
Scott A. Mitchell*
Stephen A. Vavasis**
TR 92-1267
Febwary 1992
Department of Computer Science
Cornell University
Ithaca, NY 14853-7501
*Center for Applied Mathematics, Cornell University, Ithaca, NY 14853. Supported in
part by Hughes Research Laboratories, Malibu, CA, and by DARPA under ONR
contract N00014-88K-0591, ONR contract N00014-89J-1946 and NSF grant IRI-
9006137, by a summer internship and an NSF PYl award.
** Department of Computer Science, Cornell University, Ithaca, NY 14853. Supported
by an NSF Presidential Young Investigator award.
Quality Mesh Generation in Three
Dimensions
Scott A. Mitchell*			Stephen A. Vavasist
February 7,1992
Abstract
We show how to triangulate a three dimensional polyhedral region
with holes. Our triangulation is optimal ill the following two senses.
First, our triangulation achieves the best possible aspect ratio up to
a constant. Second, for any other triangulation of the same region
into m triangles with bounded aspect ratio, our triangulation has size
n = 0(m). Such a triangulation is desired as an initial mesh for a
finite element mesh refinement algorithm. Previous three dimensional
triangulation schemes either worked only on a restricted class of in-
put, or did not guarantee well-shaped tetrahedra, or were not able to
bound the output size. We build on some of the ideas presented in
previous work by Bern, Eppstein, and Gilbert, who have shown how
to triangulate a two dimensional polyhedral region with holes, with
similar quality and optimality bounds.
1			Introduction
Triangulation of polyhedral regions is a fundamental geometric problem for
numerical analysis. In particular, if one wishes to solve an elliptic boundary
*Center for Applied Mathematics, Cornell University, Ithaca, NY 14853, Supported in
part by Hughes Research Laboratories, Malibu, CA, and by DARPA under ONR contract
NOOOl4-88K-0591, ONR contract N00014-89J-1946, and NSF grant IRI-9006137, by a
Xerox summer internship and an NSF PY! award.
tDepartment of Computer Science, Cornell University, Ithaca, NY 14853. Supported
by an NSF Presidential Young Investigator award.
value problem on a three-dimensional domain, it is necessary to discretize
the domain with a mesh. If the domain is sufficiently complicated, then the
method of finite elements is commonly used. In this method the domain is
divided into small convex polyhedral regions with a fixed number of faces
called elements. A common choice for an element in three dimensions is the
tetrahedron. Thus, a triangulation of the domain is required.
For numerical stability in the finite element method, it is necessary that
the tetrahedra have bounded aspect ratio. This means that the angle between
any adjacent pair of edges of the tetrahedron, or between any edge and a 2-
dimensional face not containing it, is bounded below by a constant. For
information about aspect ratio bounds in numerical analysis, see Babuska
and Aziz [1976].
Our algorithm generates a triangulation for a nonconvex bounded poly-
hedral domain with holes. In addition, the triangulation is optimal in two
respects. In particular, the best possible aspect ratio is achieved for the
tetrahedra, and the number of tetrahedra is within a constant factor of the
best possible for any triangulation with bounded aspect ratios.
Our work is closely based on earlier work by Bern, Eppstein and Gilbert
[1990j, who solved the corresponding problem for two-dimensional polyhedral
domains. The main complication in extending the work to three dimensions
concerns vertices of the domain. In particular, for a two-dimensional domain
a vertex is adjacent to only two edges, and there are essentially three cases
for vertices (the two edges make an acute angle, an obtuse angle, or an angle
greater than 1800). In three dimensions, however, a large number of faces
can meet at a single vertex, forming a complicated cross-section with some
convex and some concave angles. Handling vertices properly necessitates a
host of changes in the algorithm. Our running time bounds are probably
not the best possible; this is discussed further at the end of the paper. Our
optimality condition is slightly stronger than Bern et. al.'s condition. In
particular, they show that overall the number of triangles they generate is no
more than a constant above optimal. The triangles used in a particular region
of the domain, however, could be arbitrarily smaller than what is needed
for the triangulation to achieve good aspect ratio. In our triangulation, no
tetrahedron is ever more than a constant factor smaller than the tetrahedron
at the same location in the best possible triangulation.
Other authors have considered three-dimensional regions, but, to our
knowledge, no previous work has simultaneously addressed the problems of
2
optimality, aspect ratio and of complicated nonconvex regions. Indeed, as
far as we know, there is no previous algorithm to triangulate a nonconvex
three-dimensional region with guaranteed aspect ratio (regardless of the opti-
mality of the triangulation). Our technique is based on an octree partitioning
of the domain. This technique has been used before, for example, by Carey,
Sharna and Wang [1988]. Chazelle and Palios [1989] consider optimality but
are not interested in angle bounds. Chew [1989b] presents an algorithm for
twodimensional polyhedral regions, but does not address the problem of op-
timal number of triangles. Moore and Warren [1990] address the problem
of adaptive mesh generation for box-shaped regions. Mitchell [1987] and
[1988] consider the problem of adaptive mesh generation of complex regions
in regard to finite element error bounds. Because of the importance of mesh
generation, the literature on this problem is very extensive; see for example,
the conference proceedings edited by Hauser and Taylor [1986].
Recently, Bern and Eppstein [1990] have surveyed the literature on tri-
angulations, giving more details about earlier work.
Our optimality proof is based on characteristic length functions. Similar
functions have been used before by Miller and Thurston [1990], Miller and
Vavasis [1990], and Miller, Teng and Vavasis [1991]. We have not seen them
used, however, to analyze the construction of a triangulation.
The remainder of this paper is organized as follows. In Section 2 we give
more details on the aspect ratio of a tetrahedron. In Section 3 we discuss the
formation of the octree. In Section 4 we demonstrate a relation between
the boxes of the octree. In Section 5 to Section 7 we we discuss the steps
necessary to produce a triangulation given an octree. In Section 8 to Section
10 we provide the optimality proof. For the remainder of this introduction,
we describe our assumptions.
Our triangulation is denoted AocT. We assume we are given a three-
dimensional polyhedral region as follows. There is a list of vertices with
coordinates specified. There is also a list of edges and faces. The three lists
are mutually linked. The list of edges for each vertex is ordered as they
occur around the vertex. Similar orderings are assumed for the other lists.
We assume that P is connected; if not, each component could be triangulated
separately.
A face of P, or another polyhedral region, can have either zero dimensions,
called a vertex, one dimension, called an edge, or two dimensions, called
a facet. In some circumstances we want to regard P itself as the three-
3
dimensional face. Note that the containment relation induces a partial order
on faces.
We assume that the polyhedral region is nondegenerate in the following
senses. Every two-dimensional face of P has the interior of the polytope on
exactly one side. Every one-dimensional edge is incident on exactly two two-
dimensional faces. For every zero-dimensional face (vertex) v, for every small
enough open neighborhood N of v, the set (N n P) --H ?vJ has exactly one
connected component. These assumptions may be dropped in future work.
2 Aspect ratio definition
To argue about the aspect ratio of triangulations of P, we wish to define
the angle between any two faces of P, or faces in a triangulation of P, that
intersect each other non-trivially. Various concepts such as dihedral angle and
solid angle have been used in the past, but we adopt a different approach.
Interior Angle We define the interior angle or angle between two faces
f and g in the cases where f and g are a facet and an edge, two facets, or
two edges that have a common intersection v. We additionally require that
one face is not a subset of the other. We say that two rays r1, r2 with a
common endpoint v can see each other if there are points v1, v2 on r1, r2 each
at distance t > 0 from v such that the sector vviv2 lies in P.
If f and g are edges, then we measure angles in the usual way for rays,
except we require that the rays can see each other. We always measure angles
with respect to the interior of P, so that f and g could make an angle close
to 3600.
If f is an edge and g is a facet, then we consider all rays r in the plane of
g, such that r lies between the two edges of g that meet at v, and such that
f and r can see each other. Then the angle between f and g is defined to be
the minimum angle between the ray covering f and all rays r meeting the
preceding conditions. We have minimum rather than infimum because faces
are assumed to be closed.
If f and g are both facets whose intersection is a vertex, then we let both
the ray for f and the ray for g vary as in the last paragraph, and take the
minimum angle between the rays to be the angle between f and g that can
see each other. If the facets share a common edge, then the interior angle
between the facets is the angle between the rays in the half planes of f and g
4
perpendicular to the common edge. Once again, we only consider rays that
see each other.
We use the symbol "a" throughout the paper to refer to the smallest such
interior angle between any two faces of P. We note that if the interior angle
between two facets is P, then one of them has a subface edge that makes an
angle of no more than p with the other. Thus there are always either two
edges, or an edge and a facet, that have interior angle a between them.
Aspect ratio. A useful measure of the shape of a tetrahedron is its aspect
ratio, which is the ratio of the radius R of the smallest containing sphere, to
the radius r of the largest inscribed sphere. In general the smallest containing
sphere is not the sphere through the vertices of the simplex. A simplex with
bounded aspect ratio has bounds on the angles made between any pair of
faces that intersect, and not just on the angles between edges.
It can be proved that the aspect ratio of a tetrahedron is bounded above
and below by constant multiples of the reciprocal of the sharpest interior
angle (as defined above) in the tetrahedron. We do not prove this claim
because it is implicit in theorems and lemmas in upcoming sections.
3 Subdivision of P into cubes, the octree
3.1 Data structures and definitions
The main data structure we use for our algorithm is an octree. We commonly
refer to each node of the tree as a box. We associate with each box, b, a
polyhedral region of JR3 called the embedding of the box and denoted 1(b).
During the generation of the octree, 1(b) is exactly a three dimensional cube.
Later boxes are warped and triangulated, changing their geometric structure.
An octree node is either a leaf, or has eight children. The embedding
of the eight children of a node are the eight cubes obtained by dividing the
embedding of the box in half in each of the three dimensions. We say a node
is split if it is not a leaf. The process of creating the eight children of a node
is called splitting.
For clarity, we use the term box face ("box facet," "box edge," etc.) to
denote a face of a box and P face to denote a face of P. Similarly for a?
and ?b.
Duplicate. Some boxes in the octree will be duplicated into the original
5
node and several new nodes, called duplicates or duplicate boxes. We create
duplicates whenever the intersection of the box with P has more than one
connected component. Each duplicate represents the same geometric cube in
IR?, but is associated with one connected component of P n 1(b). We use the
notation P A b to denote the component of P n 1(b) assigned to a particular
box b.
Whenever we split a box b, if P A b is nonconvex a child box b' may
have more than one component (i.e., (P A b) n b' might have more than one
component even though P A b consists of one component). Whenever a box
is split, we immediately determine whether any of its children contains more
than one component of P, and if so, we make duplicates of the child box.
If a P face is incident upon P A b for a box b, we say that a box contains
the face. We maintain pointers between each box and the faces it contains.
Extended Box. For a given box b, we define the extended box of b, ex(b),
such that I(ex(b)) is the cube concentric with 1(b) and with each dimension
expanded by a factor of 5. We use the notation P A ex(b) to denote the
component of P n I(ex(b)) that contains P A b. The extended box contains
only the P faces of P A ex(b). Note that ex(b) is not a box of the octree, but
may be constructed from boxes of the octree.
Adjacent. Two boxes are called neighbors if their embeddings intersect
non-trivially, and one is not a duplicate of the other. The size of a box is
the length of an edge. We say that box b1 is adjacent to box b2 if they are
neighbors and there is a point of P common to both, i.e. if PAb1flPAb2 $
In certain cases such as when a box contains faces meeting at a reflex angle,
a box may be adjacent to two duplicates of a second box. We say b1, b2 are
balance-adjacent if they are neighbors and there is a point of P common to one
of the boxes and the extended box of the other, i.e. if (PAex(b1))fl(PAb2) # $
or if (PAex(b2))fl(PAb1) $ $. Boxes that are adjacent are balance-adjacent,
but not all balance-adjacent boxes are adjacent, see for example Figure 3.
Balance Condition. If we split a box, we immediately split other boxes
to maintain the following invariant called the balance condition: No box is
balance-adjacent to a box more than twice its size. Certain boxes containing
vertices or edges of P are called protected boxes, and are exempt from being
split by the balance condition later in the algorithm. Nonetheless, the ratio
of the size of a protected box to an adjacent box is bounded above by a
constant that depends linearly on 1/a, as we prove in the next section.
6
3.2 Creating the octree
We generate the octree by selectively splitting and duplicating nodes. The
goal of splitting and duplicating boxes is to make the boxes small enough
so that the intersection of P with the embedding of any box is as simple as
possible. However, boxes should not be made too small, as this would lead
to an excessive number of tetrahedra in the final triangulation.
The octree generation algorithm is divided into three phases, the vertex
phase, the edge phase, and the facet phase.
It is easy to state the conditions under which we duplicate b, namely,
whenever P A 1(b) has more than one component. Nonetheless, the box
duplication process is actually the most computationally complicated part of
the octree generation algorithm, because determining components of PA 1(b)
is a nontrivial computational task. The details of the duplication process are
described below.
On the other hand, for splitting, the conditions under which a box is split
are more complicated, but the computational tasks are straightforward. We
now describe the various steps of the splitting algorithm.
3.3 The vertex phase
How finely we split boxes depends on the following definition.
Vertex Cone. A vertex cone of a box b is a set of P faces, F1,F2,... , Fk,
that satisfy the following:
(a) The set is precisely one vertex and all of its superfaces. This vertex is
called the apex of the vertex cone.
(b) The apex is in b.
(c) The faces F1,F2,... , Fk are exactly the P faces incident upon PA ex(b).
Vertex Crowded. We say that a box b is vertex crowded if the following
are satisfied.
(a) There is a P vertex v in b.
(b) The superfaces of v are not the only faces of P incident upon PA ex(b).
7
Th
Th
Figure 1: Here a box is duplicated for two components, and one duplicate
box is split because it is crowded.
Equivalently, a box is vertex crowded if P A b contains a vertex that is not
the apex of a vertex cone. We say that a face G of P A ex(b) interfereswith
v E P A b if v is not a subface of G. Thus b is crowded if it contains a P
vertex v, and ex(b) contains a face that interferes with v.
We recursively split and duplicate boxes, maintaining the balance condi-
tion at each step, until every box containing a vertex contains exactly one
vertex cone. Clearly this procedure will terminate once the box sizes become
a constant factor smaller than the minimum path length in P between a
vertex and a face that does not contain that vertex.
The details of the procedure for determining whether a box b should be
split are straightforward given an enumeration of the P faces bounding P A
ex(b). Such an enumeration is obtained from the procedure for determining
whether a box should be duplicated described below.
3.3.1 The vertex centering step
The conclusion of the vertex phase is a special one-time reorganization of
the boxes so that each vertex of P is far from the boundary of the box that
contains it. The value of this property will become clear in Section 5 and
later sections.
We first make the following observation concerning the boxes at the con-
clusion of the vertex phase. Let bbe a box containing a P vertex v. We claim
that any box balance-adjacent to bhas size at least that of b. We denote the
size of bwith the notation h(b). Let b1 be balance-adjacent to b. The balance
condition ensures that h(b1) is either equal to h(b), 2h(b), or h(b)/2. Hence
we need only prove that h(b1) $ h(b)/2.
8
To prove this, suppose that h(b1) = h(b)/2. Let U1 be the parent box of
b1. Then either U1 is vertex crowded, or it is split because of the balance
condition. We may assume that b'1 is split by the balance condition, since
if it is vertex crowded we simply take k = 1 in the following discussion and
the proof holds. There is a chain of boxes b'1,. . . ,b'k such that b'k is vertex
crowded, Wj is balance-adjacent to Uj+1 (for i = 1,... , k --H 1), and such that
= 2h(b1,.+1). These are just the conditions necessary for the balance
condition to propagate from b'k back to b'1. Note that I(b'k) lies interior to
I(ex(b)), no matter how large k is, because ex(b) is 5 times larger than b in
every direction, and the sizes of the b'j'S form a decreasing geometric series.
Thus, the vertex causing b'k to be crowded is in P A ex(b), contradicting the
fact that b is not vertex crowded.
Thus, every box b with a vertex is surrounded by boxes that are either
twice the size or the same size as b. The first part of the vertex centering
step is as follows. For every box b with a vertex, we split the boxes balance-
adjacent to b that are twice as large, and then we propagate the balance
condition (which causes every other box in the octree to be split at most one
time).
We claim that after this procedure, every box balance-adjacent to a box
b with a vertex is the same size as b. This claim is clearly true before the
balance condition gets propagated. Thus, we have to show that when the
balance condition is propagated, this does not cause any box balance-adjacent
to a box b with a vertex to be split.
Consider a setup of boxes that might cause a box b1 balance-adjacent to
a box b to be split by propagation of the balance condition. Suppose b1 is
a box balance-adjacent to b, where b contains a vertex v, such that b1 is the
same size as b. Suppose elsewhere there is a box b with balance-adjacent bk
such that bk is the same size as b, and bk is balance-adjacent to a box bk?l
four times as large, which is balance-adjacent to a box bk?2 twice as large,
and so on, up to b1, which is twice as large as b2. Under these conditions
(and only these conditions) propagating the balance condition from bk would
cause b1 to split.
However, we claim that the conditions described in the previous para-
graph cannot occur. This is because such conditions would contradict the
fact that every facet contained in ex(b) has vertex v as a subface (same
argument as above).
Thus, after the first part of the vertex centering step, every box has been
9
b
Figure 2: No matter how small we split a box b', if it is outside of ex(b) then
the balance condition will not split b.
split at most once, and every box with a vertex is balance-adjacent to boxes
of the same size. Consider such a box b. Without loss of generality, b has 26
balance-adjacent boxes surrounding it of the same size, arranged in a 3 x 3 x 3
group. These boxes may have been duplicated. If necessary, we adjoin empty
boxes in the spaces of this group where no box is balance-adjacent to b.
Now, the second part of the vertex centering step is as follows. We re-
organize the octree by uniting b with some of the boxes of the 3 x 3 x 3
group, and perhaps some duplicates of these as well. We form a new box
B by uniting b with seven boxes in a 2 x 2 x 2 group, and perhaps some of
their duplicates, to obtain the result that 1(B) is a cube of side length 2h(b),
and B contains the component of P A B containing the P vertex in b. See
Figure 3 for an illustration of this procedure in two dimensions. Note that
we can always pick the correct 7 neighbors of b, so that v is distance at least
h(b)/2 from any facet of B. Since we have effectively doubled the size of b,
B may be balance-adjacent to boxes four times smaller than B. We do not
take any steps to repair this violation of the balance condition, as the size
ratio of adjacent boxes is still bounded.
From now on, B is protected, meaning that it is never split again during
the algorithm.
10
Figure 3: Uniting boxes around a vertex box in two dimensions.
3.4 The edge phase
The next phase of the octree generation algorithm focuses on edges of P, and
is analogous to the vertex phase.
We define "edge cones" and "edge crowded" and split and duplicate boxes
as for vertices. The main difference is that we take an approach different from
centering to ensure that no box face is close to an edge. The splitting and
duplicating algorithm proceeds as if the protected vertex boxes do not exist.
In particular, the extended box of any box in this phase does not extended
into a protected box, and the vertices of P are ignored.
Edge Cone. An edge cone of a box b is a set of P faces, F1, F2, F3, that
satisfy the following:
(a) The set is precisely one edge and its two superfaces. This edge is called
the apex of the edge cone.
(b) The apex is contained in b, and
(c) The faces F1, F2, F3 are exactly the P faces incident upon P A ex(b).
Edge Crowded. We say that a box b is edge crowded if it contains an
edge of P, but it and its superfacets are not the only P faces incident on
11
P A ex(b). Equivalently, b is edge crowded if any edge contained in b is not
the apex of an edge cone.
Just as in the vertex case, we recursively duplicate an unprotected box for
each edge cone it contains, and split it if it is edge crowded. We maintain the
balance condition immediately after a box is split. Recall that the protected
vertex cone boxes are not changed in this step.
Identifying edge cones is straightforward if we have a list of facets bound-
ing PA ex(b), which is obtained from the algorithm for duplication described
below.
3.4.1 Protecting edge cone boxes
We seek the same goal as the centering vertices step but for edges. The fact
that an edge cone box may be balance-adjacent to a protected vertex box
prompts us to take an approach somewhat different from the one used for
vertex cones. Here if a P edge is near a box face 0, we surround the face
with protected boxes, as in the vertex case, but do not require them to be
the same size. Unlike the vertex centering phase, we do not eliminate 0 at
this time, because 0 may be a face of a protected vertex box. Instead, after
the entire octree is generated but before we triangulate, we warp the cube
(see Section 5) into a more general shape that increases the distance between
0 and the P edge, and eliminates any small angles as well.
We begin by splitting every unprotected box containing a P edge and its
balance-adjacent boxes and enforcing the balance condition. Every box is
split at most once. This ensures by the "crowded" definition that if a box b
(other than a protected vertex box) is balance-adjacent to a box b' containing
an edge E, then b cannot be balance-adjacent to a box containing any other
edge E'.
We always protect every box containing a P edge, and every box balance-
adjacent to a box containing a P edge, regardless of whether any face of the
box is close to a P edge. These boxes are never split again.
3.5 The facet phase
The third phase of the octree generation algorithm focuses on facets of P.
We define "facet cones" and "facet crowded" as in the previous phases. We
12
split and duplicate boxes as we did for vertices. As for edges, the main com-
plications arise when centering facets, especially when a facet cone box is
balance-adjacent to a vertex or edge cone boxes. The splitting and dupli-
cating algorithm proceeds as if the protected vertex and edge boxes do not
exist. In particular, the extended box of any box does not extended into a
protected box, and the edges and vertices of P are ignored.
Facet Cone. A facet cone of an unprotected box b is a single P facet,
F1, that satisfies the following:
(a) Facet F1 is contained in b.
(b) Facet F1 is the one and only P face incident upon P A ex(b).
For the sake of consistency, we will call F1 the apex of the facet cone.
Facet crowded. We say that a box b is facet crowded if it contains a
facet that is not the only P facet incident upon P A ex(b). Equivalently, b is
facet crowded if it contains a facet that is not the apex of a facet cone.
Just as in the vertex case, we recursively split a box if it is facet crowded.
We maintain the balance condition immediately after a box is split. Recall
that the protected vertex and edge cone boxes are not changed in this step.
3.5.1 Protecting facet cone boxes
We take the same approach as for edges. We begin as in the edge case
by splitting every unprotected box containing a P facet and its balance-
adjacent boxes and enforcing the balance condition. This ensures that a box
is balance-adjacent to at most one facet cone apex. Thus, in the warping
and triangulation step to follow in Section 5, we may consider the boxes
containing a facet without reference to any other facets, except for boxes
balance-adjacent to a protected edge or vertex box.
Technically, there is no need to protect boxes because there are no more
splitting phases to consider. However, for the sake of consistent terminology,
we protect facet cone boxes and the boxes balance-adjacent to them.
3.6 The duplication process
Recall that we duplicate b whenever it is determined that P A 1(b) has more
than one component. In this section we explain how to identify the compo-
nents of PA 1(b). This same techniques allow us to determine the components
13
of PflI(ex(b)), which is necessary for determining whether boxes are crowded,
and when to propogate the balance condition.
Because we believe that the algorithm proposed here is not the most
efficient algorithm possible for the task, we skip over some of the details
involved.
The first part of the duplication algorithm is a preprocessing algorithm,
which is run before the octree generation begins. In the preprocessing al-
gorithm we first identify the separate components of ap. If P is a convex
domain (or some other fairly simple shape) then a? has a single component
and no further preprocessing is necessary. Otherwise, suppose a? has several
components. One component is the exteriorcomponent; we call all the other
components floaters. The component for each facet, edge, and vertex can be
determined in 0(n) steps with a standard graph search, where n is the total
number of facets edges and vertices of P.
Once the components are identified, we now perform a second procedure
in which we identify a tether for each floater. First, for each floater C, we
identify the point v (necessarily a vertex) with the largest x coordinate on
This point is called the base of the tether. From the base we shoot a ray in
the positive xcoordinate direction until we encounter the first point v'on ap.
This point cannot be on the same floater. Such a point exists, since the ray
must eventually pass through the exterior component. The directed segment
from v to v' is cMled the tetherfor C, and v'is called the headof the tether.
See Figure 4 for a two dimensional example of the tethering structure of a
polygon. We can carry out the search for v'in a fairly naive fashion: Simply
loop over all facets, and for each facet F carry out a point-in-polygon test for
the point where the ray v passes through the plane of F. This requires time
proportional to the number of edges of F. Thus, the total time to locate v'
for a particular tether (looping over all facets) is 0(n). Therefore, the time
to identify all tethers is 0(n2).
Note the following two properties of tethers: (1) Each tether is contained
completely in P. (2) When regarded as a digraph on the components of ap,
the tethers form an in-tree rooted at the exterior component.
Now we discuss the procedure for determining when a box should be
duplicated. We first make a list of all P facets, edges and vertices contained
in the box. A particular P facet F in the box may intersect the box in
more than one component if the facet is nonconvex. We therefore must first
determine all the components of each facet passing through the box. This
14
Figure 4: Floaters and tethers for a two dimensional polygon.
may be done by sweeping a line segment across the facet, requiring O(v log v)
steps if v is the number of vertices in the facet. The total time for identifying
all components of all facets for a box is O(n log n).
Then we perform a combinatorial graph search to identify the various
components of a? n 1(b). Let these components be called sheets. There are
two kinds of sheets: those that intersect the boundary of 1(b), and those
that are completely internal to 1(b). Note that a sheet entirely inside 1(b)
must be a floater; such a floater is called an internal floater for b. The other
sheets are called external sheets. Note that a floater that intersects 1(b) is
an internal floater if and only if it lies entirely in 1(b), otherwise it is one or
more external sheets.
Clearly any sheet is incident upon a single connected component of P n
1(b) since the sheet itself is connected. Thus, the remaining task is to deter-
mine which sheets are connected to each other in P n 1(b).
For external sheets, the following algorithm exactly determines how they
are connected. We consider the intersection of each sheet with the boundary
of 1(b) (the surface of a cube). This intersection is a collection of disjoint
simple polygons. We form the list of polygons for all sheets. We call these
polygons the "rims.? Note that these polygons break across edges of the
15
cube. We then sweep over this list of rims, forming a tree indicating which
polygons are contained in others. If these polygons were in the plane, we
could sweep with a vertical line. This would require sorting x coordinates of
the polygons, hence requiring O(n log n) steps. On the surface of the cube
we can sort by the sum of x, y, z coordinates, which means that we sweep
over the surface of the cube with a closed polygonal segment (initially the
boundary of a triangle).
In the previous paragraphs we made reference to point-in-polygon tests
and plane sweeps. These are well-known tools in the computational geometry
literature. The reader unfamiliar with these tools should consult Edelbrunner
[1987] or Preparata and Shamos [1985].
From the sweep we can build a rooted tree indicating the containment
hierarchy among the polygons. If we assume that the first point swept is
outside P, then the following simple rule determines connectivity in P: Every
polygon at an even level i (counting 0 as the top level) is connected in P to
its children at odd level i + 1. The parities are interchanged if the first point
swept is inside P.
Thus, from an O(n log n) procedure we can determine exactly which ex-
ternal sheets are in the same components of PnI(b). This leaves the problem
of determining the component of P n 1(b) for an internal floater C. This is
determined from the tether (v, v') of C. If the head v' is inside 1(b), then
the floater lies in the same component as the facet containing v' (because the
entire tether is in P). This facet might be a facet of another internal floater,
but because the tethers form an acyclic graph, eventually we will reach a
point v' on an external sheet. See Figure 5.
The other case is that the tether's head v' lies outside 1(b). In this case we
find the point v" where the tether intersects the surface of 1(b). If we knew
where in the polygon tree v" lay, then we could determine the component
containing v" (which is the same component containing the floater). The
polygon regions containing v" are easily determined if, before starting the
polygon containment sweep, we first insert v" into the list of polygons as a
singleton polygon. This should be done for all tethers. This does not change
the running time bound of O(n log n).
Thus, in time O(n log n) we can determine all the components of P n
1(b). Moreover, it is easy to see that the rims allow us to determine which
neighboring boxes are adjacent to a given duplicate box.
Finally, we need to determine the facets adjacent to P A ex(b) for testing
16
vi'
Figure 5: Determining the components of a box. Sweeping the boundary de-
termines two components, and following tethers determines that all internal
floaters belong to the same component, and there are no internal components.
crowdedness. This is done by running the component-determination algo-
rithm on ex(b), and then saving the component of P n I(ex(b)) that contains
a point in P A b. This also allows us to determine which neighboring boxes
are balance-adjacent to b.
3.7 Summary of the octree algorithm
The octree algorithm has now been completely specified. The octree struc-
ture, including adjacency, duplication, and P face containment properties
remain unchanged by the algorithms that follow.
1. We initialize by setting the octree equal to a single node representing
a large box with all of the faces of P in its interior.
2. Vertex Phase. We split and duplicate boxes until it is determined that
no box is vertex crowded. The balance condition is enforced each time
a box is split.
3. Every box b' balance-adjacent to a box b with a P vertex is split, in
the case that h(b') = 2h(b). The balance condition is enforced.
17
4. Every vertex is centered in its containing box by merging eight boxes
around a vertex. Every box containing a vertex is now protected, so is
never split again.
5. Edge Phase. We split and duplicate unprotected boxes until it is deter-
mined that no box is edge crowded. The balance condition is enforced
each time a box is split.
6.
Every box containing a P edge and its balance-adjacent boxes are split
to further separate edges from each other. The balance condition is
enforced.
7. Every box containing an edge and its balance-adjacent boxes are pro-
tected, and never split again.
8. Facet Phase. We split and duplicate boxes until it is determined that
no box is facet crowded. The balance condition is enforced each time
a box is split.
9.
Every box containing a P facet and its balance-adjacent boxes are split
to further separate facets from each other. The balance condition is
enforced.
An important feature of the algorithm is that when it terminates, every
leaf box is either empty, lying entirely inside P, or contains a vertex, edge,
or facet cone. Also recall that the running time per box is proportional to
n log n, where n is the number of faces of P.
We now consider the running time of constructing the octree. We believe
that our current running time analysis is suboptimal, so we omit some of the
details. First, there is a slight addition to the algorithm as described so far
that gives us a better bound. If a box is crowded, with 0(n) additional time
per box we determine the closest face 0 interfering with some F, where n
is the number of faces of all dimension of P. llere for closest, distance is
distance within the box. We immediately split the box nonuniformly down to
a small level around the closest point of F to 0, so that 0 does not interfere
with F in one of the descendent boxes containing F. We thus get at least
one uncrowded box every time we test if a box is crowded. The running time
of this split is dependent on the distance in the box between F and 0. This
18
time may be amortized against the boxes created by the balance condition
when balancing the octree after this split.
This allows us to claim that the total number of boxes constructed by
the octree algorithm is bounded by a constant multiple of the number of leaf
boxes that contain a point of P. From the warping and triangulating rules to
follow, each such leaf box leads to at least one tetrahedron. In particular, the
number of such leaf boxes is bounded by ?, the size of the output. Finally,
the total amount of time spent on each box is at most n log n where n is the
number of faces of the polyhedral region P. Thus, a bound on the running
time is O(?n log n). Note that this subsumes the 0(n2) preprocessing to find
tethers, since ? = ?(n). See Section 11 for more remarks on this bound.
4 Relative sizes of adjacent octree boxes
In the last section, a box may have been protected in a vertex or edge phase
of the octree generation algorithm. A protected box is not subject to the bal-
ance condition in subsequent phases, so boxes adjacent to it may be smaller
than half its size. In this section we prove that the ratio of the sizes of two
adjacent boxes cannot be smaller than a constant multiple of a. Here a is
the smallest interior angle of P as defined in Section 2. The various kinds of
boxes must be considered as separate cases.
Observe that from the centering algorithm, if the vertex was close to any
face of its containing box, that face was deleted. Thus we have the following
lemma.
Lemma 1 The distance between a vertex of P and any face of a box is at
least one ei9hth the size of the box containing the P vertex.
Immediately after boxes are split in the vertex phase, each P vertex is a
vertex cone apex in the box that contains it. That is, the only P faces in
the extended box of a box containing a P vertex are superfaces of that P
vertex. The extra splitting to further separate vertices decreases the size of
boxes (and hence extended boxes). The merging of boxes to center the P
vertices increases the size of boxes, but boxes are no bigger than they were
before the extra splitting step. Hence, even after the separating and merging,
the superfaces of a P vertex are the only P faces contained in the extended
19
box of the box containing that P vertex. So each P vertex still satisfies the
definition of a vertex cone apex.
4.1 Vertex cone boxes
Theorem 1 Let B be a box with vertex cone apex v. For any box b adjacent
to B, h(b) > k a h(B), where k is a constant.
Proof. The proof divides into a series of lemmas, depending on the phase in
which the parent of b was split. If the parent of b was split in a separating
phase, we note that any box is split at most once in a centering phase, so
that the size of b is bounded by the size of its first ancestor that was split in a
crowded box phase. Let b1 be this ancestor. So, we can ignore the separating
phases, and the proof of the theorern divides into three lemmas depending
on when b1 was split, and one corollary.
Lemma 2 If b1 was split during the vertex box stage, then h(b) > h(B)/2.
Proof. The balance condition ensures that h(b1) is at least one half of the
size of the box containing v before the centering phase begins. 1
Lemma 3 If b1 was split during the edge box stage, then h(b) > k OL
where k is a constant and a is as in Section 2.
Proof. The proof breaks up into two cases that consider why b1 was split.
We first consider the more difficult case, that when b1 was split because of
the balance condition. Let d be the edge crowded box that was split causing
the balance condition to split b1. We now divide into two subcases depending
on the position of d relative to B.
Subcase 1. If 1(d) lies entirely outside I(ex(B)) or intersects ?1(ex(B))
the balance condition gives the result: Now matter how small a box is that
intersects ?i(ex(B)), a box that intersects ?I(B) will be split to a size no
smaller than half the distance from ?I(B) to ?I(ex(B)). This is very similar
to Figure 2. If 1(d) actually lies outside of I(ex(B)), then the same bound
holds.
Subcase 2. Otherwise 1(d) lies strictly interior to I(ex(B)). Observe
then that for two boxes, one either contains the other or their intersection
is a subset of their boundaries. So, there is a box at least half the size of d
20
between 1(d) and ?I(ex(B)). If dis split down to a size one half that of this
buffer box, for each of the children d1 of d, I(ex(di)) is a subset of 1(ex(B)).
Note that d may have to be split at most twice to achieve this size, so it is
sufficient to bound the size of d1, by showing that a box d that is actually
half the size of the buffer box has a lower bound on its size. That is, we may
assume that dis actually half the size of the buffer box, if it was larger then
this implies that b is larger.
The P vertex in B can not interfere with the edge in d, so the only way
that a face can interfere with the edge contained in d is for that face to be
a facet or edge of the vertex cone for B. Let f be the edge in d and g the
closest face that interferes with f. Note that g may not be the face that
makes the sharpest angle with f at v because of alignment odities, but of
course g makes an angle p with f at least this sharpest angle.
The distance from fflex(d) to v is at least the distance from a face of 1(B)
on the shortest straightline path to v. So by Lemma 1 and trigonometry, we
have that dist(fflex(d),gflex(d)) > k p h(B). The dependence is actually
on sin(P/2). If ? is large, then the balance condition between dand B in the
vertex centering step will be the limiting factor to how large dis. Otherwise,
sin(P/2) is bounded above and below by constant multiples of p.
Since g intersected the extended box of d, and f intersected B, we also
have dist(f A ex(d),g A ex(d)) < 3$3h(d). Hence we have the result of the
lemma for d. But since the balance condition caused b1 to split after d was
split, b1 must be smaller than d, and we have h(b1) > 2h(d).
The second case is if b1 was split because it itself was crowded. If b1
intersects ?I(ex(B)) as well as ?I(B), we see that h(b) > h(B) as before.
Otherwise, b is a special case of a box d in the interior of I(ex(B)) in the
previous case, so that the proof of the last case proves the desired result. 1
Lemma 4 (Corollary) If d is an edge protected box such that 1(d) is con-
tained in I(ex(B)) and further I(d)A?I(ex(B)) = ?, then h(d) > k.P.h(B),
where P is the sharpest interior angle an edge in B makes with any superface
of v.
Proof. This was almost proved in the proof of the previous lemma. For the
completion of the proof of the corollary, we need to consider a box strictly
interior to I(ex(B)) that contains an edge, but whose parent was split because
of the balance condition and not because it was edge crowded. We note that
21
the last time its ancestor was split because it was crowded it had the desired
size bound from the proof of the theorem. If the balance condition splits
this ancestor, it must first split a balance-adjacent non-crowded box. Once
this is done, the children of this balance-adjacent box form a buffer of two
same sized non-crowded boxes between the ancestor and crowded boxes. The
balance condition can never penetrate the inner box of this buffer by splitting
crowded boxes balance-adjacent to the outer box. Also, splitting other boxes
balance-adjacent to the ancestor only once will not cause the ancestor to
split again, so we may assume that there are two layers of same sized boxes
surrounding the first generation children of the ancestor, and so the ancestor
is never split again.
Hence we have the result with a constant k half as large as for the case in
the lemma where the box was split because it was edge crowded. The proof
of the corollary for a protected box that does not actually contain an edge
but is protected because it is balance-adjacent to a box containing an edge
follows with a constant one quarter as large, from the balance condition and
the fact that it is balance-adjacent to a box where the result holds. 1
Lemma 5 If b1 was split dnring the facet box stage, then h(b) > k ? h(B).
Proof. This is the most difficult of the three lemmas comprising the proof
of Theorem 1. Many cases may be handled in the same way as in the edge
lemma. Difficulties arise when considering a box that contains two facets
that meet at an edge, and this is where the corollary to the previous lemma
is crucial. The dependence on 0 in the lemma arises from the possibility of
small angles in the following discussion: If all of the angles considered below
are large, the factor of a may be eliminated from the result. As such, the
proof concentrates on the possibility that the angles considered below are
small. Certain anomalies occur in the following discussion when the angles
under consideration are large. Except where noted, these anomalies are easily
handled, and the resulting tedious discussion omitted. In this proof there will
be a series of constants k1, k2These labels are independent of the labels
of other such constants appearing elsewhere.
If b1 was split because of the balance condition, consider the crowded
facet box d that caused us to split b1. In fact, we consider d to be the last
descendent of this box that was split because it was facet crowded.
If 1(d) does not intersect I(ex(B)) or intersects ?I(ex(B)) the same ar-
gument as in the last lemma gives the desired result.
22
So it remains to consider the case where d is contained in the strict interior
of I(ex(B)). Let f be the facet contained in d, and g the closest interfering
facet in ex(d). If f and g meet only at v and do not share an edge, then
the same argument as in the last lemma, substituting faces with the same
symbols, proves this lemma.
We are left with considering a facet crowded box d containing a facet f,
whose closest interfering facet g shares a common edge e with f. See Figure
6. We seek to obtain a lower bound on the distance between the two facets f
and g, where for distance paths are restricted to outside of protected boxes
but inside 1(ex(B)). We bound this distance by a constant times Qh(B).
Then since the box d can only be facet crowded if both f and g intersect
1(ex(d)), the box sizes are at least a constant fraction of this distance.
We first make an argument that shows we need only consider a box ad-
jacent to a box where the corollary to the last lemma holds. As one travels
parallel to e to boxes d farther from the vertex v, but still inside I(ex(B)),
the distance within 1(ex(d)) between f and g increases. So d is at least half
the size of some box adjacent to B containing f. We have "half" and not
equality because of possible odities in the alignment of the edge with the
coordinate (octree) directions leading to nonuniformity in which boxes are
protected. We now consider the first box d1 in the facet phase, such that d1
contains f and is adjacent to both B and a box Be, where Be is an edge cone
box protected because it contains e or is adjacent to a box containing e.
We now describe how badly d1 can be crowded, that is we seek to de-
scribe the distance between f and g. Because edges are well surrounded by
protected boxes, the distance from e to the boundary of Be shared with d1
is at least a constant fraction of h(Be). We recall that since Be is an edge
protected box, it is by definition not part of the extended box of d1. Hence,
the closest point of f fl ex(d1) to e is on the shared boundary of Be and d1.
So we have that the distance from ffll(ex(d1)) to efll(B) is at least kh(Be).
Let e2 be the face making sharpest interior angle p with e at v. The ratio
of the size of Be to B is bounded above by a constant times this angle as in
Lemma 4. Combining this with the distance results of the last paragraph,
the angle at v between any point of f fl 1(ex(d1)) and face e is at least a
constant factor k1 times the angle p from e to e2, as in Figure 6 left.
Let F be an imaginary facet intersecting e and e2. It may be that F is
coincident with either f or g. In Figure 6 left F is coincident with f, in
Figure 6 right F is distinct from f and g. If e2 is a facet, we take F to
23
v
c
e2
Figure 6: The distance between two superfacets of an edge inside an extended
vertex box I(ex(B)) but outside an edge box is at least koh(B).
intersect it forming an edge that makes angle p with e. Note that F passes
between g and f in I(ex(d1)).
The two right triangles at the top of Figure 6 left are similar. First, if we
assume that P is less than it/4, then the ratio of their sizes is no more than
k1/2. If p > 7r/4, then it follows from the definition of crowded that h(Be)
is at least a constant fraction of h(B), and the result of the lemma follows
immediately from the fact that e is well centered in Be and f and g make
angle at least c at e. Second, we let be the angle between F n e2 and the
half plane of g with boundary e. Note that the larger triangle's leg between
F and g is at least k2 p21 h(B). Combining these two results shows that the
smallertriangle'slegbetweenFandgisatleast? P2,.h(B).
In Figure 6 the smaller triangle's leg between f and g is completely inside
Be and makes a right angle with b n ?I(B), so that the the distance from
g A 1(ex(d1)) to F A I(ex(d1)) is at least the length of the smaller triangle's
leg. If the orientation of Be with respect to f and g is such that this leg
is not entirely inside Be, then this distance may be at most a constant k3
smaller than the leg's length, because the facet that e is well centered in
Be implies that the angle g can make with ?I(Be) is bounded below by a
constant. Hence dist(g A 1(ex(d1)), FA 1(ex(d1))) > k4 p91 h(B) for constant
24
k4 = kik22ks Similarly for f (recall that Figure 6 left is the case where F is
coincident with f, so Pj' happens to be zero). We will prove below that at
least one of Pj' or is at least a constant factor times an interior angle of
P. Assume that is at least k? for some constant k. Substituting this in
the above yields dist(g fl I(ex(d1)), F n I(ex(d1))) > k5 a
Since F passes between f and 9, we have a reverse triangle inequal-
ity of distances between sets, that is dist(f A I(ex(d1)),g A I(ex(d1))) >
dist(f A I(ex(d1)),F A I(ex(d1))) + dist(g A I(ex(d1)),F A I(ex(d1))). By
the last paragraph, one of these two distance summands must be at least
k5.a.h(B). Hence
dist(f A I(ex(d')),9 A I(ex(d1))) > k5 ?
Recall our assumption that d1 was split because g interfered with f. In
particular, f and g are contained in 1(ex(d1)), so h(d1) is at least a constant
fraction of the distance between f and g in 1(ex(d1)), that is h(d1) > kk5
a h(B). From the definition of facet crowded, as soon as box d1 and its
descendants are split down to a size that is a constant fraction of this distance,
the descendants are not facet crowded. Hence d1 `s descendants are at least
a constant fraction of this distance, and by construction so is d, b1, and b.
That is, for some constants k and k',
h(b) > k'.h(d1) >k.a.h(B).
The proof is complete except for showing that one of Pj' or is within
a constant factor of an interior angle of P. We first show that the angle
Pj' is an angle interior to P, but not necessarily an "interior angle" between
P faces as defined in Section 2. Similarly for p;. If p91 + Pj' > r/4, then
the proof is trivial. Otherwise, we may draw a cone C with apex at v, axis
along e, and F A e2 on its boundary. See Figure 6 right for a spherical cross
section of this cone, which is nearly a "top" view of Figure 6 left. Since by
assumption the P face making the smallest interior angle with e is e2, the
only P faces visible to e within this cone are f, g, and e2. The angle at v
between F A e2 and the intersection of this cone with f is an interior angle.
But this angle is just P1,. Similarly for g.
We finally show that one of Pf' or P91 is within a constant factor of an
interior angle between P faces. If e2 is a facet that meets f at an edge, then
the angle Pj' between e2 A F and f A C may be arbitrarily small, and not
25
bounded by alpha. We show that even if e2 is a facet meeting f at an edge
and also g at an edge, both angles P11 and P21 may not be small. See Figure 6
right. The proof is as follows. If one of the angles is greater than half of P,
there is nothing to prove. Otherwise they are both less than half P, but then
they are a constant fraction of the angle ? between e2 A F and e2 A f and
the angle ? between e2 A F and e2 A g, respectively. The sum of these angles
is either a true interior angle between the two P edges of f A e2 and g A e2
meeting at v, or if it is not interior, there is a P edge inside the triangular
cone formed by f, g, and e2 but outside the circular cone C, such that this
P edge makes an interior angle with f or g smaller than this sum. In either
case, this sum is at least a.
For a final note in this proof, it is surpising that the result depends on
and that there is no explicit dependence on ?, where P' = Pj1 + P91. What is
hidden is the close relationship between P' and P, which may be shown by the
following mental exercise. For simplicity assume that e2 is a subedge of f,
i.e. F = f. Fix the angle made by the facets f and g at e at a value less than
r/4. Then P' depends linearly on ?. llowever, if P is fixed, then P' depends
linearly on the angle between f and g at e. That is, provided all angles are
small compared to 7r/4, P' is bounded above and below by the product of P
and the angle between f and g at e. In this sense, P' acts as a lower bound
on P, and the explicit dependence on P is subsumed by dependence on P1. It
is important to realize that P1 is (within a constant factor) of a true interior
angle, so that the box size has dependence on a, and not a2. I
This sequence of lemmas concludes the proof of Theorem 1. I
4.2 Edge cone boxes
Theorem 2 Let B be a box with edge cone apex v. For any box b adjacent
to B,
h(b) >k.a.h(B),
where k is a constant.
Proof. We again consider when b was last split. We know that B was not
protected in the vertex centering phase, since it is an edge protected box.
If b was last split in the vertex crowded phase, regardless of whether B is
protected or not, the balance condition gives us h(b) > h(B)/2. Similarly if
b was split in the edge crowded phase.
26
So we are left to consider if b was last split in the facet crowded phase.
Let f and g be the two superfacets of e. Let I be the union of the embeddings
of all of the boxes protected because of e or its vertices. Then since e is well
centered in its protected boxes, the distance between a point of f fl 81 and
g fl 81 is at least k p h(B), where P is the angle between f and g at e. We
may now repeat the argument of Lemma 3, with facets f and g meeting at
e taking the place of edge f and face g meeting at v. That is, a rewording
of the proof that edge boxes adjacent to a vertex box have the correct size
proves that facet boxes adjacent to an edge box have the right size. 1
4.3 Facet Cone Boxes
Theorem 3 Let B be a box with facet cone apex F. For any box b adjacent
to B, h(b) > h(B)/2.
Proof. We know that B was not protected in the vertex and edge protecting
phases. if b was last split in any phase, the balance condition gives us h(b) >
h(B)/2. I
if a box is not adjacent to a protected box, then it is empty, and the
balance condition provides us with h(b) > h(B)/2. We may now summarize.
Theorem 4 Let B be any box of the octree whose embedding contains a point
of P, and b an adjacent box. Then
h(b) >k.ot.h(B),
where k is a constant.
5 Warping
All of the boxes we will now consider are leaf boxes of the octree. if a box is
adjacent to boxes that are smaller, then we replace large faces of the box with
the smaller faces of the adjacent boxes where they intersect. For example,
if a box is completely surrounded by boxes one-half its size, then each of
the six sides of the box is now four square facets with a common vertex and
shared edges, and no longer strictly speaking a cube. We warp the faces of a
27
A
Figure 7: The special case of sharing a vertex with two boxes.
box away from faces P when necessary to achieve good aspect ratio when we
triangulate. We warp box faces that are close to P-faces, where the definition
of "close" is below.
The warping consists of moving box vertices. We consider adjacent boxes
to share faces of their common boundary, so that when we warp a box face,
there is a change in all of the adjacent boxes that contain that face. The
warping rules are selected so that no inconsistencies arise.
Note that a protected vertex box (after the centering step described in
Section 3.3.1) is always adjacent to boxes the same size or smaller. This
means that, without loss of generality, we do not have to triangulate and
warp facets of a protected vertex box; instead we can triangulate and warp
the facets of the adjacent boxes.
A special note is needed in the case of protected edge boxes when the
superfacets of the edge make an angle close to 3600. If an edge box B is
adjacent to two boxes b and c, such that b and c are not adjacent to each
other, but 1(b) A 1(c) is non-empty, then there may be two duplicate copies
of a face that should be shared by one face of 1(B).
See Figure 7. In this case, we retain the two distinct copies of the same
face. Each face is warped separately for b and c. The portion of the face for
b is retained during the triangulation of P A b but ignored when considering
P A c. Similarly for c.
The purpose of warping is to make sure no P face comes too close to a
box face. Note that P vertices are already guaranteed to be well separated
from box faces because of the vertex centering procedure described in Section
28
3.3.1. Thus, we only need to warp for P edges and P faces.
The warping is done in two passes. The first warping pass is for box edges
close to P edges. We say that a box edge E is close to a P edge F if the
distance from E to F is less than h/8. Here, h is the size of the smallest
box sharing the box edge E. Note that the boxes containing E are all edge-
protected, since they are all adjacent to a box containing F. This means
that the sizes of these boxes are all within a factor of two from each other.
For every close edge E, we move it away from F as follows. We move
every point on the edge, including the endpoints, by distance h/8 in the
direction that is orthogonal to E and F, oriented away from F. (If E and F
happen to be parallel, we move the points of E in the direction orthogonal
to E and coplanar with F.)
Even though the boxes after warping have complicated shapes, we still
let the "size" of a box b be its prewarped edge length, and we continue to
use notation h(b) for this size.
The second warping pass is for box vertices close to P facets. We say that
a box vertex v is close to a P facet F if the distance between them is at most
h/16, where h is the size of the smallest box containing v. For such a vertex,
we move it distance h/1(3 away from the facet in the direction orthogonal to
the facet.
In general, because the octree algorithm separates cones, there is never an
inconsistency in these rules that could cause a vertex or edge to be close to
two P facets. There is one exception to this, namely, the case of a protected
edge box containing two facets that make a small interior or exterior angle.
If the two facets make a reflex angle, and are both close to the same vertex,
then we are in the special sharing case described before: there are actually
two copies of the vertex, and each is warped separately for the appropriate
facet. See Figure 7. If the facets make an acute angle, and are both close to
a vertex of the protected edge box, then we project the vertex orthogonally
onto the plane halfway between the facets (i.e., the plane that bisects the
edge cone).
Whenever we warp a vertex in the second warping pass, we carry its
adjacent edges with it. Thus, box edges continue to be straight segments
between box vertices during the warping. The box facets, however, are no
longer well defined because the bounding edges are no longer coplanar.
We make a new definition as follows. Let an 1 vertex be the intersection
of a P edge with a box facet, or the intersection of a P facet with a box edge,
29
or a box vertex.
In the last paragraph we referred to the intersection of a P edge with a
box facet; the box facet is no longer well defined after warping. Accordingly,
we let the 1 vertex be the intersection of the edge with the prewarped box
facet.
Lemma 6 After the warping, no two distinct I vertices x, y are closer than
min(t, h/16). Here h is the size of the smallest box containing x, y. Variable
t is defined only for the case when when x, y are on two facets making a
sharp interior angle in a protected edge box, or when one is a box vertex in
between two such facets. In this case, t = k1ha where a is the angle between
the facets, and k1 is a constant.
Proof. This proof has numerous cases, all similar. For the sake of brevity,
we prove the lemma in two cases only.
We consider the case that x and y both arise as the intersection of a P
facet F with box edges E1, E2. (This is one of the more complicated cases
of the lemma.) If E1, E2 are disjoint edges then the lemma is clear. Points
x, y before the first warping pass are separated by at least distance h if they
lie on disjoint edges. After the first warping pass x, y must be separated by
at least 3h/4. After the second warping pass, x, y must be separated by at
least 5h/8.
Next, consider the other case that E1, E2 have a common box vertex v.
(In the third case that E1, E2 are the same edge, then x, y are not distinct.)
This case has two subcases, namely, v is not considered close to F during the
second warping pass, or that v is considered close. This case is illustrated in
Figure 8.
First, observe that in either subcase, the angle /xvy is at least 600; this
angle is originally 900, and the warping steps cannot reduce it more than 300.
Now, if v is not close to F during the second warping pass, this means
that the altitude of triangle xvy (from v to the xy leg) is at least h/16, which
means by the angle bound in the last paragraph that the base of this triangle,
that is, the xy leg, is at least h/16 long.
The same argument holds if v is close to F.
Another case of the lemma is the case of a protected edge box, where x
is the intersection of a P facet F with a box edge E, and y is a box vertex.
First, we know that, as a result of the first warping pass, y has distance at
30
F			E1
Figure 8: One case in the proof of Lemma 6.
least h/8 from the P edge E in the protected edge box, where h is the size of
the protected edge box, or half that size (if the protected edge box is adjacent
to a box half its size). Now, suppose a point on F is distance less than h/16
from y after the first warping pass. Then we either move y distance h/16
from F, or distance t, where t is half the distance from z1 to z2. Here, z1, z2
are the endpoints of a path though y orthogonal to the bisector of the two
P facets in the box. See Figure 9. But since y after the first warping pass is
distance at least h/8 from E, the same holds for z1, z2. This means that the
distance between z1 and z2 is at least 2 sin(a/2) h/8, and the distance from
y after the second warping pass is at least tan(a/2) h/8. This formula has
a lower bound of k1ha. I
6 Two dimensional triangulation
After warping, we triangulate facets of boxes. This is a preliminary step to
the three dimensional triangulation algorithm of the next section. Suppose
two boxes b, b' are adjacent, and suppose their intersection is a facet S. if we
assume that h(b) > h(b'), then S is a facet of b'. Under these circumstances,
we always triangulate 5 by considering it a facet of b', not b. (If b, b' are the
same size, we break ties arbitrarily). This means that we can always assume
that when triangulating a square 5, there is no box vertex in the interior
of 5. There are four box vertices at the corners of 5, and there may be
additional box vertices at points along edges of 5.
31
F
x
pi
z2
Figure 9: Another case in the proof of Lemma 6.
We now describe the triangulation of a box facet 5. Note that after
warping, the points of 5 are not necessarily coplanar, although the points
are nearly coplanar because the warping distances are small.
The triangulation of a box facet 5 breaks down into four cases, and cases
C and D have two subcases. The cases depend on which P faces pass through
the facet. Note that we only have to triangulate the portion of 5 interior
to P, which is bounded by a closed path of line segments. Call this path of
line segments Ir. We can also think of r as a "perturbed polygon" (since the
vertices are nearly coplanar). The six subcases are illustrated in Figure 10;
the regions Ir are shaded. The points on the boundary of ir are all I vertices.
Case A, no P faces pass through the facet 5. In this case, we put in a
new point v at the center of the prewarped box facet, and we connect v to
every segment along the boundary of 5.
Case B, two P facets F, 0 and their common edge E pass through facet 5.
Let v be the I vertex where E passes through the prewarped version of facet
5. Then we connect v to all segments on the boundary of r, triangulating 7r.
Case C, one P facet F passes through (two edges of) the facet 5. Let x, y
be the two points where F passes through two edges of 5, and let v be the
midpoint of x, y. Let c be the center of 5 (before warping). Let h be the
edge length of 5 (before warping).
Subcase Cl, v is within h/4 of c. Then we connect v to every segment
along the boundary of r. This divides the region into triangles.
Subcase C2, v is further than h/4 from c. Then we connect c to every
segment along the boundary of r. This divides polygonal region :r into
32
B
v
c
C2
v
Dl
Cl
c
D2
Figure 10: The six cases for triangulating a box facet.
triangles and quadrilaterals. Then we insert a diagonal (arbitrarily) into
each quadrilateral, triangulating 7r.
Case D, two P facets F, G pass through facet 5, but no P edge. Let x, y
and x', y' be points where F, G respectively cross through the boundary of 5.
Let v, v' be the midpoints of these segments. Let c be the center of 5 (before
warping). Let h be the edge length of 5 (before warping).
Subcase Dl, v or v' is within h/4 of c. Say v is within h/4 of c. Then
we connect v to every segment along the boundary of r. This divides ir into
triangles.
Subcase D2, v and v' are both further than h/4 from c. Then we connect
c to every segment along the boundary of r. This divides r into triangles
and quadrilaterals. Then we insert a diagonal (arbitrarily) into each quadri-
lateral, triangulating ?.
We claim that this triangulation has the following property.
Lemma 7 Let 5 be a facet of a box b, triangulated with the above algorithm.
Let h be the size of 5, and let h' be the size of the smallest box adjacent to
33
S. Let T be a triangle in this triangulation. If two facets pass through S, let
t be as in the last lemma. Then the radius of the largest inscribed circle in T
is at least k2 min(h', t), where k2 is a constant.
Proof. Standard geometry tells us the radius of the largest inscribed circle
in a triangle T = Axyz is bounded below by
k3 xy min(I/zxyl, Izzyxi),
where k3 is a constant. This fact is used below.
This proof breaks down into cases. First, we take the case that T has
three vertices pqv, such that p, q are vertices on the boundary of r, and v is as
in the above cases. This is the type of triangle we get in all cases except C2
and D2. First, observe that the distance from p to q is at least min(h', t/2);
this follows from the preceding lemma. Thus, we only need a constant lower
bound on /vpq and Zvqp. But these angle bounds are immediate, because
v is near the center of the facet, and p and q are along the boundary of the
facet.
A slightly different analysis is needed for the triangle vx'y' arising in Case
Dl, where x' and y' are vertices arising from the intersection of the P facet
G and 5. We know from Lemma 6 that the distance between two I vertices,
one from F A 5 and one from GA 5, is at least t. In fact, any point of F A 5
is at distance at least kt from & A 5, from the fact that the common edge
between F and G is far from 5 after warping. To be technically correct, since
the vertices of 5 are not actually coplanar, F A 5 refers to any triangulation
of 5 intersected with F, and similarly for &. Hence the distance from v to
x'is at least kt, which gives us an edge length, and the distance from v to
G A 5 is also at least kt, which gives us angle bounds.
A different kind of triangle arises from splitting quadrilaterals in Cases
C2 and D2. For each of the quadrilaterals in these cases, we can see that the
two opposite edges e1 and e2 have length at least k4 min(h', t/2). See Figure
11. For e1, this follows from the same argument as in the last paragraph.
For e2, this follows because facet F in 5 is well separated from the center c.
Therefore, the length of e2 is at least a constant fraction of the length of e1.
Similarly, it is clear that edges e3 and e4 have length at least k4 min(h', t/2).
Similarly, we can show that the angles made by edges e3, e4 with e1, e2
are bounded below by constants. Again, this follows because c is centered,
whereas e1, e2 are away from the center.
34
e1
Figure 11: Triangulating a quadrilateral in Case C2 (closeup)
This is enough to show that the inscribed circles of the triangles ob-
tained as a result of inserting a diagonal are at least a constant multiple of
min(h',t/2). I
Theorem 5 Let B be a box. For any triangle T on a facet of B, let r be
the radius of the largest inscribed circle. Then r > k a h(B), where k is a
constant.
Proof. This is an immediate consequence of the preceding lemma and Theo-
rem 4. in particular, h' in the previous lemma is at least kah(B) by Theorem
4. Moreover, tin the previous lemma was shown to be at least k a h(B) in
the proof of Lemma 6. I
General results have been obtained for two-dimensional triangulations
with guaranteed inscribed-circle radius bounds. Bern, Edelsbrunner, Epp-
stein, Mitchell, and Tan [1991j have a result concerning optimal two-dimensional
triangulation of polygons, assuming new points cannot be introduced. Bern,
Dobkin, and Eppstein [1991] have similar results in the case that new points
can be introduced. We have not been able to incorporate these algorithms
in the present version of our triangulation algorithm because we need to in-
troduce new points, but the new points have to be introduced in a carefully
controlled fashion.
35
7 Three dimensional triangulation
We now describe how to form tetrahedra from the triangles in the last section.
We triangulate on a box by box basis. The details of how we triangulate
depends on a case analysis of what is contained in a box, and how it intersects
the box. However, the general principle is to find a central vertex, and then
form one tetrahedron or prism for each triangle in the box by taking the
convex hull of this vertex and the triangle. The organization of the argument
is very similar to the two-dimensional triangulation in the last section.
7.1 Empty box tetrahedra
We first describe how to triangulate in three dimensions a box b containing
no P faces. We take as a central vertex v the centroid of b. (The prewarping
centroid may be used). For every triangle on the surface of the box we form
a tetrahedron by taking the convex hull of v and the triangle.
7.2 Vertex box tetrahedra
For a box b containing a vertex of P, the three-dimensional triangulation is
particularly easy. We take as the central vertex v the P vertex itself. We
form tetrahedra by taking the convex hull of v each triangle on the surface
of b. The vertex centering phase ensures that the distance from the central
vertex to the triangles is at least a constant times the box size.
7.3 Edge box tetrahedra
We now describe the next case, that of triangulating in three dimensions the
protected edge boxes that contain edges. Let edge E be the edge contained
in box b. We let v be the midpoint of E n 1(b). We now connect v to every
triangle on the surface of b, creating tetrahedra.
7.4 Facet box tetrahedra
Consider a box b containing one P facet F. We introduce a central vertex
v at the centroid of the prewarped box. If v is within distance h(b)/4 to F,
then we move v to the closest point of F to v.
36
We now have three cases. The first case is if v lies in the interior of P,
but not on F. Here we form tetrahedra by taking the convex hull of every
triangle of F and v. This requires the creation of a triangulation on the
intersection of F with b, which is a polygonal region. This polygonal region
is nearly convex (it would be convex if the triangles on each the surface of the
box were coplanar). To form a triangulation of F n 1(b), we pick a centrally
located point on the F n 1(b) and connect it to all the edges. Using similar
arguments as in the previous section, it is not hard to show that the largest
inscribed radius of any triangle in this two-dimensional triangulation is at
least a constant fraction of h', where h' is the size of the smallest box adjacent
to b.
The second case is that v lies on F. We take the central vertex v, and
proceed as if it were a vertex of P. That is, for each triangle of 1(b), we
generate a tetrahedron by taking the convex hull of it and v.
The third case is if v lies outside of P. We form tetrahedra by taking
the convex hull of v and every triangle of 1(b). These tetrahedra are clipped
by F into "prisms" with nonparallel sides, and a triangular top and bottom.
Zero, one or two of the vertices of the top may also be vertices of the bottom.
If two vertices are shared, then the prism is a tetrahedron. If one vertex is
shared, the prism is split into two tetrahedra by introducing a facet between
the shared vertex, and one distinct vertex of the top and one distinct vertex
of the bottom. If no vertices are shared between the top and bottom, the
prism is split into three tetrahedra by introducing two facets. One facet is
between two of the top vertices and the opposite bottom vertex. The other
facet is between the opposite bottom vertex, one of the other bottom vertices,
and the opposite top vertex.
7.5 Two-facet box tetrahedra
Perhaps the most complicated case is the case of a box b containing two
facets F1, F2 adjacent to an edge E, but the edge itself is not in the box. Let
v be the centroid of box b (prewarped). If v is within distance h(b)/4 from
F1 or F2, then we move it to the closest point on either of these facets (break
ties arbitrarily). Say F1 is closer.
There are several possible arrangements for v, F1, F2. For example, v
could be on F1, or it could be between F1, F2, or it could be outside P, with
F1 facing away from v and F2 facing towards P.
37
For each facet Fi, (i = 1,2) we triangulate Fj fl 1(b) provided that Fi
does not contain v and faces v. This means that either one or both of Fj is
triangulated. We use the same algorithm described in the facet box case for
triangulating the facet.
Now we connect v to all the triangles on the boundary of b that are
interior to P. We also connect v to all the triangles of F1, F2.
If v is contained in F1, or if v is between F1, F2, then we have now divided
the box into tetrahedra. If v is outside P with F1 facing away, then by
connecting v to the triangles of F2 we have divided the interior of P into
tetrahedra and prisms. We now split the prisms into tetrahedra as in the
facet box case.
8 Aspect ratio of tetrahedra in AocT
Our triangulation is optimal in two respects. First, the maximum aspect
ratio among tetrahedra of our triangulation is optimal up to a constant mul-
tiple. Second, compared to all other triangulations of fixed aspect ratio, our
triangulation has the minimum number of tetrahedra up to a constant mul-
tiple. In this section and the next section we establish the first optimality
property, focusing on the optimality of the aspect ratio.
Recall that the tetrahedra we form are either the convex hull of a central
vertex and a triangle that does not contain the central vertex on the surface
of a cube or a facet of P, or else a portion of a prism that is divided into three
tetrahedra. We now want to argue about the aspect ratio of all tetrahedra
arising in the algorithm of the previous section.
There are three types of tetrahedra arising in the triangulation. A Type
A tetrahedron has one vertex v centrally located in the box, and the other
vertices on a box face. This type of tetrahedron arises in the cases of trian-
gulating a box with a vertex, a box with an edge, or a box with a facet in
which the vertex v in the last section lies near the center of the box, is in P,
and is not close to the second facet (if the box has two facets). Also some of
the tetrahedra arising from the case of two facets in a box and v on one of
the facets are of this type.
A Type B tetrahedron arises only in the two-facet case in the last section,
when the vertex v is on one facet, and the triangle it is joined to lies on
the other P facet. These tetrahedra arise in cases analogous to one of the
38
triangles in case Di of Figure 10.
A Type C tetrahedron arises from a vertex v in the last section outside
P. This happens only with boxes with one or two facets and no edges. This
causes P A b to be divided into prisms, and then each prism is divided into
tetrahedra. These tetrahedra arise in cases analogous to triangles in cases
C2 and D2 in Figure 10.
Intuitively, all three types of tetrahedra have aspect ratio at most 1/a,
because they involve a base (the base being a triangle) whose inscribed radius
is between kiah(b) and k2h(b), and the apex is well-centered over the base,
and is distance between k3ah(b) and k4h(b) from the base. Here, ..... . , k4
are constants.
The aspect ratio computations are very technical and involved. We first
prove the aspect ratio bounds for Type A tetrahedra using a construction of
an inscribed sphere. The proofs for Type B and C tetrahedra follow from
a different analysis of the same construction. A Type A tetrahedron has a
base size between kiah(b) and k2h(b), and has a height at least k5h(b).
Consider a particular box b and its central vertex v, and a Type A tetra-
hedron t3 in this box. In this section we let subscripts denote dimension, so
that 12 is one of the triangles on the surface of PA b, and t3 is the tetrahedron
formed by the convex hull of t2 and v.
Lemma 8 The smallest containing sphere, S3, of 13 has radius R3 at most
kh(b), where k is a constant.
Proof. The smallest containing sphere of t3 has radius less than the radius
of a sphere that encloses the entire box. If the box was unchanged by the
warping step, we see that a sphere with radius equal to one half a diagonal
of the box will do. If the box was warped, then it may be enclosed in a
concentric box with side lengths 5/4 times that of the original box by the
definition of "close." Thus we can take k = 5$3/4. 1
This concludes the analysis of the smallest sphere containing 13. Next, we
analyze the largest inscribed sphere of 13 by actually constructing an inscribed
sphere. Consider the triangle 12 on ?(P A b). We form a new curved triangle
u2 on the surface of a sphere as follows. Let a be the closest point of 12 to
v. It does not matter if a is a vertex of 12 or not, although in the figures of
this section a is a vertex. The intersection of 13 and the sphere centered at v
39
a
V
Figure 12: The flat (t3) and curved (u3) tetrahedra.
with radius dist(v, a) is the curved triangle u2. We let u3 be the convex hull
of u2 and v. Note that u3 is contained in t3, and has three flat facets, each
one contained in a facet of t3 that intersects v. This is illustrated in Figure
12. We will construct a sphere inscribed in u3 instead of t3.
Lemma 9 The distance from v to its base t2 is at least h(b)/8.
Proof. This is because, by assumption for the case, t3 is a Type A tetrahe-
dron. Therefore, v is centered in b. 1
Lemma 10 Let u2 be defined as in the last paragraph. Then the radius of
the largest inscribed curved circle in u2 is a constant fraction of the radius of
the largest inscribed circle in t2.
Proof. Let C be the largest inscribed circle in t2. We project this circle
onto sphere S3. This yields an ellipse, whose minor axis is the product of
the original inscribed radius and a factor depending on the angle between t2
and a tangent plane to u2. We claim that the angles between v and any edge
of t2 are bounded below by constants. This is because, before warping, the
vertices of t2 were coplanar with a box facet, and v is at least distance h/4
40
from that facet, where h is the edge length of the box. After t2 is warped,
since the warping distances are small, this can only change its inclination
with respect to the original facet by a small constant angle that is arrived at
via a case analysis (which we omit).
This same angle bound applies to the angle between t2 and the tangent
plane. I
Lemma 11 (Corollary) If v is within kah(b) of the plane of t2, then the
radius of the largest inscribed curved circle in u2 is a constant times a times
the radius of the largest inscribed circle in t2.
Proof. The key observation is that v is still within a constant of the box size
to the farthest point of the triangle. The angle between v and an edge of t2
is bounded below by a term linear in the ratio of the distances from v to the
plane of t2 and from v to the farthest point of t2. In the lemma this was a
constant, here in the corollary it is bounded by a. This will be used only for
the analysis of Type B tetrahedra. I
In a series of lemmas we derive a lower bound for the radius of the largest
sphere inscribed in t3 by the radius of the largest curved circle inscribed in
u2. Let 53 denote the largest sphere inscribed in u3, r3 its radius, and c3 its
center. Similarly, we let ?2 denote the largest circle inscribed in u2, r2 its
radius, and c2 its center.
Lemma 12 The sphere 53 is tangent to u3 at points equidistant from v.
Proof. Pick a facet f of u3 adjacent to v. Consider the triangle formed by
v, the point of tangency between 53 and f, and the center of 53. Note that
this is a right triangle. Observe that the length of the hypotenuse and the
length of one leg depends only on the position of 53 and not on on which
facet was chosen. Therefore, the lengths of the other leg (the leg between v
and the points of tangency) is the same regardless of which facets is selected.
See Figure 13. This theorem applies to arbitrary tetrahedra with any vertex
chosen for v. I
Lemma 13 There is a sphere, 53', that is tangent to u3 at exactly the same
three points as ?2 is tangent to u2.
41
s3,-center
v
Figure 13: A sphere tangent to three facets of a tetrahedron.
Proof. We call these points of tangency a, b, e. We construct 53'. Note that 53'
is not contained in u3 or t3, but will serve as a guide for finding an inscribed
sphere. We first determine uniquely the center, c', of 53'. Consider the line
vc2. Now c' is the unique point on vc2 on the opposite side of c2 as v, such
that Zvac' is right. Additionally, since ac2 is perpendicular to the boundary
of u2, ac' is also perpendicular to the facet, f, of u3 containing v and a in
the plane tangent to the vertex sphere at a. Hence, c'a is perpendicular
to the facet f. Similarly for b and e. The hypotenuse-leg theorem gives us
congruent triangles, and we see that the distances from c' to each of a, b, e
are equal. So, c' is indeed the center of a sphere with the desired properties,
53'. I
Lemma 14 (Corollary)_Any sphere tangent to the three faces of u3 con-
taining v has center on vc', and points of tangency on Wa, vb and We.
We now seek to find a sphere inscribed in u3, which is almost as large as
the largest inscribed sphere. Let d be the distance from v to the flat triangle
t2. We see that we can chose an inscribed sphere tangent to the three facets
of u3 meeting at v, with radius r, and distance from center to v equal to 1,
such that
d= 1+r%
See Figure 14 for a picture of the two dimensional case of this construction.
Let 0 be the angle between vc' and Wa. It is obvious from the figure that
42
v
Figure 14: Constructing a large inscribed sphere of u3.
= csc(0), so that
d
r1 +cs?0)
We note that the largest inscribed sphere has radius at least that of this
1			__ sin(8)
inscribed sphere, r3 > r. Note that 0 < 0 < r/2. Hence 1+csc(8) --H sin(O)+1
sin(())/2 > kO. Hence
r3 > kdO.
Also from Lemma 10 and Figure 14,
0 = r2/d.
Recall the result of Theorem 5,
r2 > kah(b).
Hence r3 > kah(b), and we have proved the following theorem for Type A
tetrahedra.
Theorem 6 (Aspect ratio of AoCT) Let B be a box. For any tetrahedron
t inside B, the aspect ratio of t is at most k/ot, where k is a trtte constant,
and ? is the sharpest interior angle of P.
This theorem is the main result for this section, which the preceding lemmas
prove for Type A tetrahedra. We now wish to establish this theorem for the
remaining types of tetrahedra, Types B and C.
We consider Type B tetrahedra briefly because the analysis is so similar
to that for Type A tetrahedra. We perform the same construction as for Type
43
A tetrahedra, defining d and 9 as before. Consider a Type B tetrahedron,
the convex hull of vertex v on P facet F and triangle T on the other P facet
& in box b. If the distance from v to G is at least kh(b), then the analysis of
Type A tetrahedra suffices. Otherwise, we have d> kah(b) as in Lemma 6.
As before, the construction gives r3 > kr2. Because d may be closer to the
plane of T, the result of Lemma 10 does not hold, and instead we use the
result of Lemma 11. This is sufficient because we have a stronger bound on
the shortest altitude of T than when considering Type A tetrahedra, namely
att(T) > kh(b).
We now consider Type C tetrahedra. These tetrahedra arise from prisms.
Since the central vertex v is well centered in the box, and no P facet is close
to it (and no P edges are in the box), the angle ? between the shortest
segment from a triangle T to v and the normal to the plane of T is no
more than a constant slightly larger than 7r/4, in contrast with the Type B
tetrahedra. We will make several constructions as for Type A tetrahedra,
only we will chose different vertices for "v". However, our choice of "v" will
be related to the box's central vertex, so that Lemma 10 will hold for Type
C tetrahedra. Hence we have r3 > kr2 > kah(b), and we now need only
specify the constructions.
Consider the various types of prisms as in Section 7.4, having a triangular
top and bottom, and one, two, or three quadrilateral sides depending on
how many vertices of the top are shared by the bottom. Prisms formed by
clipping the convex hull of v and a P triangle have zero shared top and
bottom vertices, and will be considered last. Otherwise, if two vertices of the
top of a prism are shared with the bottom, then the prism is a tetrahedron.
We take the top to be a box triangle. We now perform the same construction
as for Type A tetrahedra, with the non-shared bottom vertex as "v".
If one vertex of the top is shared with the bottom, then the prism is split
into two tetrahedra. Define the top as in the last type of prism. The bottom
is a projection of the first onto a P facet, with center of projection v far
from the top and bottom. Recalling that ? is bounded by a constant, the
bottom altitude is at least a constant fraction of the top altitude, kah(b). For
the tetrahedron that is the convex hull of the top and a non shared bottom
vertex, we chose "v" to be the bottom vertex. For the tetrahedron that is
the convex hull of the bottom and a non shared top vertex, we chose "v" to
be the top vertex.
If no top vertex is shared with the bottom, then the prism is split into
44
three tetrahedra. We let the top be the box or P triangle. As in the previous
case, the bottom triangle has altitude at least a constant fraction of the top
altitude, k?h(b). The analysis of the one shared vertex case is sufficient for
the tetrahedra with the top as a facet and the tetrahedra with the bottom
as a facet.
For the third tetrahedron, the one that is the convex hull of an edge of
the top and an edge of the bottom, we chose the bottom vertex opposite the
two top vertices as "v". The triangle of the tetrahedron opposite t in this
subcase arises from introducing a diagonal on one of the quadrilateral sides
of the prism. If both the top and bottom are subsets of P facets, or in general
if the box is a protected edge box not containing an edge, the P edge is at
least kh(b) away from t, and hence t has altitude at least kQh(b). If the top
is a subset of a box facet, then t has altitude at least kh(b). In both the case
where the top is a subset of a P face and in the case where the bottom is a
subset of a box face, the angle between the planes of the top and bottom is
between 0 and a constant slightly larger than 3r/4, by the fact that neither
plane is close to the central vertex, the central vertex is outside of P, and
both planes pass through the box. For the definition of angle between the
planes, we chose the normal to the bottom pointing towards the interior of
the prism, and the normal to the top pointing towards the exterior of the
prism. In particular, the angle between the shortest segment from t to "v"
and the normal to t is bounded by a constant, so that again Lemma 10 holds.
9 Aspect ratio tetrahedra in A*
In the last section we established an upper bound on the aspect ratio of our
triangulation. In this section we show that we are a constant factor from
optimal, by analyzing the aspect ratio of the optimal triangulation.
We wish to show that 1/a acts as a lower bound on the aspect ratio of
any triangulation of P, where a is the smallest interior angle of P. Note that
if a > r/2, then in any triangulation of P there is trivially a tetrahedron
with aspect ratio at least k/a for some constant k, since aspect ratios are
always greater than 1. Hence, we need only prove the following theorem
Theorem 7 Let P have smallest interior angle a. If a < r/2, any three
dimensional triangulation of P has a tetrahedron with aspect ratio at least
1/ sin(a).
45
Before starting the proof, we state the following simple geometric lemma,
whose proof is omitted.
Lemma 15 Let Q be a plane in JR3, and let x be a point on Q. Let r1,r2
be two rays starting from point x such that plane Q separates them (or one
or both r1, r2 may be contained in Q). Then the angle between r1 and Q is
bounded above by the angle between r1 and r2.
Proof ofTheorem 7. We first show that for any three dimensional trian-
gulation of P, there is some tetrahedron of the triangulation with angle at
most a between two edges or an edge and a facet not containing the edge.
Consider the faces F and 0 of P that intersect at p, with the angle between
F and 0 equal to a. Note that we assume that F and 0 can see each other in
the sense of Section 2. Recall that we may always chose F to be an edge. If
F makes angle a with any edge, we let that edge be 0, otherwise 0 is facet.
In either case there is some segment in the relative interior of 0 along the
ray in 0 making angle a with F. See Figure 15. We let Cj be this segment,
and note that p is contained in 0. Now, let Cm be a segment joining a point
on F and a point on Ci so that Ci, Cm, and F make a triangle and so that F
and Cm make a right angle. Let PI be the intersection of F and Cm. Let p"
be another point on 0m close to p'. Given a three dimensional triangulation
A*, we may chose 0m close enough to p and p" close enough to PI, such that
every tetrahedron of A* containing PI' also contains p and PI.
Let A* be one of these tetrahedron containing p". Then A* contains
p, p', and also has an edge, say F', that is coincident with a subsegment
of F. Without loss of generality, Cm is moved so that PI is a vertex of A*
terminating edge F'. Let 0' be the face of A; opposite PI. Then the angle
made by 0' and F' is no larger than a, because C' is on a plane separating
F' and C1 as in the last lemma.
We now show that this small angle leads to a large aspect ratio for A*
We first bound the altitude of A* at p1, which bounds the radius of the largest
sphere inscribed in A*. If 0' is a facet then we chose H = 0', otherwise we
chose H as the unique facet of A* containing edge 0' but not edge F'. That
?s, we let H be the facet of A* that intersects F1 at exactly p. Then the
angle a' between F' and the plane of H is at most a. It is clear that since
tetrahedra are convex, the radius of the largest sphere inscribed in A* is at
most half of the altitude between p' and H, where as before P' is the vertex of
46
F
p			G1
Figure 15: Finding the tetrahedron with an angle at most a.
F' that is not p. This altitude is sin(a')/IF'I, where IF' denotes the length
of edge F'. Hence we have an upper bound on the size of the largest sphere
inscribed in A* and we now find a lower bound on the size of the smallest
sphere containing A;. Since a containing sphere must have diameter at least
the longest edge, we have that the smallest containing sphere has radius at
least IF'I/2. Combining these two shows that the aspect ratio of A* is at
least 1/ sin(a'), which is at least 1/ sin(a) as desired. I
Let T be a three dimensional triangulation. Let A(T) be its aspect ratio,
defined to be be the maximum over all tetrahedron t Ei T, of the aspect ratio
of t. We may now state the theorem concerning the optimality of the aspect
ratio of our triangulation of P.
Theorem 8 (Aspect ratio optimality)
A(AocT) < kA(A*)
where k is a true constant, independent of P and a, and AocT is the trian-
gutation of P arising from our algorithm, and A* is any other triangulation
of
47
10 Optimality of the cardinality of AocT
In this section we prove that the number of tetrahedra in our triangulation
is no more than a constant factor larger than any other triangulation with
a bounded aspect ratio. The reasoning in this section is as follows. We first
prove that any triangulation with bounded aspect ratio, which we denote A*,
has certain geometric properties concerning how fast the sizes of the tetrahe-
dra can change. We then derive lower bounds on how small the tetrahedra
of our triangulation AocT can be. We then compare the two sets of bounds
to show that our triangulation is within a constant factor of optimal.
Conformal. Any triangulation A of P that we consider in this paper
must be conformal to P. This means that the following condition must hold.
Let x be a point on a P face F. Let F' be the lowest dimensional face of
A containing x. Then F' c F. It is easy to see that the triangulation we
construct has this property; vertices of P are covered by vertices of AocT,
and edges of P are covered by edges of AocT. A more general, but equivalent
definition of conformality is the following lemma, whose proof is obvious from
the definition.
Lemma 16 If a triangulation A is conformal to P, then the following con-
dition is satisfied. Let F, G be two faces of P, and x, y two points on F, C.
Let F', G' be the lowest dimensional faces of A containing x and y. Then
F' n G' c F n 0. In particular, if F, 0 are disjoint, so are F', 0'
In this section there will be a series of positive constants c1,c2,... de-
pending only on the maximum aspect ratio A of the triangulation under
consideration. For example, the maximum number of faces incident to a
vertex of A* is bounded by a constant c1.
Characteristic length function. We define fT : P ,` JR to be the
characteristic length function of a triangulation AT. This function is defined
as follows. If x is a point of P, then we define fT(x) to be the longest edge
among all tetrahedra that contain x.
The next group of lemmas focus on a triangulation A* of bounded aspect
ratio and its characteristic length function f*. We will show below that
the longest edge of any tetrahedron containing x is bounded by a constant
times the longest edge of any other tetrahedron containing x, so that any
tetrahedron containing x can be used to bound f*(x) above and below.
48
Lemma 17 Suppose x, y C P, and suppose the two tetrahedra defining f*(x)
and f*(y) are T1 and T2 (in particular, x E T1 and y E T2). If T1 and T2
share a common edge, then
C2 f*(Y) <			< C3
Proof. Let e denote the common edge. Note that because of the aspect ratio
bound of A*, the ratio of the longest to shortest side of any simplex in A*
is bounded above and below. Hence
lel < f*(x) < C4 lel.
This also holds for y. 1
Wedge. A wedge is a set of tetrahedra, T1, ..... . , Tk, that share a com-
mon vertex.
Note that the bounded aspect ratio of A* bounds the smallest interior
angle of any simplex in A*, so that the number of simplices in a wedge is
bounded by a constant c1.
Lemma 18 Given a bounded aspect ratio triangulation A*, two points x, y,
and any two faces F, G containing x, y respectively, define H = F n G. If
H # ?, then f*(x) > C5
Proof. Note that H is a face of A*. Let v be a vertex of H. We consider
the tetrahedra of a maximal wedge with common vertex v. Since no vertex
of P is degenerate, between any two tetrahedra there is a sequence of tetra-
hedra of the wedge such that each shares an edge with its successor. Since
wedges consist of a bounded number of tetrahedra, the function f* on any
tetrahedron is within a constant multiple of f* on any another by Lemma
17. This wedge includes F and G. Finally, there are tetrahedra T1, T2, such
that T1 contains x and determines f*(x) (i.e., the longest edge of T1 is equal
to f*(x)), and similarly T2 determines f*(y). Then T1 and F have a common
vertex (since they both contain x) and therefore the lengths of their edges
are within a constant of each other by the same argument. The same holds
for T2 and G. 1
49
Geodesic distance. We let distp(F, G) denote the geodesic distance in
P between two closed subsets F, G of P. In other words, this is the length of
the shortest path in P between a point of F and a point of G. This distance
is always at least as great as the Euclidean distance.
As an example of this notation, we have the following lemma. We omit
the proof of the lemma, which follows from the meaning of aspect ratio for a
tetrahedron.
Lemma 19 Let T be a tetrahedron of a bounded aspect ratio triangulation,
let x be a vertex ofT, and let F be a subface of T not containing x. Then
distp(x, F) > ?
Similarly, the following lemma follows from elementary trigonometry ap-
plied to bounded aspect ratio.
Lemma 20 Let T be a tetrahedron with bounded aspect ratio. Let F be an
edge or facet of T, and let x be a point on F whose distance from any subface
of F is at least t. Let y be a point lying on any facet G of T, such that G
does not contain F. Then the distance from y to x is at least c7t.
Proof. An illustration of the claim made by this lemma is in Figure 16 in
the case that F is a facet. If G is a subface of F, then the lemma is trivial
(we could take c7 = 1). Otherwise, let v be the point on the boundary of
F closest to x. (Note that if F is an edge, v will be one of its endpoints).
Consider the triangle Axvy. Let zFG denote the interior angle between F
and G. Then Ixvy > min(ZFG, 7r/2). But /FG is bounded below by a
constant c8 depending on the aspect ratio of the tetrahedron. Moreover, the
vx leg of this triangle has length at least t by assumption. This means that
the distance from x to y is at least t sin c8. I
We now wish to establish the following theorem about the distance be-
tween faces in a triangulation. The theorem has two cases, depending on
whether the faces have a common point or not.
Theorem 9 Let F, G be two faces of A*. Let H = F n C. Let x, y be two
points such that x ? F, y E C. Then
1. If H # ? and distp(x,H) > t, then distp(x,y) > c9t
50
Figure 16: An illustration of the conditions stated in Lemma 20.
2. If H = ?, then distp(x,y) > c10f*(x).
Proof. This proof is broken down into eight cases, depending on whether
we are in Case 1 or 2 or the theorem, and depending on whether F is a
vertex, edge, triangle or tetrahedron. The cases are labeled Al to D2. Note
that if we are in Case 1 of the theorem, we can assume that H is a proper
subset of both F and C. If H is not a proper subset, say F = H, then
distp(x, H) = 0, and the claim is trivial. Similarly, if G = H then y ? H so
distp(x, y) > distp(x, H).
In each of the upcoming cases, let p be the shortest path in P from y to
Case Al, F zs a vertex and H $ ?. This is a trivial case described in the
previous paragraphs, since H must be the same vertex as
Case A2, F is a vertex and H = ?. Let T be the first tetrahedron on
p such that T contains F. Then T must be entered on a face F' of T not
containing vertex F (since G does not contain F). Thus, from Lemma 19,
we know that distp(F, F') > cGf*(x).
Case Bi, F is an edge and H $ ?. Then G is an edge or a triangle, and
H is a vertex v common to F and G. Note that t is bounded above by f*(x)
since t cannot be longer than the length of F. Let w be the other endpoint
of F. We take two subcases: either x is within distance c7f*(w)/2 of w or
it is further than this distance. If x is within c7f*(w)/2 of w, then we apply
Case A2 to get a lower bound of ?f*(w) on the distance from y to w. This
means that the distance from y to x must be at least ?f*(w)/2, which is at
least c11t. (Note that f*(w) > f*(x) since any tetrahedron containing x will
51
also contain w.)
The other subcase is that x is distance at least ?f?(w)/2 from w and at
least t from v. Let T be the first tetrahedron on the path p that contains
edge F. Say that tetrahedron is entered at point y' on face F'. Then F' does
not contain F. By Lemma 20 the distance from y' to x is at least
c7 min(?f*(w)/2, t).
Thus, in all cases, the distance from F to G is at least c12t.
Case B2, F is an edge and H = ?. This is handled with exactly the same
argument as Case Bi, namely, we observe that either x is within ?
of an endpoint w of F, in which case the distance bound can be reduced to
A2, or else x is further than this distance from either endpoint, in which case
we use Lemma 20. Again, in this case, we get a distance bound of c12f*(x)
on the distance from x to y.
Case Cl, F is a triangle and H # ?. Let us take two subcases:
Subcase Cl (a), H is an edge of triangle F. This means that C is a second
facet containing edge H. Then we first check whether x is within distance
min(c12, 1)t/3 from a point y' on one of the other two edges E1, E2 of triangle
F. If so (say y' is on E1), then y' is at least distance 2t/3 from H, so we
use Case Bi (applied to E1 and G) to get a lower bound of 2ci2t/3 on the
distance from y' to H, and this yields a lower bound of c12t/3 on the distance
from x to H. Otherwise, x is further than distance min(c12, 1)t/3 from all
three edges of F. Then we look at the first tetrahedron T on path p that has
F as a facet. This tetrahedron must be entered by a point y' not on F. Then
Lemma 20 tells us that the distance from y' to x is at least c7 min(c12, 1)t/3.
Subcase Cl(b), H is a vertex of triangle F. Let E1, E2 be the two edges of
F containing vertex H, and let E3 be the third edge. If x is within distance
min(c12, 1)t/3 from a point y' on either E1 or E2, we can prove the theorem
as in Subcase Cl(a) by applying the argument of Case Bi to E1 and G, or to
E2 and G. Another possibility is that x is within distance c12f*(x)/2 from a
point y' E E3, in which case we can apply Case B2 to E3 and G, since they
are disjoint. This gives a lower bound of c12f*(x) on the distance from y to
G. Otherwise, the distance from x to any of the edges of F is bounded below
by either min(c12, 1)t/3 or c12f*(x)/2. Let T be the first tetrahedron on p
that has F as a facet. Then T is entered from a point y' not on F. Finally,
52
Lemma 20 gives us a lower bound of the form
c7 min(c12f*(x)/2, min(c12, 1 )t/3)
on the distance from x to y.
Cases DJ--HD2, F is a tetrahedron. These cases are not needed by the
upcoming theorems, so we omit them and leave them to the reader. I
Theorem 10 There exist constants c13 and c14 dependent on A sitch that
for all x,y E P,
< max(ci3 f*(y), ci4 distp(x, y)).
Proof. If x and y are in a common tetrahedron of A*, then the result follows
from Lemma 18.
Otherwise, let xy denote the shortest path in P from x to y. Label
the triangles of A* crossed by xy by Fi, i = 1,2,... , m. We assume in
this proof that xy intersects each Fj in a single point; the degenerate case
of exact alignment by Fj and xy can be analyzed by considering a slight
perturbation to xy. Note any two consecutive triangles of the sequence are
in a common tetrahedron. Let Xj denote Fj A xy. By Lemma 18 we know
that f*(x) < c5 f*(xi) and f*(xrn) < c5
Define
i* =maxt?:FinF?#?)
If i* = m, then x1 and Xm are in a common wedge, and f*(x1) < c5 f*(xm).
Combining these inequalities shows that f*(x) < c13 f*(y). The other case
is that i* <m, and we note that distp(xi,xj*+1) > cio.f*(xi) by Theorem 9.
Moreover, distp(x, y) > distp(xi, xj*+1). Combining these inequalities gives
the result that f*(x) < c14 distp(x,y). I
We now observe the following characterization of our triangulation AocT
Lemma 21 Let x E P be a point lying in octree box b. Then focT(x) >
k1oth(b), where k1 is a constant.
Proof. Let T be the tetrahedron of AocT lying in b containing x. Then each
edge of T on the boundary of PAb has length at least k1?h(b) as determined in
previous sections. A case analysis verifies that every tetrahedron we construct
in b has an edge on the boundary of P A b. This includes those tetrahedra
arising from prisms. I
53
We now have a lemma about geodesic distances in boxes. Note that in a
crowded box b, it is possible to have two points x, y in b such that distp(x, y)
is arbitrarily larger than h(b). For example, P might be shaped like a coil
inside b. There are two cases, however, when we can get bounds on geodesic
distances.
Lemma 22 Let b be an uncrowded box in the octree, and let x, y be two
points in P A b. Then
distp(x,y) < 2$3h(b).
Proof. Recall that in a box, b, P A b is a single component C, so x, y are
both in this component. Either there is a vertex, edge or facet cone in the
box, or b c P. In the latter case, the distance from x to y in b is the same
as the distance from x to y in P, which is bounded by $3h(b).
The other case is that b contains a cone. Let z be a point on the apex of
the cone. Then x and y are both connected to z by a straight-line path in
so distp(x,z) < $3h(b) and distp(y,z) < $3h(b). 1
We wish to establish a series of lemmas that prove that the characteristic
length function at a point of A* is bounded by a constant times the size of the
octree box containing that point in our algorithm. The first lemma involves a
construction that is useful for each of the remaining lemmas. The remaining
lemmas prove the result for a point in a protected vertex box, protected edge
box, protected facet box, and finally for an unprotected (empty) box.
Lemma 23 Let bi be a box split in the vertex, edge or facet phase, and x a
point of P face F in b1. Then there exists a point y of P face C such that
1. G does not contain F, and
2. distp(x,y) < 6$3h(b1).
Proof. If b1 is split because of the balance condition, let bi be the crowded
box that caused it to split. There may be a long chain of boxes that split
to lead to a splitting of bi; say the chain is b2...., b1. In this chain, bt is
crowded, and it causes bi?i to split, which causes b!?2 to split, ... which
causes bi to split. Without loss of generality, we can assume that all of the
boxes ..... . , bt?i are not crowded. For example, if bt', 1' < 1 is crowded, then
we can replace 1 with 1'.
54
The sizes of these boxes must be decreasing geometrically, so that h(bi) =
2?-1h(b1), because otherwise the balance condition would not be enforced
from bj+1 to bi. This means that the Euclidean distance from a point in bi to
a point in b1 is at most 2$3h(b). In other words, i(ex(b,)) c i(ex(b1)) for
all j = 1,2,... 1. Recall that the balance condition is only enforced between
balance-adjacent boxes, and neighboring boxes are balance-adjacent only if
the smaller box contains a point of the extended box of the larger box.
Hence, the single connected component PA ex(bi) contains P A ex(bj) for all
j = 1, 2,... 1. This line of reasoning is also correct in the case that b1 is itself
crowded, in which case we take 1 = 1.
Recall that xis a point on a P face Fin PAb1. We wish to find aface C
close to F and x that does not contain F. We now have two cases, depending
on why b1 was split.
Case 1. The first case is that b1 is crowded because a face interferes with
F, and this face is visible to x. In this case, we let G be the closest interfering
face to F, and y the closest point of C to x. Then we conclude 0 does not
contain F, and distp(x,y) < 4$3h(b1).
Case 2. The second case consists of one of three subcases. This first
subcase is that b1 is split because of the balance condition. The second
subcase is that b1 is crowded because some face in b1 is interfered with, but
this face is not F. The third subcase is that b1 is crowded because some face
interferes with F, but this face is not visible to x. In the second and third
subcases 1 = 1. In any subcase, since b1 is crowded, we may let J denote a P
face in b1 such that some P face in ex(bi) interferes with it. It may be that
J is F. Let H be the face interfering with J such that the geodesic distance
between H and J in P A ex(bi) is minimum. Since P A ex(b1) is a single
component containing P A b1, there is a shortest path entirely in P A ex(b1)
from x to a point z of J in bj visible to H. This path may be considered
as a chain of maximal length edges, from x to xi, x1 to x2, ... Xm to z. A
cross-sectional view of an example path is given in Figure 17.
If the smallest dimensional face containing x1 is not F, we choose 0 to
be this face, and y to be x1. Our choice of 0 as the smallest dimensional
face containing y ensures that 0 does not contain F. Since there is always a
straightline path in ex(bi) from x to x1 = y, we conclude that distp(x, y) <
3$3h(b1).
Otherwise the smallest dimensional face containing xi is F, that is xi is
in the relative interior of F. Because the path construction uses maximal
55
b
b1
Figure 17: An example path in the proof of Lemma 23. Here all P vertices
represent edges, and edges represent facets, except for the P vertex J =
Note b3 is crowded and b1 is split by the balance condition.
length edges, the only way this can happen is if x1 is actually z. Hence
there is a straightline path on F from x to z, and J is contained in F. See
Figure 18. We choose G to be H, the face interfering with J, and y to be
the closest point of il to z. Since G does not contain J, and J in contained
in F, it follows that G does not contain F. Since there exist straightline
paths from x E b1 to z ? bi to y E b1, with h(b?) < h(b1), we conclude that
distp(x,y) < 6?h(b1). I
Lemma 24 For any vertex v of P in a box b, f*(v) < c15h(b).
Proof. Since the protecting and merging steps for a protected vertex box
change the size of any box by at most a factor of two, we can assume b is a
box at the end of a vertex phase. Let b1 be the parent of b.
We apply Lemma 23 with x = v in box b1. For any vertex, if another
face does not contain it, the two are disjoint. Hence the point y on P face G
is disjoint from F = v and from the construction distp(v, y) < 4$3h(b). We
apply Theorem 9 to conclude f*(v) ? 4cV03h(b) I
Lemma 25 (Corollary) For any point p in a box bv protected for vertex v,
< c16h(bv)
Proof. In a protected box Euclidean distance to the apex is the same as
geodesic distance, so distp(p, v) < $3h(bv). From Theorem 10 and Lemma
24, f*(p) < max(c13. f*(v), c14 distp(p, v)) < max(c13c15, $3ci4)h(b). I
56
bv
x
v
y
Figure 18: An illustration of an instance of Case 2 of the proof of Lemma
23 leading to the case in the proof of Lemma 26 where F' and G' have a
common vertex v.
Lemma 26 For anypointp ofP edge F in a box b, f*(p) < c17h(b).
Proof. As in the vertex case, Lemma 24, since the splitting steps within a
protecting step subdivide any box at most once, we may assume that the
parent b1 of b was split in the vertex or edge phase. Assume that the box
containing pat the end of the octree algorithm is not vertex protected, since
otherwise Lemma 25 applies.
We apply the Lemma 23 with x = p. Unlike the vertex case, Lemma
24, we may have that G and F intersect nontrivially. To complete the proof
we need to consider the particular triangulation A*. Let F' and G' be the
smallest dimensional faces in A* such that x E F' c Fand y E G' c G. If F'
and G' are disjoint, we can apply Theorem 10 as in Lemma 24 and the proof
is complete. Otherwise, by conformality, F'A G' is a P vertex v. Let bv be
the box containing v at the end of the vertex phase. See Figure 18. Since v is
centered in bv, and by assumption p is outside bv, distp(p, v) > h(bv)/8. But
h(bv) > c15f*(v) by Lemma 24. Hence by Theorem 9 Case 1, distp(p,y) >
8?q5f*(v). Recalling from Lemma 23 that distp(p, y) < 6$3h(bi) completes
the proof. I
Lemma 27 (Corollary) For any pointp in a box bF protectedfor edge F,
? c18h(bF).
Proof. The proof is identical to the proof of Lemma 25. I
57
Lemma 28 For any point p of P facet F in a box b, f*(p) < ci9h(b).
Proof. As in the previous lemmas, we may assume that the parent b1 of b
was split in the vertex or edge phase, and not a protecting step. Assume
that the box containing p at the end of the edge phase, is not vertex or edge
protected, since otherwise Lemma 25 or Lemma 27 applies.
We apply Lemma 23 with x = p. As in Lemma 26, F and G may not
be disjoint, so we must consider the particular triangulation A*. Let F'
and 0' be the smallest dimensional faces in A*, such that x ? F' c F and
y ? 0' c 0.
If F' and 0' are disjoint, we can apply Theorem 10. Otherwise, by
conformality, F' n 0' is a P vertex or a subset of a P edge. Let v be
the closest point of F' A 0' to p. Let bv be the box containing v at the
end of the edge phase. Since v is well centered in protected boxes, and
p is not in a protected box by assumption, distp(p, v) > h(bv)/8. De-
pending on whether bv is vertex or edge protected, we rely on Lemma 25
or Lemma 27 to show h(bv) > min(c16, ci8)f*(v). Hence, by Theorem 9
Case 1, distp(p, y) > 8? min(c16, ci8)f*(v). Recalling from Lemma 23 that
distp(p, y) < 6$3h(b1) completes the proof. I
Lemma 29 (Corollary) For any point p of a box bF protected for facet F
< c20h(b?).
Proof. The proof is identical to the proof of Lemma 25. I
Lemma 30 For any point p in an unprotected box b (containing no P faces),
< c21h(b).
Proof. As in the previous lemmas, we may assume that the parent b1 of b
was split in the vertex, edge, or facet phase. Either b1 is split by the balance
condition, or b1 is actually crowded, but no P face is in the octant of its child
b. By imitating the path construction of Lemma 23, with the only difference
being that x is not on any P face, we may identify as 0 the first face we
meet on a path from p to a face in the crowded box b1. That is, we may
show that distp(p,y) < 3$3h(b1), and h(b?) ? h(b1), for some point y on
face 0 in box b,, at the end of the octree generation algorithm. Depending
on whether b? is a vertex, edge or facet protected box, we apply one of
the previous three corollaries, Lemma 25, Lemma 27 or Lemma 29 to show
58
max(ci6,ci8,c2o)f*(v) < h(b?). We apply Theorem 10 to p, so that f*(p) <
max(c13f*(y), ci4distp(p, y)) < max(ci3/ min(c16, cis, c2o), 3Vm3ci4)h(b). 1
We may now summarize the results of the previous lemmas.
Theorem 11 Let x E P lie in a box b. Then f*(x) < c22h(b), where c22 is a
constant depending on the aspect ratio bound of A*.
The following result combines the bounds derived for focT and f*.
Theorem 12 For all polytopes P and constants A, there exists a constant
c', dependent only on A and a, such that for any z E P and any triangulation
A* of P with maximum aspect ratio A,
focT(z) >
Proof. This follows directly from Lemma 21 and Theorem 11 I
Theorem 13 For all polytopes P and constants A, there exists a constant
c", dependent only on A and a, such that
I AocT I < c"I A* ,
where I A* denotes the number of tetrahedra of A* and similarly for I AocT I
Proof. Let AT be a triangulation with bounded aspect ratio. We derive
some inequalities that depend on aspect ratio; these inequalities apply to
either A* or AocT. By the definition of fT and bounded aspect ratio A,
there are constants c23 and c24 depending only on A such that
c23f?(x)3 < vd(T?) < c??f?(x)3
for all x in the interior of Ti and for any tetrahedron Tj.
Thus
c23 ? ITi fTk)3 dx < c24
for every tetrahedron. Summing this chain of inequalities over all tetrahedra
in the triangulation, we get
c23o+IATI<Jpf?(lx)3dx=?J?tf?(lx)3dx<c24IATk
59
This holds for all triangulations, so in particular for AocT and A*. By
Theorem 12
C23 AocT < IP ?oc??x?3 dx < f? ?????x??3 dx < c??/(&)? A*
11 Conclusions
An important open question is the running time. We demonstrated a run-
ning time bound of O(?n log n) for our algorithm, which is inferior to Bern,
Epstein and Gilbert's running time of O(n log n + ?) for two dimensional re-
gions. However, we have recently found out by private communication with
those authors that this bound does not actually hold for the algorithm that
they propose for two-dimensional nonconvex polygons. Therefore, it is not
clear what is the best possible running time for triangulation of nonconvex
polyhedral regions in either two or three dimensions.
It would be interesting to generalize our algorithm to work for higher
dimensional regions. Also, the non-degeneracy conditions on the input region
could be relaxed.
Another open question concerns optimizing the tetrahedra for properties
other than aspect ratio. For example, is there an algorithm for triangulating a
three-dimensional nonconvex polyhedral region to maximize inscribed sphere
radius of the tetrahedra?
Finally, we have plans to implement this algorithm. The two-dimensional
analog of this algorithm has been implemented by the first author in Ct+.
Acknowledgement
Thanks to Paul Chew of Cornell for discussing various triangulations with
us. Thanks to Marshall Bern, David Eppstein and John Gilbert of Xerox for
help with their earlier work. Thanks to Joe Mitchell of Cornell for helping
with the plane sweep algorithm.
60
References
I. Babuska and A. K. Aziz [1976], On the angle condition in the finite element
method, SlAM J. Num. AnaL 13:214--H226.
M. Bern, D. Dobkin, and D. Eppstein [1991], Triangulating polygons without
large angles, preprint.
M. Bern, H. Edelsbrunner, D. Eppstein, 5. Mitchell, and T. 5. Tan [1991],
Edge-insertion for optimal triangulations, preprint.
M. Bern, D. Eppstein, J. Gilbert [1990 Provably Good Mesh Generation,
Xerox Palo Alto Research, Proc. 1990 Sympos:um on Foun-
dations of Computer Science.
M. Bern and D. Eppstein [1991], Mesh Generation and Optimal Triangula-
tion, preprint.
G. Carey, M. Sharna, and K. Wang [1988], A Class of Data Structures for
2-D and 3-D adaptive mesh refinement, Int. J. Numer. Meth. in Enyr.
26:2607--H2622.
B. Chazelle, L. Palios [1989 , Triangulating a Nonconvex Polytope, Proc. 5th
ACM Symposium on omputational Geometry. 393--H399.
L. P. Chew [1989a], Constrained Delaunay Triangulation, Algorithmica 4:97--H
108.
L. P. Chew [1989b], Guaranteed-Quality Triangular Meshes, C. 5. Cornell,
TR 89--H983.
H. Edelsbrunner [1987], Algorithms in combinatorial geometry, Springer Ver-
lag, Berlin.
J. Hauser and C. Taylor [1986], Numerical Grid Generation in Computational
Fluid Dynamics (Proceedings), Pineridge Press, Swansea, U.K.
G. L. Miller, S.-H. Ten , and 5. A. Vavasis [1991], A Unified Geometric
Approach to Graph eparators, Proc. 32nd Symp. Foundations of Comp.
Sci.
G. L. Miller and W. Thurston [1990], Separators in Two and Three Dimen-
sions, Proc. 22nd Symp. on Theory of Computing, 300--H309.
G. L. Miller and 5. A. Vavasis [1990], Density Graphs and Separators, De-
partment of Computer Science, Cornell, Tech Report 90-1169, and to
appear, Symposium on Discrete Algorithms.
W. F. Mitchell [1987], A Comparison of Adaptive Refinement Techniques
for Elliptic Problems, Department of Computer Science, University of
Illinois at Urbana-Champaign, Report No. UIUCDCS-R-87-1375.
61
W. F. Mitchell (1988], Unified Multilevel Adaptive Finite Element Methods
for Elliptic Problems, Department of Computer Science, University of
Illinois at Urbana-Champaign, Report No. UIUCDCS-R-88- 1446.
D. Moore, J. Warren [1990j, Multidimensional Adaptive Mesh Generation,
Department of Computer Science, Rice University, Rice COMP TR 90-
106.
F. Preparata and M. Shamos [1985], Computational Geometry: An Introduc-
tion, Springer Verlag, New York.
62
