OUTPUT-SENSITIVE CONSTRUCTION OF CONVEX HULLSByTIMOTHY MOON-YEW CHANB.A., Rice University, 1992A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFDOCTOR OF PHILOSOPHYinTHE FACULTY OF GRADUATE STUDIESDEPARTMENT OF COMPUTER SCIENCEWe accept this thesis as conformingto the required standardTHE UNIVERSITY OF BRITISH COLUMBIAOctober, 1995© Timothy Moon-Yew Chan, 1995In presenting this thesis in partial fulfilment of the requirements foran advanceddegree at the University of British Columbia, I agree that the Library shall make itfreely available for reference and study. I further agree that permission for extensivecopying of this thesis for scholarly purposes may be ‘granted by the head of mydepartment or by his or her representatives. It is understoodthat copying orpublication of this thesis for financial gain shall not be allowed withoót my writtenpermission;(Signature)Department of Comput-r cneThe University of British ColumbiaVancouver, CanadaDate October 4, 1995DE-6 (2)88)AbstractThe construction of the convex hull of a finite point set in a low-dimensional Euclideanspace is a fundamental problem in computational geometry. This thesis investigatesefficient algorithms for the convex hull problem, where complexity is measured as afunction of both the size of the input point set and the size of the output polytope.Two new, simple, optimal, output-sensitive algorithms are presented in two dimensions and a simple, optimal, output-sensitive algorithm is presented in three dimensions.In four dimensions, we give the first output-sensitive algorithm that is within a polylogarithmic factor of optimal. In higher fixed dimensions, we obtain an algorithm thatis optimal for sufficiently small output sizes and is faster than previous methods forsublinear output sizes; this result is further improved in even dimensions.Although the focus of the thesis is on the convex hull problem, applications of ourtechniques to many related problems in computational geometry are also explored, including the computation of Voronoi diagrams, extreme points, convex layers, levels inarrangements, and envelopes of line segments, as well as problems relating to ray shooting and linear programming.11Table of ContentsAbstract iiTable of Contents iiiList of Tables vList of Figures viAcknowledgements vii1 Introduction 11.1 The Convex Hull Problem 11.2 The Number of Faces of the Convex Hull . . . 41.3 Output-Sensitive Algorithms 61.4 Previous Convex Hull Algorithms 81.5 Results in This Thesis 112 Two- and Three-Dimensional Convex Hulls 142.1 A Simplified “Ultimate Planar Convex Hull Algorithm” 142.1.1 The prune-and-divide algorithm in the plane 152.1.2 Analysis of the prune-and-divide algorithm in the plane . . 182.2 An Optimal Convex Hull Algorithm in Two and Three Dimensions . . . 202.2.1 The group-and-wrap algorithm in the plane 212.2.2 The group-and-wrap algorithm in three dimensions 232.2.3 Refinements of the group-and-wrap method 272.3 Application: Lower Envelopes of Line Segments in the Plane . . . . 291113 Four-Dimensional Convex Hulls 323.1 Preliminaries on the Divide-and-Conquer Construction of Convex Hulls 333.1.1 The upper hull 333.1.2 Facets and their duals 343.1.3 Cuttings for divide and conquer 363.2 A Prune-and-Divide Convex Hull Algorithm in Four Dimensions . 363.2.1 Primal dividing 383.2.2 Dual pruning 413.2.3 Converting from dual to primal 423.2.4 Specializing for d = 4 433.2.5 The prune-and-divide algorithm in four dimensions . . . 463.2.6 Analysis of the prune-and-divide algorithm in four dimensions . 483.3 Application: Three-Dimensional Voronoi Diagrams 514 Higher-Dimensional Convex Hulls 544.1 Ray Shooting Queries 554.2 Linear Programming Queries 624.3 A Convex Hull Algorithm and an Extrema Algorithmin Any Fixed Dimension 664.4 Application: Convex Layers and Depths 704.5 Further Applications: Levels in Arrangements andLinear Programming with Violations 744.6 The Prune-and-Divide Convex Hull Algorithm in Even Dimensions . 804.7 Appendix: Using Randomization in Linear Programming Queries . . 845 Conclusion 87Bibliography 91ivList of Tables4.1 Known data structures for ray shooting queries in polytopes and linearprogramming queries 575.2 Summary of output-sensitive results for the convex hull problem 88vList of Figures1.1 The boundary of the convex hull of a planar point set 22.2 Pairing and pruning points in the plane 162.3 Wrapping a set of rn/mi convex polygons of size m 222.4 The lower envelope of a set of line segments 303.5 (a) The upper hull of a point set in E3 and (b) the vertical projection ofits facets 343.6 A simple region S of the point set in Figure 3.5 383.7 (a) The region S and (b) the intersection of its interior with the interiorof the simple region S from Figure 3.6 433.8 Tracing a boundary component (d = 3) 453.9 The Voronoi diagram of a planar point set . . 524.10 The convex layers of a planar point set 714.11 The boundary of the 1-level in an arrangement of lines 75viAcknowledgementsI wish to thank my research supervisor, Jack Snoeyink, for providing many helpfulcomments on various drafts of this work, sharing many enlightening discussions, andintroducing me to the exciting subject of computational geometry. His guidance, encouragement, and friendship are invaluable to me. I am also much indebted to DavidKirkpatrick and Nick Pippenger for their advice and support. I would also like to thankmy University Examiners, Maria Klawe and Maurice Queyranne.To the entire Computer Science Department, I would like to extend my gratitude forproviding a wonderful environment for my graduate study at the University of BritishColumbia. The financial support of the Killam Trusts and the Natural Sciences andEngineering Research Council of Canada is gratefully acknowledged.Finally, I am deeply thankful to my parents for their love, patience, and encouragement, and for all that they have done for my sake; this work is dedicated to them.viiChapter 1IntroductionComputational geometry [Ede87, Mu193, O’Ro94, P585] studies the design and analysisof algorithms for solving geometric problems. One central problem that has receivedconsiderable attention in the area is the problem of constructing convex hulls. Theimportance of the problem stems not only from its many applications (such as to patternrecognition, statistics, and image processing) but also from its usefulness as a tool forsolving a variety of problems in computational geometry. The concept of convex hulls iswell-studied in mathematics and is appealing both mathematically and intuitively: givena point set, its convex hull is simply the smallest convex set that encloses it.In the following section, we define the convex hull problem in more precise terms.In Section 1.3, we describe the type of algorithms seeked in this thesis, namely, output-sensitive algorithms. Section 1.4 gives a brief review of previous convex hull algorithms.Finally, in Section 1.5, we outline the results obtained in this thesis.1.1 The Convex Hull ProblemLet P = {pi,. . . , p} be a set of n points in d-dimensional Euclidean space Ed, whered 2 is a small fixed constant. An assumption we make throughout the thesis is that theinput points are in general position. This means that no d+ 1 points of P lie in a commonhyperplane; in certain places, we also require that no d points of P lie in a vertical hyperplane. (In this thesis, the terminology “vertical,” “above/below,” “upward/downward,”Chapter 1. Introduction 2Figure 1.1: The boundary of the convex hull of a planar point set.and “highest/lowest” are always with respect to the last coordinate.) There are generalperturbation methods [EM9O, EC92] to cope with point sets not in general position. Theidea behind these methods is to eliminate degenerate configurations by applying an arbitrarily small perturbation to the input; the perturbation is done only conceptually, soeach primitive operation on the perturbed input has to be simulated.Recall that a set S ç Ed is convex if for every p, q 5, the line segment iscontained in S. The convex hull of F, denoted by conv(P), is the “smallest” convex setcontaining F, that is, the intersection of all convex sets containing P. Equivalently, it canbe defined as the intersection of all halfspaces containing P. Alternatively, it is the setof all convex combinations of F: conv(P) = {D cp2 : > c = 1, a, . . . , > O}.An intersection of a finite set of (closed) halfspaces is called a polytope. By a wellknown fact [Grü67, MS71], conv(P) is a polytope. For example, if d = 2, it is a convexpolygon (see Figure 1.1). We can represent the boundary of the polygon by its sequenceof edges in, say, counterclockwise order.In higher dimensions, we can describe the boundary of the polytope by its faces:Chapter 1. Introduction 3Suppose that P E’ is a polytope with a non-empty interior. If h is a hyperplane thatintersects the boundary of P but not its interior, then h fl P is a face of P. A face isa i-face (0 j < d) if it has dimension j—that is, if it is contained in some j-flat butnot in a (j — 1)-flat. A (d — 1)-face is called a facet, a (d — 2)-face is called a ridge, a1-face is called an edge, and a 0-face is called a vertex. The faces of P, together withthe empty set 0 and P itself, form a lattice under inclusion, and the union of the facetsis the boundary of P. Thus, we can represent the boundary of P by the Hasse diagramcorresponding to this lattice of faces.The vertices of conv(P) belong to the point set P and are also called the extremepoints of P. The following statements are equivalent given a point p F: (i) p isextreme in P; (ii) p conv(P— {p}); (iii) conv(P— {p}) conv(P); and (iv) there existsa vector e E’ such that p > q for all q e P (q p). In (iv), we say that p isextreme/maximal along direction (or minimal along —c). The convex hull of P is thesame as the convex hull of the extreme points.“Constructing the convex hull of F” then means producing a complete representationof the boundary of conv(P)—the ordered sequence of edges if d 2, or the facial latticestructure if d > 2. This construction problem is the main topic of this thesis.A closely related and equally important problem is the computation of an intersectionof a set of halfspaces H = {h1, . . . , h,}. Denote this intersection by fl H. The well-knownlinear programming problem seeks a point in fl H that is maximal along a given direction.In contrast, the halfspace intersection problem asks for a complete representation of fl H.(In the context of linear programming, each halfspace in H is called a (linear) constraint,a point in fl H is called a feasible solution, and the intersection fl H is called the regionof feasibility.)To solve the halfspace intersection problem, we can first find a point in the interiorof fl H (if the intersection is non-empty), using linear programming techniques [Meg84].Chapter 1. Introduction 4By translation, we can move this point to the origin o so that each halfspace h is ofthe form {x e Ed : x < 1} for some vector tj e Ed. Then there is a one-to-onecorrespondence [Ede87] between the j-faces of fl H and the (d — j — 1)-faces of theconvex hull conv({i,. . ., }). This shows that computing intersections of halfspacesand computing the convex hulls are in fact equivalent problems.Notice that we have just applied a form of duality (or polarity) when we map ahalfspace {x e Ed : • x < 1} to a point e Er’. Duality [CGL85, Ede87] is extremelyimportant in computational geometry as it allows one to transform a problem involvingpoints to a problem involving halfspaces/hyperplanes and vice versa. Sometimes we maygain more insight into the geometry of a problem by examining the problem in both itsprimal and dual setting.For a more in-depth exposition of the concepts discussed so far, we refer the readerto the standard computational geometry textbooks [Ede87, Mu193, O’Ro94, PS85].1.2 The Number of Faces of the Convex HullGiven a set P of n points in Ed, how many faces can conv(P) have? For d 2, thenumber is clearly at most 2n, since conv(P) is a polygon with at most n vertices andat most n edges. For d = 3, Euler’s formula implies that the number of edges andfacets relate linearly to the number of vertices, so the total number of faces is also 0(n).However, for d = 4, 5, the number can be quadratic in n, and for d> 6, it can be as highas e(nld/2i), as the following theorem shows.Theorem 1.2.1 Let d be a fixed constant and n be any number.(i) Every n-point set P ç Ed in general position has at most 0(n[°/2J) faces in itsconvex hull.Chapter 1. Introduction 5(ii) There exists an n-point set P C Ed in general position that has (nL”/2i) faces inits convex hull.Proof: We first prove the upper bound in (i). Since the number of faces of the convexhull is at most 2d times the number of facets, it suffices to bound the number of facets inconv(P). Since the number of j-faces is at most (‘4) = O(n’) (as a j-face is incidentto j + 1 vertices), it suffices to bound the number of facets in terms of the number of([d/2J — 1)-faces in conv(P). If we dualize points to halfspaces as described in Section 1.1,then our task is to bound the number of vertices in terms of the number of Id/21 -facesin an intersection fl H of n halfspaces. This is done by charging vertices to Ed/21 -facesin fl H as follows.Consider a vertex v of fl H. By the general position assumption, there are preciselyd edges incident to v. Orient the edges in an upward direction. Then we can find eitherEd/21 edges oriented towards v, or Fd/21 edges oriented away from v. In the former case,these Ed/21 edges define a Ed/21-face a that has v as its highest vertex. In the lattercase, the Ed/21 edges define a Ed/21-face a that has v as its lowest vertex. In eithercases, we charge v to a. Since each face has a unique highest/lowest vertex (assumingno degeneracies), at most two vertices are charged to a Ed/21-face. Thus, the number ofvertices is no more than twice the number of Ed/21 -faces, and (i) is proved.To prove the lower bound in (ii), we choose a point set P on the moment curve:P = {M(t1),. . . , M(t,j}, where t1,. . . , t, are n distinct real numbers and M(t) =(t,t2,.. . , td) é EE for any t. Observe that a hyperplane can pass through at most dof the points in P, since a function that is linear in t, t2,. . . , is a polynomial in tof degree < d and thus has at most d distinct zeroes. Therefore, the points in P arein general position. We now show that there are at least ([d22j)= (nLd/2i) faces inconv(P). In fact, we prove a stronger property (the [d/2J -neighborly property) aboutChapter 1. Introduction 6the polytope conv(P): every [d/2J-subset Q of P defines a ([d/2J — 1)-face of conv(P).The proof of this property is not difficult. Given a subset Q P with IQ [d/2j,we need to show that there is a hyperplane h such that all points of Q lie on h and allpoints of P— Q lie strictly on one side of h. This follows if we define the functionF(t)= fi (t-t)2M(t1)Qand observe that F is linear in t, t2,. . . , td, F(t) = 0 for all M(t) e Q, and F(t) > 0for all other ti’s. UFor the point set P used in the above lower-bound argument, conv(P) is called acyclic polytope. Using a tighter upper-bound analysis, McMullen [McM7O, MS71] in factshowed that the cyclic polytopes (or more generally, the [d/2J-neighborly polytopes)attain the maximum number of j-faces among all d-dimensional, n-vertex polytopes, forany 0 <j <d. This result is the celebrated Upper Bound Theorem in polytope theory.1.3 Output-Sensitive AlgorithmsSince the number of faces of the convex hull can be as large as e(nLd/2i)by Theorem 1.2.1,any algorithm for constructing the convex hull needs to spend at least (nLd/2i) time in theworst case just to write out a representation of the hull. Hence, 2(nL’/i) is automaticallya lower bound for the convex hull problem.However, for point sets that occur in practice, the convex hull has usually a smallnumber of faces. For example, if the points in P are chosen independently at random froma uniform distribution on a convex r-gon in E2, then the expected number of faces of theconvex hull is O(rlogn) [RS63]. For points chosen uniformly and independently from theinterior of a d-dimensional hypercube, the expected number of vertices of the convex hullChapter 1. Introduction 7is 0(logd_l n) [BK+78]. For points from the interior of a d-dimensional ball, the expectednumber of faces is O(n(d_1)/(1)) [Ray7O], which is far from the pessimistic 0 (n1d/2J)bound. Even when the points are chosen on the surface of the unit paraboloid so thatall points are extreme (which is the case in applications to Voronoi diagrams [Aur9l]),the expected number of faces in the convex hull is still 0(n) if the points are lifted froma uniformly distributed point set in the interior of a (d — 1)-dimensional ball [Dwy9l].In analyzing the performance of convex hull algorithms, it is therefore desirable totake into account the number of faces of the convex hull. From now on, this number willbe denoted f (although in a few occasions f is also used to denote a facet; this shouldbe clear from the context). Our goal is to develop efficient algorithms for constructingconvex hulls, where the measure of efficiency is asymptotic worst-case running time as afunction of both n (the input size) and f (the output size). In general, algorithms withcomplexity measured in terms of not only the input size, but also the output size, aresaid to be output-size sensitive, or output-sensitive for short. Ideally, one should need tospend less work when the output size is small.Besides the number of faces f, a related parameter one could use to measure thecomplexity of convex hull algorithms is the number of hull vertices or extreme points.This is denoted by h (when it does not represent a hyperplane or halfspace; again, thisshould be clear from the context). While the number f ranges from 0(1) to 0(nL°/2i),the number h ranges from 0(1) to n. For d < 3, f and h are asymptotically equivalent,so for historical reasons, we will use h instead of f. For d 4, we will return to using f.Remark: As stated in the beginning of Section 1.1, we frequently need to assume thatthe input points are in general position for our algorithms to work properly; when thisassumption does not hold, we have to rely on perturbation methods [EM9O, EC92].However, perturbing the input points may cause the output size to increase, so we needChapter 1. Introduction 8to redefine f to be the maximum number of faces of conv(P’) over all “perturbations” F’of P. For d < 3, we can simply redefine h to be the number of input points on theboundary of the convex hull (since these are the only points that can become a vertexafter perturbation). These redefinitions are important only when there are a large numberof degeneracies. For our two- and three-dimensional convex hull algorithms, we are ableto handle degenerate cases directly, so redefining h is unnecessary.1.4 Previous Convex Hull AlgorithmsThe convex hull problem has had a long history going back to the beginning of computational geometry and has been an intensively studied subject even up to the presentday. Early papers dealt primarily with the planar case d = 2. The first O(nlogn)-timealgorithm in E2 was given by Graham [Gra72] and is commonly known as Graham ‘s scan.In terms of n alone, Graham’s algorithm is (asymptotically) optimal, since an i2(n log n)lower bound for the convex hull problem can be obtained by a reduction from sorting.Using more complex arguments [Ben83, Yao8l], the (nlogn) lower bound applies alsoto the weaker problem of identifying the vertices of convex hull (the extreme points).However, all of these lower bound arguments assume that the number of hull vertices his at least a fraction of n, so it is conceivable that there is an algorithm that beats Graham’s scan if h is substantially smaller than n. This was indeed shown by Jarvis [Jar73],who gave a simple O(nh)-time algorithm, dubbed Jarvis ‘s march. This algorithm is thusoutput-sensitive.After the publication of Graham’s and Jarvis’s algorithm, a number of different convexChapter 1. Introduction 9hull algorithms (as well as variants on previous algorithms) were proposed in E2. We mention here two divide-and-conquer algorithms [PS85]—”MergeHull” and “QuickHull”—modeled after the sorting algorithms MergeSort and QuickSort. The former divide-and-conquer algorithm, due to Preparata and Hong [PH77], runs in O(nlogn) time. Thelatter algorithm, discovered independently by several researchers around the late 1970s,runs in O(nh) time in the worst case, but is usually faster in practice (as its name suggests). At the time, no asymptotic improvement to the original bounds by Graham andJarvis was known, so the true complexity, in terms of n and h, of the convex hull problem in E2 remained unresolved. Finally, in 1986 Kirkpatrick and Seidel [KS86] gave thedefinitive answer in a paper entitled “The Ultimate Planar Convex Hull Algorithm?”by giving an output-sensitive Ofri log h)-time algorithm and providing a matching lowerbound. It would appear that Kirkpatrick and Seidel’s optimal algorithm is thus the“ultimate” convex hull algorithm for d = 2—or is it?The convex hull problem in its three-dimensional setting (d = 3) has also been studiedintensively. Preparata and Hong [PH77] presented the first O(nlogn)-time algorithmin E3, based on divide-and-conquer. The first output-sensitive algorithm is the giftwrapping method of Chand and Kapur [CK7O], which works in arbitrary dimensions andis a generalization of Jarvis’s march in two dimensions (although historically, Chand andKapur’s method appeared before Jarvis’s). As analyzed by Swart [Swa85], the methodruns in O(nh) time.For a long time the gift-wrapping method was the only output-sensitive algorithmknown for three-dimensional convex hulls Then, in their seminal work on randomizationtechniques in computational geometry [CS89], Clarkson and Shor gave an algorithm withan optimal O(n log h) expected running time, where the expectation is based solely on therandom choices made by the algorithm [Mu193]. Afterwards, Edelsbrunner and Shi [ES91]Chapter 1. Introduction 10proposed a deterministic algorithm with an O(n log2 h) running time, by following theparadigm of Kirkpatrick and Seidel. An optimal 0 (n log h)-time deterministic algorithmwas finally obtained when Chazelle and Matouek [CM95] applied recently-developedderandomization techniques to Clarkson and Shor’s convex hull algorithm. (A differentrandomized 0(n log h) algorithm was recently reported by Clarkson [C1a94] and, according to his paper, can also be derandomized.) The resulting algorithm is not very practical,since derandomization techniques tend to be complicated, using tools like the methodof conditional probabilities. The problem of finding a simple optimal output-sensitivealgorithm for computing convex hulls in E3 thus remained.For dimensions d > 4, much attention was directed to devising efficient worst-caseconvex hull algorithms. The gift-wrapping method by Chand and Kapur [CK7O] wasshown to run in O(n[d/2H1)time [Swa85] in the worst case. Seidel [Sei8l] improved thistime bound to O(nFd/21), using a different approach called the beneath-beyond method.In a later paper [Sei86J, Seidel exploited a shelling order to obtain a second algorithmwith an O(nL4/2J log n) worst-case running time. The randomized incremental construction technique of Clarkson and Shor [CS89] yields an algorithm with expected 0(nL/2i)time, which is optimal for worst-case output. Another randomized algorithm with thesame complexity was given by Seidel [Sei9l]. A deterministic O(n[d/2i )-time algorithmwas finally obtained by Chazelle [Cha93bJ; his method is based on Clarkson and Shor’srandomized solution combined with new ideas on derandomization.Although Chazelle’s algorithm has settled the complexity of the convex hull problemin arbitrary fixed dimension for worst-case output, the output-sensitive complexity is farfrom resolved. Among the algorithms we have discussed, only the gift-wrapping methodand Seidel’s second deterministic algorithm are output-sensitive. For an f-face output(recall that the number f is 2(1) and e(nLd/2i)), Swart proved an 0(nf)-time boundChapter 1. Introduction 11for the gift-wrapping method, and Seidel proved an O(n2 + f log n)-time bound for hisalgorithm. Note the substantial improvement of these bounds over Chazelle’s ofriLd/2i)bound for small values of f.In terms of n and f, Q(n log f + f) is the oniy lower bound known (the first term isa consequence of Kirkpatrick and Seidel’s two-dimensional lower bound). It is suspectedthat there is an algorithm with running time close to 0 (i-i log f + f), but finding oneappeared difficult and remained an outstanding problem, even for d = 4. There weresome improvements to the upper bound: Matouek [Mat93] showed that the running timeof Seidel’s method [Sei86] can be brought down to 0(n22/(Ld/2i+1)+e + f logn) if datastructures for linear programming queries are used. (Here and throughout this thesis,> 0 denotes an arbitrarily small but fixed constant.) The n factor can even be reducedto log°’ n using Matouek’s static structures. Further improvements seem to requirenew ideas.Remark: In this brief look at previous convex hull algorithms, we have not mentionedresults (e.g. [BS78, Dwy9l]) that assume a certain probability distribution on the input,since we are not focussing on “average-case” complexity. Moreover, we have discussedonly sequential algorithms; see [AGR94] for a recent work on parallel convex hull algorithms in a fixed dimension. We have not examined the dynamic maintenance problem [OvL8l], nor considered the construction of convex hull of geometric objects otherthan point sets.1.5 Results in This ThesisIn this thesis, we obtain new methods for the output-sensitive construction of convexhulls. An outline of our results is as follows. In the planar case, we discover two simpleChapter 1. Introduction 12algorithms that compute the convex hull in optimal O(nlogh) time. The first one canbe viewed as a simplification of Kirkpatrick and Seidel’s “ultimate” algorithm and thesecond one uses a grouping idea to speed up Jarvis’s march. Considering the long historyand the fundamental nature of the planar convex hull problem, it is surprising thatthese two simple algorithms have been left undetected. The second algorithm is alsointeresting from a practical viewpoint, as it does not require a linear-time subroutine forfinding medians, unlike Kirkpatrick and Seidel’s original solution. More importantly, thisalgorithm can be extended to yield an optimal O(n log h)-time method for constructingconvex hulls in E3. This 3-d algorithm does not require the complex derandomizationtools used by Chazelle and Matouek’s optimal algorithm, and uses relatively simpledata structures, namely, Dobkin-Kirkpatrick hierarchies [DK83, DK9O]. These resultsare described in Chapter 2.Then in Chapter 3, we present what we regard as the central part of the thesis: anoutput-sensitive convex hull algorithm in E4 that runs in O((n + f) log2 f) time. Thisis the first algorithm in four dimensions that is within logarithmic factors of optimalover the whole range of output sizes f. This result has an important consequence asit immediately leads to an output-sensitive algorithm for computing three-dimensionalVoronoi diagrams, which have numerous applications [Aur9l, 0BS92]. Although our 4-dconvex hull algorithm is the most intricate algorithm we study, its guiding principle isnevertheless the same simple divide-and-conquer strategy used by our first planar convexhull algorithm.In Chapter 4, we enter into higher-dimensional space. We observe that the giftwrapping method can be improved using the data structures for ray shooting queries inpolytopes developed by Agarwal and Matouek [AM93] and refined by Matouek andSchwarzkopf [MS93]. Together with the grouping idea from our second planar convexhull algorithm, this ‘implies an O(n log f + (nf)11/(Ld/2J+1) log°’ n) time bound for theChapter 1. Introduction 13construction of convex hulls in Er’. If f = O(nh/1d/2i / logK n) for a sufficiently large K,then the O(ri log f) term dominates; hence, our method is optimal for small output sizes.Furthermore, our time bound improves Matouek’s previous O(n22/(Ld/2H1) 1og°’ n +flogn) bound for sublinear output size f, i.e., for f = O(n/1og’n). We manage toimprove this result further in even dimensions by combining these higher-dimensionaltechniques with our 4-d divide-and-conquer algorithm. The running time obtained isO((n + (f)l_l/[d/2]+fnl_211d/) log° n).Our techniques are not limited to the convex hull problem. We apply these to manyproblems in computational geometry throughout the thesis. These applications includeVoronoi diagrams (as mentioned above), envelopes of line segments, enumeration ofextreme points, convex layers and depths of point sets, levels in arrangements of hyperplanes, and linear programming with few violated constraints. We hope that these“digressions” demonstrate the important role that convex hulls have in computationalgeometry.Chapter 5 summarizes our work and concludes with open problems and remarks ondirections for further research.Remark: Most of the results of this thesis have been presented in conference papers, andtheir full versions have been submitted for publication in journals. See [CSY95b] for thesimplification of Kirkpatrick and Seidel’s algorithm and its extension to four dimensions.A dual version of the 4-d algorithm in the halfspace-intersection setting is describedin [CSY95a]. Most of our higher-dimensional results appear in [Cha95b]; specializationto 2-d and 3-d can be found in [Cha95a].Chapter 2Two- and Three-Dimensional Convex HullsIn this chapter, we present two O(n log h)-time convex hull algorithms in the plane E2,the second of which is also extended to E3 with the same optimal time complexity.Although algorithms with this time complexity were known, our methods are simplerthan the previous and also illustrate the basic ideas to be used in the rest of this thesis forconstructing convex hulls in higher dimensions. In particular, the first planar convex hullalgorithm we present, which can be considered as a simplification of the original optimalmethod of Kirkpatrick and Seidel [KS86], serves as the basis of the four-dimensionaloutput-sensitive algorithm in the next chapter.2.1 A Simplified “Ultimate Planar Convex Hull Algorithm”This section describes a simple O(n log h) convex hull algorithm in the plane. Since ouralgorithm can be viewed as a simplification of Kirkpatrick and Seidel’s planar convexhull algorithm, we first sketch here their method for comparison.Given an n-point set P C E2, we want to construct the convex hull of P. It sufficesto compute just the upper hull of F, consisting of the sequence of hull edges that havean upward normal vector. Then the lower hull can be computed in a similar manner byreflection and the convex hull can be obtained by joining these two hulls.Kirkpatrick and Seidel’s algorithm constructs the upper hull of P as follows: (i) find apoint p E P with the median x-coordinate, (ii) compute the edge of the upper hull14Chapter 2. Two- and Three-Dimensional Convex Hulls 15that intersects the vertical line through p (x[p1] < x[p’] < x[p2]), and (iii) recursivelycompute the upper hull of all points left of (and including) Pi and the upper hull of allpoints right of (and including) P2To find the edge ij (the bridge) that intersects a given vertical line in step (ii),Kirkpatrick and Seidel used a prune-and-search procedure, similar to the prune-and-search linear programming algorithm of Dyer [Dye84] and Megiddo [Meg83b, Meg84].(In fact, bridge-finding can be formulated as a linear program in dual space.) Here is ahigh-level description what is involved in this prune-and-search procedure: First, pointsare paired, the slope of the line through each pair is calculated, and the median slope mis computed. Then the upper-hull vertex p with a supporting line of slope m is found.A comparison involving Pm and the given vertical line is then performed, which allowsone point in half of the pairs be pruned. This step eliminates 1/4 of the points and theprocedure is repeated.This ends our brief sketch of Kirkpatrick and Seidel’s algorithm. As a summary,we can say that their algorithm has two levels: the lower level is a prune-and-searchprocedure, and on top of that is a divide-and-conquer method. Our main observation isthat we can get a simpler algorithm if we combine these two levels into one, i.e., if we usepruning directly for divide-and-conquer rather than for searching. As a result, we canskip the step that computes the point p with median x-coordinate and avoid actuallysearching for the bridge at each recursive step.2.1.1 The prune-and-divide algorithm in the planeWe now give the details of our simplified planar convex hull algorithm. As in Kirkpatrickand Seidel’s algorithm, only the upper hull of the given n-point set P ç E2 is computed.We first pair the points of P arbitrarily and calculate the slope of the line through eachpair. We then find the median slope m and compute the upper-hull vertex Pm that has aChapter 2. Two- and Three-Dimensional Convex Hulls 16slope=mRmedian slope mFigure 2.2: Pairing and pruning points in the plane. Points marked L belong to P,points marked R belong to Pr, and points marked X belong to neither sets.supporting line of slope m; this vertex can be computed by taking the maximum along aprojection of P parallel to m. The x-coordinate of p is then used to divide P into twoparts: P, which contains p and all points to its left, and Pr, which contains pr-,., and allpoints to its right.Now, if a pair has slope less than m, then the right point in the pair cannot participatein the upper hull of P and thus can be pruned from P. Similarly, if a pair has slopegreater than m, then its left point in the pair cannot participate in the upper hull of Prand can be pruned from Pr. Since half (n/4) of the pairs have slope less than the medianm and half have slope greater than m, pruning ensures that P and Pr each containat most 3n/4 points. We then recursively compute the upper hull of P and Pr. SeeFigure 2.2 for an example.Lupper hull,,Chapter 2. Two- and Three-Dimensional Convex Hulls 17The pseudocode of the algorithm is given below. For convenience, we assume thatthe leftmost and rightmost points p and Pr of P have been identified and we let n be thecardinality of the set P = P—{p, Pr } instead. In the interest of practical efficiency, line 1has been added to the algorithm; it does not affect asymptotic worst-case performance.Algorithm DivideHull2d(P, p, Pr)[Given n-point set P C E2 and points P,Pr e E such that x[pj] <x[p] <X[Pr] forall p P, return the sequence of edges of the upper hull of P = P U {p,pr}. I1. discard points from P that lie below pp2. if P = 0 then return (p)if P’ = {p} then return (hp,3. arbitrarily choose [n/2j disjoint pairs {{s1,t1},. . . , {Sn/2j, }} from F’and order each pair so that x[s] <x[t24. let m = (y[t] — y[Sj]) / (x[t] — x[s]), i = 1,..., [n/2Jand m median of (mi,.. . , m[n/2j)5. let p = point in P that maximizes y[pm] — m X[Pm]6. letP’={pEP:x{p]<x[pmj}—{tj:mjm}P={pEP’:x[p]>x[pm]}—{sj:mj>m}7. if p Pr then return DivideHu112d(Pi,p,pr)if p Pe then return DivideHu112d(P,p,pr)otherwise return the concatenation ofDividellull2d(P7, p, p) and DivideHull2d(P, pm, Pr)Remark: It is not difficult to modify Dividellull2d() to work for point sets P not ingeneral position. When there are more than one point in P that maximize ypm]—mx[pmjin line 5, we simply pick the leftmost one. When two points in a pair share the samex-coordinate, we can eliminate the bottom one.Chapter 2. Two- and Three-Dimensional Convex Hulls 182.1.2 Analysis of the prune-and-divide algorithm in the planeLet T(n, h) be the running time of algorithm Dividellull2d() on a point set with n + 2points (i.e., n points excluding p and Pr) and h +1 upper-hull vertices (i.e., h upper-hulledges). By noting that median-finding (line 4) can be done in linear time, we obtain thefollowing recurrence for T(n, h), where c denotes a constant:c ifn<1T(n h) <T(ne, h) + cn if Ti> 2 and h,- 0— T(n,-, h) + cn if n> 2 and h 0T(ri, h) + T(n,-, h,-) + cn if n> 2 and h, h,- 1for some 0 <flj, n,- < f3n/41 and h, h,- > 0 with n + n,- <n and h + hr = h.Using the concavity of the logarithm, one can then prove that T(n, h) = O(n log h) byinduction. Here, we observe an alternative proof that is perhaps simpler as it avoids theuse of induction. The proof is more general and provides better insight into recurrencesof this kind by examining their recursion trees.Let T be a rooted tree in which each node v is assigned a cost c(zi’) e [0, oo). We saythat the cost function c is a-fading for a constant a e (0, 1) if c(i) < a c(v) for everynode 1 and its parent v. As part of the analysis of their 3-d output-sensitive convexhull algorithm, Edelsbrunner and Shi [ES91, Lemma 3.1] proved that the total cost insuch a tree is asymptotically bounded by the per-level cost times the logarithm of thenumber of nodes. Their proof uses a path compression operation that transforms T into abalanced tree. We give a simple, short proof of their result that avoids path compressionaltogether; we then improve the bound to depend on the number of leaves rather thanthe number of nodes.Lemma 2.1.1 In a recursion tree T with m nodes and £ leaves and an a-fading costfunction c, if the sum of the costs at each level is bounded by C, then the sum of the costsof all nodes in T is (i) at most C(1og1,m+2) and (ii) at most C(log1/£+ 1 + 1/(1 — a)).Chapter 2. Two- and Three-Dimensional Convex Hulls 19Proof: Number the levels of the tree 0, 1, 2, . . . with the root at level zero. Letk = [log11mj. The sum of the costs at levels 0, 1,. . . , k is bounded by C(k + 1)C(log11 rn+ 1). Furthermore, by the o-fading property, each node on a level greater thank has cost bounded by Cak <C/rn; hence, the sum of the costs at level k + 1, k + 2,is bounded by C. Part (i) follows.To prove part (ii), we choose k = [logia £] instead. As before, the sum of the costsat levels 0, 1,. . . , k is bounded by C(k + 1) < C(log11£ + 1). Thus, we just have toaccount for the costs of nodes at levels greater than k. Note that each node belongs tosome root-to-leaf path in T. By the o-fading property, the sum of the costs at levelsk + 1, k + 2,... along such a path is bounded by__CCa4 +Cac+2+...1—ce (1—cSince there are root-to-leaf paths in total, the sum of the costs at levels k + 1, k + 2,...is bounded by C/(1 — ce). Part (ii) follows. 0With Lemma 2.1.1, it is now easy to show that the running time of algorithmDividellull2d() is 0(nlogh). Consider the recursion tree generated by the calls toDividellull2dO. It ‘is clear that the sum of the costs at each level of the tree is boundedby cn and that the cost function satisfies the (3/4)-fading property. Since the number of leaves is at most h (as a new edge is discovered at every leaf), Lemma 2.1.1(u)immediately implies that the total cost of the algorithm is bounded by cn log413 h + 0(n).The storage requirement of the algorithm is clearly linear. We have thus shown:Theorem 2.1.2 Algorithm Dividellull2d() computes the (h + 1)-vertev upper hull ofan (n + 2)-point set P ç E2 in 0(nlogh) time and 0(n) space.Remarks:1. Compared to the algorithm by Kirkpatrick and Seidel, Dividellull2d() is fasterChapter 2. Two- and Three-Dimensional Convex Hulls 20by a constant factor: if finding the median of n numbers takes bn time, then in the worstcase, Kirkpatrick and Seidel’s algorithm spends 3bn log2 h+ 0(n) time and our algorithmspends bn log413 h + 0(n) 1.2bn log2 h + 0(n) time in median-finding (which is themost costly operation in both algorithms).2. Besides viewing DivideHull2d() as a simplification of Kirkpatrick and Seidel’salgorithm, one can also view Dividellull2dQ as a variant of the “QuickHull” algorithm [PS85], since QuickHull recursively uses an extreme vertex to divide a convex hullinto two. But as pruning is not done in QuickHull, its worst-case complexity can bee (nh). We learned recently that Wenger [Wen94] has proposed a randomized version ofQuickHull that performs pruning. His algorithm, with an 0(n log h) expected runningtime, is similar to ours, except that finding the median slope in line 4 is replaced byrandomly selecting a slope. In the next section, we give another 0(nlogh) algorithmthat avoids median-finding but is deterministic.2.2 An Optimal Convex Hull Algorithm in Two and Three DimensionsWe now give a different 0(n log h)-time convex hull algorithm in the plane, along withits extension in three dimensions. In the worst case, this 2-d algorithm is faster than themethod from the previous section since this algorithm does not perform median-findingoperations. The extension in 3-d is also simpler than the previous optimal derandomization method of Chazelle and Matouek [CM95]. Our idea here is to improve Jarvis’smarch and the gift-wrapping method by applying a common grouping trick. This grouping idea can be applied to other problems besides the construction of convex hulls, andwe will consider one example later in Section 2.3. Some variants of the method for convexhulls are discussed in Section 2.2.3.Chapter 2. Two- and Three-Dimensional Convex Hulls 212.2.1 The group-and-wrap algorithm in the planeLet P C E2 be a set of n 3 points. The algorithm Jarvis’s march [Jar73, O’Ro94, PS85]computes the h vertices of the convex hull one at a time, in counterclockwise (ccw) order,by a sequence of h wrapping steps, as follows: if Pk—1 and Pk are the previous two verticescomputed, then the next vertex Pk+1 is set to be the point p e P that maximizes theangle LPk_1PkP with p Pk. One wrapping step can obviously be done in 0(n) time byscanning all n points; with an appropriate initialization the method constructs the entireconvex hull in 0(nh) time.We observe that a wrapping step can be done faster if we preprocess the points.Choose a parameter rn between 1 and n and partition P into rn/mi groups each ofsize at most rn. Compute the convex hull of each group in 0(m log rn) time by, say,Graham’s scan [Gra72]. This gives us rn/mi possibly overlapping convex polygons eachwith at most rn vertices, after a preprocessing time of 0((mlogm)) = 0(nlogm).Now, a wrapping step can be done by scanning all En/mi polygons and computingtangents or supporting lines of the polygons through the current vertex Pk, as shownin Figure 2.3. Since tangent finding takes logarithmic time for a convex polygon bybinary or Fibonacci search [CD87, PS85j (the dual problem is to intersect a convexpolygon with a ray), the time required for a wrapping step is then 0( log m). Ash wrapping steps are needed to compute the hull, the total time of the algorithm becomes0(nlogm + h( logm)) = 0(n(1 + h/rn) logrn).The following is the pseudocode of the algorithm just described. The procedurealways runs within 0(n(1 + H/rn) log m) time and successfully returns the list of edgesof conv(P) in ccw order when H> h.Chapter 2. Two- and Three-Dimensional Convex Hulls 22Algorithm GroupHull2d(P, m, H), where P C E2, 3 m n, and H 11. partition P into subsets F1,. .. ,En/rn1 each of size at most m2. for i = 1,..., rn/mi do3. compute conv(P) by Graham’s scan and store its vertices in an arrayin ccw order4. P0 (0, —cc)5. joi +— the rightmost point of P6. fork=1,...,Hdo7. for i = 1,..., rn/mi do8. compute the point qj E P that maximizes7pklpkqj (qi pk)by performing a binary search on the vertices of conv(P)9. Pk+1 — the point q from {q1, . . . ,qn/m]} that maximizes Lpk_lpkq10. ifpk+1 =p then return..., Pk—lPk, PkP1)11. return incompleteFigure 2.3: Wrapping a set of rn/mi convex polygons of size m.Chapter 2. Two- and Three-Dimensional Convex Hulls 23By choosing m = H, the complexity of the algorithm is then O(n(1 + H/rn) log m) =O(ri log H). Since the value of h is not known in advance, we use a sequence of H’s to“guess” its value as shown below (the same strategy is used in Chazelle and Matouek’salgorithm [CM95]):Algorithm Groupllull2d(P), where P ç E21. fort=1,2,... do2. L +— Groupllull2d(P, m, H), where m = H = min{22t,n}3. if L incomplete then return LThe procedure stops with the list of hull edges as soon as the value of H in the for-loopreaches or exceeds h. The number of iterations in the loop is log log hi (using base-2logarithms), and the t-th iteration takes 0(nlogH) = 0(n2t) time. Therefore, the totalrunning time of the algorithm is n2t) = O(n2El0l0I1) = 0(n log h). Thestorage requirement is clearly linear.Theorem 2.2.1 Groupllull2d() computes the h-vertex convex hull of an n-point set P CE2 in O(nlogh) time using 0(n) space.Remark: We can handle point sets that are not in general position as follows: whenthere are more than one point q that maximize the angle tp1pkq in line 9 ofGroupHull2d(P, rn, H), pick the point q that is farthest from pk; use the same rule tobreak ties in line 8.2.2.2 The group-and-wrap algorithm in three dimensionsLet P C E3 be a set of n 4 points. Assuming that the points are in general position,we know (by Euler’s formula) that there are precisely 2h — 4 facets (triangular faces)of the convex hull if there are h hull vertices. It suffices to construct these 2h — 4 hullChapter 2. Two- and Three-Dimensional Convex Hulls 24facets; with the aid of a dictionary, we can easily generate the list of vertices and edges,together with the facial lattice structure of conv(P), in additional 0(hlogh) time.The higher-dimensional analogue of Jarvis’s march is Chand and Kapur’s giftwrapping method [CK7O, PS85, Swa85], which computes the hull facets one at a timeas follows: from a given facet f, we generate its three adjacent facets f3 by performinga wrapping step about each of the three edges e3 of f (j = 1, 2,3). Here, a wrappingstep about e3 is to compute a point Pj e P that maximizes the angle between f andconv(e3 U {p3}) with P3 e3. Since such a step can be done in 0(n) time, we can findthe facets adjacent to f in 0(n) time. Assuming an initial facet f is given (which canbe found in two wrapping steps), a breadth-first or depth-first search can then generateall facets of the convex hull. Using a dictionary to detect duplication, we can ensurethat each facet is processed once. This implies that the algorithm performs 3(2h — 4)wrapping steps and thus runs in 0(nh) time.We can use the same grouping idea from the previous subsection to improve the timecomplexity to optimal 0(n log h) while maintaining linear space. The calls to Graham’sscan (line 3 of GroupHull2d(P, m, H)) are now replaced by calls to Preparata and Hong’sthree-dimensional convex hull algorithm [PH77], which has the same complexity. Toextend line 8 to three dimensions, we need to calculate tangents or supporting planes of3-dimensional polytopes through a given line. In order to obtain the same running time asin the previous subsection, we need a method to perform each of these tangent operationsin logarithmic time. A data structure that can do precisely this is the polyhedral hierarchyof Dobkin and Kirkpatrick [DK83, DK9O].A polyhedral hierarchy can be defined as a monotone sequence of 3-dimensional polytopes P1 C P2 C ... C P with the property that each connected component (or cap) of— Pk is of constant complexity. Each Pk is called a level of the hierarchy. If P1 hasChapter 2. Two- and Three-Dimensional Convex Hulls 25constant size and P = 7’, then the sequence is said to be an inner hierarchical representation of P. Similarly, if P1 = P and P has constant size, then it is an outer hierarchicalrepresentation of P. For any 3-dimensional polytope P with m vertices, Dobkin andKirkpatrick showed that an inner/outer hierarchical representation of 0(log m) levelsexists and can be computed in 0(m) time. The polytope at each level is not explicitlystored in the representation; instead, pointers between two adjacent levels are providedso that one can easily traverse up or down the hierarchy.The hierarchical representation provides a very useful data structure for manipulatingwith polytopes in 3-space. For instance, by “walking up” the inner hierarchy, we can findthe tangent of a given polytope passing through a given line in 0(log m) time. (Thebinary-search solution to the tangent-finding problem for convex polygons can be seenas an implicit use of the hierarchy.) By “walking down” the outer hierarchy, we can alsosolve the dual problem of intersecting a polytope with a ray in 0(logm) time. There aremany other applications of the Dobkin-Kirkpatrick hierarchy; for example, we mentiona polylogarithmic algorithm by Eppstein [Epp9l] for detecting whether three polytopesin E3 have a common intersection.We can now give the pseudocode of our three-dimensional convex hull algorithm. Byusing the Dobkin-Kirkpatrick hierarchical representation to store the polytopes in line 3(which require only linear-time preprocessing), we can perform line 11 in logarithmictime for each polytope conv(Pj. The analysis of this algorithm is thus identical to thatof the two-dimensional algorithm.Chapter 2. Two- and Three-Dimensional Convex Hulls 26Algorithm GroupHull3d(P, m, H), where P C E3, 4 < m < n, and H 11. partition P into subsets Fi,. . . , .P1n/m each of size at most m2. for i = 1,..., rn/mi do3. compute conv(P) by Preparata and Hong’s algorithm and store it ina Dobkin-Kirkpatrick hierarchy4. F, Q {fo}, where fo is some initial facet of conv(P)5. fork=1,...,2H—4do6. ifQ=OthenreturnF7. picksomefeQandsetQ—Q—{f}8. let e3 be the edges of f (j = 1,2,3)9. forj=1,2,3do10. for i = 1,..., rn/mi do11. compute the point qj e P that maximizes the angle between f andconv(e3 U {qi}) by searching the hierarchy of conv(P)12. P2 — the point q from {qi,.. . ,q1/]} that maximizes the angle betweenf and conv(e3 U {q}) (q ç’ e)13. f3 — conv(e3 U {p3})14. if f3 ‘ F then15. F-FU{f}, Q÷—QU{f}16. return incompleteWe can use a queue or a stack to implement Q and a dictionary to implement F. Asthere are only 0(h) dictionary operations, they can be carried out in 0(hlogh) time. Infact, more clever implementations of the gift-wrapping method via a shelling order [Sei86]replace the need for dictionaries with just a priority queue.As before, we choose the group size m = H and guess the value of h with a sequenceof H’s:Algorithm GroupHull3d(P), where P ç E31. fort=1,2,... do2. L — Groupllull3d(P,m,H), where m = H = min{22t,n}3. if L incomplete then return LChapter 2. Two- and Three-Dimensional Convex Hulls 27Theorem 2.2.2 GroupHull3d() computes the h-vertex convex hull of an n-point set P çE3 in 0(nlogh) time using 0(n) space.Remark: We can handle point sets that are not in general position as follows: In line 8,let e3 = J6 with a3 and b3 oriented in a counterclockwise order around f. When thereare more than one point q that maximize the angle between f and conv(e3U{q}) in line 12of GroupHull3d(P, m, H), pick the point q that maximizes the angle Lbjaq; and if thereare more than one q that achieve this maximum, pick the one farthest from a3. Use thesame rule to break ties in line 11. For degenerate point set, it is easier to keep track ofedges rather than facets, since facets can be convex polygons rather than triangles. So,make F and Q sets of edges instead and in line 15, add the oriented edges bja andto F and Q. Although we may not have a complete description of the facet incident tothese two edges, we know the equation of the plane containing the facet; this equation issufficient to perform wrapping about these edges.2.2.3 Refinements of the group-and-wrap methodIn this subsection, we suggest ideas on possible improvements that may make algorithmsGroupllull2d() and Groupllull3d() run even faster in practice.Idea 1. First, points found to be in the interior of conv(P) in line 3 ofGroupllull2d(P, m, H) or GroupHull3d(P, in, H) can be eliminated from further consideration. This may potentially save work during future iterations of the algorithm,although it does not affect the worst-case complexity.Idea 2. In Groupllull2d(P) and Groupllull3d(P), we choose the group size m = Hso as to balance the 0(nlogm) preprocessing cost and the O(H(logm)) cost for theChapter 2. Two- and Three-Dimensional Convex Hulls 280(H) wrapping steps. Alternatively, we can choose m = min{H log H, n} or reverselyset H = rn/log m. This choice of m does not affect the first cost except in the lower-order terms, but it reduces the second cost from 0(ri log H) to 0(n) and thus results ina smaller constant factor overall.Idea 3. With Idea 2, the dominant cost of algorithm Groupllull2d(P, m, H) lies inthe preprocessing, i.e., the computation of the convex hulls of the groups in line 3. Toreduce this cost, we may consider reusing hulls computed from the previous iterationand merging them as the group size is increased. Suppose rn’ is the previous groupsize. Since the convex hull of two convex polygons can be computed in linear time(the dual problem is to intersect two convex polygons), we can compute the convexhull of Frn/rn’l convex m’-gons in 0(m log(m/rn’)) time by the standard “MergeHull”divide-and-conquer algorithm [PS85]. Thus, the rn/mi hulls in line 3 can be constructedin 0(nlog(m/rn’)) rather than 0(nlogm) time. The same can be said for the three-dimensional case, but merging two 3-dimensional polytopes, though possible in lineartime [Cha92], is more expensive.Idea 4. In GroupHull2d(P), we use the sequence of group size m= 22, t = 1,2,..., toguess h. The improvements from Ideas 2 and 3 in fact permit us to choose slower growing sequences and still retain optimal 0 (n log h) complexity. For example, one possiblesequence is simply m = 2, t = 2, 3,..., which corresponds to doubling the group sizeafter each iteration. Note that a coarser sequence approximates h less well while a densersequence requires more iterations. We may try to optimize the worst-case constant factor and lower-order terms using sequences with different growth rates. We suggest thesequence m = 22, t = 2,3,...Chapter 2. Two- and Three-Dimensional Convex Hulls 29Idea 5. E. Weizi has observed that the binary search in line 8 of algorithmGroupllull2d(P, m, H) can be replaced by a simpler linear search without changing thetime complexity of the algorithm. The following monotonicity property provides the justification: during the course of the algorithm, the variable qj in line 8 can ollly advancein the ccw direction along conv(P) for each fixed i. As a result, the h-vertex convexhull of p convex polygons with a total of n vertices can be computed in O(n + hp) timeby gift-wrapping; the two-polygon (p = 2) version of the algorithm is in fact the dual ofan intersection algorithm by O’Rourke et al. [OC+82j (see also [O’Ro94, PS85]). Thetotal cost of Groupllull2d(P, m, H) can then be reduced to O(nlogm + H(n/m)) time,which is a log m factor saving in the second term. Although the overall constant factor is unaffected by the saving if Idea 2 is employed (as the first term is the dominantone), the linear search is easier to implement. There does not seem to be an analogoussimplification in three dimensions.2.3 Application: Lower Envelopes of Line Segments in the PlaneIn the previous section, we have presented optimal output-sensitive convex hull algorithms in both E2 and E3. Besides simplicity, the approach of the previous section hasthe advantage that it is applicable to a variety of other problems. This section gives oneapplication that illustrates this idea particularly well.Consider the problem of computing the lower envelope £(S) of a set S of n line segments in the plane, which we define as the boundary of U8ES where a is the trapezoidformed by extending the segment s vertically to +00. See Figure 2.4 for an example. (Convex hulls correspond to lower envelopes of lines in the dual.) Let h be theoutput size, i.e., the number of edges in the envelope; it is known that h is at mostO(nc(n)) [HS86]. Hershberger [Her89] has given a worst-case optimal algorithm thatChapter 2. Two- and Three-Dimensional Convex Hulls 30Figure 2.4: The lower envelope of a set of line segments. (Shown in dotted lines.)computes lower envelopes in Ofri log n) time. We now describe how his algorithm can bemade output-sensitive with our technique.First, observe that we can trace the h edges in £(S) from left to right by performing h ray shooting operations, where a ray shooting operation is: given a ray p emanating from a point on or beneath £(S), find the first trapezoid a (s e S) that pcrosses. As such an operation can be done in 0(n) time, this gives us a naïve O(nh)method, like Jarvis’s march. To improve the running time, partition S into [n/rn] groupseach of at most m segments and compute the lower envelope of each group by Hershberger’s algorithm; this takes O(nlogm) time in total. Using knowil data structures suchas [CE+91, GH+87, HeS93], we can perform ray shooting under each of these [n/mienvelopes in 0(log m) time after O(mc(m)) preprocessing (the ray shooting methodscan be simplified in our case since envelopes are monotone). This implies that the h rayshooting operations on £(S) can be done in 0(h( log rn)) time. Choosing an appropriate group size rn and guessing the output size h give us an optimal output-sensitiveChapter 2. Two- and Three-Dimensional Convex Hulls 31O(n log h) algorithm for computing the lower envelope.Theorem 2.3.1 The h-edge lower envelope of a set of n line segments in E2 can beconstructed in O(nlogh) time.Other applications of the same approach, including higher-dimensional ones, can befound in Chapter 4. In many cases, our grouping technique, combined with appropriatedata structures, reduces 0 (n log n) terms in the running time to 0 (n log h), where hrepresents output size.Chapter 3Four-Dimensional Convex HullsIn the previous chapter, we have presented some optimal output-sensitive algorithms forconstructing convex hulls in two and three dimensions. Although there is not a drasticdifference between the 0(n log h) performance of these algorithms and the performanceof the 0(n log n) algorithms, in higher-dimensional space output-sensitivity can make areal difference.This chapter considers the four-dimensional space E4, where we give a near-optimaloutput-sensitive convex hull algorithm with an 0((n + f) log2 f) running time. The algorithm is therefore quite efficient for the whole range of output sizes f from 0(1) to0(n2). For example, when f = 0(n), the algorithm runs in O(nlog2n) time, whichis a significant improvement over the 0(n2) running time of a worst-case optimal algorithm {Sei8l]. (The previous output-sensitive method by Seidel [Sei86], combined withMatouek’s improvement [Mat93], achieves O(n4/3log°1n) time in this case.)The basic strategy behind our algorithm is divide-and-conquer. In order to obtainan output-sensitive method, the subproblems we solve cannot have asymptotically morefaces than the original polytope. Therefore, we make each subproblem compute somerestricted portion of the original polytope. For each of the subproblems defined, we showthat sufficiently many input points can be removed without changing the subproblem.As in quicksort-like recursions, the merge step is trivial once we have devised a partitioning scheme for dividing a problem into subproblems. We have already seen this type32Chapter 3. Four-Dimensional Convex Hulls 33of recursion before, namely, in the planar algorithm of Section 2.1; in fact, our 4-d algorithm is an extension of this 2-d algorithm. Further extension of the method to higherdimensions will be given in the next chapter in Section 4.6.Our result on 4-d convex hulls has an important application, namely, to the computation of Voronoi diagrams in 3-space. We will examine this application in Section 3.3.3.1 Preliminaries on the Divide-and-Conquer Construction of Convex HullsBefore we give our four-dimensional convex hull algorithm, we first set up notation,introduce key concepts concerning the divide-and-conquer computation of convex hulls,and then describe a tool that we need.3.1.1 The upper hullFor the most part, our 4-d algorithm is described in an arbitrary constant dimension d.In this setting, one can then tell when the four-dimensionality of the problem is usedso that extensions to higher dimensions can be described more easily later. Figures aredrawn for the case d = 3 only.Given a set P c Ed of n points in general position, our goal is compute the facialstructure of the convex hull conv(P). Following the approach of Section 2.1, we focus ourattention only on the upper hull, consisting of faces of the convex hull with an upwardnormal vector (see Figure 3.5(a)). The upper hull of P can be thought of as the boundedfaces of conv(P U {(O,.. . , 0, —oo)}). Once we have a method for computing the upperhull of P, we can also compute the lower hull of P in a similar manner by reflection andjoin the two hulls to form the convex hull of P.Notation. Let F(P), R(P), and V(P) be the set of all facets, ridges, and vertices(respectively) of the upper hull of P.Chapter 3. Four-Dimensional Convex Hulls 34/X2 (a)Figure 3.5: (a) The upper hull of a point set in E3 and (b) the vertical projection of itsfacets.To simplify representational issues, we require our algorithm to output only the setF(P) of all facets of the upper hull of P. From this set, we can then generate allfaces and build the complete lattice structure of the faces (the Hasse diagram) using adictionary in O(F(P) log F(P)) time; this additional cost will be absorbed in the costof the algorithm. Our algorithm for computing F(P) is based on divide-and-conquer: tocompute all the facets in F(P), we partition F(P) into suitable subsets and recursivelycompute these subsets of facets.3.1.2 Facets and their dualsThe following provides a simple characterization of the set of facets F(P). First by thenon-degeneracy assumption, F(P) consists only of (d 1)-dimensional simplices withvertices all from P. Let f be such a simplex and let h(f) denote the unique hyperplanex3......(b)Chapter 3. Four-Dimensional Convex Hulls 35containing f. Then f e F(P) if all points of P lie on or below h(f).We now introduce some useful notation used throughout the chapter.Notation. Let 4. denote the vertical projection operator: p. = (Xi,. . . , Xd_i) if p =(xi,. . . , x), and P4. = {p. : p e P} for any P c Ed. Given P ç Ed andS C Ed_i, let P15 = {p e P : pJ. E S} be the restriction ofF to S. Let ms denotethe interior of S and 8S denote the boundary of S.Observe that the vertical projection of the facets in F(P) forms a collection of (d — 1)-dimensional simplices in Ed_i that have disjoint interiors, that is, mt (ft) n mt (f’--) = 0for any two distinct facets f, f’ e F(P). In fact, the vertical projection of all faces ofthe upper hull forms a simplicial complex. For example, if d = 3, then {f: f e F(P)}forms a triangulation in the plane, as shown in Figure 3.5(b). Thus, one possible divide-and-conquer approach is to use these vertical projections to partition F(P).An equally natural approach is to use the vertical projections of the facets’ duals topartition F(P). For each f e F(P), we can define a point fD E Ed via the standardduality transformation of [Ede87]: if the hyperplane h(f) is given by { (x1,. . . , Xd) : Xd =21x + ... +2d_1Xd_— d}, then we letfD= (‘,.. . ,d). The projection of thedual fD4 = (ci,.. . , ed—i) is geometrically just the gradient of h(f) (ignoring the scalarmultiple 2); for example, if d = 2, then this is just the slope.To allow us to speak about the two divide-and-conquer approaches more succinctly,we make the following definitions:Definition. Given sets 5, C E_i, let F5(P) = {f e F(P) : f1. C S} be the primalrestriction of F(P) to S and F’(P) = {f e F(P) : z} be the dualrestriction of F(P) to L.Before we describe our divide-and-conquer algorithm, we first need to introduce ageneral tool for geometric divide-and-conquer known as the (1/r) -cutting.Chapter 3. Four-Dimensional Convex Hulls 363.1.3 Cuttings for divide and conquerLet H be a set of ri hyperplanes. A cutting in Ed is a covering of Ed with closed(possibly unbounded) simplices with disjoint interiors; the size of the cutting is thenumber of simplices. A cutting is a (1/r)-cutting if any simplex of = intersects atmost n/r hyperplanes of H. The (1/r)-cutting and its relatives have recently become apopular tool in the design of divide-and-conquer algorithms for many geometric problems;see [AGR94, Cha93b] for examples in the context of convex hulls.In [CF9O], Chazelle and Friedman showed that a (1/r)-cutting of size O(rd) existsfor any finite set of hyperplanes in Ed. Since then, many researchers, notably Matouek [Mat9lb, Mat9lc] and Chazelle ECha93a], have looked into the problem of howsuch a cutting can be constructed deterministically. The following theorem is one of theresults that have been established, using derandomization techniques. We remark thatimprovements to this theorem were known, but since we need it only for the special casewhen r is constant, it is more than sufficient for our purposes.Theorem 3.1.1 ([Mat9la, Theorem 6.1]) Given n hyperplanes in a fixed dimensiond, a (1/r)-cutting of size O(r’) can be computed in O(nrd) time.3.2 A Prune-and-Divide Convex Hull Algorithm in Four DimensionsWe are now ready to provide the details of our four-dimensional convex hull algorithm,which is an extension of the planar algorithm Dividellull2d() from Section 2.1. First,here is a high-level description how the extension works.Recall that in line 4 of algorithm Dividellull2dQ, the median of a set of n/2 numbersis computed. Since the median can be thought of as a one-dimensional (1/2)-cutting fromSection 3.1.3, we extend this step to d dimensions by computing the (1/2)-cutting of a setChapter 3. Four-Dimensional Convex Hulls 37of n/2 hyperplanes in Ed_i. In line 5, a vertex p, used for dividing the hull, is computedby taking the maximum of a set of numbers formed by projecting the input points alonga direction of slope m. Of course, the maximum of a set of numbers can be interpretedas the upper hull of a one-dimensional point set. In d dimensions, Pm then becomes acollection of ridges computed by projecting the input points along certain directions andtaking the upper hulls of the resulting (d — 1)-dimensional point sets. For d = 4, theseupper hulls are 3-dimensional and are therefore of linear size.In the 2-d algorithm, Pm divides the upper hull of P into two parts: the portion ofthe hull to the left Of pm and the portion to the right Of p. Observe that the left hull isalso the portion with slope less than m, and similarly the right hull is the portion withslope greater than m. We have thus used Pm to partition the upper hull in two ways:(i) by x-coordinate and (ii) by slope. The restriction of the upper hull with x-coordinateinside a given interval is just primal restriction, in the terminology of Section 3.1.2, andthe restriction of the upper hull with slope within a given interval is just dual restriction.In our extension to 4-d, we adopt the same strategy of using both primal and dualrestrictions to partition the upper hull.In the planar case, dividing the point set by x-coordinate ensures that the two sub-problems do not share any input points except for the vertex Pm, and pruning by slopeensures that each of the two subproblems has at most 3/4 of the input points. In the samemanner for E4, dividing by primal restrictions controls the sum of the sizes of the subproblems and pruning by dual restrictions guarantees that no subproblem receives morethan a fixed fraction of the input. Subproblems can now share more than one input point,but we argue that the number of points shared is proportional to the size of the output.The analysis then follows from an application of the general cost lemma (Lemma 2.1.1)from Section 2.1.2: primal dividing bounds the per-level cost of the recursion tree anddual pruning ensures the of-fading property.Chapter 3. Four-Dimensional Convex Hulls 38Figure 3.6: A simple region S of the point set in Figure 3.5.3.2.1 Primal dividingOur convex hull algorithm computes the primal restrictions of F(P) to certain regionsin Ed_i recursively. The regions are not arbitrary but are of a special form that we callsimple regions.Definition. A set S C E’ is a simple region of P if it is the vertical projection of aunion of facets in F(P).Figure 3.6 shows an example of a simple region. A simple region S may be disconnected. There may even exist a non-empty open (d — 1)-ball centered on OS, of an arbitrarily small radius, whose intersection with mt S is not homeomorphic to any (d— 1)-ball.This intersection cannot be empty however, as S is a union of full-dimensional simplices;in particular, this rules out “spikes,” e.g., (d — 2)-simplices that are attached to theboundary.The following lemma lists some useful properties concerning simple regions and primalrestrictions. Part (a) is an identity that follows from definition and is important forChapter 3. Four-Dimensional Convex Hulls 39proving other parts of the lemma. Parts (b) and (c) discuss when poillts can be removedwithout changing the primal restriction. Parts (d) and (e) provide bounds on the numberof vertices and ridges restricted to a simple region. Finally, (f) and (g) describe propertiesof the boundary of a simple region.Lemma 3.2.1 Let S be a simple region of P. The following statements are true:(a) UfEFS(P) f = S. (The projection of the facets of the primal restriction to S coversthe region S.)(b) If Q C P contains all vertices of the facets inF8(P), then S is a simple region ofQ and Fs(P) =F8(Q). (Points that do not contribute to facets inF8(P) can beremoved from P without affectingF8(P).)(c) S is a simple region of the restricted point set P1s and the restricted facets F8 (Pis) =F8(P).(d) V(P8) < dF8(P). (The number of facets in the primal restriction gives a boundfor the number of vertices.)(e) {r E R(P) : r4. c S}j < dF8(P). (The number of facets in the primal restrictiongives a bound for the number of ridges.)(f) 3S is the vertical projection of a union of ridges in R(P). Thus, we can representOS as a set of at most d F8(P) ridges.(g) Ps = {v : v is a vertex of some ridge r in OS}, and IFas < dFs(P). (Thenumber of vertices in the boundary OS is bounded.)Proof: Recall that {int (ft) : f E F(P)} are disjoint. Then (a) is immediate from thedefinition of the primal restriction F8 (F).Chapter 3. Four-Dimensional Convex Hulls 40To prove (b), first note that Fs(P) C Fs(Q) follows directly from the hypothesis.This implies that S is a simple region of Q. Now, (a) says that UfEFS(P) f4. = S =UfEFs(Q) f4.. We therefore must have equality: Fs(P) = Fs(Q).Statement (c) is a direct consequence of (b).To prove (d), let p be a point in V(F1s). Since the projection p4. E 5, by (a) we havep4. e f4. for some facet f Fs(P). Since p e V(P15), p must be a vertex of f. Then (d)follows as each facet has d vertices incident on it (by the non-degeneracy assumption).To prove (e), observe that each ridge r with r4. C S is incident on some facet inFs (F), by (a). Then (e) follows as each facet has d ridges incident on it.The first part of (f) is immediate from the definition of a simple region. The cardinalitybound is just a consequence of (e).The first part of (g) follows from (f) and the non-degeneracy assumption. In particular, this implies that P105 C V(P) and consequently, Pos V(P15). So the second partfollows from (d). UTo compute the primal restriction F5(P) for a simple region S of F, our divide-and-conquer algorithm first subdivides S into smaller simple regions {S} with disjointinteriors and then recursively computes F5 (F) for each of the Si’s. In computingF5 (P) we may consider only those input points that belong to P5 = F1mts1 Uby Lemma 3.2.1(c). Since mt S are disjoint, the points shared between subproblemsare points restricted to the boundary of the Si’s and we can bound the size of theseboundaries in terms of the number of output facets by Lemma 3.2.1(f,g).Note that when d = 2, the boundary of a connected simple region consists of just twopoints. In higher dimensions, the boundary becomes more complex and its manipulationdemands more care.Chapter 3. Four-Dimensional Convex Hulls 413.2.2 Dual pruningTo find a good strategy for subdividing a simple region, we switch to dual space. Weshow how to partition Ed_i into a constant number of simplices such that in computingthe dual restriction of F(P) to each of the simplices, a fraction of the points of P can bepruned.Lemma 3.2.2 In O(P) time, one can find closed simplices L,. . . , z.k C Ed_i andP1,. . . , P C P such that (i) uL1 = E”1 and {int z.} are disjoint, (ii) F(P) =F’(P), and (iii) PZ < a P, for all i = 1,. . . , k. Here, k and 0 < a < 1 are bothconstants depending only on d (assuming that P exceeds a certain constant).Proof: A proof of this lemma in the dual setting can be found in Edelsbrunner’s exposition [Ede87J of Megiddo’s linear programming algorithm [Meg84]. The algorithm isbased on the prune-and-search paradigm, and this lemma represents its “prune step.” InMegiddo’s original approach, the constant k is quite large and a is very close to 1. Weobserve an alternative solution using results on cuttings.Let = n. First, form the set H of dual hyperplanes by mapping each pointp=(pi,...,pd) ofPtothehyperp1ane{(1,...,d) :d=2pii+...+2pd_id_i—pd}.Then each facet f of the upper hull of P corresponds to a vertex of the lower envelopeof H. (The lower envelope of H consists of faces of the polytope P = { eis below every hyperplane of H}.) Furthermore, if L C Ed_i, then a facet f in thedual restriction F’(P) corresponds to a vertex of the lower envelope that has verticalprojection in /i.Arbitrarily pair the n hyperplanes in H, compute the intersection of each pair, andvertically project these intersections. This gives us n/2 hyperplanes in Ed_i. Computea constant-sized (1/2)-cutting {} of these (d — 1)-dimensional hyperplanes by Theorem 3.1.1. Consider a simplex /. from the cutting. Half of the n/2 pairs have anChapter 3. Four-Dimensional Convex Hulls 42intersection whose vertical projection lies completely outside For such a pair, one ofthe two hyperplanes cannot participate in the restriction of the lower envelope of H to/. and is thus redundant. Therefore, in computing the dual restriction F(P), n/4 ofthe points can be pruned. This proves the lemma with ci = 3/4. 03.2.3 Converting from dual to primalIn this subsection we show how to convert the partitioning {z} of the dual space obtainedfrom Lemma 3.2.2 into a partitioning in the primal space. Specifically, we show that forany simplex / ç Ed_i, we can define a region S = S(P) Ed_i such that theprimal restriction of F(P) to S is the same as the dual restriction of F(P) to z (i.e.,F8(P) = F’(P)).We start with the case when L is just a halfspace. By a transformation of coordinates, we can make L the halfspace {(, . . . , d_1) : < O}. Define a projectionEd_Ed_i that sends a point (x1,. . . , Xd) to (x2,. . . , Xd). Consider the upperhull of the (d — 1)-dimensional point set ir(P): its facets are projection of ridges inthe upper hull of F, that is, F(ir(P)) C {r(r) : r e R(P)}. Let “boundary” B bethe union of all r4. with 7r(r) e F(ir(P)). Then B is monotone in the first coordinate: given any . . , _) e Ed_i, there is at most one point (a,.. . , ewith 2 = , ..., = We define S to be the region right of B, that is,= {(,.. . : (, . , d—1) e B for some < }, as in Figure 3.7(a).The following lemma can now be established for a halfspace / by verifying definitions.Lemma 3.2.3 F5(P) = F(P).Chapter 3. Four-Dimensional Convex Hulls 43> x1(a) (b)Figure 3.7: (a) The region S and (b) the intersection of its interior with the interior ofthe simple region S from Figure 3.6.Remark: In the dual setting, as described in the proof of Lemma 3.2.2, the (d — 1)-dimensional upper hull of the projected point set 7rx(P) corresponds to the (d — 1)-dimensional intersection of the lower envelope of H with Thus, ridges appearing inB correspond to edges of the lower envelope of H that intersect 8L.We now extend the definition of S. to the case when L is a simplex rather thana halfspace: Since Li is a simplex, write / as an intersection of a constant number ofhalfspaces {ö}. Then define S to be a region with interior fl3 mt Sc,. It is not difficultto see that Lemma 3.2.3 holds for simplices z as well.3.2.4 Specializing for d = 4We now show that the region S as defined in the previous subsection satisfies some nicecomputational properties if d = 4. We first coilsider the case in which L is a halfspace.Chapter 3. Four-Dimensional Convex Hulls 44For d = 4, the projected point set ir(P) is 3-dimensional. We can compute the facetsof the upper hull F(7r(P)), and thus, the boundary B, in O(P log V(P)) time byTheorem 2.2.2. This permits computations involving the region S to be done efficiently,such as deciding if a point lies in the interior of S.Lemma 3.2.4 Suppose that d = 4. Then the restricted point set P mt s. can be computedin O(PlogV(P)) time using O(P) space.Proof: Compute the facets of the 3-dimensional upper hull F(’ir(P)) and store {ir(r):ir (r) e F(r (F)) }, which is a set of 0 ( V(P) ) triangles in E2 with disjoint interiors, in aplanar point location structure [EGS86, Kir83, Pre9O, ST86]; this takes O(P logV(P))time. For each p E F, we can then test whether p4. E mt S in logarithmic time by findinga facet ir(r) of F(r(P)) with ir(p). e ir(r)j. and then determining which side of rj.the point p.j.. lies on. CAnother operation on the region S that we need is that of intersecting S witha simple region (see Figure 3.7(b)). To ensure that the resulting region is simple, weintersect their interiors only. We represent a simple region S by its boundary, which is aset ofO(F8(P)) ridges by Lemma 3.2.1(f). We assume that each ridge r of 8S is givenan orientation to indicate which side of r. the region S lies on.Lemma 3.2.5 Suppose that d = 4. Given the boundary t9S for a simple region S ofF, one can construct the boundary aS’ for a new simple region S’ of F with mt 5’ =intS fl intS in 0((P + Fs(P))logV(P)) time using O(F + F8(F)) space.Proof: Call a subset of S a subregion if it is the closure of a connected component ofEd— (8S U Ba). A subregion is a simple region, so we can define the new simple region5’ to be the union of all subregions contained in S. The boundary OS’ is made up ofChapter 3. Four-Dimensional Convex Hulls 45boundary components, which are connected components of the boundary of subregions.To decide whether a given boundary component B contributes to OS’, take a point q nearB but inside the region bounded by B (this requires an examination of the orientationof a ridge in B), and then test if q E mt S using the point location method from theprevious lemma. Therefore, to compute OS’, it suffices to produce all the boundarycomponents. This is done using depth-first search as follows.We first record the ridges of the boundaries OS and B in a dictionary; as the(d — 1)-dimensional upper hull F(ir(P)) and the boundary B can be computed inO(P log jV(P) ) time for d = 4, this takes O((P + IFs(P) ) log jV(P) ) time. We maketwo copies of a ridge to represent the two “sides” of a ridge and assign different orientations to them. We then generate all the (d — 3)-subfaces of these ridges, and for eachsuch (d — 3)-face u, we create the list of ridges incident to a in sorted order and storea in a dictionary for (d — 3)-faces. The ordering of these ridges is based on the anglesmade by their vertical projections with a fixed hyperplane through a.j. in Ed—i.Then, given an (oriented) ridge in a boundary component B, we can identify its d — 1adjacent ridges in B in constant time by following pointers. (Here, two oriented ridgesFigure 3.8: Tracing a boundary component (d = 3).Chapter 3. Four-Dimensional Convex Hulls 46r1 and r2 are adjacent in B if there is a common (d — 3)-subface a incident on bothridges and there is no other ridge r’ in B that a is incident on, such that r’. lies withinthe angle range defined by r1j. and r2. around at.) Using a depth-first search to visitthe adjacent ridges recursively, we can then trace all ridges that belong to the sameboundary component B, as indicated in Figure 3.8. All boundary components can thenbe generated by ensuring that all ridges in OS are visited. The time required by thedepth-first search is proportional to the number of ridges in OS and B and is thereforeonly O(F8(P) + V(P)j).Note: if both sides of a ridge appear in the boundary OS’, then the ridge can beremoved. DIt remains to extend the above lemmas to the case when L is a simplex rather thana halfspace. Recall that we have defined S such that mt S = fl3 mt S, if L is writtenas an intersection of a constant number of halfspaces {S3}. By applying Lemmas 3.2.4and 3.2.5 to each halfspace 3j individually, we see that the lemmas are also true for thesimplex L.3.2.5 The prune-and-divide algorithm in four dimensionsWe now have all the pieces needed for an output-sensitive convex hull algorithm inE4. Let Dual—Partition(P) represent a dual partitioning {(P, z)}1 obtained fromLemma 3.2.2. Let Restrict—Interior(P, z) represent the restricted point set P1intas computed by Lemma 3.2.4 and let Restrict—Boundary(P, B, ) be the boundary ofthe simple region 5’ returned in Lemma 3.2.5 for B = OS. The following provides anoutline of our recursive algorithm.Chapter 3. Four-Dimensional Convex Hulls 47Algorithm DivideHull4d(P, B)[Given P = mts and B = 8S for a simple region S of a point set P E4 whereB 0 is represented as a set of (oriented) ridges, return the set of facets F5(F).1. P+—PU{v:visavertexofsomeridgerinB}2. if P < n for a constant n0 then returnF5(P) in constant time3. {(P, L)}L1 +— Dual-Part ition(P) by computing a (1/2)-cutting (Lemma 3.2.2)4. fori=1,...,kdo5. P’ — P fl P fl Restrict—Interior(P, by computing a 3-d upper hulland performing 2-d point location (Lemma 3.2.4)6. B — Restrict-Boundary(P, B, by computing a 3-d upper hull andperforming depth-first search on the boundary ridges (Lemma 3.2.5)7. return U {DivideHull4d(F, B) : B O}We first argue that the algorithm indeed computes the primal restriction F5(F). Inthe first line of the algorithm, we reset P to the point set P1ms U PI8 = Pis, accordingto Lemma 3.2.1(g); the justification is provided by Lemma 3.2.1(c): F5(P) =F5(P15).Line 2 provides the base case. Line 3 gives us a constant number of simplices {z}with disjoint interiors, covering Ed_i; for each /j, we are also given a subset P of F, ofcardinality at most P, such that F’(P) = F(P).Let S denote the simple region with interior mt S fl mt S. Since F5.(F) FZXi (F)by Lemma 3.2.3, we know that the Si’s have disjoint interiors and that their union isS. Furthermore, asF5(P) CF5.(P) = F(P) = F(P), all facets in Fs(F) havevertices from Pj, which implies thatF5(P) =F5(P) by Lemma 3.2.1(b).In line 5 we set P’ = Pfl F’ fl its = Piints and in line 6 we let B be theboundary of S. Then line 7 returns lJ Fs(F) = IJF5(P) =F5(P), as claimed.Having argued that DivideHull4d() correctly computes the primal restrictionF5(P),we can use the algorithm to compute the set F(P) of all facets of the upper hull. Theinitial simple region S we use is just the convex hull of P, which can be computedusing the 3-dimensional algorithm from Section 2.2. Thus, by letting F’ = {p e FChapter 3. Four-Dimensional Convex Hulls 48pj. is not a vertex of So} and B = OS0, a call to Dividellull4d(P, B) then returns F(P),as desired.3.2.6 Analysis of the prune-and-divide algorithm in four dimensionsWe now analyze the running time of the algorithm. We do so by counting the cost of therecursion tree produced by the calls to Dividellull4dQ, in a way similar to our analysisof Dividellull2d() in Section 2.1.2. Let n be the number of input points and f be thenumber of facets of the upper hull. Let P and S denote the input point set and thesimple region associated with a node v of the recursion tree. Let n = andf =By Lemmas 3.2.2, 3.2.4, and 3.2.5, the non-recursive part of the algorithm (lines 1—6)requires O((P + = O((P + fV)logfM) time at node v, since= V(Ps,) < dFs(P) by Lemma 3.2.1(d). To get the total running time,we just have to sum this cost over all nodes in the recursion tree.We first analyze the cost contributed by the O(P logf) term. By Lemma 3.2.2(iii),this cost is n-fading, so we can apply Lemma 2.1.1. To sum the costs on a given levelof the tree, we write = intS, + P < n1, + df by Lemma 3.2.1(g). Sincethe Sn’s have disjoint interiors over all nodes ii of one level, we have n, < n and, f1, = f for each level of the recursion tree. This gives us an 0((n + f) log f) boundon the cost-per-level. The tree has 0(f) leaves, as each leaf discovers at least one facet(note thatF8(P) = 0 only if 8S = 0 by definition of a simple region). Lemma 2.1.1(u)says that the total contribution is O((n + f) log2 f).Next we analyze the cost contributed by the O(f log f) term. This cost may not beof-fading, so we cannot apply Lemma 2.1.1. But since f = f and the recursion treehas depth at most logj, n by Lemma 3.2.2(iii), we can bound the sum of these costs by0(flogflogn), which never dominates 0((n + f)log2f).Chapter 3. Four-Dimensional Convex Hulls 49We conclude that the total running time of the algorithm is O((n + f) log2 f). Totalspace is O(n + f) as long as we free up the space used to store the boundary B beforewe make the recursive calls in line 7.Theorem 3.2.6 Algorithm Dividellull4d() computes the f-face upper hull of an n-pointset P C E4 in O((n + f) log2 f) time and O(n + f) space.Remarks:1. Algorithm DivideHull4d() can be considered as a primal-based divide-and-conquer algorithm since it recursively computes the primal restriction of F(P) to a simpleregion. Alternatively, one may consider an algorithm that computes the dual restrictionof F(P) to a simplex recursively. This dual-based approach is perhaps less complex sincesimplices are easier to handle than boundaries of simple regions. However, the problemwith this approach is that the dual analogue of Lemma 3.2.1(b) is not true in general:one can construct a point set P ç E3 and a triangle C E with F(P) F’(Q) forQ = {v : v is a vertex of some facet f F’(P)}.2. Although the cutting techniques used for dual pruning in our algorithm havebeen well-studied, our strategy for primal dividing appears new. This strategy providesa simple way to guarantee that the total problem size at any level of the recursion isO(n+f); it would be difficult to obtain such a bound using the existing cutting techniquesalone. Previously, primal dividing was used only in two and three dimensions, notably inthe algorithms of Kirkpatrick and Seidel [KS86] and Edelsbrunner and Shi [ES91], andour algorithm can be regarded as an extension of these approaches. In fact, the threedimensional version of our algorithm simplifies Edelsbrunner and Shi’s algorithm in thesame way as our two-dimensional algorithm Dividellull2d() simplifies Kirkpatrick andSeidel’s (see the first remark after Theorem 2.1.2). The “contour”-based approach usedin a recent parallel 3-d convex hull algorithm by Amato, Goodrich, and Ramos [AGR94]Chapter 3. Four-Dimensional Convex Hulls 50can also be interpreted as a form of primal dividing. There, contours play a role similarto the lower-dimensional upper hulls of 7r(P) in Section 3.2.3 and are used to ensurethat the total problem size at any level of the recursion remains 0(n); but since theydescribe their method in the dual setting, its geometry is less apparent in some places.3. Dividellull4d() can return not only F(P) but also a point location structure forthe set {fJ. : f 6 F(P)} of tetrahedra in E3. We simply maintain the recursion treeand store the planar point location structures from Lemma 3.2.4 at every node; thisrequires O((n + f) log f) space. Then we can find a facet of which the vertical projectioncontains our given query point by just following a path down the recursion tree. Sincethe tree has depth at most log11 n, the query time is O(log2n). For small output size f,we can further reduce the space and query time bound to O(f log f) and 0(log2f) byfirst calling DivideHull4d() to identify the vertices V(P) and then building the pointlocation structure for V(P) instead of P (as F(V(P)) = F(P)). We thus achieve thesame performance as Goodrich and Tamassia’s 3-d point location structure [GT91].4. Some practical issues. With no additional work, Dividellull4d() can return theincidence structure between facets and ridges; this fact can be used to reduce the numberof dictionary operations needed. Moreover, with an appropiate choice of coordinatesystem, it is not necessary to compute the upper and lower hulls separately; we chooseto do so here merely because vertical projections are easier to visualize. In Section 3.2.3,a transformation of coordinates is used to define the projection ir for an arbitraryhalfspace ; this is used only to simplify presentation and such a transformation needsnot be carried out explicitly. Finally, we should mention that degeneracies may occur inthe projected point set ir(P) even though the point set P is itself non-degenerate; insuch a situation, we may wish to apply a perturbation to L.5. In the next chapter, we describe one way that algorithm DivideHull4d() can beextended to higher dimensions. Since this extension requires higher-dimensional dataChapter 3. Four-Dimensional Convex Hulls 51structures for ray shooting and linear programming queries, we postpone its discussionuntil we come to Section 4.6. See that section for the difficulties that arise when thedimension is beyond 4.3.3 Application: Three-Dimensional Voronoi DiagramsIn this section we discuss one important application of our four-dimensional convex hullalgorithm, namely the output-sensitive computation of Voronoi diagrams in three dimensions. Let P be a set of n point sites in Ed. The Voronoi region of a site p E P consistsof all points q E Ed such that q is closer to p than to any other point in P (with respectto Euclidean distance). The Voronoi diagram of P is the collection of all Voronoi regions(see Figure 3.9). The Voronoi diagram is an extremely useful structure in computationalgeometry, since it is a powerful technique for dealing with problems related to proximity,such as answering nearest neighbor queries. Furthermore, the dual of the Voronoi diagram is the Delaunay triangulation, another fundamental geometric structure with manyapplications. See the survey by Aurenhammer [Aur9l] or the book by Okabe, Boots, andSugihara [0BS92] for more on these and other applications.We first describe how an algorithm for constructing convex hulls automatically yieldsan algorithm for constructing Voronoi diagrams in one dimension lower by “lifting” thesites. The connection between convex hulls and Voronoi diagrams is well known: it wasfirst noted by Brown [Bro8O] and further developed by Edelsbrunrier and Seidel [ES86].We also describe how to get an output-sensitive algorithm for computing a portion of theVoronoi diagram clipped to a polytope.For a given site p = (p’,. . . ,pd) E F, we can use a “lifting map” [Ede87, ES86]to define a halfspace p’ in E1: p = {(x1,. . . ,xj) : Xd+1 2p1x + ... +2Pdxd —p.p}. It is well known that the Voronoi regions are just the vertical projection ofChapter 3. Four-Dimensional Convex Hulls 52Figure 3.9: The Voronoi diagram of a planar point set. Bold lines indicate edges of thediagram clipped to the polygon W.wChapter 3. Four-Dimensional Convex Hulls 53the facets of the polytope fl2,pp*. Thus, the computation of a Voronoi diagram inEd is reduced to the computation of an intersection of halfspaces in E’. Since thecomputing an intersection of halfspaces is equivalent to the computing convex hulls byduality (Section 1.1), Theorem 3.2.6 has the following consequence:Theorem 3.3.1 The Voronoi diagram of 71 point sites in E3 can be computed inO((n + f) log2 n) time and O(n + f) space, where f is the size of the Voronoi diagram((n) = f = Q(n2)).In certain applications, only the portion of a Voronoi diagram lying in a given areais needed, and the size of this portion may be much smaller than the size of the entireVoronoi diagram. What we want is then an output-sensitive algorithm to compute theVoronoi diagram of P clipped to a region W, defined simply as the collection of all non-empty intersections of the Voronoi regions with W (see Figure 3.9).Suppose that W is a k-dimensional polytope (fl F) fl F, where F is a set of m halfspacesand F is a k-fiat in Ed. Lift each halfspace -y e F to a vertical halfspace denoted by ‘y andlift the k-fiat to a vertical (k + 1)-flat denoted by F*. Then the clipped Voronoi diagramis just the vertical projection of the facets of the polytope flpEpp* fl fly* fl F*.Thus, the clipped Voronoi diagram can be computed by constructing the intersection ofthe halfspaces {p* n F* : p P} U {7* fl F* : -y e F} inside the (k + 1)-flat F*. Ifk = 3, we can use Theorem 3.2.6 to compute the intersection of these n + m halfspacesof dimension k + 1.Theorem 3.3.2 Let d > 4 be a constant. The Voronoi diagram of n point sites inEd clipped to a 8-dimensional polytope defined by m halfspaces can be computed inO((n + m + f) log2 f) time and O(n + m + f) space, where f is the size of the clippedVoronoi diagram ((1) = f = 0(n2)).Chapter 4Higher-Dimensional Convex HullsIn this chapter, we consider the convex hull problem in Ed where the dimension d isany fixed constant. We begin by generalizing the two- and three-dimensional convex hullalgorithms in Section 2.2. Recall that these algorithms are based on the gift-wrappingmethod. Since one can interpret a wrapping step as performing a ray shooting query indual space, we recast the grouping technique of Section 2.2 in a more general setting interms of ray shooting queries (Section 4.1). With a little more work, the same techniquecan actually be applied to linear programming queries as well (Section 4.2). Data structures for answering ray shooting and linear programming queries are then our main toolsin higher-dimensional space.Using known data structures for ray shooting queries in polytopes by Agarwal, Matouek and Schwarzkopf [AM93, MS93J, we obtain an O(nlogf +(nf)’/(L”/2H4log°1n)-time convex hull algorithm in Ed. Using known data structuresfor linear programming queries by Matouek [Mat93], we also obtain an O(n log°1 h +(nh)l_h/([d/2Hl)log°’ n)-time algorithm for identifying the extreme points (i.e., the vertices of the convex hull) of a set of n points in Ed. (Recall that f = O(nLd/2i) is thenumber of hull faces and h < n is the number of hull vertices/extreme points.) Theseoutput-sensitive results are described in Section 4.3.Ray shooting and linear programming queries turn out to have numerous applications,and we examine how our techniques lead to improved results in some of these applicationsin Sections 4.4 and 4.5. Finally, using ray shooting and linear programming queries in our54Chapter 4. Higher-Dimensional Convex Hulls 554-d divide-and-conquer convex hull algorithm from Chapter 3, we obtain an extensionof the algorithm that computes d-dimensional convex hulls in O((ri + (nf)’h/E0/2l+fnl2IFd/1)1og°’ n) time for any even d> 4 in Section 4.6.4.1 Ray Shooting QueriesWe first investigate the problem of ray shooting in a convex polytope. Given a collection H of n (closed) halfspaces in Ed, where each halfspace contains a known point, say,the origin o, a ray shooting qnery is to determine the first bounding hyperplane h of fl Hthat is crossed by a query ray originating from fl H (a ray crosses a hyperplane h if itintersects h but is not contained in h).In two dimensions, the ray shooting problem can be solved as follows: first computethe polygon fl H and store its vertices in an array in counterclockwise order; then a querycan be done by a simple binary search. Observe that computing the intersection fl H isequivalent to computing a convex hull in the dual space, and thus takes O(n log n) timeby Graham’s scan for example [Gra72]; and the binary search takes 0(log n) time. Hence,this method requires 0(nlogn) preprocessing time, 0(n) space, and 0(logn) query time.The same preprocessing time, space, and query time can be obtained in three dimensions: in the preprocessing, compute the polytope fl H by the dual of Preparata andHong’s convex hull algorithm [PH77] and construct its Dobkin-Kirkpatrick hierarchicalrepresentation [DK83, DK9O]; then a query can be answered in logarithmic time (seeSection 2.2.2).Our first observation is that a preprocessing time/query time tradeoff is possible usingthe grouping idea from Section 2.2. Using this observation, we can perform q queries in0(nlogq) time rather than 0(nlogn) time for small q’s.Chapter 4. Higher-Dimensional Convex Hulls 56Lemma 4.1.1 There is a (static) data structure for ray shooting in a polytope definedby a set H of n halfspaces in E orE3 with 0(nlogm) preprocessing time, 0(n) space,and 0((n/m) logm) query time, where m is a parameter between 1 and n.Proof: Partition H into rn/mi subsets (“groups”) H1,. .. ,Hn/mi, each of size atmost m and build the above structures for each H. The total preprocessing time is0((mlogm)) = 0(nlogm), and the space complexity remains 0(n). Now, ray shooting is a decomposable problem [BS8O], i.e., the answer to a query on H’ U H” can becomputed from the answers to the queries on H’ and H” in constant time. Therefore,a query on H can be computed directly by querying on each H, taking 0((n/m) log m)time. ElCorollary 4.1.2 An (online) sequence of q ray shooting queries in a polytope defined bya set H of n halfspaces in E2 or E3 can be performed in 0(nlogq + qlogn) time and0(n) space.Proof: By Lemma 4.1.1, the total time needed to answer q queries is 0(nlogm +q( logm)), where 1 < m < n. Choose m = q when q < n and choose m = n whenq>n. ElFor d-dimensional polytopes with d > 3, Agarwal and Matouek [AM93] were thefirst to obtain efficient ray shooting data structures. Subsequently, Matouek andSchwarzkopf [MS93] proposed a simpler and slightly faster approach, with the resultsshown in Table 4.1.Structure 1 in the table originates from a data structure by Matouek [Mat92] for theproblem of preprocessing an n-point set P ç E” for halfspace range queries. (A halfspacerange query on P is to report all points of P that lie in a given query halfspace.)Chapter 4. Higher-Dimensional Convex Hulls 57preprocessing update time ray shooting linear programmingStructures time, space (amortized) query time [M593] query time [Mat93]1 nlogn, n N/A l_h/Lo/2jlog°’n nl/Ld/2i1og°n2 m1og°’ n N/Am1/I2i ml/[Z/2J1og2’logn3 n/i 1og°’ n N/A log n log’ j1’ nlogn, n log2 n 1_h/Ld/2j nl_l/Ld/2i+E2’ m’ m’/nmh/[d/2i mh/Ld/21log2dnlogn3’ nLdh/2i+E n[d/2Jl+ logn log jTable 4.1: Known data structures for ray shooting queries in polytopes and linear programming queries. For Structures 2 and 2’, m is a parameter between n and n[d/2i. Thebounds are all deterministic, with the big-Oh notation omitted.A collection {(P1,si), . . . , (Pk, Lk)} is called a simplicial partition of size k if(i) U P = F, (ii) the Pj’s are disjoint, and (iii) L\ is a simplex containing P foreach i. The key behind Matouek’s structure is the following Partition Theorem:given 1 < r < n/2, there exists a simplicial partition {(P, z)} of size 0(r), withn/r I P < 2n/r, such that any hyperplane with fewer than n/r points of P on oneside crosses at most O(r’’/L°/2i+ logr) of the simplices A method for constructingsuch a simplicial partition is described in [Mat92]. This theorem suggests a data structure for storing the point set P called the partition tree: the simplices zX are stored atthe root of the tree, and a subtree is generated in a recursive fashion for each of the pointsets Ps’s.Matouek and Schwarzkopf [MS93] observed that the partition tree can be used to answer ray shooting queries on the polytope fl H if we dualize the halfspaces in H to pointsand augment each node of the partition tree with a (1/r)-net [HW87, Mat9lc]. Choosingr = n7 for a suitable constant 7> 0, they then obtained a linear-space structure that canChapter 4. Higher-Dimensional Convex Hulls 58be built in O(n log n) time and can answer ray shooting queries in O(nl_h/L/2i iog°’ n)time.Structure 3 in Table 4.1 is obtained in a different manner. Here, the key concept isa shallow cutting. A 0-level shallow cutting for H is a covering of fl H with simpliceswith disjoint interiors; the size of is the number of simplices. Furthermore, is a0-level shallow (1/r)-cutting if any simplex of intersects at most n/r of the boundinghyperplanes of H. In [Mat92], it is shown that a 0-level shallow (1/r)-cutting of sizeO(r[d/2i) exists.Following an earlier randomized algorithm of Clarkson {C1a87] for the polytope membership problem, Matouek and Schwarzkopf defined the following tree-like structure: a0-level shallow (1/r)-cutting for H is stored at the root of the tree, and for each simplex Lin the cutting, a subtree is generated recursively for the set of halfspaces that do not completely contain zS. If one sets r = n7, such a structure can be built in Q(nl.d/2i log°’ n)time. It was shown by Matouek and Schwarzkopf [MS93] that a ray shooting query canbe answered in logarithmic time with this structure.Finally, a combination of the linear-space approach used in Structure 1 and the large-space approach used in Structure 3 yields a continuous tradeoff between preprocessingand query time: with O(m logO(1) n) preprocessing, one can answer a ray shooting queryin O((n/m1/[d/2)logn) time, for any n m < Ld/2J This is named Structure 2 inTable 4.1.We observe here that the grouping scheme used in Lemma 4.1.1 can in fact be usedto obtain further preprocessing time/query time tradeoffs for Structure 1.Lemma 4.1.3 There is a (static) data structure for ray shooting in a polytope definedby a set H of n halfspaces in Ed (d> 3) with O(nlogm) preprocessing time, 0(n) space,and O((n/m’/1-’/2i)1og°’ m) query time, where m is a parameter between 1 and n.Chapter 4. Higher-Dimensional Convex Hulls 59Proof: By partitioning H into rn/mi groups as in Lemma 4.1.1 and using Structure 1to store each group, the preprocessing time becomes O((m log m)) O(n log m) andquery time becomes O((ml_h/1d/2ilog°’ m)) = 0(mh//2J 1og°’’ m). E]Corollary 4.1.4 A sequence of q ray shooting queries in a polytope defined by a set Hof n halfspaces in Ed (d> 3) can be performed in O(nlogq + (nq)1_h/(L/2i+1) log°’ n +qiogn) time (using O(n + (nq)l_h/(Ld/2H)1og° n) space).Proof:CASE I. q < nh/Ld/2i/logK n, where K is a sufficiently large constant. UseLemma 4.1.3’s modification of Structure 1 with m = (q1og’ q)Ld/2] (1 <m < n). Thenthe running time iso (nlogm + 1d/2j logO(1) m) = O(nlogq).CASE II. nh/Ld/2j < q < n’’2. Use Structure 2 with m = (nq)’’/(1/2H4)(n mnL°V2i). Then the running time iso (mlog0 n + 1d,2j logn) = O((nq)’’2’1og°’ n).CASE III. q> nld/2J 1og’ n. Use Structure 3. Then the running time is0 (nLd/2J log°’ n + qlogn) = 0(qlogn).Remark: In some applications, the number of queries q may not be known in advance.In that case, the parameter m cannot be set directly. This problem can be avoided bybreaking the q queries into k clusters of q, . . . , q, queries, where q1,q2,... is a knownsequence and q + . .. + <q q1 + . . . + q,. For example, in Case I of the proof ofChapter 4. Higher-Dimensional Convex Hulls 60Corollary 4.1.4, if we choose the sequence q = 22t (i = 1, 2, . . .), then the total runningtime is 0(Z nlogq)= 0(ZPoo1 n2) = 0(nlogq), as before. (Logarithms are inbase 2.) Similarly, in Case II, we see that the complexity remains unchanged by choosingthe sequence qj = nu/Ld/2j2z (i = 1,2,. . .). This method of guessing the value of q usingan increasing sequence is analogous to the method used in Section 2.2 for guessing theoutput size.We now discuss dynamic ray shooting in polytopes, where halfspaces may be insertedor deleted. In two dimensions, a data structure by Overmars and van Leeuwen [OvL8l]has 0 (n log n) preprocessing time, 0(n) space, 0 (log2 n) update time, and 0 (log n) querytime; the data structure uses a balanced binary search tree and concatenable queueoperations to maintain the dual convex hull. It is straightforward to extend Lemma 4.1.1to get a method with 0 (n log rn) preprocessing time, 0(n) space, 0 ((n/rn) log2 rn) updatetime, and 0((n/m) logrn) query time (1 <rn < n). We can then obtain a dynamic planarversion of Corollary 4.1.2:Lemma 4.1.5 A sequence of q ray shooting queries in a polygon defined by a dynamicset H of at most n halfplanes in E2, and q insertions/deletions on H can be performedin 0(nlogq+qlog2ri) time and 0(n) space.Proof: By the above method, the total time needed to perform q queries and updates is0(nlogm + q(log2m)), where 1 <rn n. Choose rn = qlogq when q < n/logn andchoose m = n otherwise.In higher dimensions, Matouek and Schwarzkopf [MS93] have provided dynamic versions of their data structures, as shown in the bottom half of Table 4.1. Structure 1’ usesthe dynamic partition trees from Agarwal and Matouek [AM91] and is similar to thestatic Structure 1, except that (1/r)-nets are replaced by more robust simplicial partitionsChapter 4. Higher-Dimensional Convex Hulls 61and the parameter r is set to a large enough constant rather than the number n7. Thiscauses the polylogarithmic factor in the query time to increase to an nE factor. Structure 3’ is based on Structure 3 and again uses dynamization techniques from [AM91].Higher-level shallow cuttings are used, and the parameter r is also set to a constant; as aresult, space and preprocessing time is now O(m) instead of O(m1og°’ n). Finally,Structure 2’ combines the approaches in the other two dynamic structures to yield acontinuous tradeoff.We can apply our grouping scheme to get a preprocessing time/query time tradeoff forStructure 1’. This modification is easily shown to have 0(nlogm) preprocessing time,0(n) space, 0((n/m) log2 m) amortized update time, and 0(n/mh/14/2J_6) query time(1 <m < n). As a consequence, we get the following:Lemma 4.1.6 A sequence of q ray shooting queries in a polytope defined by a dynamicset H of at most n halfspaces in Ed (d> 2) can be performed in(i) 0(nlogq + (nq)’_h/(Ld/2H)+qnl_/(Ld/H1)) time, if the number of insertions/deletions is 0(q);(ii) 0(n log2 n+ (nq) 1/([d/2]+1)+6 + q log n) time, if the number of insertions/deletionsis 0(n).Proof: For (i), we consider three cases:CASE I. q nh/L/2j_E. Use the above modification of Structure 1’ withm””2= q(1 <m <n). Then the running time is0 (nlogm+ qlog2m + 1/L2J_E) 0(nlogq).CASE II. nh/Ld/2J < q < n. Use Structure 2’ with m = (nq)l_h/(Ld/2H) (n mn [d/2J). Then the running time is0 (mi+6 +qm+ 1d/2j logn) = 0((nq)’112’).Chapter 4. Higher-Dimensional Convex Hulls 62CASE III. q n. Use Structure 2’ with m =n22/(Ld/i+l) (n < m < n[/21). Thenthe running time iso (m1 + + 1/[d/2J lOn) =O(qnh_2/i+1)+2E).For (ii), we perform a similar analysis. We use Structure 1’ when q <Structure 2’ with m = (nq)’’/(L’/2i when nh/1d/2j < q < Ld/2i and Structure 3’ whenq > nLd/2j+6. E4.2 Linear Programming QueriesIn this section, we apply similar techniques to those of the previous section to answerlinear programming queries. Given a collection H of n halfspaces in Ed, each containingthe origin o, a linear programming query is to determine the vertex v of the polytopefl H that maximizes • v for a query vector e Er’.We begin by extending the grouping technique of Lemmas 4.1.1 and 4.1.3 to handlelinear programming with a small number of queries. This is not trivial because linearprogramming, unlike ray shooting, is not a decomposable problem.Lemma 4.2.1 There is a dynamic data structure for linear programming queries ona set H of n halfplanes in E2 with O(nlogm) preprocessing time, 0(n) space, and0((n/m) log2 m) update and query time, where m is a parameter between 1 andProof: We consider the static case first. Partition H into rn/mi groups H1, . . . , HEn/mi,each of size at most m, compute the convex polygon H fl H for each i, and storeeach of them in an ordered array. The total preprocessing time is then 0 ( (m log m)) =0(nlog m), while space is linear. Reichling [Rei88a] showed that in 0(k log2 m) time, onecan detect whether the intersection of k convex m-gons is empty, and if not, report theChapter 4. Higher-Dimensional Convex Hulls 63point in the intersection that is extreme in a given direction ; his method is based onMegiddo’s prune-and-search technique. Using Reichling’s algorithm on the k = rn/mipolygons.H En/mi, we can answer a linear programming query in O((n/m) log2 m)time.The dynamic part can be proved using Overmars and van Leeuwen’s data structure {OvL8l] to store each of the Hi’s, which requires O((n/m) log2 m) update time;Reichling’s time bound still applies. 0As a result of this lemma, q linear programming queries in the plane can be answeredin 0 (n log q) time for q < n/log n.Corollary 4.2.2 A sequence of q linear programming queries and q insertions/deletionson a dynamic set H of at most n halfplanes in E2 can be performed in O(nlogq+qlog2n)time and 0(n) space.In higher dimensions, Matouek [Mat93] obtained data structures for linear programming queries achieving the complexities shown in Table 4.1. His approach is as follows:He first showed how a multidimensional parametric search technique can reduce the problem of answering linear programming queries to that of answering halfspace-emptinessqueries with witness. (In the dual setting, a halfspace-emptiness query on H is to determine whether a given query point p belongs to fl H, and if not, provide a witnesshalfspace h e H that does not contain p.) Then he obtained data structures for thehalfspace-emptiness problem using techniques completely analogous to those used in theray shooting problem. (In fact, the halfspace-emptiness is a special case of ray shooting.)We now show how to obtain a preprocessing time/query time tradeoff for Structure 1in the case of linear programming queries. The bounds we get are similar to thoseChapter 4. Higher-Dimensional Convex Hulls 64obtained in Lemma 4.1.3, except for an extra polylogarithmic factor in n in the querytime; this causes an additional O(nloglogn) term in the overall time bound.Lemma 4.2.3 There is a (static) data structure for linear programming queries on a setH of n halfspaces in Ed (d > 2) with 0(nlogm) preprocessing time, 0(n) space, and0((n/m’IL°V2J)log°’ m logd n) query time, where m is a parameter between 1 and n.Proof: In this proof, we assume that the reader is familiar with Matouek’s technique [Mat93].We consider the halfspace-emptiness queries first. Partition H into In/mi groupsH1,. . . , Hnm1, each of size at most m. For each of the Hi’s, we build a datastructure [Mat93] with 0(mlogm) preprocessing time and 0(m) space, so that eachhalfspace-emptiness query on H, can be answered in 0(logm) parallel steps using0(m1_h/1d/2Jlog°’ m) processors. The total preprocessing time is then 0((m logm)) =0(n log m) and the space requirement remains linear. Since the halfspace-emptinessproblem is decomposable, a query on H can be performed in r(n, m) = 0(log m) parallelsteps using r(n, m) = 0((m1/L’/2Jlog°’ m)) = 0(m1/’d/2j log° m) processors; orsequentially, in t(n, m) = 0((ml_h/Ld/ilog°’ m)) = 0(m1/d/2j log°’ m) time.Matouek has shown that any data structure for halfspace-emptiness queries (satisfying some reasonable conditions) can be used to answer linear programming queriesby a multidimensional version of Megiddo’s parametric search method [Meg83a]. Theresulting query time is given by 0(t(n,m)T(n,m)dlogdK(n,m)), which, in our case, is0(mh//2j log°’m logdn). ElCorollary 4.2.4 A sequence of q linear programming queries on a set H of n halfspacesin Ed (d > 2) can be performed in 0(nloglogn + nlogq + (nq)1_/(Ld/2H1)logO n +q1ogd n) time.Chapter 4. Higher-Dimensional Convex Hulls 65Proof: The proof is as in Corollary 4.1.4, except that for Case J (q < i/ [d/2] / logK n)we use Lemma 4.2.3 with m = (qlog’ n)Ld/2i (1 <m < n). The running time for Case Inow becomeso (nlogm + 1d/2J log°1m logd n) = O(nloglogn + nlogq).Remark: Again, the complexity remains the same even if the value of q is not known inadvance. (Use the sequence qj = (logn)2 (i = 1, 2,...) for Case I.)Since Matouek’s data structures can be made dynamic (see Table 4.1), the followinganalogue of Lemma 4.1.6 is straightforward:Lemma 4.2.5 A sequence of q linear programming queries on a dynamic set H of atmost n halfspaces in E° (d> 2) can be performed in(i) O(nloglogn + nlogq + (nq)l_h/(Ld/H4+ qn 2/(L0V2H4)) time, if the numberof insertions/deletions is 0(q);(ii) O(nlog2n + (nq)l_h/(Ld/2H) + qlog4n) time, if the number of insertions/deletions is 0(n).Finally, we observe that for the semidynamic case, where there are no deletions,Lemma 4.2.5(u) may be improved somewhat.Lemma 4.2.6 A sequence of q linear programming queries and n insertions onan initially empty set of halfspaces in E’ can be performed in O(nlog2n +(nq)11/([0/2H4)log°’ n + qlog°’ n) time.Proof: As in the proof of Lemma 4.2.3, we consider the halfspace-emptiness problemfirst. Since, this problem is decomposable, the techniques by Bentley and Saxe [BS8O]Chapter 4. Higher-Dimensional Convex Hulls 66may be applied to convert a static structure to a semidynamic one (which increases building time and query time by a logarithmic factor). We then apply Matouek’s parametricsearch to use this structure for answering linear programming queries. The resulting timebound is only a polylogarithmic factor increase on the static bound in Corollary 4.2.4. 0In Section 4.7, we show how to eliminate the log n factors in Lemma 4.2.3, and thus remove the nloglogn term from both Corollary 4.2.4 and Lemma 4.2.5(i), if randomizationis allowed.4.3 A Convex Hull Algorithm and an Extrema Algorithmin Any Fixed DimensionIn this section, we apply the results from Sections 4.1 and 4.2 to obtain output-sensitivealgorithms for the problem of constructing convex hulls and the related problem of enumerating extreme points, in any fixed dimension.We first show that the f-face convex hull of an n-point set in Ed can be constructed byperforming 0(f) ray shooting queries in a d-dimensional polytope defined by n halfspaces.The algorithm that we use is just the well-known gift-wrapping method [CK7O, PS85,5wa85] dualized, since a wrapping step (recall Section 2.2) corresponds to shooting a rayin the dual polytope. If the ray shooting queries are performed directly by scanning thehalfspaces, then we get an 0(nf)-time bound. We observe that this can be improvedusing the data structures from Section 4.1.Theorem 4.3.1 The convex hull of a set P of n points in Ed can be constructed in0(nlogf + (nf)1_u/(Ld/2H1og°’ n) time (and 0(n + (nf)1_h/(Ld/2H1)log°’ n) space),where f is the number of hull faces.Chapter 4. Higher-Dimensional Convex Hulls 67Proof: In the dual setting, our problem becomes computing an intersection of a set Hof n halfspaces in Ed (assumed to be in general position), each containing the origin o.It suffices to compute the vertices of the intersection fl H, from which one can easilygenerate the complete facial lattice structure of fl H in 0(f log f) time by a dictionary.First, an initial vertex v0 can be found by performing d ray shooting queries in fl H,since shooting a ray from o gives a point in a (d — 1)-face, and shooting a ray from apoint in a j-face inside its affine hull (i.e., inside the j-fiat containing the face) gives apoint in a (j — 1)-face (1 <j < d). Furthermore, given a vertex v, the vertices adjacentto v in the 1-skeleton (the graph formed by the vertices and edges of fl H) can be foundby performing d ray shooting queries: if h1, . . . , hd are the hyperplanes defining v, thenshoot a ray from v along each of the d lines formed by intersecting d — 1 hyperplanesfrom {h1,. ..,hd}.Since the 1-skeleton is connected, we can use a depth-first search or a breadth-firstsearch to visit all vertices of fl H; we can ensure that each vertex is visited only once byusing a dictionary to detect replication (as in the 3-d algorithm in Section 2.2.2). Thisshows that the vertices of fl H can be computed by performing 0(f) ray shooting queriesin fl H. The theorem then follows by applying Corollaries 4.1.2 and 4.1.4 (recall thatf=O(nld/2J)). DNext, we consider the problem of finding the h extreme points of an n-point set P inEd. Although in two and three dimensions this problem has the same complexity as theconvex hull problem, the extreme point problem in higher dimensions is generally a lessdemanding problem, since we are required to output only the vertices of the convex hull,not the entire facial lattice structure.Determining whether a given point is extreme is reducible to solving a linear programon a set of n dual halfspaces. By testing all n points in F, we can find all the extremeChapter 4. Higher-Dimensional Convex Hulls 68points in n linear programming queries. Since Megiddo [Meg84j showed that a linearprogram can be solved in linear time for any constant d, the extreme point problemcan be solved in quadratic time, as is well known. Using the data structures fromSection 4.2 (Corollary 4.2.4), the n linear programming queries can in fact be solvedin O(n2_2/(Ld/2i+1) 1og°’ n) time, which implies a subquadratic bound for enumeratingextreme points, as Matouek [Mat93] has observed. We now show that the bound canbe improved to an output-sensitive one (i.e., one that depends on both n and h).As we have noted, n linear programming queries on n halfspaces are sufficient to solvethe extreme point problem. We observe that the number of queries or the number ofhalfspaces can be reduced if h is small. Specifically, we give a simple algorithm thatfinds the extreme points using h queries on n halfspaces together with n queries onh halfspaces. Using Megiddo’s linear programming algorithm, this leads to an O(nh)time extrema algorithm which was also discovered recently by Clarkson [C1a94] andOttmann et al. [0SS95]. We observe that the time bound can be further improved usingthe results from Section 4.2.Theorem 4.3.2 The h extreme points of a set P of n points in E° can be identified inO(nlog’ h + (rih)l_h/(Ld/Hl)) time or in O(nlog°’ h + (nh)’’/(L’1/2H )1og°’ n)time.Proof: Without loss of generality, assume that the origin o is in the interior of conv(P).Consider the following incremental algorithm, which is essentially the same as Clarkson’salgorithm [C1a94] and the algorithm by Ottmann et al. [0SS95}:Chapter 4. Higher-Dimensional Convex Hulls 69Algorithm Extrema(P)1. Q-O2. for each p e P (in any order) do3. ifp Q and p ‘ conv(Q U {o}) then4. if p is an extreme point of P then5. Q—QU{p}6. else find the facet f of conv(P) that intersects ray7. let v be a vertex of f that is not in Q8. Q÷-QU{v}9. return QObserve that v must exist in line 7, because otherwise all vertices of f would be in Q;since p e conv(f U {o}), this would imply that p E conv(Q U {o}): a contradiction withline 3. It is then clear that the algorithm correctly returns the set of extreme points ofP.We now analyze the cost of the algorithm. Note that line 3 can be done by solvinga linear program on Q in the dual and lines 4 and 6 can be done by solving a linearprogram on P in the dual. (Line 7 takes constant time since each facet has d verticesby the general position assumption.) Observe that although line 3 is executed n times,lines 4—8 are executed only h times since each execution adds a new point to Q. Thus,the algorithm requires h linear programming queries on P, a static set of size n, and nlinear programming queries on Q, a semidynamic set of size at most h.By Corollary 4.2.4, the h queries on P can be done in O(nloglogn + nlogh +(nh)1_h/(Ld/2H4log°’ n) time. By Lemma 4.2.5(u), the n queries and h insertions on Qcan be done in O((nh)l_h/(Ld/2H)+ nlog’ h) time. The total running time is thenO(nloglogn + nlog’1h + (nh)l_h/(Ld/2H4)+E).Notice that when h < n for a constant o < (1/[d/2J)2,the number of hull faces isso we can compute the entire convex hull in optimal O(nlogh) time andChapter 4. Higher-Dimensional Convex Hulls 700(n) space by Theorem 4.3.1. This allows us to remove the 0(nloglogn) term in thetime bound. The first part of the theorem is thus proved, and the second part followssimilarly, using Lemma 4.2.6 instead of Lemma 4.2.5(u) for Q. DTheorem 4.3.2 has an interesting corollary: the h-vertex convex hull of an n-point setin d dimensions can be constructed in 0(n1og°’ h + hLd/2i) time. In terms of n and h,this bound is within a polylogarithmic factor of optimal in the worst case for the convexhull problem, since (n log h+h1d/2J) is a lower bound. This bound is not output-sensitivethough, since the output size f can range anywhere from (h) to e(hLd/2i).Corollary 4.3.3 The convex hull of a set P of n points in Ed can be constructed in0(n1og°’ h + hL”2i) time, where h is the number of hull vertices.Proof: Compute the extreme points by Theorem 4.3.2 and then construct the convexhull of these h points by Chazelle’s algorithm [Cha93b] in 0(hL”12i) time (note that whenh = l(n1I’/2J), we have hL°121 = 2((nh)l_h/(Ld/H4))). D4.4 Application: Convex Layers and DepthsWe now discuss an application of the output-sensitive convex hull and extrema algorithmsfrom the previous section. Let P C Ed be a n-point set. The convex layers (or “onion”layers) of P are defined iteratively as follows: layer 1 is convex hull of F, and if layer i isnon-empty, then layer i + 1 is defined as the convex hull of the points of P that are notvertices of the previous layers i,. . . , i (see Figure 4.10). If P denotes the vertices of the ith layer, then P+’ can be expressed recursively as the extreme points of P— (P1U. . .UPZ).Every point p e P is a vertex of exactly one layer; that is, there is an unique i for whichp E P. This index i is called the depth of p.Chapter 4. Higher-Dimensional Convex Hulls 71Two problems come to mind: (i) computing the facial structure of each of the convexlayers, and (ii) computing just the vertices of the convex layers. The convex layers problemhere refers to the first problem. We call the second problem the depth problem, becauseit is equivalent to finding the depths of all the points. As Preparata and Shamos [PS85]described, one motivation of these problems is from robust estimation in the area ofstatistics.For d = 2, Chazelle [Cha85] gave an optimal O(nlogn)-time algorithm for the convexlayers problem (and thus, the depth problem as well). For d = 3, a near-optimal O(n’jtime method was provided by Agarwal and Matouek [AM93]. For d> 4, the problemshave been less well studied. As the total number of faces over all convex layers can rangefrom e (n) to 0 (n Ld/2i), the output size of the convex layers problem can be huge inhigher dimensions. In applications where we only need to know which point lies in whichlayer, the depth problem is therefore of more interest.Since an O(n3)-time solution to the depth problem in any constant dimension d canbe easily derived using linear programming, Edelsbrunner [Ede87, Problem 10.3(c)] askedFigure 4.10: The convex layers of a planar point set.Chapter 4. Higher-Dimensional Convex Hulls 72whether there is an o(n3)method for d> 4. By repeatedly computing extreme points andremoving points, one can solve the depth problem in 0(n2) linear programming queriesand 0(n) deletions on n halfspaces, which by Lemma 4.2.5(u), take 0(n3_3/(Ld/2i+1))time. However, as Ottmann et al. [OSS95j have pointed out recently, a simpler andbetter solution is to apply an 0(nh)-time extrema algorithm (like the one in Section 4.3)to compute the vertices of each layer. Since the sum of the number of extreme pointsover all layers is just n, the vertices of all layers are identified in only 0(n2) time.Here, we show how a hybrid of the convex hull and extrema algorithms from Section 4.3 leads to an improved subquadratic solution to the depth problem for any d; forexample, the time bound obtained is 0(n8/5)for d = 4 or d = 5. We then show howthis result implies an output-sensitive algorithm for the convex layers problem.Theorem 4.4.1 The depth of all points in a set P of n points in Ed can be computedin 0(n2jtime, where = 2/([d/2J + 1).Proof: We iteratively compute the vertices of the i-th layer (i = 1, 2,...) as follows. Weuse the convex hull algorithm in Theorem 4.3.1 to construct the i-th layer, but as soonas more than n vertices are discovered in the layer, we stop the computation and switchto the extrema algorithm in Theorem 4.3.2 to compute the vertices of the layer. We thenremove the vertices of the i-tb layer from P and proceed to the (i+ 1)-st layer. At the end,we will have the depths of every point in P. For the calls to the convex hull algorithm, wewill use a dynamic ray shooting data structure instead of a static one so that structuresdon’t have to be rebuilt as points are removed from P after each iteration; for the callsto the extrema algorithm, however, we will leave the data structures unchanged.Let h, denote the number of vertices of the i-tb layer ( h = n). We first analyze thecost of the calls to the convex hull algorithm in Theorem 4.3.1, which involve a number ofray shooting queries and n deletions on a dynamic set of at most n halfspaces; the numberChapter 4. Higher-Dimensional Convex Hulls 73of queries is proportional to the number of facets discovered. Since we stop the computation in a layer when n vertices are found, we make at most O(min{h’2,/3Ld/2i })queries for the i-th layer. The total number of queries is then asymptotically boundedby+[d/2J nd/2j) ( h) + (Ld/2i-1) (h<nI3 h>ni3 h<nI h>n3< -i3(Ld/2i1) h, < nH3(Ld/2j_l)By Lemma 2.6(u), we see that the cost of these queries is o((n2(Ld/2i_1))1—h/(Id/2i+1))= O(n2°) by our choice of 3 (where c is an appropriate constant).Next, we analyze the cost of the calls to the extrema algorithm in Theorem 4.3.2.Note that the extrema algorithm is called only for the layers i with h > n, and thenumber of hi’s with h > is at most n’ (since h = n). Ignoring logarithmicfactors, the cost is then(n + (nh)1_h/(Ld/2i+1)) < ( 1) + nl_l/([d/2]+l) ( h_h//2i+1))/ \ 1—1/(d/2j+1) / \ 1/([d/2j+1)< n (n1) + 1—1/(ld/2]+l) ( > h ) ( 1\h>n1 J \h2>n13< ,2—/3 + 1_1/(Ld/2]+1) (nl_l/(Ld/2i+l)) (fll_/3)1/(Ld/2J+1) = O(n2),by Holder’s inequality. Therefore, the entire method runs in O(n2)time. 0Corollary 4.4.2 The convex layers of a set P of n points in Ed can be constructed inO(n23 + flogn) time, where /3 = 2/([d/2J + 1) and f is now the total number offaces in all convex layers.Proof: Let P be the vertices of layer i (i.e., the points of depth i) and let h and f bethe number of vertices and faces of the layer (Z h = n, > f2 = f). We first computeChapter 4. Higher-Dimensional Convex Hulls 74P for all i in O(n213) time by Theorem 4.4.1 and then construct the convex hull ofeach P using Seidel’s algorithm [Sei86j with Matouek’s improvement [Mat93]. The totaltime needed is O(n2+(h21/2]+1)++ f log hi)) = O(n2’+f logn). 0Remarks:1. A worst-case optimal convex layers algorithm for d 4 is not difficult to get:just use an O(nh)-time extrema algorithm to compute the vertices of each layer and useChazelle’s convex hull algorithm [Cha93b] to construct the layers; then the running timeis 0(n2 + Ld/2i)2. For a more direct output-sensitive convex layers algorithm, we can simply dothe following: iteratively use the convex hull algorithm in Theorem 4.3.1 to constructthe i-th layer (i = 1, 2,.. .) and delete points from P that are vertices of a layer aftereach iteration. This method is the same as the one by Agarwal and Matouek [AM93](proceedings version) for the three-dimensional case. It requires 0(f) ray shooting queriesand ii deletions, and by Lemma 4.1.6(u), takes 0((nf)1_h/(Ld/2H1)) time, which issuperior to the bound in Corollary 4.4.2 only when the output size f is near linear (recall(n) = f = O(nLd/2J)).4.5 Further Applications: Levels in Arrangements andLinear Programming with ViolationsThis section describes some more applications of the results in Sections 4.1 and 4.2.Problems considered here include the output-sensitive construction of an interesting geometric structure known as the k-level [Ede87], and also a variant of the familiar linearprogramming problem.We first consider the construction of a k-level in a hyperplane arrangement. Given aChapter 4. Higher-Dimensional Convex Hulls 75set H of n hyperplanes in Ed, the k-level in the arrangement A(H) is defined as the set ofall points in Ed that have at most k hyperplanes of H above it (0 < k <ri). The 0-level(the upper envelope) is just the dual of a convex hull. Figure 4.11 shows an example ofa k-level with k = 1. Let f denote the size (combinatorial complexity) of the k-level.In the plane, an output-sensitive algorithm for constructing the k-level was givenby Edelsbrunner and Welzl [EW86]. We improve its running time from 0 (n log n +flog2n) to O(nlog f + flog2 n). (Alternatively, Cole, Sharir, and Yap [CSY87] gave analgorithm that runs in 0(n log n+(n+f) log2 k) time.) In higher dimensions, Agarwal andMatouek [AM93] proposed a depth-first search method based on ray shooting queries,which runs in O(nlogn+f1)time ford = 3 (actually they state a weaker O((n+f)nE)bound). We improve this to 0(nlog f +f16) and analyze the running time in higherdimensions. Tight worst-case bounds on the size f of a k-level is currently not knowneven for d = 2; thus, achieving output-sensitivity is particularly important here.Theorem 4.5.1 A k-level in an arrangement A(H) of n hyperplanes in E’ can be constructed inFigure 4.11: The boundary of the 1-level in an arrangement of lines. (Shown in bold.)Chapter 4. Higher-Dimensional Convex Hulls 76(i) O(nlog f + flog2n) time, if d = 2;(ii) 0(n1ogf+f’) time, if d= 3;(iii) 0(nlogf + (nf)l_h/(Ld/2J)+fn2/(L(/J+l)+6) time, if d 4;where f is the output size.Proof: The depth-first search algorithm by Agarwal and Matouek [AM93] (proceedingsversion) constructs the k-level using 0(f) polytope ray shooting queries and 0(f) insertions/deletions on two dynamic sets of at most n halfspaces. (In two dimensions, theiralgorithm is the same as Edelsbrunner and Welzl’s [EW86].) Hence, the theorem followsfrom Lemmas 4.1.5 and 4.1.6(i). DNext, we consider the following problem: given a set H of n hyperplanes in Ed, adirection and a number 0 < k <n, find a point in the k-level of A(H) that is minimalalong ; in other words, find a minimal point that lies on or above all but at most k ofthe hyperplanes in H. This is the feasible case of the linear programming problem with atmost k violated constraints. If k = 0, it is just an ordinary linear programming problemand can be solved in 0(n) time by Megiddo’s algorithm [Meg84]. We are interested inthe case when k > 0 is a small integer; that is, we allow only a few violations of theconstraints.For this k-violation problem, Matouek [Mat94] has devised a method that not onlyfinds the minimum but also enumerates all 0(kd) local minima in the (< k)-levels. Hismethod is a depth-first search procedure that uses linear programming queries. It runs in0(n log n + k2 log2 n) time if d = 2 and 0(nlog ii + k) time if d = 3; the running timeis 0(nlogn) if d 4 and k is sufficiently small. We show how the 0(nlogn) terms inthese bounds can be reduced to 0(nlogk) in two dimensions or to 0(nloglogn+nlogk)Chapter 4. Higher-Dimensional Convex Hulls 77in higher dimensions. As we will see in Section 4.7, the extra O(nloglogn) term caneven be removed using randomization.Theorem 4.5.2 The linear programming problem on n constraints in Ed with at mostk violations for the feasible case can be solved in(i) O(n log k + k2 log2 n) time, if d = 2;(ii) O(nloglogn + nlogk +k3) time, if d = 3;(iii) O(nloglogn+nlogk) time, ifd 4 andkd <nl/Ld/2J_E.Proof: The depth-first search algorithm by Matouek [Mat94] solves this problem usingO(kd) linear programming/membership queries and O(kd) insertions/deletions on twodynamic sets of at most n halfspaces. (A membership query is just a special case ofa ray shooting query.) Hence, part (i) of the theorem follows from Lemma 4.1.5 andCorollary 4.2.2, and parts (ii) and (iii) follow from Lemmas 4.1.6(i) and 4.2.5(i). DRemark: As Roos and Widmayer [RW94] observed, the planar feasible linear programming problem with at most k violated constraints can be solved in O(n log2 n) time, ifone uses an O(nlogn) slope selection algorithm [CS+89] and performs binary search.Roos and Widmayer also studied a related but different planar problem, which in ourterminology corresponds to finding the y-maximum point on the boundary of the k-level.They give a slightly improved O(nlogn + klog2 k)-time method for this problem; however, this method does not seem to apply to our k-violation linear programming problem.The techniques here may also be applicable to the infeasible case of linear programming with k violated constraints (which has applications to separation and transversalChapter 4. Higher-Dimensional Convex Hulls 78problems [ERvK93]), or to the smallest k-enclosing circle problem; see Matouek’s paper [Mat94].As we now show, Matouek’s approach on linear programming with violations can infact be used to improve an output-sensitive algorithm by Mulmuley [Mul90] for constructing (< k)-levels——i.e., the i-level for all i = 0, 1,. . . , k—of a non-redundant arrangementA(H) of n hyperplanes in Ed. Here, A(H) is non-redundant if for every h E H, theupper envelope of H — {h} coincides with the upper envelope of H.Mulmuley’s algorithm can be regarded as a generalization of Seidel’s output-sensitiveconvex hull algorithm [Sei86]. While Seidel’s algorithm constructs convex hulls in O(n2+f log n) time, Mulmuley’s algorithm constructs ( k)-levels in O(n2kd_l + flog n) time.In both cases, f denotes output size. Matouek [Mat93] had already shown how to reducethe first term of the running time of Seidel’s algorithm to O(n2_2/(Ld/H1)) (actually,O(n2_/(1d/1+l) log°1n)). Here, we show how to reduce the first term of the runningtime of Mulmuley’s algorithm to O(n2_21(Ld/21 )+Ekd_l). Using the lifting map fromSection 3.3, this result can be used in the construction of Voronoi diagrams of order0, 1,.. . . , k in one dimension lower; see [ES86, Mu190].Theorem 4.5.3 We can compute i-levels in a non-redundant arrangement A(H) of nhyperplanes in Ed for all i = 0,1,. . . , k inO(n22/(Ld/1+l)k_l + flogn) time, wheref is the output size.Proof: Let L,(H) denote the boundary of the i-level in A(H) and let f be its size= f). For each h e H, let Hh = {hn h’ : h’ E H — {h}}, which is a set of(d — 1)-dimensional hyperplanes.Mulmuley [Mu190] gave an algorithm which constructs the facial structure of L (H)in O((f + f_) log n) time, given the following information:1. the local minima in L (H) — L_1(H) (along some predefined direction),Chapter 4. Higher-Dimensional Convex Hulls 792. the local minima in L(Hh) — L_l(Hh) for each h E H, and3. the facial structure ofL1(H).Matouek [Mat94] has shown that the local minima of all i-levels in A(H) (i =0, 1,. . . , k) can be enumerated by performing O(kd) linear programming/membershipqueries and O(kd) insertions/deletions on two dynamic sets of at most n halfspaces.Similarly, the local minima of all i-levels in A(Hh) (i = 0 1,. . . , k) can be computed usingO(k’1) linear programming/membership queries and O(kd_l) insertions/deletions, foreach h E H. Observe that we do not need separate structures to store each Hh as the datastructures from Section 4.2 can perform linear programming queries restricted to any jfiats [Mat93]. The total number of queries and updates is then O(kd+rjkd_l) = O(nkd_l).By Lemmas 4.1.6(i) and 4.2.5(i), this takesO(n2_/(Ld/H1)k(_l) time.Thus, items 1 and 2, for all i = 0, 1,.. . , k, can be computed in O(n2_2/(Ld/2i+1)kd_1)time. Now, Mulmuley’s algorithm can be used to construct the facial structure ofL0(H),L1(H),. . . , Lk(H) incrementally, in additional O(flogn) time. DRemark: In the worst case, the total size f of the (< k)-levels can be as large ase(nLd/2ikld1) and this bound is tight [CS89]. Currently, algorithms for constructingthe ( k)-levels that are optimal for worst-case output are known in dimension 2 byEverett et al. [ERvK93] and in dimensions > 4 by Mulmuley [Mu191]. The latter resultof Mulmuley is randomized.Before we close this section, we remark that further applications of our ideas arepossible. For example, Theorem 4.3.1 can be extended to compute the intersection ofa convex hull with a j-fiat in an output-sensitive manner; in the dual, this correspondsto computing projections (shadows) of an intersection of halfspaces. More generally, wecan obtain output-sensitive bounds for computing skeletons in a halfspace intersection, orChapter 4. Higher-Dimensional Convex Hulls 80with the known methods for ray shooting in a collection of hyperplanes {AM93], skeletonsin a hyperplane arrangement; see [Ede87, Chapter 9]. With suitable data structures, thisapplies to arrangements of different objects as well, such as line segments in the plane.4.6 The Prune-and-Divide Convex Hull Algorithm in Even DimensionsIn this section, we return to the convex hull problem. We show that Theorem 4.3.1 canbe further improved in even dimensions. The method is an extension of the divide-and-conquer algorithm Dividellull4d() given in Chapter 3. To describe the extension, let usfirst recall the terminology and notation from Chapter 3 and particularly, Sections 3.2.3and 3.2.4. It suffices to provide higher-dimensional analogues of Lemmas 3.2.4 and 3.2.5from Section 3.2.4, since other parts of the algorithm work in any fixed dimension. Asin Section 3.2.4, we may assume that LI is a halfspace.Recall that in the algorithm in Chapter 3, the facets F(P) of the upper hull of apoint set P C Ed are constructed by recursively computing the primal restrictionF8(P)to various simple regions S of P. To divide a simple region into smaller simple regions,(d — 1)-dimensional upper hulls of the projection ir(P) Ed_I for certain halfspaces zare used.The main difficulty that arises when d > 4 is that we do not have control on thesize of these (d — 1)-dimensional upper hulls, i.e., the number of facets in F(ir,(P)).For d < 4, the size is O(V(P)), which we can bound by O(F8(P)) if P.4. C S byLemma 3.2.1(d). For d > 4, we can only bound the size by O(F(P)) and this can bemuch larger than the actual output size jF5(P), since there may exist facets in F(P)with vertical projection outside S even if P4. C S. To overcome this difficulty, we do notconstruct all the hull facets in F(ir(P)) but instead apply the results in Sections 4.1and 4.2 to perform queries on F(ir(P)).Chapter 4. Higher-Dimensional Convex Hulls 81Lemma 4.6.1 Suppose that d> 4. Then the restricted point set can be computedin O((jPI + (P v(p)Dl_l/rd/2) log°’ P) time.Proof: As in the proof of Lemma 3.2.4, we can test whether p4. E mt S for a given pointp e P by finding a facet ir1(r) of F(r(P)) with ir(r).I. e r(r). and then determiningwhich side of r4. the point p4. lies on. This facet can be found by performing a linearprogramming query on a polytope P defined by P dual halfspaces in Ed_i correspondingto the P points in the projected point set ir(P) C Ed_i. As we need P queries foreach point p e P, the cost is O(P2—21Fd/21 log0 P) by Corollary 4.1.4.To further reduce the running time, we first identify the vertices ir(v) of V(ir(P))(v e V(P)); using the extrema algorithm from Section 4.3, this requires at most O((P +(P v(p))1_i/rd/21) 1og°’ V(P)) time by Theorem 4.3.2 (remember that the point setir(P) is just of dimension d — 1). Because F(ir(P)) = F(V(r,(P))), we only needthe halfspaces corresponding to these < V(P) vertices to define our polytope P. Thecost of the P queries is then no more than O((P + (P {V(P)I)i_uIFd/21) log0 P)according to Corollary 4.1.4. DLemma 4.6.2 Suppose that d > 4. Given 88 for a simple region S of F, one canconstruct 8S’ for a new simple regionS’ ofF with jutS’ = intS fl intS in O((P +F5(P) + (P jF5(P) Fd/21) log0 P) time.Proof: As in the proof of Lemma 3.2.5, it suffices to compute all the boundary components. Then we can select which boundary components contribute to the boundary OS’—that is, which boundary components correspond to a subregion inside S—by the techniques of the previous lemma. This requires at most O(IFs(P)D linear programming queries on P halfspaces in Ed_i, and thus can be done withinO((P + F5(P) + (PjF5(P)I)1_uIFd/21) log0 P) time.Chapter 4. Higher-Dimensional Convex Hulls 82The boundary computation is again based on depth-first search, but now we cannotgenerate all the ridges of the boundary B in advance as we cannot afford to computeall the facets in F(ir(P)); rather, a ridge is generated when it is needed.As before, we store the ridges in OS (but not B) and their (d — 3)-subfaces in adictionary and sort these ridges around each (d — 3)-face u. Suppose B is a boundarycomponent and we are given a ridge r in B (with its orientation). We first describe howwe can generate the d — 1 ridges that are adjacent to r in B.These adjacent ridges can be classified into two types: (i) ones that are in OS, and(ii) ones that are in the boundary B. We deal with the adjacent ridges that are in OSfirst. A ridge adjacent to r must share a common (d — 3)-subface, so let us consider one(d — 3)-subface o of r. Look up the dictionary to see if a is a (d — 3)-face of OS. If so,by performing a binary search on the list of ridges that a is incident on, we can identifya candidate ridge in OS that r may be adjacent to. Repeating this procedure for every(d — 3)-subface a of r, we get all the ridges in OS that r may be adjacent to.Next we deal with the adjacent ridges that are in the boundary B. Again we considera (d — 3)-subface a of r. Determine whether ir(a) is a ridge of R(ir(P)); this test canbe reduced to a linear programming query on a (d — 1)-dimensional polytope definedby P dual halfspaces. If the test is true, the linear programming query can be usedto get a facet ir(r’) of F(7r(P)) (r’ e R(P)) that r(u) is incident on. Then a rayshooting query in dual space can be used to find (if it exists) the other facet ir (r”) ofF(ir(P)) (r” e R(P)) that ir(a) is also incident on. These ridges r’ and TI’ in B givetwo possible candidates for the adjacent ridges of r. We repeat this procedure for every(d — 3)-subface a of r.We now have a list of possible candidates for ridges that may be adjacent to r inB. By performing some local tests, we can deduce which of these ridges are actuallyadjacent. As in the proof of Lemma 3.2.5, we can then trace the complete boundaryChapter 4. Higher-Dimensional Convex Hulls 83component B by visiting the adjacent ridges recursively in a depth-first manner. Tocompute all boundary components, we ensure that all ridges in OS are visited.To evaluate total time needed by this computation, observe that the number of ridgesvisited by the depth-first search is O(IFs(P)) by Lemma 3.2.1(e) since we only generateridges r with r. C S. The work is then dominated by O(jF5(P)) linear programmingand ray shooting queries in Ed_i, which, by Corollaries 4.1.4 and 4.2.4, require no morethan O((P + F5(P) + (P Fs(P))l_u1’’2l 1og0 PD time. 0To get a convex hull algorithm in Ed, we just have to replace Lemmas 3.2.4 and 3.2.5by Lemmas 4.6.1 and Lemma 4.6.2 in the algorithm outline for Dividellull4d() fromSection 3.2.5. We follow the same notation from Section 3.2.6 to analyze the runningtime.By Lemmas 4.6.1 and 4.6.2, the non-recursive part of the algorithm now takesO((PI + FsjP) + (V(P) + Fs(P)D)l_h/rdI21 logO(1) P) = O((P + f,, +(PV f)u/rd/21) log° n) time at node 11, if we recall that V(P11) = V(PII is)djFs(P)j by Lemma 3.2.1(d).To sum this cost, we recall, from Section 3.2.6, that , n, < n and >, f,, = fover every level of the recursion tree. Since P = P lint s + P n, + df, byLemma3.2.1(g), we also have , PV <n + df. Using Holder’s inequality, we obtain thefollowing cost-per-level bound, ignoring polylogarithmic factors:(jpj + f + (IPVI f)1_i/rd/21) = o ( + f + i-2/Fd/21 (jp 1/Ed/2 fi_1/rd/21))= O(ri + f + 1—2/Ed/21( + f)d/21f11’1V21)= O(n + (nf)’21 +fi-2/fd/1)Summing over all O(1og1 n) levels, we get O((n + (nf)i_114/2+fnl_2/Ed11)log° n)as the total running time.Chapter 4. Higher-Dimensional Convex Hulls 84Theorem 4.6.3 Let d > 4 be a constant. The convex hull of an n-point set in Ed canbe computed in O((n + (nf)1_h/Ed/2+fnl_211d/) log°’ n) time.For odd dimensions d, this method is of no use since it is more complicated and notbetter than the convex hull algorithm in Section 4.3; however, for even dimensions d,we do obtain improvement over previous results for a certain range of f. For example,when f = e(n), our previous method (and also Matouek’s method [Mat93]) achievesO(n2_/(Ld/H4)log°1n) time; the method here achieves O(n22/rd/11og°’ n) time. Ingeneral, if the output size is linear or sublinear, Theorem 4.6.3 provides the best upperbound currently known for the convex hull problem, ignoring polylogarithmic factors.4.7 Appendix: Using Randomization in Linear Programming QueriesIn this appendix, we describe how randomization can be used to improve the results fromSection 4.2.Consider the following linear programming problem: given k preprocessed convexpolytopes Hi,. . . , ilk C Ed, each defined by m halfspaces containing the origin o, computethe vertex v of fl fl. . . n fI, that maximizes . v for a given E Ed. Suppose that linear-space static structures (Structure 1) from Table 4.1 are used to store these polytopes. Asis demonstrated in the proof of Lemma 4.2.3, a direct application of Matouek’s multidimensional parametric search technique would yield an O(k ml/L0/2j 1og°’ m logd k)time solution. Here we describe how the loge’ k factor can be eliminated by using Sharirand Welzl’s randomized algorithm for generalized linear programming [SW92] (which isbased on Seidel’s linear programming algorithm [Sei9l]). This in turn improves the querytime in Lemma 4.2.3.We first observe that the problem of finding an extremum in a non-empty intersectionof k convex objects in E’ belongs to the class of LP-type problems of combinatorialChapter 4. Higher-Dimensional Convex Hulls 85dimension d as defined by Sharir and Weizi [SW92]. Sharir and Weizi presented a simplerandomized algorithm for solving LP-type problems of fixed combinatorial dimension thatrequires an expected number of 0(k) primitive operations. The primitive operations, inour case, are: (i) to test whether a given point lies inside one of the objects (violationtests), and (ii) to find the extremum in an intersection of d + 1 of the objects (basiscomputations).For our application, the objects are polytopes. A violation test is simply a membershipquery and costs O(m1_h/Ld/2Jlog°1m) time. A basis computation involves solving ourlinear programming problem on d + 1 of the polytopes; since the number of polytopesis now constant, we can apply our previous method, via Matouek’s parametric search,to solve this problem in 0(m’’/L”/2Jlog°’ m) time. Because 0(k) violation tests andbasis computations are expected to be performed by Sharir and Welzl’s algorithm (theexpected number of basis computations is actually only 0(logd k) [We191j), we obtain arandomized 0(k ml_l/Id!2]log°1m)-time solution to our linear programming problemon k polytopes.The above method carries through if the polytopes are stored in linear-space dynamicstructures (Structure 1’ from Table 4.1); we simply replace the log°’ m factors withm6. With slightly more effort, we can even remove the assumption that the polytopesall contain the origin; the method can detect whether k preprocessed polytopes have acommon intersection.Note that in the two-dimensional case both violation tests and basis computations canbe performed in O(log m) time. Thus, Sharir and Welzl’s algorithm achieves expectedO(k log m) time, which is an improvement over the previous O(k log2 m) algorithm byReichling [Rei88aj, as used in our proof of Lemma 4.2.1. It is also interesting to compare the techniques here with those used in the previous deterministic and randomizedmethods by Reichling [Rei88b] and Eppstein [Epp9l] for the three-dimensional problem.Chapter 4. Higher-Dimensional Convex Hulls 86The (expected) query time in Lemma 4.2.3 can now be improved toO((n/mL/2i)logO m) since it uses k = En/mi. As a consequence, the O(nloglogn)term in Corollary 4.2.4 can be eliminated; the same is true for the dynamic case(Lemma 4.2.5(i)), which leads to corresponding improvements in Theorem 4.5.2.Chapter 5. Conclusion 88rdimension running time references2 O(nlogh) [KS86]3 O(nlogh) [CM95]4 O(rif) [CK7O, Swa85]O(n4/3log°1n + f log n) [Mat93, Sei86]d> 4 O(nf) [CK7O, Swa85}O(n2_21(Ld/2H4) 1og°’ n + flogn) [Mat93, Sei86]2 O(nlogh) Theorems 2.1.2 and 2.2.13 O(n log h) Theorem 2.2.24 O((n + f) log2 f) Theorem 3.2.6O(nlogf + (nf)2/3iog° n) Theorem 4.3.1d> 4 O(nlogf + (nf) 1/([d/2j+1) iog°’ n) Theorem 4.3.1O((ri + (nf)1_111d/2+fn1_2/E/2l)1og°’ n) Theorem 4.6.3Table 5.2: Summary of output-sensitive results for the convex hull problem. The tophalf of the table shows the best running time achieved by previous algorithms, and thesecond half shows the running time achieved by the algorithms of this thesis.Chapter 5. Conclusion 89In trying to answer the question of whether e(nlogf + f) is the true complexityof the convex hull problem, we may consider some special cases, for instance, whenoutput size is very small or very large. For sufficiently small values of f, we do indeedattain optimal 0(nlogf) time in arbitrary dimension (Theorem 4.3.1). A theoreticallyinteresting question is whether 0(f) time can be achieved for sufficiently large valuesof f. The method by Seidel [Sei86} achieved 0(f log n) time for f super-quadratic.A result not listed in Table 5.2 that is worth mentioning here is that we can constructthe h-vertex convex hull in Ed in O(nlog°’ h + hL”/2J) time (Corollary 4.3.3). Thus,there is a near-optimal convex hull algorithm if the output polytope is of the worst kind,i.e., if f = e (h1d/2i). Unfortunately, polytopes of this kind are rare in practice, so thisresult is more of a theoretical nature. (Unlike Chazelle’s method [Cha93b], this is at leastsensitive to the number of hull vertices.)Note that all the algorithms listed in Table 5.2 are deterministic. Can randomizationlead to faster or simpler output-sensitive algorithms for higher-dimensional convex hulls?Most of our new algorithms in Table 5.2, especially the four- and higher-dimensionalalgorithms, assume general position in the input. To what extent can this assumptionbe removed?This thesis has explored new techniques for solving the convex hull problem. Alongthe way, we have also touched on applications to many related problems. For instance,we have obtained output-sensitive algorithms for computing Voronoi diagrams (Theorems 3.3.1 and 3.3.2), extreme points (Theorem 4.3.2), convex layers (Corollary 4.4.2),k-levels (Theorems 4.5.1 and 4.5.3), and even envelopes of line segments (Theorem 2.3.1).We have also obtained improved algorithms for non-constructive problems such as computing depths in a point set (Theorem 4.4.1) and linear programming with few violatedconstraints (Theorem 4.5.2). Theorem 4.4.1 on depths provides a nice example on howChapter 5. Conclusion 90output-sensitive algorithms can be used to produce better worst-case algorithms for certain problems not requiring output-sensitivity.Our discussion on these related problems raises more interesting questions. For one, isour extrema algorithm (Theorem 4.3.2) close to optimal? A weaker question is: in termsof n alone, is e(n2_2/(Ld/2H1)) near the true complexity of the extreme point problem?We do not know of even an Q(n’+E) lower bound on the problem. Obtaining nontriviallower bounds for problems of this type may be a topic of future research. Another openquestion is: can the depth problem be solved in close toO(22/(Ldh/i+1)) time, like theextreme point problem? One can also consider improving our results on convex layersand k-levels; for example, can we construct four-dimensional k-levels in O((n + f)+)time, as we can for convex hulls? Our result on linear programming with violations mayalso be improved; for example, can one get an O(n log k)-time algorithm for all 0 k <nin the planar (feasible) case?Bibliography[AM91] P. K. Agarwal and J. Matouek. Dynamic half-space reporting and its applications. Technical Report CS-1991-43, Duke University, 1991.[AM93] P. K. Agarwal and J. Matouek. Ray shooting and parametric search. SIAMJournal on Computing, 22:764—806, 1993. Also in Proceedings of the 2thAnnual ACM Symposium on Theory of Computing, pages 517—526, 1992.[AGR94] N. M. Amato, M. T. Goodrich, and E. A. Ramos. Parallel algorithms forhigher-dimensional convex hulls. In Proceedings of the 35th Annual IEEESymposium on Foundations of Computer Science, pages 683—694, 1994.[Aur9l] F. Aurenhammer. Voronoi diagrams: a survey of a fundamental data structure. ACM Computing Surveys, 23:345—405, 1991.[Ben83] M. Ben-Or. Lower bounds for algebraic computation trees. In Proceedings ofthe 15th Annual Symposium on Theory of Computing, pages 80—86, 1983.[BK+78] J. L. Bentley, H. T. Kung, M. Schkolnick and C. D. Thompson. On theaverage number of maxima in a set of vectors. Journal of the Association forComputing Machinery, 25:536—543, 1978.[BS8O] J. L. Bentley and J. Saxe. Decomposable searching problems I: static-to-dynamic transformations. Journal of Algorithms, 1:301—358, 1980.[BS78] J. L. Bentley and M. I. Shamos. Divide-and-conquer for linear expected time.Information Processing Letters, 7:87—91, 1978.[Bro8O] K. Q. Brown. Geometric transforms for fast geometric algorithms. PhD thesis,Carnegie—Mellon University, Pittsburg, Penn., 1980.[Cha95a] T. M. Chan. Optimal output-sensitive convex hull algorithms in two and threedimensions. Submitted to Discrete é4 Computational Geometry.{Cha95b] T. M. Chan. Output-sensitive results on convex hulls, extreme points, andrelated problems. Submitted to Discrete é4 Computational Geometry. Also inProceedings of the 11th Annual ACM Symposium on Computational Geometry,pages 10—19, 1995.91Bibliography 92[CSY95a] T. M. Chan, J. Snoeyink, and C.-K. Yap. Output-sensitive construction ofpolytopes in four dimensions and clipped Voronoi diagrams in three. In Proceedings of the Sixth Annual ACM-SIAM Symposium on Discrete Algorithms,pages 282—291, 1995.[CSY95b] T. M. Chan, J. Snoeyink, and C.-K. Yap. Primal dividing and dual pruning: output-sensitive construction of 4-d polytopes and 3-d Voronoi diagrams.Submitted to Discrete é4 Computational Geometry.[CK7O] D. R. Chand and S. S. Kapur. An algorithm for convex polytopes. Journalof the Association for Computing Machinery, 17:78—86, 1970.[Cha85] B. Chazelle. An optimal algorithm for computing convex layers. IEEE Transactions on Information Theory, IT-31:509—517, 1985.[Cha92] B. Chazelle. An optimal algorithm for intersecting three-dimensional convexpolyhedra. SIAM Journal on Computing, 21:671—696, 1992.[Cha93a] B. Chazelle. Cutting hyperplanes for divide-and-conquer. Discrete é4 Computational Geometry, 9:145—158, 1993.{Cha93b] B. Chazelle. An optimal convex hull algorithm for point sets in any fixeddimension. Discrete é4 Computational Geometry, 10:377—409, 1993.[CD87] B. Chazelle and D. P. Dobkin. Intersection of convex objects in two and threedimensions. Journal of the Association for Computing Machinery, 34:1—27,1987.[CE+91] B. Chazelle, H. Edelsbrunner, M. Grigni, L. Guibas, J. Hershberger, M. Sharir,and J. Snoeyink. Ray shooting in polygons using geodesic triangulations. InProceedings of the 18th International Colloquium on Automata, Languages,and Programming, Lecture Notes in Computer Science, volume 510, Springer-Verlag, pages 661—673, 1991.[CF9O] B. Chazelle and J. Friedman. A deterministic view of random sampling andits use in geometry. Combinatorica, 10(3):229—249, 1990.[CGL85] B. Chazelle, L. J. Guibas, and D. T. Lee. The power of geometric duality.BIT, 25:76—90, 1985.[CM95] B. Chazelle and J. Matouek. Derandomizing an output-sensitive convex hullalgorithm in three dimensions. Computational Geometry: Theory and Applications, 5:27—32, 1995.Bibliography 93[C1a87] K. L. Clarkson. New applications of random sampling in computational geometry. Discrete é1 Computational Geometry, 2:195—222, 1987.[C1a94] K. L. Clarkson. More output-sensitive geometric algorithms. In Proceedingsof the 35th Annual IEEE Symposium on Foundations of Computer Science,pages 695—702, 1994.[CS89] K. L. Clarkson and P. W. Shor. Applications of random sampling in computational geometry, II. Discrete é4 Computational Geometry, 4:387—421, 1989.[CS+89] R. Cole, J. Salowe, W. Steiger, and E. Szemerédi. An optimal-time algorithmfor slope selection. SIAM Journal on Computing, 18:792-810, 1989.[CSY87] R. Cole, M. Sharir, and C.-K. Yap. On k-hulls and related problems. SIAMJournal on Computing, 16:61—77, 1987.[DK83] D. P. Dobkin and D. G. Kirkpatrick. Fast detection of polyhedral intersection.Theoretical Computer Science, 27:241—253, 1983.[DK9O] D. P. Dobkin and D. G. Kirkpatrick. Determining the separation of preprocessed polyhedra: a unified approach. In Proceedings of the 17th InternationalColloquium on Automata, Languages, and Programming, Lecture Notes inComputer Science, volume 443, Springer-Verlag, pages 440—413, 1990.[Dwy9l] R. A. Dwyer. Higher-dimensional Voronoi diagrams in linear expected time.Discrete é.4 Computational Geometry, 6:343—367, 1991.[Dye84] M. E. Dyer. Linear time algorithms for two- and three-variable linear programs. SIAM Journal on Computing, 13(1):31—45, 1984.[Ede87] H. Edelsbrunner. Algorithms in Combinatorial Geometry. Springer-Verlag,Berlin, 1987.[EGS86] H. Edelsbrunner, L. Guibas, and J. Stolfi. Optimal point location in a monotone subdivision. SIAM Journal on Computing, 15:317—340, 1986.[EM9O] H. Edelsbrunner and E. P. Mücke. Simulation of simplicity: a technique tocope with degenerate cases in geometric algorithms. ACM Transactions onGraphics, 9(1):66—104, 1990.[ES86} H. Edelsbrunner and R. Seidel. Voronoi diagrams and arrangements. Discreteé4 Computational Geometry, 1:25—44, 1986.[ES91] H. Edelsbrunner and W. Shi. An O(nlog2h) time algorithm for the threedimensional convex hull problem. SIAM Journal on Computing, 20:259—277,1991.Bibliography 94[EW86] H. Edelsbrunner and E. Weizi. Constructing belts in two-dimensional arrangements with applications. SIAM Journal on Computing, 15:271—284, 1986.[EC92] I. Emiris and J. Canny. An efficient approach to removing geometric degeneracies. In Proceedings of the Eighth Annual ACM Symposium on ComputationalGeometry, pages 74—82, 1992.[Epp9l] D. Eppstein. Dynamic three-dimensional linear programming. In Proceedingsof the 32nd IEEE Symposium on Foundations of Computer Science, pages488—494, 1991.[ERvK93] H. Everett, J.-M. Robert, and M. van Kreveld. An optimal algorithm forcomputing (< k)-levels, with applications to separation and transversal problems. In Proceedings of the Ninth Annual ACM Symposium on ComputationalGeometry, pages 38—46, 1993.[GT91J M. Goodrich and R. Tamassia. Dynamic trees and dynamic point location.In Proceedings of the 23rd Annual ACM Symposium on Theory of Computing,pages 523—533, 1991.[Gra72] R. L. Graham. An efficient algorithm for determining the convex hull of afinite planar set. Information Processing Letters, 1:132—133, 1972.[Grü67] B. Grünbaum. Convex Polytopes. John Wiley & Sons, London, 1967.[GH+87] L. J. Guibas, J. Hershberger, D. Leven, M. Sharir, and R. E. Tarjan. Linear-time algorithms for visibility and shortest path problems inside triangulatedsimple polygons. Algorithmica, 2:209—233, 1987.[HS86] S. Hart and M. Sharir. Nonlinearity of Davenport-Schinzel sequences and ofgeneralized path compression schemes. Combinatorica, 6:151—177, 1986.[HW87] D. Haussler and E. Welzl. Epsilon-nets and simplex range queries. DiscreteComputational Geometry, 2:127—151, 1987.[Her89] J. Hershberger. Finding the upper envelope of n line segments in O(n log n)time. Information Processing Letters, 33:169—174, 1989.[HeS93] J. Hershberger and S. Sun. A pedestrian approach to ray shooting: shoot aray, take a walk. In Proceedings of the Fourth Annual A CM-SIAM Symposiumon Discrete Algorithms, pages 54—63, 1993.[Jar73] R. A. Jarvis. On the identification of the convex hull of a finite set of pointsin the plane. Information Processing Letters, 2:18—21, 1973.Bibliography 95[Kir83] D. Kirkpatrick. Optimal search in planar subdivisions. SIAM Journal onComputing, 12:28—35, 1983.[KS86] D. Kirkpatrick and R. Seidel. The ultimate planar convex hull algorithm?SIAM Journal on Computing, 15:287—299, 1986.[Mat9la] J. Matouek. Approximations and optimal geometric divide-and-conquer. InProceedings of the p3rd Annual ACM Symposium on Theory of Computing,pages 505—511, 1991.[Mat9lb] J. Matouek. Cutting hyperplane arrangements. Discrete é4 ComputationalGeometry, 6:385—406, 1991.[Mat9lc] J. Matouek. Efficient partition trees. In Proceedings of the Seventh AnnualACM Symposium on Computational Geometry, pages 1—9, 1991.[Mat92] J. Matouek. Reporting points in halfspaces. Computational Geometry: Theory and Applications, 2:169—186, 1992.[Mat93] J. Matouek. Linear optimization queries. Journal of Algorithms, 14:432—448,1993. Also with 0. Schwarzkopf in Proceedings of the Eighth Annual ACMSymposium on Computational Geometry, pages 16—25, 1992.[Mat94] J. Matouek. On geometric optimization with few violated constraints. InProceedings of the Tenth Annual ACM Symposium on Computational Geometry, pages 312—321, 1994.[MS93] J. Matouek and 0. Schwarzkopf. On ray shooting in convex polytopes. Discrete é4 Computational Geometry, 10(2):215—232, 1993.[McM7O] P. McMullen. The maximal number of faces of a convex polytope. Mathematika, 17:179—184, 1970.[MS71] P. McMullen and G. C. Shephard. Convex Polytopes and the Upper BoundConjecture. Cambridge University Press, 1971.[Meg83a] N. Megiddo. Applying parallel computation algorithms in the design of serialalgorithms. Journal of the Association for Computing Machinery, 30:852—865,1983.[Meg83b] N. Megiddo. Linear time algorithm for linear programming in R3 and relatedproblems. SIAM Journal on Computing, 12:759—776, 1983.[Meg84] N. Megiddo. Linear programming in linear time when the dimension is fixed.Journal of the Association for Computing Machinery, 31:114—127, 1984.Bibliography 96[Mu190] K. Mulmuley. Output sensitive construction of levels and Voronoi diagramsin Rd of order ito k. In Proceedings of the 22rd Annual ACM Symposium onTheory of Computing, pages 322—330, 1990.[Mu191] K. Mulmuley. On levels in arrangements and Voronoi diagrams. Discrete é4Computational Geometry, 6:307—338, 1991.[Mu193] K. Mulmuley. Computational Geometry: An Introduction Through Randomized Algorithms. Prentice-Hall, Englewood Cliffs, N.J., 1993.[0BS92] A. Okabe, B. Boots, and K. Sugihara. Spatial Tessellations: Concepts andApplications of Voronoi Diagrams. John Wiley & Sons, Chichester, England,1992.[O’Ro94] J. O’Rourke. Computational Geometry in C. Cambridge University Press,1994.[OC+82] J. O’Rourke, C.-B. Chien, T. Olson, and D. Naddor. A new linear algorithmfor intersecting convex polygons. Computer Graphics and Image Processing,19:384—391, 1982.[0SS95] T. Ottmann, S. Schuierer, and S. Soundaralakshmi. Enumerating extremepoints in higher dimensions. To appear in Proceedings of the 12th Symposiumon Theoretical Aspects of Computer Science, 1995.[OvL8i] M. H. Overmars and J. van Leeuwen. Maintenance of configurations in theplane. Journal of Computer and System Sciences, 23:166—204, 1981.[Pre9O] F. P. Preparata. Planar point location revisited. International Journal ofFoundations of Computer Science, i(i):7i—86, 1990.[PH77] F. P. Preparata and S. J. Hong. Convex hulls of finite sets of points in twoand three dimensions. Communications of the Association for ComputingMachinery, 20:87—93, 1977.[PS85] F. P. Preparata and M. I. Shamos. Computational Geometry: An Introduction. Springer-Verlag, New York, 1985.[Ray7O] H. Raynaud. Sur l’enveloppe convexe des nuages des points aléatoires dansRTh, I. Journal of Applied Probability, 7:35—48, 1970.[Rei88a] M. Reichling. On the detection of a common intersection of k convex objectsin the plane. Information Processing Letters, 29:25—29, 1988.Bibliography 97[Rei88b] M. Reichling. On the detection of a common intersection of k convex poiyhedra. In Computational Geometry and its Applications, Lecture Notes inComputer Science, volume 333, Springer-Verlag, pages 180—186, 1988.[RS63] A. Rényi and R. Sulanke. Uber die konvexe Hülle von n zufällig gerwähtenPunkten I. Z. Wahrsch. Verw. Gebiete, 2:75—84, 1963.[RW94] T. Roos and P. Widmayer. k-Violation linear programming. InformationProcessing Letters, 52:109—114, 1994.[ST86] N. Sarnak and R. E. Tarjan. Planar point location using persistent searchtrees. Communications of the Association for Computing Machinery, 29:669—679, 1986.[Sei8l] R. Seidel. A convex hull algorithm optimal for point sets in even dimensions. Technical Report 81-14, Department of Computer Science, Universityof British Columbia, Vancouver, B.C., 1981.[Sei86] R. Seidel. Constructing higher-dimensional convex hulls at logarithmic costper face. In Proceedings of the 18th Annual ACM Symposium on Theory ofComputing, pages 404—413, 1986.[Sei9l] R. Seidel. Small-dimensional linear programming and convex hulls made easy.Discrete Computational Geometry, 6:423—434, 1991.[SW92] M. Sharir and E. Welzl. A combinatorial bound for linear programming andrelated problems. In Proceedings of the Ninth Symposium on Theoretical Aspects of Computer Science, Lecture Notes in Computer Science, volume 577,Springer-Verlag, pages 569—579, 1992.[Swa85] G. F. Swart. Finding the convex hull facet by facet. Journal of Algorithms,6:17—48, 1985.[We191] E. Welzl. Smallest enclosing disks (balls and ellipsoids). In H. Maurer (Ed.),New Results and New Trends in Computer Science, Lecture Notes in Computer Science, volume 555, Springer-Verlag, pages 359—370, 1991.[Wen94] R. Wenger. Randomized quick hull. Manuscript, 1994.[Yao8l] A. C. Yao. A lower bound to finding convex hulls. Journal of the Associationfor Computing Machinery, 28:780—787, 1981.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Output-sensitive construction of convex hulls
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Output-sensitive construction of convex hulls Chan, Timothy Moon-Yew 1995
pdf
Page Metadata
Item Metadata
Title | Output-sensitive construction of convex hulls |
Creator |
Chan, Timothy Moon-Yew |
Date Issued | 1995 |
Description | The construction of the convex hull of a finite point set in a low-dimensional Euclidean space is a fundamental problem in computational geometry. This thesis investigates efficient algorithms for the convex hull problem, where complexity is measured as a function of both the size of the input point set and the size of the output polytope. Two new, simple, optimal, output-sensitive algorithms are presented in two dimensions and a simple, optimal, output-sensitive algorithm is presented in three dimensions. In four dimensions, we give the first output-sensitive algorithm that is within a poly logarithmic factor of optimal. In higher fixed dimensions, we obtain an algorithm that is optimal for sufficiently small output sizes and is faster than previous methods for sublinear output sizes; this result is further improved in even dimensions. Although the focus of the thesis is on the convex hull problem, applications of our techniques to many related problems in computational geometry are also explored, including the computation of Voronoi diagrams, extreme points, convex layers, levels in arrangements, and envelopes of line segments, as well as problems relating to ray shooting and linear programming. |
Extent | 1699842 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2009-04-15 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0051478 |
URI | http://hdl.handle.net/2429/7188 |
Degree |
Doctor of Philosophy - PhD |
Program |
Computer Science |
Affiliation |
Science, Faculty of Computer Science, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 1995-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_1995-059324.pdf [ 1.62MB ]
- Metadata
- JSON: 831-1.0051478.json
- JSON-LD: 831-1.0051478-ld.json
- RDF/XML (Pretty): 831-1.0051478-rdf.xml
- RDF/JSON: 831-1.0051478-rdf.json
- Turtle: 831-1.0051478-turtle.txt
- N-Triples: 831-1.0051478-rdf-ntriples.txt
- Original Record: 831-1.0051478-source.json
- Full Text
- 831-1.0051478-fulltext.txt
- Citation
- 831-1.0051478.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0051478/manifest