Delaunay Triangulations with Predictions

Algorithms with predictions have received considerable attention in recent years. Given a problem instance and a prediction for the solution, we would like to devise an algorithm that solves the problem in such a way that the algorithm performs “better” when the prediction is “close” to the correct solution. “Better” may be in terms of running time, or approximation factor (for optimization problems), or competitive ratio (for on-line problems). Different measures of “closeness” may be used depending on the problem, or alternatively we may assume that predictions obey some probabilistic models. The paradigm has been applied to problems in diverse areas, such as scheduling, classical graph problems, resource allocation, mechanism design, facility location, and graph searching to name a few; see the survey in [52] and the regularly updated online list of papers in [1].

In contrast, work on algorithms with predictions in computational geometry have been few and far between. The goal of this paper is to promote this direction of research. Primarily, we concentrate on the use of predictions to improve the running time of geometric algorithms.

As many fundamental problems in classical computational geometry may be viewed as two- or higher-dimensional generalizations of standard one-dimensional problems like sorting and searching, we start our discussion with the 1D sorting problem. Here, algorithms with predictions are related to the historically well-explored topic of adaptive sorting algorithms [32] – if we are given a sequence that is already close to sorted, we can sort faster. There are different ways to measure distance to sortedness. For example:

• $\mathrm{■}$

  Consider $D_{\text{inv}}$, the number of inversions in the given sequence (i.e., the number of out-of-order pairs). In terms of this parameter, good old insertion-sort runs in $O(n+D_{\text{inv}})$ time and already beats $O(n\log n)$ algorithms when $D_{\text{inv}}=O(n)$. Better algorithms are known running in $\mathrm{\Theta }(n\log (1+D_{\text{inv}}/n))$ time, which is optimal in comparison-based models.

• $\mathrm{■}$

  Or consider $D_{\text{rem}}$, the minimum number of elements to remove (the “wrong” elements) so that the remaining list is sorted. Here, the optimal time bound is known to be $\mathrm{\Theta }(n+D_{\text{rem}}\log n)$.

• $\mathrm{■}$

  Or consider $D_{\text{runs}}$, the number of runs (i.e., the number of increasing contiguous subsequence that the input can be partitioned into). Here, the optimal time bound is known to be $\mathrm{\Theta }(n\log D_{\text{runs}})$. Unfortunately, even when $D_{\text{runs}}$ is as low as $n^{0.1}$, it is not possible to beat $n\log n$ under this distance measure.

The survey in [32] gives an extensive background and a (long) list of other parameters that have been considered in the algorithms literature. Adaptive sorting has continued to be studied (some random recent examples include [7, 53, 4]), and some recent papers [6, 47] examine newer models for sorting with predictions.

Can we obtain similar results for geometric problems in 2D? In this paper, we pick arguably one of the most basic problems in computational geometry to start this study: namely, the computation of the Delaunay triangulation, denoted $\text{DT}(P)$, of a set $P$ of $n$ points in the plane [27]. One way to define $\text{DT}(P)$ is via the empty circle property: three points $p,p^{\prime },p^{\prime \prime }\in P$ form a Delaunay triangle iff the circumcircle through $p,p^{\prime },p^{\prime \prime }$ does not have any points of $P$ strictly inside. (Throughout this paper, we assume for simplicity that the input is in general position, e.g., no 4 input points are co-circular.)

Delaunay triangulations are combinatorially equivalent to Voronoi diagrams (under the Euclidean metric) by duality. (It is a little more convenient to work with Delaunay triangulations than Voronoi diagrams when specifying the prediction $G$ or defining closeness, and so all our results will be stated in terms of the former.) Delaunay triangulations and Voronoi diagrams are among the most ubiquitous geometric structures, with countless applications, and lie at the very core of computational geometry.

In the prediction setting, we are given not just the input point set $P$, but also a triangulation $G$ of $P$ (i.e., a decomposition of the convex hull of $P$ into nonoverlapping triangles, using all the points of $P$ as vertices, but no Steiner points). The goal is to devise an algorithm to compute $\text{DT}(P)$ faster, when the “distance” between $G$ and $\text{DT}(P)$ is small. In the traditional setting, several standard textbook algorithms are known for computing the Delaunay triangulation in $O(n\log n)$ time, e.g., via divide-and-conquer, plane sweep, and randomized incremental construction [8, 27, 5, 34], which is optimal in the worst case. Without closeness, an arbitrary triangulation $G$ of $P$ still provides some extra information about the input, but not enough to beat $n\log n$: it is not difficult to see that $\mathrm{\Omega }(n\log n)$ remains a lower bound by reduction from 1D sorting (see Figure 1). We want to do better than $n\log n$ when $G$ is close. As we will see, designing faster adaptive algorithms for 2D Delaunay triangulations is technically more difficult, and so more interesting, than for 1D sorting.

