Minimum Convex Partitions and Maximum Empty Polytopes

Let $S$ be a set of $n$ points in $\mathbb{R}^d$. A Steiner convex partition is a tiling of ${\rm conv}(S)$ with empty convex bodies. For every integer $d$, we show that $S$ admits a Steiner convex partition with at most $\lceil (n-1)/d\rceil$ tiles. This bound is the best possible for points in general position in the plane, and it is best possible apart from constant factors in every fixed dimension $d\geq 3$. We also give the first constant-factor approximation algorithm for computing a minimum Steiner convex partition of a planar point set in general position. Establishing a tight lower bound for the maximum volume of a tile in a Steiner convex partition of any $n$ points in the unit cube is equivalent to a famous problem of Danzer and Rogers. It is conjectured that the volume of the largest tile is $\omega(1/n)$. Here we give a $(1-\varepsilon)$-approximation algorithm for computing the maximum volume of an empty convex body amidst $n$ given points in the $d$-dimensional unit box $[0,1]^d$.


Introduction
Let S be a set of n ≥ d + 1 points in R d , d ≥ 2. A convex body C is empty if its interior is disjoint from S. A convex partition of S is a partition of the convex hull conv(S) into empty convex bodies (called tiles) such that the vertices of the tiles are in S. In a Steiner convex partition of S the vertices of the tiles are arbitrary: they can be points in S or Steiner points. For instance, any triangulation of S is a convex partitions of S, where the convex bodies are simplices, and so conv(S) can be always partitioned into O(n ⌊d/2⌋ ) empty convex tiles [13].
In this paper, we study the minimum number of tiles that a Steiner convex partition of every n points in R d admits, and the maximum volume of a single tile for a given point set. The research is motivated by a longstanding open problem by Danzer and Rogers [2,6,9,18,33]: What is the maximum volume of an empty convex body C ⊂ [0, 1] d that can be found amidst any set S ⊂ [0, 1] d of n points in a unit cube? The current best bounds are Ω(1/n) and O(log n/n), respectively (for a fixed d). The lower bound comes from decomposing the unit cube by n parallel hyperplanes, each containing at least one point, into at most n + 1 empty convex bodies. The upper bound is Theorem 1. For n ≥ d + 1, we have g d (n) ≤ n−1 d . For d = 2, this bound is the best possible, that is, g 2 (n) = ⌈(n − 1)/2⌉; and for every fixed d ≥ 2, we have g d (n) ≥ Ω(n).
We say that a set of points in R d is in general position if every k-dimensional affine subspace contains at most k + 1 points for 0 ≤ k < d. We show that in the plane every Steiner convex partition for n points in general position, i of which lie in the interior of the convex hull, has at least Ω(i) tiles. This leads to a simple constant-factor approximation algorithm.
Theorem 2. Given a set S of n points in general position in the plane, a ratio 3 approximation of a minimum Steiner convex partition of S can be computed in O(n log n) time.
The average volume of a tile in a Steiner convex partition of n points in the unit cube [0, 1] d is an obvious lower bound for the maximum possible volume of a tile, and for the maximum volume of any empty convex body C ⊂ [0, 1] d . The lower bound g d (n) ≥ Ω(n) in Theorem 1 shows that the average volume of a tile is O(1/n) in some instances, where the constant of proportionality depends only on the dimension. This implies that a simple "averaging" argument is not a viable avenue for finding a solution to the problem of Danzer and Rogers.
Maximum empty polytope among n points in a unit cube. In the second part of the paper, we consider the following problem: Given a set of n points in rectangular box B in R d , find a maximum-volume empty convex body C ⊂ B. Since the ratio between volumes is invariant under affine transformations, we may assume without loss of generality that B = [0, 1] d . We therefore have the problem of computing a maximum volume empty convex body C ⊂ [0, 1] d for a set of n points in [0, 1] d . It can be argued that the maximum volume empty convex body is a polytope, however, the (number and) location of its vertices is unknown and this represents the main difficulty. For d = 2 there is a polynomial-time exact algorithm (see Section 5) while for d ≥ 3 we are not aware of any exact algorithm. Thus the problem of finding faster approximations naturally suggests itself.
There exist exact algorithms for some related problems. Eppstein et al. [16] find the maximum area empty convex k-gon with vertices among n points in O(kn 3 ) time, if it exists. As a byproduct, a maximum area empty convex polygon with vertices among n given points can be computed exactly in O(n 4 ) time with their dynamic programming algorithm.
By John's ellipsoid theorem [30], the maximum volume empty ellipsoid in [0, 1] d gives a 1/d dapproximation. Here we present a (1 − ε)-approximation for the maximum volume empty convex body C opt by first guessing a good approximation of the bounding hyperrectangle of C opt of minimum volume, and then finding a sufficiently close approximation of C opt inside it. We obtain the following two approximation algorithms. The planar algorithm runs in near-quadratic time in n, however, the running time degrades with the dimension.
As far as the problem of Danzer and Rogers is concerned, one need not consider convex sets-it suffices to consider triangles-and for triangles the problems considered are much simpler. Specifically, the largest area triangle amidst n points in the unit square U yields a constant-factor approximation of the largest area convex body (polygon) amidst the same n points. That is, the asymptotic dependencies on n of the areas of the largest empty triangle and convex body are the same. Even in higher dimensions, if d is fixed, the asymptotic dependencies on n of the volumes of the largest empty simplex and convex body are the same. For d = 2 there is a polynomial-time exact algorithm for computing the largest empty triangle contained in U (see Section 5) while for d ≥ 3 we are not aware of any exact algorithm.
Related work. Decomposing polygonal domains into convex sub-polygons has been also studied extensively. We refer to the article by Keil [25] for a survey of results up to the year 2000. For instance, when the polygon may contain holes, obtaining a minimum convex partition is NPhard, regardless of whether Steiner points are allowed. For polygons without holes, Chazelle and Dobkin [11] obtained an O(n + r 3 ) time algorithm for the problem of decomposing a polygon with n vertices, r of which are reflex, into convex parts, with Steiner points permitted. Keil [25] notes that although there are an infinite number of possible locations for the Steiner points, a dynamic programming approach is amenable to obtain an exact (optimal) solution; see also [26,35].
Fevens et al. [17] designed a polynomial time algorithm for computing a minimum convex partition for a given set of n points in the plane if the points are arranged on a constant number of convex layers. The problem of minimizing the total Euclidean length of the edges of a convex partition has been also considered. Grantson and Levcopoulos [19], and Spillner [36] proved that the shortest convex partition and Steiner convex partition problems are fixed parameter tractable, where the parameter is the number of points of P lying in the interior of conv(P ). Dumitrescu and Tóth [15] proved that every set of n points in R 2 admits a Steiner convex partition which is at most O(log n/ log log n) times longer than the minimum spanning tree, and this bound cannot be improved. Without Steiner points, the best upper bound for the ratio of the minimum length of a convex partition and the length of a minimum spanning tree (MST) is O(n) [27].
A largest area convex polygon contained in a given (non-convex) polygon with n vertices can be found by the algorithm of Chang and Yap [10] in O(n 7 ) time. The problem is known as the potatopeeling problem. On the other hand, a largest area triangle contained in a simple polygon with n vertices, can be found by the algorithm of Melissaratos and Souvaine [31] in O(n 4 ) time. Hall-Holt et al. [21] compute a constant approximation in time O(n log n). The same authors show how to compute a (1 − ε)-approximation of the largest fat triangle inside a simple polygon (if it exists) in time O(n). Given a triangulated polygon (with possible holes) with n vertices, Aronov et al. [4] compute the largest area convex polygon respecting the triangulation edges in O(n 2 ) time.
For finding a maximum volume empty axis-parallel box amidst n points in [0, 1] d , Backer and Keil [5] reported an algorithm with worst-case running time of O(n d log d−2 n). An empty axisaligned box whose volume is at least (1−ε) of the maximum can be computed in O 8ed ε 2 d n log d n time by the algorithm of Dumitrescu and Jiang [14].
Lawrence and Morris [29] studied the minimum integer k d (n) such that the complement R d \S of any n-element set S ⊂ R d , not all in a hyperplane, can be covered by k d (n) convex sets. They prove k d (n) ≥ Ω(log n/d log log n). It is known that covering the complement of n uniformly distributed points in [0, 1] d requires at least Ω(n/d log n) convex sets, which follows from the upper bound in the problem of Danzer and Rogers.

