This algorithm works for a triangulation, which is no loss of generality: a plane graph can be triangulated in O(n log n) time (see text).
|
The idea: suppose a vertex v of the triangulation is deleted, as are its incident edges.
The result is a face that is the union of the triangles incident to v.
Q: the original set of triangles Q, If a point is in one of Q', it can be found in Q by looking through each of the triangles of Q. If v has degree K, this takes O(K) time. |
|
Since the retriangulated graph is planar, and has 1 fewer vertex, it has fewer faces and edges.
The approach:
find many vertices, all of low degree, to delete. Retriangulate, and recursively build a data structure for the coarser graph.
If a fixed fraction of the vertices, say n/10, can be removed each time, then the recursion depth is O(log n), as is the query time.
One more condition: the set of vertices must be independent,
with no incident edge in common.
Otherwise, the faces to be retriangulated and searched will be not be O(K)
when degree K vertices are removed.
The goal:
Given a planar triangulation, find a large set of independent vertices, all of low degree.
Since a planar graph has at most 3n-6 edges, the average degree of a vertex is at most 6.
As with Markov's inequality,
half the nodes must have degree less than twice average,or in general,
(1-z)n of the nodes must have degree <= 6/z.
Pick z=2/3: n/3 nodes have degree <= 9.
Now pick greedily from such nodes: place a "mark" on all nodes of degree > 9, and repeat the following until no unmarked node remains:
pick any unmarked node; put the node in independent set I; mark the node and its neighbors.
The nodes in I are independent;
each node picked results in a graph with 10 fewer unmarked nodes,
and the number of unmarked nodes is at least n/3, so
|I| > n/30.
(29/30)k n <1,or
k=log30/29n = O(log n).
Preprocessing:
each independent set takes linear time to find, and the number of nodes decreases geometrically, so O(n) time is needed overall.
Space: also O(n).
But: the constants are large!
Better constants are possible, including use of the four-color theorem.
Suppose a planar subdivision has faces that are x-monotone:
the intersection of a vertical line with a face is an interval. Define the << relation for faces in a monotone subdivision as follows:
Since (the interiors of) faces do not meet, and the faces are monotone, this relation is well-defined. |
|
Lemma: The << relation is acyclic.
Proof.
Suppose
F0 << F1 <<...Fk-1 << Fk = F0is a shortest cycle. Suppose Fi is the face whose rightmost point is leftmost. Since Fi and Fi+1 meet a common vertical line, and Fi and Fi-1 meet a common vertical line, by choice of Fi there is a vertical line meeting all three. But this implies Fi-1 << Fi+1, contradicting the assumption that the cycle was shortest, QED.
For any i with 1 <= i < m and any vertical line L,
L intersect union1 <= j <=i Fjis below
L intersect unioni < j <=m Fj,and these two sets meet at a common point. Let si be the union of all such points (for the given i).
si is called a separator.
Then si is connected, and is the union of bounding edges of the faces, and, roughly speaking,
Fj << si for j<=i, and
Fj >> si for j>i.
For a given point p, to locate p with respect to si:
use binary search to find the intersection of si with the vertical line through p.
Check if p is above or below that intersection.
To search:
do binary search on the separators: check locate p with respect to sn/2,
then s3n/4 or sn/4, etc.
In this way,
find the separators si and si+1 with si << p << si+1.
Then the edge above p on si+1 bounds the face containing p.
The query time is O(log n)2;
Although separators are not disjoint, the storage can be made to be O(n).
The constants are not bad.
With "fractional cascading", O(log n) query time is possible.
From a known point, walk through triangulation to query point.
For fat triangles or uniformly distributed points, O(sqrt(n)) time.
|
For a set S of n sites (points) in the plane,
A point p will have 1 or more nearest sites: For p a point in the plane, define
|
|
| VR(B) = |
{p | NN(p)=B}, the set of points with the same nearest sites B For most B, VR(B) is empty For B with one element, VR(B) is a region in the plane |
|---|
| Voronoi diagram VD(S): | the partition of the plane {VR(B) | B a subset of S} |
|---|
Suppose S={s1, s2}.
Then VD(S) is has three regions:
|
|
| empty disk: | a disk that contains no sites in its interior. | |
|---|---|---|
| Nearest Neighbor disk: |
for a point p, the disk centered at p bounded by
a circle through NN(p). That is, an empty disk with a site on its bounding circle. |
As p moves away from its nearest site s,
the NN disk expands but keeps its bounding circle on s.
For p moving in VR(s1, s2), the bounding circle must always contain s1 and s2.
VR(B) is nonempty:
iff there is a point p with NN(p)=B;
iff there is an empty disk with the sites of B on its boundary.
a·b is axbx + aybySuppose p is in VR(S), so s1 and s2 are equidistant from p
Euclidean norm ||a|| =sqrt(a2), where a2=a·a.
a·(b+c) = (a·b) + (a·c):
dot product distributes over addition.
a·b == 0 iff a and b are perpendicular.
||p-s1||2 == ||p-s2||2
(p-s1)·(p-s1) == (p-s2)·(p-s2) ==>
p2 + s12 - 2 p·s1 == p2 + s22 - 2 p·s2, ==>
s12 - 2 ·s1 == s22 - 2 p·s2 ==>
p·(s2 - s1) == (s22 - s12)/2
This is the equation of a line, perpendicular to s2 - s1.
Note that (s2+s1)/2, the midpoint of the segment s1--s2, satisfies the equation.
A point in VR(s1)==VR({s1}) satisfies
p·(s2 - s1) < (s22 - s12)/2
|
S={s1, s2, s3}
Here
|
|
Any p in VR(S) solves
p·(si - s1) == (si2 - s12)/2, for i = 2 or 3,two linear equations in two unknowns.
This has a unique solution unless s2 -s1 = z(s3-s1) for some z;
this means that the points of S are collinear.
p·(t - s) < (t2 - s2)/2,that is, p is in a halfplane bounded by the perpendicular bisector of s and t.
Any point p that satisfies
p·(si - s) < (si2 - s2)/2 for si in S, si!=sis in VR(s): p has s uniquely closest among all sites.
That is, VR(s) is an intersection of halfplanes; if bounded, a convex polygon.
|
The Delaunay subdivision DG(S) of S has vertices == the sites |
|
There is a duality correspondence between the parts of VD(S) and DG(S), where
vertices <--> faces
edges <--> edges
faces <--> vertices
the mapping is incidence preserving:
if vertex a is dual to face a', and a is an endpoint of edge e, and e is dual to e',
then e' is a bounding edge
of a'.
More formally,
For B a subset of S with VR(B) nonempty,
| cell(B): the "relative interior" of the convex hull of B |
|---|
| |B|=1 : | the site in B |
|---|---|
| |B|=2 : | the line segment between the sites, not including the sites |
| |B|=3 : | the interior of the triangle with endpoints in B, assuming not collinear. |
| |B|>3 : | a convex polygon, the interior of the convex hull of B, assuming not collinear. |
The subvision DG(S) can be defined as {cell(B) | VR(B) is not empty}.
| "face of" : |
If a is an endpoint of line segment b, or a is a bounding edge of the convex polygon b, say that a is a face of b. |
|---|
The incidence preserving condition amounts to:
cell(B) is a face of cell(B') iff VR(B') is a face of VR(B).
For example,
a site a is an endpoint of a Delaunay edge {a,b} iff VR(a,b) is a face of VR(a).
An unbounded edge of the Voronoi diagram is dual to an edge of the convex hull of S: one Voronoi vertex is "at infinity" == the corresponding "NN disk" is a halfplane.
VR(B) nonempty for |B|>3 implies there is some point equidistant to all sites of B, so the points of B are cocircular.
General position, for these diagrams: no four sites are cocircular, as well as no three sites collinear
If all sites are on a common line, the Voronoi edges are parallel lines, and the convex hull of the sites is 1-dimensional.
Theorem: DG(S) is is a plane graph.
That is, the edges of DG(S) don't cross.
A more general fact, that implies this theorem, is:
Theorem: If B and B' are subsets of S with B!=B', then cell(B) does not meet cell(B').
Note that this is true even if B is a subset of B': for example, cell(a,b) is a line segment that does not include its endpoints a=cell(a) and b=cell(b).
Lemma: Suppose D and D' are closed disks whose interiors meet, and line L passes through the intersection points of their bounding circles. Then D-D' and D'-D are separated by L.Proof: "left to the reader"
Proof of Theorem: Suppose B and B' are two subsets of S with nonempty Voronoi regions. Since VR(B) is not empty there is an empty disk D with B on its boundary, and an empty disk D' with B' on its boundary. Since D and D' are empty, every site of B is in D-D', unless the site is one of the intersection points of the bounding circles of D and D'. Since the analogous claim is true for B', all sites of B and B' are either separated by L, or are one of the intersection points. Hence the relative interiors of the convex hulls of B and B' are separated by L, and disjoint, QED
Bounds follow from the planarity of DG(S):
Theorem: The number of Voronoi vertices (and Delaunay triangles) is no more than 2n-4 and the number of Voronoi or Delaunay edges is no more than 3n-6.
Proof: From Euler's relation, as described in Lecture 2, applied to DG(S). QED
Call a triangulation locally Delaunay if for any triangles abc and abd, d is not inside the circumcircle of abc.
Clearly, any Delaunay triangulation is locally Delaunay. In fact, the converse holds:
Theorem:A locally Delaunay triangulation is Delaunay.
First,
Lemma: Given two disks D and D', with line L through the intersection points of their bounding circles, letIn short, D contains D' on "its" side of L, and vice versa.Then
H be the open halfplane bounded by L and containing D-D', and H' be the open halfplane bounded by L and containing D'-D.
D' intersect H is a subset of D intersect H, and
D intersect H' is a subset of D' intersect H.
Proof: also "left to the reader".
Proof of the theorem: Given a Delaunay triangle abc and a site d, it must be shown that d is not inside the circumcircle of abc.
Pick a point p inside abc so that the line segment dp contains no sites. Suppose dp crosses edge bc. The triangles meeting dp can be written
Ti == aiai+1ai+2, for i=1 to m-2, for some m,where
a1=a, a2=b, a3=c, and am=d.The edges crossed by dp can be written
ei == aiai+1, for i=2 to m-1.
By the first "disk" lemma and the locally Delaunay condition, the line through edge ei separates ai-1 from ai+2, and by the above "disk" lemma, the intersection of the Delaunay disk of Ti with dp is farther from d than the intersection of the Delaunay disk of Ti+1 with dp.
Since the Delaunay disk of Tm-2 does not contain d in its interior, none of the disks of Ti do, for i=1..m-2. Hence the disk of T1=abc does not. QED
|
Suppose a triangulation is not locally Delaunay,
so some triangles abc and abd have site d inside the circumcircle of abc.
The quadrilateral acbd is convex, and has two triangulations: T = {abc, abd}, with edge ab, the original one, and |
|
Facts:
T' is Delaunay (since the other one isn't)
Angle increase: The minimum angle of T' (of the six angles) is larger than the minimum angle of T.
| Edge flip: |
The transformation of the triangulation
deleting {abc, abd} and inserting {cda, cdb},that is, swapping cd for ab. |
|---|
Note: if {a,b,c,d} are cocircular, both triangulations are Delaunay and have the same smallest angle.
...for computing a Delaunay triangulation:
Find any triangulation of S; repeat find non-locally-Delaunay pair of triangles flip edge until the triangulation is locally Delaunay
Minimum angle of a triangulation: the smallest angle of any triangle in the triangulation
Theorem: The edge flipping algorithm terminates, and the Delaunay triangulation maximizes the minimum angle
Proof: Given a triangulation T, let
| Angle(T) = | (a1, a2,...) be the sequence of measures of angles of triangles of T, sorted into nondecreasing order. |
|---|
Suppose T is not Delaunay, and so not locally Delaunay, and so there is a triangulation T' resulting from an edge flip.
Then Angle(T') is lexicographically larger than Angle(T).
Hence if T is not Delaunay, Angle(T') is larger than Angle(T),
so the Delaunay triangulation is lexicographically maximum.
Since the number of triangulations is finite,
and each triangulation generated is distinct,
the algorithm terminates.
QED
The Delaunay triangulation is generally easier to compute. Some algorithms
Given sets L and R, S= L union R, separated by a vertical line, and their triangulations DG(L) and DG(R), compute DG(S).
There are RR edges, LL edges, and cross edges.
All new Delaunay edges in the merge are cross edges: if an RR edge is Delaunay now, it was before.
Some triangles of DG(L) are deleted:
those with circumcircles containing points of R.
...and symmetrically for DG(R).
The cross edges are ordered along the split line, and consecutive edges share a vertex.
The cross edges are built from bottom to top, starting from a cross edge that is a convex hull edge.
Conceptually, maintain an empty disk on the current cross edge {a,b}, pushing it up, but keeping a and b on the bounding circle.
The first site this "rising bubble" hits gives the next cross edge.
That site is on an edge incident to a or b. To find it,
walk around
the edges incident to a, deleting those bounding triangles conflicting
with b. Find the first edge {a,a'} not deleted.
Similarly walk around the edges incident to b, finding edge {b,b'}.
Either b or b' is hit first by the "rising bubble", yielding
the next cross edge {a,b'} or {a',b}.
Since O(n) work is done in the merge, O(n log n) is needed overall.
Dwyer has shown that O(n log log n) is needed for points uniform in the square, using bucketing.
There are at least four viewpoints of the sweep algorithm.
When viewed as a Delaunay triangulation algorithm:
Sweep a vertical line from left to right.
For a point p on the sweep line,
the tangent empty circle at p:
is empty: has no sites inside;
is tangent to the sweep and left of it;
contains a site (on the circle, not inside)
S(p): the sites on the tangent empty circle.
The points on the sweep form equivalence classes based on S(p).
Intervals of points have the same S(p) with |S(p)|=1, separated by points with |S(p)|=2, corresponding to Delaunay edges.
| site events: | The sweep meets a site | |
| circle events: | the sweep contains a point with S(p)=3 yields a Delaunay circle corresponds to the shrinking down of an interval |
Suppose the circle event has S(p)={a,b,c},
the previous interval had label {a}, bounded by points with
labels {a,b} and {a,c}
New label {b,c} has possible events with its neighbors on the sweep: schedule these.
At a site event, suppose site a is met by the sweep,
in an interval labeled by {b}, and bounded by {b,c} and {b,d}
Then a new interval, labeled by {a} and bounded by {a,b} on both sides,
will be seen.
Check for circle events for {a,b} and {b,c}, and for {a,b} and {b,d}
Some scheduled circle events will not actually occur: for example,
if a site event splits an interval, that interval cannot shrink to a point.
if a circle event shrinks an interval to a point, that interval cannot expand to shrink a neighboring interval to a point.
Such "false alarm" events can be detected and removed when real events occur, and there are a bounded number of them.