There are a number of natural parameters to measure distance to the Delaunay triangulation, as we define in Table 1. Some of these definitions (e.g., $D,D_{\text{flip}},D_{\text{cross}}$) may be used to compare any two triangulations, while others (e.g., $D_{\text{vio}},D_{\text{local}}$) are specific to the Delaunay triangulation. For all these parameters, the value is 0 if and only if $G=\text{DT}(P)$. Some relationships between these parameters are obvious, e.g., $D_{\text{cross}}\le O(n{d}_{\text{cross}})$ and $D_{\text{vio}}\le O(n{d}_{\text{vio}})$, and it is known that

The inequality $D_{\text{local}}\le D$ holds because all non-locally Delaunay edges can’t be actual Delaunay edges, while the last inequality $D_{\text{flip}}\le D_{\text{cross}}$ is nontrivial and was proved by Hanke et al. [36].

There are examples when $D_{\text{local}}$ is a constant and $D$ is $\mathrm{\Omega }(n)$. In a very vague sense, $D_{\text{local}}$ is more similar to $D_{\text{runs}}$ (the number of runs) for sorting, while $D$ is more similar to $D_{\text{rem}}$ and $D_{\text{flip}}$ and $D_{\text{cross}}$ are more similar to $D_{\text{inv}}$. Indeed, it is not difficult to prove an $\mathrm{\Omega }(n\log D_{\text{local}})$ lower bound for computing the Delaunay triangulation in terms of $D_{\text{local}}$, by reduction from sorting with $\mathrm{\Theta }(D_{\text{local}})$ runs, or equivalently, merging $\mathrm{\Theta }(D_{\text{local}})$ sorted lists (see Figure 2). So, $D_{\text{local}}$ may not be an ideal parameter for the problem (even when $D_{\text{local}}$ is as small as $n^{0.1}$, we can’t hope to beat $n\log n$).

Arguably the simplest and most natural of all these parameters is $D$, the number of “wrong” edges in $G$ that are not in $\text{DT}(P)$ (which is the same as the number of edges in $\text{DT}(P)$ that are in $G$, since all triangulations of the same point set have the same number of edges). To initiate the study of algorithms with predictions for the Delaunay triangulation, we thus begin with the following simple question:

Question 1: what is the best 2D Delaunay triangulation algorithm with running time analyzed as a function of $n$ and the parameter $D$?

This simple question already turns out to be quite challenging. An upper bound on time complexity in terms of $D$ would automatically imply an upper bound in terms of $D_{\text{flip}}$ and $D_{\text{cross}}$ by the above inequalities. We can ask a similar question for other parameters such as $d_{\text{cross}}$ and $d_{\text{vio}}$.

Although Question 1 provides one way to model predictions, often errors are probabilistic in nature. Thus, for various problems, algorithms with predictions have sometimes been studied under probabilistic models. To define a simple probabilistic model of prediction for the Delaunay triangulation, one could just consider tossing an independent coin to determine whether an edge of $\text{DT}(P)$ should be placed in $G$ (similar to a model studied by Cohen-Addad et al. [26] for the max-cut problem). But there is one subtlety for our problem: the input $G$ is supposed to be a triangulation, containing wrong edges (but not containing self-crossings). We propose the following 2-stage model: first, each edge of $\text{DT}(P)$ is placed in a set $R$ independently with a given probability $\rho$; second, $R$ is completed to a triangulation of $P$, but we make no assumption about how this second step is done, i.e., $R$ can be triangulated adversarially (this makes results in this model a bit more general).

Question 2: under the above probabilistic model, what is the best 2D Delaunay triangulation algorithm with running time analyzed as a function of $n$ and $\rho$?

In this model, the expected number of “wrong” edges $D$ is upper-bounded by $(1-\rho )m$, where $m\le 3n-6$ is the number of edges in any triangulation of $P$. The hope is that the problem may be solved faster under the probabilistic model than by applying a worst-case algorithm in terms of the parameter $D$.