Combinatorial bounds
In this section we prove Theorem 1. We start with the upper bound. The following simple algorithm returns a Steiner convex partition with at most ⌈(n − 1)/d⌉ tiles for any n points in R d .
Algorithm A1: Step 1. Compute the convex hull R ← conv(S) of S. Let A ⊆ S be the set of hull vertices, and let B = S \ A denote the remaining points.
Step 2. Compute conv(B), and let H be the supporting hyperplane of an arbitrary (d − 1)dimensional face of conv(B  It is obvious that the algorithm generates a Steiner convex partition of S. An illustration of Algorithm A1 on a small planar example appears in Figure 1 (right). Let h and i denote the number of hull and interior points of S, respectively, so that n = h + i. Each hyperplane used by the algorithm removes d interior points of S (with the possible exception of the last round if i is not a multiple of d). Hence the number of convex tiles is 1 + ⌈i/d⌉, and we have 1 + ⌈i/d⌉ = ⌈(i + d)/d⌉ ≤ ⌈(n − 1)/d⌉, as required.
Lower bound in the plane. A matching lower bound in the plane is given by the following construction. For n ≥ 3, let S = A ∪ B, where A is a set of 3 non-collinear points in the plane, and B is a set of n − 3 points that form a regular (n − 3)-gon in the interior of conv(A), so that conv(S) = conv(A) is a triangle. If n = 3, then conv(S) is an empty triangle, and g 2 (S) = 1. If 4 ≤ n ≤ 5, S is not in convex position, and so g 2 (S) ≥ 2. Suppose now that n ≥ 6.
Consider an arbitrary convex partition of S. Let o be a point in the interior of conv(B) such that the lines os, s ∈ S, do not contain any edges of the tiles. Refer to Figure 2 (left). For each point s ∈ B, choose a reference point r(s) ∈ R 2 on the ray − → os in conv(A) \ conv(B) sufficiently close to point s, and lying in the interior of a tile. Note that the convex tile containing o cannot contain any reference points. We claim that any tile contains at most 2 reference points. This immediately implies Suppose, to the contrary, that a tile τ contains 3 reference points r 1 , r 2 , r 3 , corresponding to the points s 1 , s 2 , s 3 . Refer to Figure  would contain all points s 1 , s 2 , s 3 in its interior. Hence conv{o, s 1 , s 2 , s 3 } is a quadrilateral, and conv{o, r 1 , r 2 , r 3 } is also a quadrilateral, since the reference points are sufficiently close to the corresponding points in B. We may assume w.l.o.g. that vertices of conv{o, s 1 , s 2 , s 3 } are o, s 1 , s 2 , s 3 in counterclockwise order. Then s 2 lies in the interior of conv{o, r 1 , r 2 , r 3 }. We conclude that every tile τ contains at most 2 reference points, as required.
Lower bounds for d ≥ 3. A similar construction works in for any d ≥ 2, but the lower bound no longer matches the upper bound g d (n) ≤ ⌈(n − 1)/d⌉ for d ≥ 3.
Recall that a Horton set [24] is a set S of n points in the plane such that the convex hull of any 7 points is non-empty. Valtr [37] generalized Horton sets to R d . For every d ∈ N, there exists a minimal integer h(d) with the property that for every n ∈ N there is a set S of n points in general position in R d such that the convex hull of any h(d) + 1 points in S is non-empty. It is known that h(2) = 6, and Valtr proved that h(3) ≤ 22, and in general that h(d) ≤ 2 d−1 (N (d − 1) + 1), where N (k) is the product of the first k primes.
We construct a set S of n ≥ d + 1 points in R d as follows. Let S = A ∪ B, where A is a set of d + 1 points in general position in R d , and B is a generalized Horton set of n − (d + 1) points in the interior of conv(A), such that the interior of any h(d) + 1 points from B contains some point in B.
Consider an arbitrary Steiner convex partition of S. Every point b ∈ B is in the interior of conv(S), and so it lies on the boundary of at least 2 convex tiles. For each b ∈ B, place two reference points in the interiors of 2 distinct tiles incident to b. Every tile contains at most h(d) reference points. Indeed, if a tile contains h(d) + 1 reference points, then it is incident to h(d) + 1 points in B, and some point of B lies in the interior of the convex hull of these points, a contradiction.
We have 2(n − d − 1) reference points, and every tile contains at most h(d) of them. So the number of tiles is at least 3 Approximating the minimum Steiner convex partition in R 2 In this section we prove Theorem 2 by showing that our simple-minded algorithm A1 from Section 2 achieves a constant-factor approximation in the plane if the points in S are in general position.
Approximation ratio. Recall that algorithm A1 computes a Steiner convex partition of conv(S) into at most 1 + ⌈i/2⌉ parts, where i stands for the number of interior points of S.
If i = 0, the algorithm computes an optimal partition, i.e., ALG = OPT = 1. Assume now that i ≥ 1. Consider an optimal Steiner convex partition Π of S with OPT tiles. We construct a planar multigraph G = (V, E) as follows. The faces of G are the convex tiles and the exterior of conv(S) (the outer face). The vertices V are the points in the plane incident to at least 3 faces (counting the outer face as well). Since i ≥ 1, G is non-empty and we have |V | ≥ 2. Each edge in E is a Jordan arc on the common boundary of two faces. An edge between two bounded faces is a straight line segment, and so it contains at most two interior points of S. An edge between the outer face and a bounded face is a convex arc, containing hull points from S. Double edges are possible if two vertices of the outer face are connected by a straight line edge and a curve edge along the boundary-in this case these two parallel edges bound a convex face. No loops are possible in G. Since Π is a convex partition, G is connected.
Let v, e, and f , respectively, denote the number of vertices, edges, and bounded (convex) faces of G; in particular, f = OPT. By Euler's formula for planar multigraphs, we have v − e + f = 1, that is, f = e − v + 1. By construction, each vertex of G is incident to at least 3 edges, and every edge is incident to two vertices. Therefore, 3v ≤ 2e, or v ≤ 2e/3. Consequently, f = e − v + 1 ≥ e − 2e/3 + 1 = e/3 + 1. Since S is in general position, each straight-line edge of G contains at most 2 interior points from S. Curve edges along the boundary do not contain interior points. Hence each edge in E is incident to at most two interior points in S, thus i ≤ 2e. Substituting this into the previous inequality on f yields OPT = f ≥ e/3 + 1 ≥ i/6 + 1. Comparing this lower bound with the upper bound ALG ≤ ⌈i/2⌉ + 1, we conclude that and the approximation ratio of 3 follows.
Tightness of the approximation ratio. We first show that the above ratio 3 is tight for Algorithm A1. We construct a planar point set S as follows. Consider a large (say, hexagonal) section of a hexagonal lattice. Place Steiner vertices at the lattice points, and place two points in S on each lattice edge. Slightly perturb the lattice, and add a few more points in S near the boundary, and a few more Steiner points, so as to obtain a Steiner convex partition of S with no three points collinear. Denote by v, e, and f , the elements of the planar multigraph G as before.
Since we consider a large lattice section, we have v, e, f → ∞. We write a ∼ b, whenever a/b → 1.
Time analysis. It is easy to show that Algorithm A1 runs in O(n log n) time for a set S of n points in the plane. We employ the semi-dynamic (delete only) convex hull data structure of Hershberger and Suri [23]. This data structure supports point deletion in O(log n) time, and uses O(n) space and O(n log n) preprocessing time. We maintain the boundary of a convex polygon R in a binary search tree, a set B ⊂ S of points lying in the interior of R, and the convex hull conv(B) with the above semi-dynamic data structure [23]. Let k = ⌈ (n − 4)/3⌉. Place (k + 1) 2 Steiner points (a, b, a 2 − b 2 ) on the saddle surface z = x 2 − y 2 for pairs of integers (a, b) ∈ Z 2 , −⌊k/2⌋ ≤ a, b ≤ ⌈k/2⌉. The four points {(x, y, x 2 − y 2 ) : x ∈ {a, a + 1}, y ∈ {b, b + 1}} form a parallelogram for every (a, b) ∈ Z 2 , −⌊k/2⌋ ≤ a, b ≤ ⌈k/2⌉ − 1. Refer to Figure 2 (right). These parallelograms form a terrain over the region {(x, y) : −⌊k/2⌋ ≤ x, y ≤ ⌈k/2⌉}. Note that no two parallelograms are coplanar. Subdivide the space below this terrain by vertical planes x = a, −⌊k/2⌋ ≤ a ≤ ⌈k/2⌉. Similarly, subdivide the space above this terrain by planes y = b, −⌊k/2⌋ ≤ b ≤ ⌈k/2⌉. We obtain 2k interior-disjoint convex regions, k above and k below the terrain, such that the common boundary of a region above and a region below is a parallelogram of the terrain. The points in R 3 that do not lie above or below the terrain can be covered by 4 convex wedges.
Enclose the terrain in a sufficiently large tetrahedron T . Clip the 2k convex regions and the 4 wedges into the interior of T . These 2k + 4 convex bodies tile T . Choose 3 noncollinear points of S in each of the k 2 parallelograms, such that no 4 points are coplanar and no 2 are collinear with vertices of T . Let the point set S be the set of 4 vertices of the large tetrahedron T and the 3k 2 points selected from the parallelograms.

Approximating the maximum empty convex body
Let S be a set of points in the unit cube [0, 1] d ⊆ R d . Our task is to approximate the largest convex body C ⊆ [0, 1] d that contains no points of S in its interior. Let C opt = C opt (S) denote this body, and let vol opt (S) denote its volume.

Approximation by the discrete hull
In the following, assume that m > 0 is some integer, and consider the grid point set Let S ⊆ [0, 1] d be a point set, and let C opt be the corresponding largest empty convex body in [0, 1] d . Given a grid G(m), we call conv(C opt ∩ G(m)) the discrete hull of C opt [22]. We need the following easy lemma. Proof. Consider a point p ∈ C \ D. Consider the set of 2d points X = {p ± 2(d/m)e i }, where e i is the unit vector having one in the ith coordinate, and 0 everywhere else. If one of the points of X is outside C, then the distance from p to the boundary of C is at most 2d/m. Otherwise, the cube p + [−2, 2] d /m is contained in the "diamond" conv(X), which is in turn contained in C. But then, the grid points of the grid cell of G(m) containing p are in C, and p can not be outside D.
We reached a contradiction.
It follows that all the points in the corridor C \ D are at distance at most 2(d/m) from the boundary of C. The volume of the boundary of C is bounded by the volume of the boundary of the unit cube, namely 2d. As such, the volume of this corridor is vol Lemma 1 implies that if vol opt (S) ≥ ρ, in order to obtain a (1 − ε)-approximation, we can concentrate our search on convex polytopes that have their vertices at grid points in G(m), where m = O(1/(ερ)). If ρ is a constant, then the maximum volume empty lattice polytope in G(m) with m = O(1/ε) is an (1 − ε)-approximation for C opt . However, for arbitrary vol opt (S) ≥ Ω(1/n), a much finer grid would be necessary to achieve this approximation.

An initial brute force approach
In this section we present approximation algorithms relying on Lemma 1 alone, approximating the maximum volume empty polytope by a lattice polytope in a sufficiently fine lattice (grid). We shall refine our technique in Subsections 4.3 and 4.4.
For the plane, we take advantage of the existence of an efficient solution for a related search problem in the spirit of those considered by Eppstein et al. [16]. Proof. The algorithm works by dynamic programming. First, we compute for all triangles with vertices from S whether they contain a forbidden point inside them; trivially this can be done in O(n 3 m) time. We then build a directed graph G on the allowable triangles, connecting two triangles ∆ and ∆ ′ if they share their left endpoint, are interior disjoint, share an edge, and their union forms a convex quadrilateral. We orient the edge from the triangle that is most counterclockwise (around the common vertex) to the other triangle. All edges are oriented "upwards", so G is a directed acyclic graph (DAG). Observe that G has O(n 3 ) vertices (allowable triangles) and the maximum out-degree in G is bounded from above by n.
The weight of a vertex corresponding to a triangle is equal to its area. Clearly, a convex polygon corresponds to a path in G, namely the triangulation of the polygon from its leftmost vertex, and its weight is the area of the polygon. Finding the maximum weight path can be done in linear time in the size of the DAG; see e.g., [12,Section 4.7]. G has O(n 3 ) vertices and O(n 4 ) edges, and as such the overall running time is O(n 3 m + n 4 ).
The combination of Lemmas 1 and 2 readily yields an approximation algorithm for the plane, whose running time depends on vol(C opt ). In dimensions d ≥ 3, we are not aware of any analogue of the dynamic programming algorithm in Lemma 2. Instead, we fall back on a brute force approach that enumerates all feasible subsets of a sufficiently fine grid. Andrews [3] proved that a convex lattice polytope of volume V has O(V (d−1)/(d+1) ) vertices. Hence a convex lattice polytope in G(m) has O(m d(d−1)/(d+1) ) vertices. By the well-known inequality Bárany and Vershik [8] proved that there are exp O(m d(d−1)/(d+1) ) convex lattice polytopes in G(m). If the polytopes can also be enumerated in this time (as in the planar case [7]), then the runtime in Lemma 4 reduces accordingly.

A better approximation in the plane
If C opt is long and skinny (e.g., ρ is close to 1/n), then the uniform grid G(m) we used in Lemmas 3 and 4 is unsuitable for finding a (1 − ε)-approximation efficiently. Instead, we employ a rotated and stretched grid (an affine copy of G(m)) that has similar orientation and aspect ratio as C opt . This overcomes one of the main difficulties in obtaining a good approximation. Since we do not know the shape and orientation of C opt , we guess these parameters via the minimum area rectangle containing C opt . Proof. The intuitive idea is now to guess a rectangle R that contains C opt such that its area, vol(C opt ), is at least a constant fraction of the area of R, namely vol(R); and then apply Lemma 3 to the rectangle R (as the unit square) to get the desired approximation.
Let B 0 be the minimum area rectangle (or arbitrary orientation) that contains C opt ; see Figure 3 (left). We guess an approximate copy of B 0 . In particular, we guess the lengths of the two sides of B 0 (up to a factor of 2) and the orientation of B 0 (up to an angle of O(1/n)), and then try to position a scaled copy of the guessed rectangle so that that it fully contains C opt . Assume for convenience that n ≥ 10. Denote by a and b the lengths of the two sides of rectangle B 0 , where a ≤ b. It is clear that b ≤ √ 2, the diameter of the unit square. We now show that a ≥ 1/n. From Theorem 1, we have vol( n , for n ≥ 10. From the other direction, vol(C opt ) ≤ √ 2 width(C opt ) ≤ √ 2a. Comparing these bounds yields a ≥ 1/n. It follows that if we want to guess the side lengths of B 0 up to a factor of two, we need to consider only O(log 2 n) possibilities. Indeed, we consider the canonical lengths ℓ i = 2 i /n, for i = 0, . . . , ⌈log 2 n⌉. Let B 1 be a minimum area rectangle with canonical side lengths and the same orientation as B 0 , so that B 0 ⊆ B 1 .
The orientation of a rectangle is given by the angle between one side and the x-axis. We approximate the orientation of B 0 by the n canonical orientations α = iπ/(5n), for i = 0, 1, . . . , 5n− 1. Let B 2 be a congruent copy of B 1 rotated clockwise to the nearest canonical orientation about the center of B 1 . We show that B 1 ⊆ 2B 2 , i.e., a scaled copy of B 2 contains B 1 . Let k ≥ 1 be the minimum scale factor such that B 1 ⊆ kB 2 . Refer to Figure 3 (right). Denote by o the common center of B 1 and B 2 , let x be a vertex of B 1 on the boundary of kB 2 , and let y be the corresponding vertex of kB 2 . Clearly, sin(∠xoy) ≤ π/(5n) since we rotate by at most π/(5n). We also have sin(∠oyx) ≥ (1/n)/ √ 2 > π/(5n), since sin(∠oyx) is the ratio between a side of kB 2 and its diameter (and this ratio is the same for B 0 ). The law of sines yields |ox| ≥ |xy|; and we have |ox| + |xy| ≥ |oy| by the triangle inequality. If follows that |oy| < 2|ox|, and so k ≤ 2 suffices.
So far we have guessed the canonical side lengths and orientation of B 2 , however, we do not know its location in the plane. If a translated copy Consider an arbitrary tiling of the plane with translates of B 2 . By a packing argument, only O(1/ρ) translates intersect the unit square [0, 1] 2 . One of these translates, say B 2 + v, intersects C opt , and hence the rectangle R = 3B 2 + v contains C opt .
We can apply Lemma 3 to the rectangle R (as the unit square) to get the desired approximation. Specifically, let T : R 2 → R 2 be an affine transformation that maps R into the unit square Since ρ = Ω(1/n), the above expression is bounded by O ε −8 n 2 log 2 n , as claimed.
Remark. If ρ = Ω(1) the running time of this planar algorithm is near linear in n.
By doing an exponentially decreasing search for ρ, the running time increases only by a constant factor (this is a geometrically decreasing series, hence the term with the last value of ρ dominates the whole running time). We summarize our result for the plane in the following.

A better approximation in higher dimensions
Given a set S ⊆ [0, 1] d of n points and a parameter ε > 0, we compute an empty convex body C ⊆ [0, 1] 2 such that vol(C) ≥ (1 − ε)vol opt (S). Similarly to the algorithm in Subsection 4.3, we guess a hyperrectangle R that contains C opt such that vol(C opt ) is at least a constant fraction of vol(R); and then apply Lemma 4 to R (as the hypercube) to obtain the desired approximation.
Consider a hyperrectangle B 0 of minimum volume (and arbitrary orientation) that contains C opt . The d edges incident to a vertex of a hyperrectangle B are pairwise orthogonal. We call these d directions the axes of B; and the orientation of B is the set of its axes.
We next enumerate all possible discretized hyperrectangles of volume Ω(1/n), guessing the lengths of their axes, their orientations, and their locations as follows: Guess the length of every axis up to a factor of 2. Since the minimum length of an axis in our case is Ω(1/n) and the maximum is √ d, the number of possible lengths to be considered is O log d n . Let B 1 be a hyperrectangle of minimum volume with canonical side lengths and the same orientation as B 0 such that B 0 ⊆ B 1 .
We can discretize the orientation of a hyperrectangle as follows. We spread a dense set of points on the sphere of directions, with angular distance O(1/n) between any point on the sphere and its closest point in the chosen set. O(n d−1 ) points suffice for this purpose. We try each point as the direction of the first axis of the hyperrectangle, and then generate the directions of the remaining axes analogously in the orthogonal hyperplane for the chosen direction. Overall, this generates Successively replace each axis of B 1 by an approximate axis that makes an angle at most α = 1/(cn) with its corresponding axis, where c = c(d) is a constant depending on d. Let B 2 be a congruent copy of B 1 obtained in this way. If c = c(d) is sufficiently small, then B 1 ⊆ 2B 2 .
Consider a tiling of R d with translates of B 2 . Note that only O(1/vol(C opt )) = O(n) translates intersect the unit cube [0, 1] d . One of these translates B 2 + v intersects C opt , and then the hyperrectangle R = 3B 2 + v contains C opt . Since C opt (S) takes a constant fraction of the volume of R, we can deploy Lemma 4 in this case, and get the desired (1 − ε)-approximation in exp O ε −d(d−1)/(d+1) log ε −1 time. Putting everything together, we obtain the following. Remark. Consider a set S of n points in R d . The approximation algorithm we have presented can be modified to approximate the largest empty tile, i.e., the largest empty convex body contained in conv(S), rather than [0, 1] d . The running time is slightly worse, since we need to take the boundary of conv(S) into account. We omit the details.

Conclusions
In this section we briefly outline two exact algorithms for finding the largest area empty convex polygon and the largest area empty triangle amidst n points in the unit square. At the end we list a few open problems.
Largest area convex polygon. Let S ⊂ U = [0, 1] 2 , where |S| = n. Let T be the set of four vertices of U . Observe that the boundary of an optimal convex body, C opt , contains at least two points from S ∪ T . By convexity, the midpoint of one of these O(n 2 ) segments lies in C opt . For each such midpoint m, create a weakly simple polygon P m by connecting each point p ∈ S to the boundary of the square along the ray mp. The polygon P m has O(n) vertices and is empty of points from S in its interior. Then apply the algorithm of Chang and Yap [10] for the potatopeeling problem (mentioned in Section 1) in these O(n 2 ) weakly simple polygons. The algorithm computes a largest area empty convex polygon contained in a given (non-convex) polygon with n vertices in O(n 7 ) time. Finally, return the largest convex polygon obtained this way. The overall running time is O(n 9 ). The running time can be reduced to O(n 8 log n) as follows. Instead of considering the O(n 2 ) midpoints, compute a set P of O(n log n) points so that every convex set of area at least 2/(n + 4) contains at least one of these points. In particular, C opt contains a point from P . The set P can be computed by starting with a O(n) × O(n) grid, and then computing an ε-net for it, where ε = O(1/n), using discrepancy [30]. The running time of this deterministic procedure is roughly O(n 2 ), and the running time of the overall algorithm improves to O(n 7 · n log n) = O(n 8 log n).
Largest area empty triangle. The same reduction can be used for finding largest area empty triangle contained in U , resulting in O(n 2 ) weakly simple polygons P m . Then the algorithm of Melissaratos and Souvaine [31] for finding a largest area triangle contained in a polygon is applied to each of these O(n 2 ) polygons. The algorithm computes such a triangle in O(n 4 ) time, given a polygon with n vertices. Finally, return the largest convex polygon obtained this way. The overall running time is O(n 6 ). Via the ε-net approach (from the previous paragraph) the running time of the algorithm improves to O(n 4 · n log n) = O(n 5 log n).

Open questions.
Interesting questions remain open regarding the structure of optimal Steiner convex partitions and the computational complexity of computing such partitions. Other questions relate to the problem of finding the largest empty convex body in the presence of points. We list some of them: (1) Is there a polynomial-time algorithm for computing a minimum Steiner convex partition of a given set of n points in R d ? Is there one for points in the plane?
(2) Is there a constant-factor approximation algorithm for the minimum Steiner convex partition of an arbitrary point set in R d (without the general position restriction)? Is there one for points in the plane?
(3) For d > 2, the running time of our approximation algorithm for the maximum empty polytope has a factor of the form n O(d 2 ) . It seems natural to conjecture that this term can be reduced to n O(d) . Another issue of interest is extending Lemma 2 to higher dimensions for a faster overall algorithm.
(4) Given n points in [0, 1] d , the problem of finding the largest convex body in [0, 1] d that contains up to k (outlier) points naturally suggests itself and appears to be also quite challenging.