We could consider other probabilistic models too. For example, one variant is to pick a random edge in $\text{DT}(P)$, perform a flip if its two incident triangles form a convex quadrilateral (or do nothing otherwise), and repeat for $(1-\rho )m$ steps. In this model, the flip distance $D_{\text{flip}}$ is guaranteed to be at most $(1-\rho )m$; in contrast, the above model allows for potentially larger flip distances and so is in some sense stronger.

We obtain a number of results addressing Questions 1 and 2 and related questions:

1. 1.

   We present several different algorithms for computing $\text{DT}(P)$ from $G$ in time sensitive to the parameter $D$. One deterministic algorithm runs in $O(n+D{\log }^{3}n)$ time. A randomized algorithm runs in $O((n+D\log n){\log }^{\ast }n)$ time, where ${\log }^{\ast }$ is the (very slow-growing) iterated logarithmic function. Our final randomized algorithm eliminates the ${\log }^{\ast }n$ factor and runs in $O(n+D\log n)$ expected time, which is optimal since an $\mathrm{\Omega }(\max \{n,D\log D\})=\mathrm{\Omega }(n+D\log n)$ lower bound follows by reduction from sorting $D$ numbers (we can just embed the configuration in Figure 1 on $\mathrm{\Theta }(D)$ points into a larger Delaunay triangulation). Thus, we have completely answered Question 1 (modulo the use of randomization). All of our algorithms do not require knowing the value of $D$ a priori.

   It is interesting to note that our optimal randomized algorithm runs in linear time even when the number $D$ of “wrong” edges is as large as $n/\log n$.

   An immediate consequence is that the problem can be solved in $O(n+D_{\text{flip}}\log n)$ expected time – we get linear running time as long as the flip distance is at most $n/\log n$.

2. 2.

   On the probabilistic model in Question 2, we obtain another algorithm for computing $\text{DT}(P)$ from $G$ in $O(n\log \log n+n\log (1/\rho ))$ time with high probability (the algorithm is always correct). For any constant $\rho$, the running time is thus almost linear except for the small $\log \log n$ factor. This is notable when compared against our $O(n+D\log n)$ algorithm, since in this probabilistic setting, $D$ is expected to be near $(1-\rho )m$. For example, if we set $\rho =0.1$ (or even $\rho =1/\log n$), we can still substantially beat $n\log n$ even though a majority of the edges are wrong!

   We can also get a similar result for the random flip model mentioned earlier.

3. 3.

   We also obtain results sensitive to other parameters. For example, in terms of $d_{\text{vio}}$, we describe a deterministic algorithm for computing $\text{DT}(P)$ from $G$ in $O(n{\log }^{\ast }n+n\log d_{\text{vio}})$ time. Note that there are examples where $d_{\text{vio}}=1$ but $D=\mathrm{\Omega }(n)$ (i.e., a fraction of the edges are wrong!), so our $D$-sensitive algorithms can’t be used directly. The parameter $d_{\text{vio}}$ is related to a notion called order-$k$ Delaunay triangulation¹¹1The name is somewhat misleading as it is not equivalent to the order-$k$ Voronoi diagram in the dual. The set of all candidate triangles for order-$k$ Delaunay triangulations does correspond to the vertices of the order-$k$ Voronoi diagram, but an order-$k$ Delaunay triangulation is not unique, comprising of a subset of the candidate triangles. in the literature [35] (with motivation very different from ours): a triangulation $G$ is an “order-$k$ Delaunay triangulation” precisely when $d_{\text{vio}}\le k$. Our result implies the following interesting consequence: given an order-$O(1)$ Delaunay triangulation, we can recover the actual Delaunay triangulation in almost linear ($O(n{\log }^{\ast }n)$) time!

4. 4.

   We similarly describe an $O(n{\log }^{\ast }n+n\log d_{\text{cross}})$-time deterministic algorithm for computing $\text{DT}(P)$ from $G$. Again, there are examples where $d_{\text{cross}}=1$ but $D=\mathrm{\Omega }(n)$. Although the parameter $d_{\text{cross}}$ may not appear to be the most natural at first, it turns out to be the key that allows us to obtain items 2 and 3 above in an elegant, unified way: it is not difficult to show that $d_{\text{cross}}=O((1/\rho )\log n)$ with high probability in the probabilistic model, and we are also able to show that $d_{\text{cross}}\le O(d_{\text{vio}}^2)$ always holds.

5. 5.

   Although we have focused on the Delaunay triangulation problem, we can obtain similar results for a few related problems. For example, consider the problem of computing the Euclidean minimum spanning tree $\text{EMST}(P)$ of a set $P$ of $n$ points in 2D. Given a spanning tree $T$, we can compute $\text{EMST}(P)$ from $T$ in $O((n+D_{\text{emst}}\log n){\log }^{\ast }n)$ expected time, where $D_{\text{emst}}$ is the number of “wrong” edges in $T$ that are not in $\text{EMST}(P)$. In fact, this result, and the idea of using a spanning tree in the first place, will be instrumental to the eventual design of our $O(n+D\log n)$-time randomized algorithm for Delaunay triangulations.

   We get similar results for EMST in the probabilistic model or in terms of other parameters, and also for computing other related geometric graphs such as the relative neighborhood graph and the Gabriel graph in the plane [38]. We also get similar results for the problem of computing the convex hull of $n$ points in 3D that are in convex position (by a standard lifting map, Delaunay triangulations in 2D maps to convex hulls in 3D) – the convex position assumption is necessary since one could show an $\mathrm{\Omega }(n\log n)$ lower bound for checking whether $n/2$ points are all inside the convex hull of another $n/2$ points, even when we are given the correct convex hull.

We are not aware of any prior work in computational geometry that explicitly addresses the paradigm of algorithms with predictions, except for a recent paper by Cabello and Giannopoulos [10] on an on-line geometric search problem (where predictions are used to obtain competitive ratio bounds). Implicitly related are adaptive data structures, and these have been studied in computational geometry; as one example, see Iacono and Langerman’s work [37] on “proximate” planar point location (where a query can be answered faster when given a finger near the answer).

Our work is related to, and may be regarded as a generalizaton of, work on the verification of geometric structures: how to check whether a given answer is indeed correct, using a simpler or faster algorithm than an algorithm for computing the answer. This line of work has received much attention from computational geometers in the 1990s [29, 50, 51], motivated by algorithm engineering. For the Delaunay triangulation problem, it is easy to verify whether $G$ is equal to $\text{DT}(G)$ in linear time, by checking that all edges are locally Delaunay. Verifying correctness corresponds to the case $D=0$, and our results provides more powerful extensions: we can correct errors (as many as $n/\log n$ wrong edges) in linear time. Viewing our results in the lens of correcting errors in the output provides further motivation, even for readers not necessarily interested in learning-augmented algorithms.

Some recent work such as [31] studies geometric algorithms in models that tolerate errors in the computation process itself (e.g., noisy or faulty comparisons), but our work seems orthogonal to this research direction. On the other hand, there has been considerable prior work in computational geometry on models that capture errors or imprecisions in the input. For example, in one result by Löffler and Snoeyink [46], an imprecise input point is modelled as a range, specifically, a unit disk (whose center represents its estimated position); it is shown that $n$ disjoint unit disks in the plane can be preprocessed in $O(n\log n)$ time, so that given $n$ actual points lying in these $n$ disks, their Delaunay triangulation can be computed in $O(n)$ time. Similar results have been obtained for different types of regions and different problems. Our model provides an alternative way to address imprecise input: in the preprocessing phase, we can just compute the Delaunay triangulation $G$ of the estimated points; when the actual points $P$ are revealed, we can correct $G$ to obtain $\text{DT}(P)$ using our algorithms.²²2One technical issue is that $G$ may have self-crossings when re-embedded into the actual point set $P$, but at least one of our algorithms works even when $G$ has self-crossings (see Corollary 6). An advantage of our model is its simplicity: we don’t need to be given any geometric range or distribution for each input point; rather, accuracy is encapsulated into a single parameter $D$.

Even without relevance to algorithms with paradigms or coping with errors in the input/output, Question 1 is appealing purely as a classical computational geometry problem. In the literature, a number of papers have explored conditions under which Delaunay triangulations can be computed in linear or almost linear time. For example, Devillers [28] showed that given the Euclidean minimum spanning tree of a 2D point set, one can compute the Delaunay triangulation in $O(n{\log }^{\ast }n)$ expected time (this was later improved to deterministic linear time [22]); two Delaunay triangulations $\text{DT}(P_1)$ and $\text{DT}(P_2)$ can be merged into $\text{DT}(P_1\cup P_2)$ in linear time [40]; conversely, a Delaunay triangulation $\text{DT}(P_1\cup P_2)$ can be split into $\text{DT}(P_1)$ and $\text{DT}(P_2)$ in linear time (by Chazelle et al. [15]); the Delaunay triangulation of the vertices of a convex polygon can be computed in linear time (by Aggarwal et al. [2, 21]); the Delaunay triangulation for 2D points with integer coordinates can be computed in the same time complexity as sorting in the word RAM model (by Buchin and Mulzer [9]); etc. A powerful work by Löffler and Mulzer [45] unifies many of these results by showing equivalence of 2D Delaunay triangulations with quadtree construction and related problems. Although some of these tools will be useful later, it isn’t obvious how they can used to resolve Questions 1 and 2.

Although we view our main contributions to be conceptual (introducing a new research direction in computational geometry), our algorithms are also technically interesting in our opinion. The solutions turn out to be a lot more challenging than 1D adaptive sorting algorithms or the standard textbook Delaunay triangulation algorithms:

• $\mathrm{■}$

  Our deterministic $O(n+D{\log }^{3}n)$-time algorithm (given in Section 2) uses planar graph separators or $r$-divisions in the given graph $G$ (for classical computational geometry problems, Clarkson–Shor-style sampling-based divide-and-conquer [23, 24, 19] is more commonly used than planar graph separators).

• $\mathrm{■}$

  Our randomized $O((n+D\log n){\log }^{\ast }n)$-time algorithm (given in Section 3) is obtained via an adaptation of the abovementioned randomized $O(n{\log }^{\ast }n)$-time algorithm for computing the Delaunay triangulation from the EMST, by Devillers [28] (which in turn was based on the randomized $O(n{\log }^{\ast }n)$-time algorithms by Seidel [54] and Clarkson, Tarjan, and Van Wyck [25] for the polygon triangulation problem). The adaptation is more or less straightforward, but our contribution is in realizing that these randomized $O(n{\log }^{\ast }n)$-time algorithms are relevant in the first place.

• $\mathrm{■}$

  Our final randomized $O(n+D\log n)$-time algorithm (given in Section 4) eliminates the extra ${\log }^{\ast }n$ factor by incorporating a more complicated divide-and-conquer involving planar graph separators or $r$-divisions, this time applied to actual Delaunay triangulations, not the input graph $G$. (Curiously, ideas from certain output-sensitive convex hull algorithms [30, 12] turn out to be useful.)

• $\mathrm{■}$

  For our deterministic $O(n{\log }^{\ast }n+n\log d_{\text{cross}})$-time algorithm (given in Section 5), which as mentioned is the key to our result in the probabilistic model as well as our $d_{\text{vio}}$-sensitive result, we return to using $r$-divisions in the given graph $G$. The way the ${\log }^{\ast }n$ factor comes up is interesting – it arises not from the use of known randomized $O(n{\log }^{\ast }n)$ algorithms like Seidel’s, but from an original recursion. Although our $D$-sensitive algorithm does not directly yield a $d_{\text{cross}}$-sensitive result, it will turn out to be an important ingredient in the final algorithm.

Many interesting questions remain from our work (which we see as a positive, as our goal is to inspire further work in this research direction):

• $\mathrm{■}$

  Could the extra ${\log }^{\ast }n$ factors be eliminated in some of our results?

• $\mathrm{■}$

  Could the extra $\log \log n$ factor in our result in the probabilistic model be removed?

• $\mathrm{■}$

  Could we get a time bound in terms of $D_{\text{flip}}$ better than $O(n+D_{\text{flip}}\log n)$? For example, is there a Delaunay triangulation algorithm that runs in $O(n+D_{\text{flip}})$ time?

• $\mathrm{■}$

  Similarly, is there a Delaunay triangulation algorithm that runs in $O(n+D_{\text{cross}})$ time? Or, better still, $O(n\log (1+D_{\text{cross}}/n))$ time? Or $O(n+D_{\text{vio}})$ or $O(n\log (1+D_{\text{vio}}/n))$ time?

More generally, we hope that our work will prompt other models of prediction to be studied for the Delaunay triangulation problem, or similar models for other fundamental problems in computational geometry and beyond.