Four Soviets Walk the Dog: Improved Bounds for Computing the Fréchet Distance
Four Soviets Walk the Dog: Improved Bounds for Computing the Fréchet Distance
Kevin Buchin 0 1 2 3
Maike Buchin 0 1 2 3
Wouter Meulemans 0 1 2 3
Wolfgang Mulzer 0 1 2 3
0 Department of Mathematics and Computer Science , TU Eindhoven, PO Box 513, 5600 MB Eindhoven , The Netherlands
1 A preliminary version appeared as K. Buchin , M. Buchin, W. Meulemans, and W. Mulzer. Four Soviets Walk the Dogwith an Application to Alt's Conjecture. Proc. 25th SODA, pp. 13991413, 2014
2 Institut für Informatik, Freie Universität Berlin , 14195 Berlin , Germany
3 Mathematik, Ruhr Universität Bochum , 44780 Bochum , Germany
Given two polygonal curves in the plane, there are many ways to define a notion of similarity between them. One popular measure is the Fréchet distance. Since it was proposed by Alt and Godau in 1992, many variants and extensions have been studied. Nonetheless, even more than 20 years later, the original O (n2 log n) algorithm by Alt and Godau for computing the Fréchet distance remains the state of the art (here, n denotes the number of edges on each curve). This has led Helmut Alt to conjecture that the associated decision problem is 3SUMhard. In recent work, Agarwal et al. show how to break the quadratic barrier for the discrete version of the Fréchet distance, where one considers sequences of points instead of polygonal curves. Building on their work, we give a randomized algorithm to compute the Fréchet distance between two polygonal curves in time O (n2√log n(log log n)3/2) on a pointer machine and in

time O(n2(log log n)2) on a word RAM. Furthermore, we show that there exists an
algebraic decision tree for the decision problem of depth O(n2−ε), for some ε > 0.
We believe that this reveals an intriguing new aspect of this wellstudied problem.
Finally, we show how to obtain the first subquadratic algorithm for computing the
weak Fréchet distance on a word RAM.
1 Introduction
Shape matching is a fundamental problem in computational geometry, computer
vision, and image processing. A simple version can be stated as follows: given a
database D of shapes (or images) and a query shape S, find the shape in D that most
resembles S. However, before we can solve this problem, we first need to address an
issue: what does it mean for two shapes to be “similar”? In the mathematical literature,
on can find many different notions of distance between two sets, a prominent
example being the Hausdorff distance. Informally, the Hausdorff distance is defined as the
maximal distance between two elements when every element of one set is mapped
to the closest element in the other. It has the advantage of being simple to describe
and easy to compute for discrete sets. In the context of shape matching, however, the
Hausdorff distance often turns out to be unsatisfactory: it does not take the continuity
of the shapes into account. There are well known examples where the distance fails to
capture the similarity of shapes as perceived by human observers [6].
In order to address this issue, Alt and Godau introduced the Fréchet distance into
the computational geometry literature [8,45]. They argued that the Fréchet distance is
better suited as a similarity measure, and they described an O(n2 log n) time algorithm
to compute it on a real RAM or pointer machine.1 Since Alt and Godau’s seminal
paper, there has been a wealth of research in various directions, such as extensions
to higher dimensions [7,23,26,28,33,46], approximation algorithms [9,10,37], the
geodesic and the homotopic Fréchet distance [29,34,38,48], and much more [2,22,25,
35,36,51,54,55]. Most known approximation algorithms make further assumptions
on the curves, and only an O(n2)time approximation algorithm is known for arbitrary
polygonal curves [24]. The Fréchet distance and its variants, such as dynamic
timewarping [13], have found various applications, with recent work particularly focusing
on geographic applications such as mapmatching tracking data [15,63] and moving
objects analysis [19,20,47].
Despite the large amount of published research, the original algorithm by Alt and
Godau has not been improved, and the quadratic barrier on the running time of the
associated decision problem remains unbroken. If we cannot improve on a quadratic
bound for a geometric problem despite many efforts, a possible culprit may be the
underlying 3SUMhardness [44]. This situation induced Helmut Alt to make the
following conjecture.2
1 For a brief overview of the different computational models in this paper, refer to Appendix.
2 Personal communication 2012, see also [6].
Conjecture 1.1 (Alt’s conjecture) Let P, Q be two polygonal curves in the plane.
Then it is 3SUMhard to decide whether the Fréchet distance between P and Q is at
most 1.
1.1 Our Contribution
We address the question by Agarwal et al. and show how to extend their approach
to the Fréchet distance between two polygonal curves. Our algorithm requires total
expected time O(n2√log n(log log n)3/2). This is the first algorithm with a running
time of o(n2 log n) and constitutes the first improvement for the general case since the
original paper by Alt and Godau [8]. To achieve this running time, we give the first
subquadratic algorithm for the decision problem of the Fréchet distance. We emphasize
that these algorithms run on a real RAM/pointer machine and do not require any
bitmanipulation tricks. Therefore, our results are more in the line of Chan’s recent
subcubictime algorithms for allpairsshortest paths [30,31] or recent
subquadratictime algorithms for minplus convolution [16] than the subquadratictime algorithms
for 3SUM due to Baran et al. [12]. If we relax the model to allow constant time
tablelookups, the running time can be improved to be almost quadratic, up to O(log log n)
factors. As in Agarwal et al., our results are achieved by first giving a faster algorithm
for the decision version, and then performing an appropriate search over the critical
values to solve the optimization problem.
In addition, we show that nonuniformly, the Fréchet distance can be computed in
subquadratic time. More precisely, we prove that the decision version of the problem
can be solved by an algebraic decision tree [11] of depth O(n2−ε), for some fixed ε > 0.
It is, however, not clear how to implement this decision tree in subquadratic time, which
hints at a discrepancy between the decision tree and the uniform complexity of the
Fréchet problem.
Finally, we consider the weak Fréchet distance, where we are allowed to walk
backwards along the curves. In this case, our framework allows us to achieve a subquadratic
algorithm on the work RAM. Refer to Table 1 for a comprehensive summary of our
results.
3 It is well known that the four Russians are not actually Russian, so we refer to them as four Soviets in
Table 1 Summary of results
O n2√log n(log log n)3/2
We distinguish the continuous and the discrete Fréchet distance in the decision (D) and the computation (C)
version. We also consider the weak continuous Fréchet distance. The computational models are the pointer
machine (PM), the word RAM (WRAM) or the algebraic decision trees (DT). The old bounds are due to
Alt and Godau [8] (continuous and weak) and Eiter and Mannila [39] (discrete)
1.2 Recent Developments
Recently, Ben Avraham et al. [14] presented a subquadratic algorithm for the discrete
Fréchet distance with shortcuts that runs in O (n4/3 log3 n) time. This running time
resembles, at least superficially, our result on algebraic computation trees for the
general discrete Fréchet distance.
When we initially announced our results, we believed that they provided strong
evidence that Alt’s conjecture is false. Indeed, for a long time it was conjectured
that no subquadratic decision tree exists for 3SUM [57] and an (n2) lower bound
is known in a restricted linear decision tree model [4, 40]. However, in a recent–
and in our opinion quite astonishing–result, Grønlund and Pettie showed that if we
allow only slightly more powerful algebraic decision trees than in the previous lower
bounds, one can decide 3SUM nonuniformly in O (n2−ε) steps, for some fixed ε > 0
[52]. They also show that this leads to a general subquadratic algorithm for 3SUM, a
situation very similar to the Fréchet distance as described in the present paper. Thus,
despite some interesting developments, the status of Alt’s conjecture remains as open
as before. However, we can now see that there exists a wide variety of efficiently
solvable problems such as (in addition to 3SUM and the Fréchet distance) Sorting
X + Y [41], Min Plus Convolution [16], or finding the Delaunay triangulation
for a point set that has been sorted in two orthogonal directions [27], for which there
seems to be a noticeable gap between the decision tree complexity and the uniform
complexity.
In our initial announcement, we also asked whether, besides 3SUMhardness, there
may be other reasons to believe that the quadratic running time for the Fréchet distance
cannot be improved. Karl Bringmann provided an interesting answer to this question
by showing that any algorithm for the Fréchet distance with running time O(n2−ε),
for some fixed ε > 0, would violate the strong exponential time hypothesis (SETH)
[17]. These results were later refined and improved to show that the lower bound holds
in basically all settings (with the notable exception of the onedimensional continuous
Fréchet distance, which is still unresolved) [18]. We believe that these developments
show that the Fréchet distance still holds many interesting aspects to be discovered
and remains an intriguing object of further study.
2 Preliminaries and Basic Definitions
Let P and Q be two polygonal curves in the plane, defined by their vertices
p0, p1, . . . , pn and q0, q1, . . . , qn . Depending on the context, we interpret P and
Q either as sequences of n and n edges, or as continuous functions P : [0, n] → R2
and Q : [0, n] → R2. In the latter case, we have P(i + λ) = (1 − λ) pi + λpi+1 for
i = 0, . . . , n − 1 and λ ∈ [0, 1], and similarly for Q. Let be the set of all continuous
and nondecreasing functions α : [0, 1] → [0, n] with α(0) = 0 and α(1) = n. The
Fréchet distance between P and Q is defined as
where · denotes the Euclidean distance.
The classic approach to computing dF ( P, Q) uses the freespace diagram
FSD( P, Q). It is defined as
FSD( P, Q) := {(x , y) ∈ [0, n] × [0, n]  P(x ) − Q(y) ≤ 1}.
In other words, FSD( P, Q) is the subset of the joint parameter space for P and Q
where the corresponding points on the curves have distance at most 1, see Fig. 1.
The structure of FSD( P, Q) is easy to describe. Let R := [0, n] × [0, n] be the
ground set. We subdivide R into n2 cells C (i, j ) = [i, i + 1] × [ j, j + 1], for i,
j = 0, . . . , n − 1. The cell C (i, j ) corresponds to the edge pair ei+1 and f j+1, where
ei+1 is the (i + 1)th edge of P and f j+1 is the ( j + 1)th edge of Q. Then the set
F (i, j ) := FSD( P, Q) ∩ C (i, j ) represents all pairs of points on ei+1 × f j+1 with
distance at most 1. Elementary geometry shows that F (i, j ) is the intersection of
C (i, j ) with an ellipse [8]. In particular, the set F (i, j ) is convex, and the intersection
of FSD( P, Q) with the boundary of C (i, j ) consists of four (possibly empty) intervals,
one on each side of ∂C (i, j ). We call these intervals the doors of C (i, j ) in FSD( P, Q).
A door is said to be closed if the interval is empty, and open otherwise.
A path π in FSD( P, Q) is bimonotone if it is both x  and ymonotone, i.e., every
vertical and every horizontal line intersects π in at most one connected component. Alt
and Godau observed that it suffices to decide whether there exists a bimonotone path
from (0, 0) to (n, n) inside FSD( P, Q). We define the reachable region reach( P, Q)
Fig. 1 Two polygonal curves P and Q, together with their associated freespace diagram. The reachable
region reach(P, Q) is shown in blue. For example, the white area in C(2, 1), denoted F(2, 1), corresponds
to all points on the third edge of P and the second edge of Q that have distance at most 1. As (5, 5) ∈
reach(P, Q), we have dF (P, Q) ≤ 1
as the set of points in FSD( P, Q) that are reachable from (0, 0) on a bimonotone
path. Then, dF ( P, Q) ≤ 1 if and only if (n, n) ∈ reach( P, Q), see Fig. 1. It is not
necessary to compute all of reach( P, Q): since FSD( P, Q) is convex inside each
cell, we only need the intersections reach( P, Q) ∩ ∂C (i, j ). The sets defined by
reach( P, Q) ∩ ∂C (i, j ) are subintervals of the doors of the freespace diagram, and
they are defined by endpoints of doors in the freespace diagram in the same row or
column. We call the intersection of a door with reach( P, Q) a reachdoor. The
reachdoors can be found in O(n2) time through a simple breadthfirsttraversal of the cells
[8]. In the next sections, we show how to obtain the crucial information, i.e., whether
(n, n) ∈ reach( P, Q), in o(n2) time instead.
2.1 Basic Approach and Intuition
In our algorithm for the decision problem, we basically want to compute reach( P, Q).
But instead of propagating the reachability information cell by cell, we always group
τ × τ cells (with 1 τ n) into an elementary box of cells. When processing a
box, we can assume that we know which parts of the left and the bottom boundary
of the box are reachable. That is, we know the reachdoors on the bottom and left
boundary, and we need to compute the reachdoors on the top and right boundary of
the elementary box. These reachdoors are determined by the combinatorial structure
of the box. More specifically, suppose we know for every row and column the order
of the door endpoints (including for the reachdoors on the left and bottom boundary).
Then, we can deduce which of these door boundaries determine the reachdoors on
the top and right boundary. We call the sequence of these orders, the (full) signature
of the box.
The total number of possible signatures is bounded by an expression in terms of τ .
Thus, if we pick τ sufficiently small compared to n, we can precompute for all possible
signatures the reachdoors on the top and right boundary, and build a data structure to
query these quickly (Sect. 3). Since the reachdoors on the bottom and left boundary are
required to make the signature, we initially have only incomplete signatures. In Sect. 4,
we describe how to compute these efficiently. The incomplete signatures are then used
to preprocess the data structure such that we can quickly find the full signature once
we know the reachdoors of an elementary box. After building and preprocessing the
data structure, it is possible to determine dF ( P, Q) ≤ 1 efficiently by traversing the
freespace diagram elementary box by elementary box, as explained in Sect. 5.
3 Building a Lookup Table
3.1 Preprocessing an Elementary Box Before it considers the input, our algorithm builds a lookup table. As mentioned above, the purpose of this table is to speed up the computation of small parts of the freespace diagram.
Let τ ∈ N be a parameter.4 The elementary box is a subdivision of [0, τ ]2 into
τ columns and rows, thus τ 2 cells.5 For i, j = 0, . . . , τ − 1, we denote the cell
[i, i + 1] × [ j, j + 1] with D(i, j ). We denote the left side of the boundary ∂ D(i, j )
by l(i, j ) and the bottom side by b(i, j ). Note that l(i, j ) coincides with the right side
of ∂ D(i − 1, j ) and b(i, j ) with the top of ∂ D(i, j − 1). Thus, we write l(τ, j ) for the
right side of D(τ − 1, j ) and b(i, τ ) for the top side of D(i, τ − 1). Figure 2 shows
the elementary box.
The doororder σ rj for a row j is a permutation of {s0, t0, . . . , sτ , tτ }, having 2τ + 2
elements. For i = 1, . . . , τ , the element si represents the lower endpoint of the door
on l(i, j ), and ti represents the upper endpoint. The elements s0 and t0 are an
exception: they describe the reachdoor on the boundary l(0, j ) (i.e., its intersection with
reach( P, Q)). The doororder σ rj represents the combinatorial order of these
endpoints, as projected onto a vertical line, i.e., they are sorted into their vertical order.
Some doororders may encode the same combinatorial structure. In particular, when
door i is closed, the exact position of si and ti in a doororder is irrelevant, as long as
ti comes before si . For a closed door i (i > 0), we assign si to the upper endpoint of
l(i, j ) and ti to the lower endpoint. The values of s0 and t0 are defined by the
reachdoor and their relative order is thus a result of computation. We break ties between si
and ti by placing si before ti , and any other ties are resolved by index. A doororder
σic is defined analogously for a column i . We write x <ic y if x comes before y in σic,
and x <rj y if x comes before y in σ rj . An incomplete doororder is a doororder in
which s0 and t0 are omitted (i.e. the intersection of reach( P, Q) with the door is still
unknown); see Fig. 3.
4 A preview for the impatient reader: we later set τ = (√log n/ log log n).
5 For now, the elementary box is a combinatorial concept. In the next section, we overlay these boxes on
the freespace diagram to obtain “concrete” elementary boxes.
Fig. 2 The elementary box. The
cell D(2, 1) is shown white. Its
boundaries—l(2, 1), l(3, 1),
b(2, 1), b(2, 2)—are indicated
t3
s3
Fig. 3 The doororder of a row (the vertical order of the points) encodes the combinatorial structure of the
doors. The doororder for the row in the figure is s1s3s4t5t3t0s2t4s0s5t1t2. Note that s0 and t0 represent
the reachdoor, which is empty in this case. These are omitted in the incomplete doororder
We can now define the (full) signature of the elementary box as the aggregation of
the doororders of its rows and columns. Therefore, a signature = (σ1c, . . . , στc, σ1r ,
. . . , στr ) consists of 2τ doororders: one doororder σic for each column i and one
doororder σ rj for each row j of the elementary box. Similarly, an incomplete signature is
the aggregation of incomplete doororders.
For a given signature, we define the combinatorial reachability structure of the
elementary box as follows. For each column i and for each row j , the combinatorial
reachability structure indicates which door boundaries in the respective column or row
define the reachdoor of b(i, τ ) or l(τ, j ).
Lemma 3.1 Let be a signature for the elementary box. Then we can determine the
combinatorial reachability structure of the box in total time O(τ 2).
Proof We use dynamic programming, very similar to the algorithm by Alt and Godau
[8]. For each vertical edge l(i, j ) we define a variable l(i, j ), and for each horizontal
edge b(i, j ) we define a variable b(i, j ). The l(i, j ) are pairs of the form (su , tv),
representing the reachdoor reach( P, Q) ∩ l(i, j ). If this reachdoor is closed, then
tv <rj su holds. If the reachdoor is open, then it is bounded by the lower endpoint of
the door on l(u, j ) and by the upper endpoint of the door on l(v, j ). (Note that in this
case we have v = i .) Once again s0 and t0 are special and represent the reachdoor on
l(0, j ). The variables b(i, j ) are defined analogously.
Fig. 4 The three cases for the recursive definition of l(i, j ). If the lower boundary is reachable, we can
reach the whole right door (left). If neither the lower nor the left boundary is reachable, the right door is
not reachable either (middle). Otherwise, the lower boundary is the maximum of l(i − 1, j ) and the lower
boundary of the right door (right)
Now we can compute l(i, j ) and b(i, j ) recursively as follows: first, we set
l(0, j ) = b(i, 0) = (s0, t0), for i, j = 0, . . . , τ − 1.
Next, we describe how to find l(i, j ) given l(i − 1, j ) and b(i − 1, j ), see Fig. 4.
Case 1: Suppose b(i −1, j ) is open. This means that b(i −1, j ) intersects reach( P, Q),
so reach( P, Q) ∩ l(i, j ) is limited only by the door on l(i, j ), and we can set l(i, j ) :=
(si , ti ).
Case 2: If both b(i − 1, j ) and l(i − 1, j ) are closed, it is impossible to reach l(i, j )
and thus we set l(i, j ) := l(i − 1, j ).
Case 3: If b(i − 1, j ) is closed and l(i − 1, j ) is open, we may be able to reach l(i, j )
via l(i − 1, j ). Let su be the lower endpoint of l(i − 1, j ). We need to pass l(i, j )
above su and si and below ti , and therefore set l(i, j ) := (max(su , si ), ti ), where the
r
maximum is taken according to the order < j .
The recursion for the variable b(i, j ) is defined similarly. We can implement the
recursion in time O (τ 2) for any given signature, for example by traversing the
elementary box column by column, while processing each column from bottom to top.
There are at most ((2τ + 2)!)2τ = τ O(τ 2) distinct signatures for the elementary
box. We choose τ = λ√log n/ log log n for a sufficiently small constant λ > 0, so
that this number becomes o(n). Thus, during the preprocessing stage we have time
to enumerate all possible signatures and determine the corresponding combinatorial
reachability structure inside the elementary box. This information is then stored in an
appropriate data structure.
3.2 Building the Data Structure
Before we describe this data structure, we first explain how the doororders are
represented. This depends on the computational model. By our choice of τ , there are o(n)
distinct doororders. On the word RAM, we represent each doororder and incomplete
doororder by an integer between 1 and (2τ )!. This fits into a word of log n bits. On the
pointer machine, we create a record for each doororder and incomplete doororder;
we represent an order by a pointer to the corresponding record.
The data structure has two stages. In the first stage, we assume we know the
incomplete doororder for each row and for each column of the elementary box,6 and we
wish to determine the incomplete signature. In the second stage we have obtained the
reachdoors for the left and bottom sides of the elementary box, and we are looking
for the full signature. The details of our method depend on the computational model.
One way uses table lookup and requires the word RAM; the other way works on the
pointer machine, but is a bit more involved.
3.2.1 Word RAM
We organize the lookup table as a large tree T . In the first stage, each level of T
corresponds to a row or column of the elementary box. Thus, there are 2τ levels. Each
node has (2τ )! children, representing the possible incomplete doororders for the next
row or column. Since we represent doororders by positive integers, each node of T
may store an array for its children; we can choose the appropriate child for a given
incomplete doororder in constant time. Thus, determining the incomplete signature
for an elementary box requires O(τ ) steps on a word RAM.
For the second stage, we again use a tree structure. Now the tree has O(τ ) layers,
each with O(log τ ) levels. Again, each layer corresponds to a row or column of the
elementary box. The levels inside each layer then implement a balanced binary search
tree that allows us to locate the endpoints of the reachdoor within the incomplete
signature. Since there are 2τ endpoints, this requires O(log τ ) levels. Thus, it takes
O(τ log τ ) time to find the full signature of a given elementary box.
3.2.2 Pointer machine
Unlike in the word RAM model, we are not allowed to store a lookup table on every
level of the tree T , and there is no way to quickly find the appropriate child for a given
doororder. Instead, we must rely on batch processing to achieve a reasonable running
time.
Thus, suppose that during the first stage we want to find the incomplete signatures
for a set B of m elementary boxes, where again for each box in B we know the
incomplete doororder for each row and each column. Recall that we represent the
doororder by a pointer to the corresponding record. With each such record, we store
a queue of elementary boxes that is empty initially.
We now simultaneously propagate the boxes in B through T , proceeding level by
level. In the first level, all of B is assigned to the root of T . Then, we go through the
nodes of one level of T , from left to right. Let v be the current node of T . We consider
each elementary box b assigned to v. We determine the next incomplete doororder
6 In the next section, we describe how to determine the incomplete doororders efficiently.
for b, and we append b to the queue for this incomplete doororder—the queue is
addressed through the corresponding record, so all elementary boxes with the same
next incomplete doororder end up in the same queue. Next, we go through the nodes
of the next level, again from left to right. Let v be the current node. The node v
corresponds to a next incomplete doororder σ that extends the known signature of
its parents. We consider the queue stored at the record for σ . By construction, the
elementary boxes that should be assigned to v appear consecutively at the beginning
of this queue. We remove these boxes from the queue and assign them to v . After this,
all the queues are empty, and we can continue by propagating the boxes to the next
level. During this procedure, we traverse each node of T a constant number of times,
and in each level of the T we consider all the boxes in B. Since T has o(n) nodes, the
total running time is O(n + mτ ).
For the second stage, the data structure works just as in the word RAM case,
because no table lookup is necessary. Again, we need O(τ log τ ) steps to process one
box. After the second stage, we obtain the combinatorial reachability structure of the
box in constant time since we precomputed this information for each box (Lemma 3.1).
Thus, we have shown the following lemma, independently of the computational model.
Lemma 3.2 For τ = λ√log n/ log log n, with a sufficiently small constant λ > 0, we
can construct in o(n) time a data structure of size o(n) such that
• given a set of m elementary boxes where the incomplete doororders are known,
we can find the incomplete signature of each box in total time O(n + mτ );
• given the incomplete signature and the reachdoors on the bottom and left boundary
of an elementary box, we can find the full signature in O(τ log τ ) time;
• given the full signature of an elementary box, we can find the combinatorial
reachability structure of the box in constant time.
4 Preprocessing a Given Input
Next, we perform a second preprocessing phase that considers the input curves P
and Q. Our eventual goal is to compute the intersection of reach( P, Q) with the cell
boundaries, taking advantage of the data structure from Sect. 3. For this, we aggregate
the cells of FSD( P, Q) into (concrete) elementary boxes consisting of τ × τ cells.
There are n2/τ 2 such boxes. We may avoid rounding issues by either duplicating
vertices or handling a small part of FSD( P, Q) without lookup tables.
The goal is to determine the signature for each elementary box S. At this point, this
is not quite possible yet, since the signature depends on the intersection of reach( P, Q)
with the lower and left boundary of S. Nonetheless, we can find the incomplete
signature, in which the positions of s0, t0 (the reachdoor) in the (incomplete) doororders
σir , σ jc are still to be determined.
We aggregate the columns of FSD( P, Q) into vertical strips, each
corresponding to a single column of elementary boxes (i.e., τ consecutive columns of cells in
FSD( P, Q)). See Fig. 5.
Let A be such a strip. It corresponds to a subcurve P of P with τ edges. The
following lemma implies that we can build a data structure for A such that, given any
Fig. 5 FSD(P, Q) is subdivided into n2/τ 2 elementary boxes of size τ × τ . The freespace diagram
is subdivided into n2/τ 2 elementary boxes of size τ × τ . A strip is a column of elementary boxes: it
corresponds to a subcurve of P with τ edges
Fig. 6 By using the arrangement A defined by unit circles centered at vertices of P , we can determine
the incomplete doororder of each segment s on Q. This is done by locating the dual point of s in the dual
arrangement B. The dual arrangement also contains pseudolines to determine when s leaves a circle of A
Lemma 4.1 Given a subcurve P with τ edges, we can compute in O (τ 6) time a data
structure that requires O (τ 6) space and that allows us to determine the incomplete
doororder of any line segment on Q in time O (log τ ).
Proof Consider the arrangement A of unit circles whose centers are the vertices of
P (see Fig. 6). The incomplete doororder of a line segment s is determined by the
intersections of s with the arcs of A (and for a circle not intersecting s by whether
s lies inside or outside of the circle). Let s be the line spanned by line segment s.
Suppose we wiggle s . The order of intersections of s and the arcs of A changes only
when s moves over a vertex of A or if s leaves or enters a circle.
We use the standard duality transform that maps a line : y = ax + b to the
point ∗ : (a, −b), and vice versa. Consider a unit circle C in A with center (cx , cy ).
Elementary geometry shows that the set of all lines that are tangent to C from above
dualizes to the curve ta∗(C ) : y = cx x − cy − √1 + x 2. Similarly, the lines that are
tangent to C from below dualize to the curve tb∗(C ) : y = cx x − cy + √1 + x 2. Define
C ∗ := {ta∗(C ), tb∗(C )  C ∈ A}. Since any pair of distinct circles C1, C2 has at most
four common tangents, one for each choice of above/below C1 and above/below C2,
it follows that any two curves in C ∗ intersect at most once.
Let V be the set of vertices in A, and let V ∗ be the lines dual to the points in V
(note that V  = O (τ 2)). Since for any vertex v ∈ V and any circle C ∈ A there are
at most two tangents through v on C , each line in V ∗ intersects each curve in C ∗ at
most once. Thus, the arrangement B of the curves in V ∗ ∪ C ∗ is an arrangement of
pseudolines with complexity O (τ 4). Furthermore, it can be constructed in the same
expected time, together with a point location structure that finds the containing cell in
B of any given point in time O (log τ ) [60, Chap. 6.6.1].
Now consider a line segment s and the supporting line s . As observed in the
first paragraph, the combinatorial structure of the intersection between s and A is
completely determined by the cell of B that contains the dual point s∗. Thus, for every
cell f (s) ∈ B, we construct a list L f (s) that represents the combinatorial structure of
s ∩ A. There are O (τ 4) such lists, each having size O (τ ). We can compute L f (s) by
traversing the zone of s in A. Since circles intersect at most twice and since a line
intersects any circle at most twice, the zone has complexity O (τ 2α(τ )), where α(·)
denotes the inverse Ackermann function [60, Thm. 5.11]. Since O (τ 2α(τ )) ⊂ O (τ 2),
we can compute all lists in O (τ 6) time.
Given the list L f (s), the incomplete doororder of s is determined by the position
of the endpoints of s in L f (s). There are O (τ 2) possible ways for this, and we build
a table T f (s) that represents them. For each entry in T f (s), we store a representative
for the corresponding incomplete doororder. As described in the previous section,
the representative is a positive integer in the word RAM model and a pointer to the
appropriate record on a pointer machine.
The total size of the data structure is O (τ 6) and it can be constructed in the same
time. A query works as follows: given s, we can compute s∗ in constant time. Then
we use the point location structure of B to find f (s) in O (log τ ) time. Using binary
search on T f (s) (or an appropriate tree structure in the case of a point machine), we
can then determine the position of the endpoints of s in the list L f (s) in O (log τ ) time.
This bound holds both on the word RAM and on the pointer machine.
Lemma 4.2 Given the data structure of Lemma 3.2, the incomplete signature for each
elementary box can be determined in time O (nτ 5 + n2(log τ )/τ ).
Proof By building and using the data structure from Lemma 4.1, we determine the
incomplete doororder for each row in each vertical τ strip in total time proportional
to
We repeat the procedure with the horizontal strips. Now we know for each
elementary box in FSD( P, Q) the incomplete doororder for each row and each column. We
use the data structure of Lemma 3.2 to combine these. As there are n2/τ 2 boxes, the
number of steps is O(n2/τ + n) = O(n2/τ ). Hence, the incomplete signature for
each elementary box is found in O(nτ 5 + n2(log τ )/τ ) steps.
5 Solving the Decision Problem
With the data structures and preprocessing from the previous sections, we have all
ingredients in place to determine whether dF ( P, Q) ≤ 1. We know for each elementary
box its incomplete signature and we have a data structure to derive its full signature
(and with it, the combinatorial reachability structure) when its reachdoors are known.
What remains to be shown is that we can efficiently process the freespace diagram to
determine whether (n, n) ∈ reach( P, Q). This is captured in the following lemma.
Lemma 5.1 If the incomplete signature for each elementary box is known, we can
determine whether (n, n) ∈ reach( P, Q) in time O(n2(log τ )/τ ).
Proof We go through all elementary boxes of FSD( P, Q), processing them one
column at a time, going from bottom to top in each column. Initially, we know the full
signature for the box S in the lower left corner of FSD( P, Q). We use the signature
to determine the intersections of reach( P, Q) with the upper and right boundary of
S. There is a subtlety here: the signature gives us only the combinatorial reachability
structure, and we need to map the resulting si , t j back to the corresponding vertices on
the curves. On the word RAM, this can be done easily through table lookups. On the
pointer machine, we use representative records for the si , ti elements and use O(τ )
time before processing the box to store a pointer from each representative record to
the appropriate vertices on P and Q.
We proceed similarly for the other boxes. By the choice of the processing order
of the elementary boxes we always know the incoming reachdoors on the bottom
and left boundary when processing a box. Given the incoming reachdoors, we can
determine the full signature and find the structure of the outgoing reachdoors in
total time O(τ log τ ), using Lemma 3.2. Again, we need O(τ ) additional time on
the pointer machine to establish the mapping from the abstract si , ti elements to the
concrete vertices of P and Q. In total, we spend O(τ log τ ) time per box. Thus, it
takes time O(n2(log τ )/τ ) to process all boxes, as claimed.
As a result, we obtain the following theorem for the pointer machine (and, by
extension, for the real RAM model). For the word RAM model, we can obtain an even
faster algorithm (see Sect. 6).
Theorem 5.2 There is an algorithm that solves the decision version of the Fréchet
problem in O(n2(log log n)3/2/√log n) time on a pointer machine.
Proof Set τ = λ√log n/ log log n, for a sufficiently small constant λ > 0. The
theorem follows by applying Lemmas 3.2, 4.2, and 5.1 in sequence.
6 Improved Bound on Word RAM
We now explain how the running time of our algorithm can be improved if our
computational model allows for constant time tablelookup. We use the same τ as above (up
Fig. 7 A cluster consists of
τ × τ elementary boxes, thus of
τ 2 × τ 2 cells. A row R and its
corresponding R for the central
elementary box are indicated
to a constant factor). However, we change a number of things. “Signatures” are
represented differently and the data structure to obtain combinatorial reachability structures
is changed accordingly. Furthermore, we aggregate elementary boxes into clusters and
determine “incomplete doororders” for multiple boxes at the same time. Finally, we
walk the freespace diagram based on the clusters to decide dF ( P, Q) ≤ 1.
6.1 Clusters and Extended Signatures
We introduce a second level of aggregation in the freespace diagram (see Fig. 7): a
cluster is a collection of τ ×τ elementary boxes, that is, τ 2×τ 2 cells in FSD( P, Q). Let
R be a row of cells in FSD( P, Q) of a certain cluster. As before, the row R corresponds
to an edge e on Q and a subcurve P of P with τ 2 edges. We associate with R an
ordered set Z = e0, z0, z1, z1, z2, z2, . . . , zk , zk , e1 with 2 · k + 3 elements. Here k
is the number of intersections of e with the unit circles centered at the τ 2 vertices of P
(all but the very first). Hence, k is bounded by 2τ 2 and Z  is bounded by 4τ 2 + 3. The
order of Z indicates the order of these intersections with e directed along Q. Elements
e0 and e1 represent the endpoints of e and take a special role. In particular, these are
used to represent closed doors and snap open doors to the edge e. The elements zi
are placeholders for the positions of the endpoints of the reachdoors: z0 represents a
possible reachdoor endpoint between e0 and z1, the element z1 represents an endpoint
between z1 and z2, etc.
Consider a row R of an elementary box inside the row R of a cluster, corresponding
to an edge e of Q. The doorindex of R is an ordered set s0, t0, . . . , sτ , tτ of size
2τ + 2. Similar to a doororder, elements s0 and t0 represent the reachdoor at the
leftmost boundary of R ; the elements si and ti (1 ≤ i ≤ τ ) represent the door at the
right boundary of the i th cell in R . However, instead of rearranging the set to indicate
relative positions, the elements si and ti simply refer to elements in Z . If the door
is open, they refer to the corresponding intersections with e (possibly snapped to e0
or e1). If the door is closed, si is set to e1 and ti is set to e0. The elements s0 and t0
are special, representing the reachdoor, and they refer to one of the elements zi . An
incomplete doorindex is a doorindex without s0 and t0. The advantage of a
doorindex over a doororder is that the reachdoor is always at the start. Hence, completing
Fig. 8 (left) Every field represents the incomplete doorindex of a row in an elementary box. (center) The
fields are grouped into words per row in a cluster. (right) Transposition yields the desired organization,
where a word represents the incomplete doorindex of the rows in an elementary box
an incomplete doorindex to a full doorindex can be done in constant time. Since a
doorindex has size 2τ + 2, the number of possible doorindices for R is τ O(τ ).
We define the doorindices for the columns analogously. We concatenate the
doorindices for the rows and the columns to obtain the indexed signature for an elementary
box. Similarly, we define the incomplete indexed signature. The total number of
possible indexed signatures remains τ O(τ 2).
For each possible incomplete indexed signature we build a lookup table T as
follows: the input is a word with 4τ fields of O(log τ ) bits each. Each field stores the
positions in Z of the endpoints of the ingoing reachdoors for the elementary box: 2τ
fields for the left side, 2τ fields for the lower side. The output consists of a word that
represents the indices for the elements in Z that represent the outgoing reachdoors for
the upper and right boundary of the box. Thus, the input of T is a word of O(τ log τ )
bits, and T has size τ O(τ ). Hence, for all incomplete indexed signatures combined,
the size is τ O(τ 2) = o(n) by our choice of τ .
6.2 Preprocessing a Given Input
During the preprocessing for a given input P, Q, we use superstrips consisting of τ
strips. That is, a superstrip is a column of clusters and consists of τ 2 columns of the
freespace diagram. Lemma 4.1 still holds, albeit with a larger constant c in place of
6. The data structure gets as input a query edge e, and it returns in O(log τ ) time a
word that contains τ fields. Each field represents the incomplete doorindex for e in the
corresponding elementary box and thus consists of O(τ log τ ) bits. Hence, the word
size is O(τ 2 log τ ) = O(log n) by our choice of τ . Thus, the total time for building
a data structure for each superstrip and for processing all rows is O(n/τ 2 (τ c +
n log τ )) = O(n2(log τ )/τ 2). We now have parts of the incomplete indexed signature
for each elementary box packed into different words. To obtain the incomplete indexed
signature, we need to rearrange the information such that the incomplete doorindices
of the rows in one elementary box are in a single word. This corresponds to computing
a transpose of a matrix, as is illustrated in Fig. 8. For this, we need the following
lemma, which can be found—in slightly different form—in Thorup [62, Lem. 9].
Lemma 6.1 Let X be a sequence of τ words that contain τ fields each, so that X can
be interpreted as a τ × τ matrix. Then we can compute in time O(τ log τ ) on a word
RAM a sequence Y of τ words with τ fields each that represents the transpose of X .
Proof The algorithm is recursive and solves a more general problem: let X be a
sequence of a words that represents a sequence M of b different a × a matrices, such
that the i th word in X contains the fields of the i th row of each matrix in M from
left to right. Compute a sequence of words Y that represents the sequence M of the
transposed matrices in M .
The recursion works as follows: if a = 1, there is nothing to be done. Otherwise, we
split X into the sequence X1 of the first a/2 words and the sequence X2 of the remaining
words. X1 and X2 now represent a sequence of 2b (a/2) × (a/2) matrices, which
we transpose recursively. After the recursion, we put the (a/2) × (a/2) submatrices
back together in the obvious way. To finish, we need to transpose the offdiagonal
submatrices. This can be done simultaneously for all matrices in time O(a), by using
appropriate bitoperations (or table lookup). Hence, the running time obeys a recursion
of the form T (a) = 2T (a/2) + O(a), giving T (a) = O(a log a), as desired.
By applying the lemma to the words that represent τ consecutive rows in a
superstrip, we obtain the incomplete doorindices of the rows for each elementary box. This
takes total time proportional to
We repeat this procedure for the horizontal superstrips. By using an appropriate lookup
table to combine the incomplete doorindices of the rows and columns, we obtain the
incomplete indexed signature for each elementary box in total time O(n2(log τ )/τ 2).
6.3 The Actual Computation
We traverse the freespace diagram cluster by cluster (recall that a cluster consists
of τ × τ elementary boxes). The clusters are processed column by column from left
to right, and inside each column from bottom to top. Before processing a cluster,
we walk along the left and lower boundary of the cluster to determine the incoming
reachdoors. This is done by performing a binary search for each box on the boundary,
and determining the appropriate elements zi which correspond to the incoming
reachdoors. Using this information, we assemble the appropriate words that represent the
incoming information for each elementary box. Since there are n2/τ 4 clusters, this
step requires time O((n2/τ 4)τ 2 log τ ) = O(n2(log τ )/τ 2). We then process the
elementary boxes inside the cluster, in a similar fashion. Now, however, we can process
each elementary box in constant time through a single table lookup, so the total time
is O(n2/τ 2). Hence, the total running time of our algorithm is O(n2(log τ )/τ 2). By
our choice of τ = λ√log n/ log log n for a sufficiently small λ > 0, we obtain the
following theorem.
Theorem 6.2 The decision version of the Fréchet problem can be solved in
O(n2(log log n)2/ log n) time on a word RAM.
7 Computing the Fréchet Distance
The optimization version of the Fréchet problem, i.e., computing the Fréchet distance,
can be done in O(n2 log n) time using parametric search with the decision version
as a subroutine [8]. We showed that the decision problem can be solved in o(n2)
time. However, this does not directly yield a faster algorithm for the optimization
problem: if the running time of the decision problem is T (n), parametric search gives an
O((T (n)+n2) log n) time algorithm [8]. There is an alternative randomized algorithm
by Raichel and HarPeled [49]. Their algorithm also needs O((T (n) + n2) log n) time,
but below we adapt it to obtain the following lemma.
Lemma 7.1 The Fréchet distance of two polygonal curves with n vertices each can
be computed by a randomized algorithm in O(n22α(n) + T (n) log n) expected time,
where T (n) is the time for the decision problem.
Before we prove the lemma, we recall that possible values of the Fréchet distance
are limited to a certain set of critical values [8]:
1. the distance between a vertex of one curve and a vertex of the other curve (vertex–
vertex);
2. the distance between a vertex of one curve and an edge of the other curve (vertex–
edge); and
3. for two vertices of one curve and an edge of the other curve, the distance between
one of the vertices and the intersection of e with the bisector of the two vertices
(if this intersection exists) (vertex–vertex–edge).
If we also include vertex–vertex–edge tuples with no intersection, we can sample
a critical value uniformly at random in constant time. The algorithm now works as
follows (see HarPeled and Raichel [49] for more details): first, we sample a set S
of K = 4n2 critical values uniformly at random. Next, we find a , b ∈ S such that
the Fréchet distance lies between a and b and such that [a , b ] contains no other
value from S. In the original algorithm this is done by sorting S and performing a
binary search using the decision version. Using medianfinding instead, this step can
be done in O(K + T (n) log K ) time. Alternatively, the running time of this step could
be reduced by picking a smaller K . However, this does not improve the final bound,
since it is dominated by a O(n22α(n)) term. The interval [a , b ] with high probability
contains only a small number of the remaining critical values. More precisely, for
K = 4n2 the probability that [a , b ] has more than 2cn ln n critical values is at most
1/nc [49, Lem. 6.2].
The remainder of the algorithm proceeds as follows: first, we find all critical values
of type vertex–vertex and vertex–edge that lie inside the interval [a , b ]. This can be
done in O(n2) time by checking all vertex–vertex and vertex–edge pairs. Among these
values, we again use medianfinding to determine the interval [a, b] ⊆ [a , b ] that
contains the Fréchet distance in O(K + T (n) log K ) time. It remains to determine
the critical values corresponding to vertex–vertex–edge tuples that lie in [a, b].
For this, take an edge e of P and the vertices of Q. Conceptually, we start with
circles of radius a around the vertices of Q, and we increase the radii until b. During
this process, we observe the evolution of the intersection points between the circle
arcs and e. Because all vertex–vertex and vertex–edge events have been eliminated,
each circle intersects e in either 0 or 2 points, and this does not change throughout the
process. A critical value of vertex–vertex–edge type corresponds to the event that two
different circles intersect e in the same point, i.e., that two intersection points meet
while growing the circles. Two intersection points can meet at most once, and when
they do, they exchange their order along e.
This suggests the following algorithm: let Aa be the arrangement of circles with
radius a around the vertices of Q, and let Ab be the concentric arrangement of circles
with radius b. We determine the ordered sequence Ia of the intersection points of
the circles in Aa with e, and we number them in their order along e. Next, we find
the ordered sequence of intersection points Ib between e and the circles in Ab. We
assign to each point in Ib the number of the corresponding intersection points in Ia .
Since Ia  = Ib, this gives a permutation of {1, . . . , Ia }, Two intersection points
change their order from Ia to Ib exactly if there is a vertex–vertex–edge event in [a, b],
so these events correspond to the inversions of the resulting permutation. Given that
there are k such inversions, we can find them in time O(Ia  + k) using insertion sort.
Thus, the overall running time to find the critical events in [a, b], ignoring the time
for computing Ia and Ib, is O(n2 + K ).
It remains to show that we can quickly find Ia and Ib. We describe the algorithm for
Ia . First, compute the arrangement Aa of circles with radius a around the vertices of
Q. This takes O(n2) time [32]. To find the intersection order, traverse in Aa the zone
of the line spanned by e. The time for the traversal is bounded by the complexity of
the zone. Since the circles pairwise intersect at most twice and intersects each circle
only twice, the complexity of the zone is O(n2α(n)) [60, Thm. 5.11]. Summing over
all edges e, this adds a total of O(n22α(n)) to the running time. To find Ib, we proceed
similarly with Ab. Thus the overall time is O(T (n) log(n) + n22α(n) + K ). The event
K > 8n ln n has probability less than 1/n4, and we always have K = O(n3). Thus,
this case adds o(1) to the expected running time. Given K ≤ 8n ln n, the running
time is O(n log n). Lemma 7.1 follows. Theorem 7.2 now results from Lemma 7.1,
Theorem 5.2, and Theorem 6.2.
Theorem 7.2 The Fréchet distance of two polygonal curves with n edges each can be
computed by a randomized algorithm in time O(n2√log n(log log n)3/2) on a pointer
machine and in time O(n2(log log n)2) on a word RAM.
8 Discrete Fréchet Distance on the Pointer Machine
As mentioned in the introduction, Agarwal et al. [1] give a subquadratic algorithm
for finding the discrete Fréchet distance between two point sequences, using the word
RAM. In this section, we explain how their algorithm for the decision version of
the problem can be adapted to the pointer machine. This shows that, at least for the
decision version, the speedup does not come from bitmanipulation tricks but from
a deeper understanding of the underlying geometric structure. Our presentation is
slightly different from Agarwal et al. [1], in order to allow for a clearer comparison
with our continuous algorithm.
We recall the problem definition: we are given two sequences P = p1, p2, . . . , pn
and Q = q1, q2, . . . , qn of n points in the plane. For δ > 0, we define a directed
graph Gδ with vertex set P × Q. In Gδ, there is an edge between two vertices ( pi , q j ),
( pi , q j+1) if and only if both d( pi , q j ) ≤ δ and d( pi , q j+1) ≤ δ. The condition is
similar for an edge between vertices ( pi , q j ) and ( pi+1, q j ), and vertices ( pi , q j ) and
Fig. 9 The discrete Fréchet distance: (left) two point sequences P (disks) and Q (crosses) with 5 points
each; (middle) the associated freespace matrix F (white = 1, gray = 0); (right) the resulting reachability
matrix M. Since M55 = 1, the discrete Fréchet distance is at most 1
( pi+1, q j+1). There are no further edges in Gδ. The discrete Fréchet distance between
P and Q is the smallest δ for which Gδ has a path from ( p1, q1) to ( pn, qn). In the
decision version of the problem, we are given δ > 0, and we need to decide whether
there is a path from ( p1, q1) to ( pn, qn) in Gδ.
We now describe a subquadratic pointer machine algorithm for the decision version.
Thus, let point sequences P, Q be given, and suppose without loss of generality that
δ = 1. The discrete analogue of the freespace diagram is an n × n Boolean matrix
F where Fi j = 1, if d( pi , q j ) ≤ 1, and Fi j = 0, otherwise, for i, j = 1, . . . , n. We
call F the freespace matrix. Similarly, the discrete analogue of the reachable region
is an n × n Boolean matrix M that is defined recursively as follows: M11 = F11,
and for i, j = 1, . . . , n, (i, j ) = 1, we have Mi j = 1 if and only if Fi j = 1 and
at least one of Mi−1, j , Mi, j−1 or Mi−1, j−1 equals 1 (we set Mi,0 = M0, j = 0, for
i, j = 1, . . . n). Then the discrete Fréchet distance between P and Q is at most 1 if
and only if Mnn = 1. We call M the reachability matrix; see Fig. 9.
Adapting the method of Agarwal et al. [1], we show how to use preprocessing and
table lookup in order to decide whether Mnn = 1 in o(n2) steps on a pointer machine.
Let τ = λ log n, for a suitable constant λ > 0. We subdivide the rows of M into
k = O(n/τ ) strips, each consisting of τ consecutive rows: the first strip L1 consists
of rows 1, . . . , τ , the second strip L2 consists of rows τ + 1, . . . , 2τ , and so on. Each
strip Li , i = 1, . . . , k, corresponds to a contiguous subsequence Pi of τ points on
P. Let Ai be the arrangement of disks obtained by drawing a unit disk around each
vertex in Pi . The arrangement Ai has O(τ 2) faces.
Next, let ρ = λ log n/ log log n, with λ > 0 as above. We subdivide each strip Li ,
i = 1, . . . , k into l = O(n/ρ) elementary boxes, each consisting of ρ consecutive
columns in Li . We label the elementary boxes as Bi j , for i = 1, . . . , k and j = 1, . . . , l.
As above, an elementary box Bi j has corresponding contiguous subsequences Pi of
τ vertices on P and Q j of ρ vertices on Q. Now, the incomplete signature of an
elementary box Bi j consists of (i) the index i of the strip that contains it; and (ii)
the sequence f1, f2, . . . , fρ of faces in the disk arrangement Ai that contain the ρ
vertices of Q j , in that order. The full signature of an elementary box Bi j consists of
its incomplete signature plus a sequence of ρ + τ bits, that represent the entries in the
reach matrix M directly above and to the left of Bi j . We call these bits the reach bits.
2
2
3 4
1 3
Fig. 10 (left) We subdivide the reachability matrix M into strips of τ rows. Each strip is subdivided into
elementary boxes of ρ columns. Each elementary box corresponds to a subsequence of length τ on P and
a subsequence of length ρ on Q. (right) The incomplete signature σ of an elementary box Bi j consists of
the index of the containing strip and the face sequence for Q j in Ai . The full signature σ f additionally
contains ρ + τ reach bits that represent the bits in M directly above and to the left of Bi j
As in the continuous case, the information in the full signature suffices to determine
how the reachability information propagates through the elementary box; see Fig. 10.
The preprocessing phase proceeds as follows: first, we enumerate all possible
incomplete signatures. For this, we need to compute all strips Li and the corresponding
disk arrangements Ai , for i = 1, . . . , τ . Furthermore, we also compute a suitable point
location structure for each Ai . Since there are O (n/τ ) strips, each of which consists
of τ rows, this takes time O ((n/τ ) · τ 2 log τ ) = O (nτ log τ ) = O (n log n log log n).
For each strip Li , since Ai has O (τ 2) faces, the number of possible face sequences
f1, . . . , fρ is τ O(ρ) ≤ n1/3, by our choice of τ and ρ and for λ small enough. Thus,
there are O (n4/3/ log n) incomplete signatures, and they can be enumerated in the
same time. Now, for each incomplete signature σ = (i, f1, . . . , fρ ) we build a
lookuptable that encodes for each possible setting of the reach bits the resulting reach
bits at the bottom and the right boundary of the elementary box. There are 2ρ+τ ≤ n1/3
possible settings of the reach bits, by our choice of τ and ρ and for λ small enough.
We enumerate all of them and organize them as a complete binary tree of depth ρ + τ .
For each setting of the reach bits, we use the information of the incomplete signature
to determine the result through a straightforward dynamic programming algorithm
[1, 18, 39] in O (τ · ρ) = O (log2 n) time, and we store the result as a linked list of
length ρ + τ − 1 at the leaf for the corresponding reach bits; see Fig. 11. Thus, the
total time for this part of the preprocessing phase is O (n5/3 log n).
Next, we determine for each elementary box Bi j its incomplete signature. For this,
we use the point location structure for Ai to determine for each vertex in Q j the face
of Ai that contains it. There are O (n2/τρ) elementary boxes, each Q j has ρ vertices,
and one point location query takes O (log τ ) time, so the total time for this step is
O ((n2/τρ) · ρ · log τ ) = O ((n2/ log n) log log n). Using this information, we can
store with each elementary box a pointer to the lookup table for the corresponding
incomplete signature.
1
0
1
0
1
0
Fig. 11 For each incomplete signature, we create a lookup table organized as a complete binary tree. Each
leaf corresponds to a setting of the reach bits for the elementary box. In the leaves, we store a linked list of
length ρ + τ − 1 that represents the contents of the reachability matrix at the bottom and at the right of the
elementary box
Finally, we can now use the lookup tables to propagate the reachability information
through M , one elementary box at a time, as in the continuous case. The time to
process one elementary box is O (ρ +τ ), because we need to traverse the corresponding
lookup table to find the reach bits for the adjacent boxes. Thus, the total running time
is O ((n2/τρ) · (ρ + τ )) = O (n2/ρ) = O ((n2/ log n) log log n). Thus, we get the
following pointer machine version of the result by Agarwal et al. [1]
Theorem 8.1 There is an algorithm that solves the decision version of the discrete
Fréchet problem in O ((n2/ log n) log log n) time on a pointer machine.
Remark Agarwal et al. [1] further describe how to get a faster algorithm for the
decision version by aggregating the elementary boxes into larger clusters, similar to
the method given in Sect. 6. This improved algorithm finally leads to a subquadratic
algorithm for computing the discrete Fréchet distance. Unfortunately, as in Sect. 6, it
seems that this improvement crucially relies on constant time table lookup, so it does
not directly translate to the pointer machine.
9 Decision Trees
Our results also have implications for the decisiontree complexity of the Fréchet
problem. Since in that model we account only for comparisons between input elements,
the preprocessing comes for free, and hence the size of the elementary boxes can be
increased. Before we consider the continuous Fréchet problem, we first note that a
similar result can be obtained easily for the discrete Fréchet problem.
Proof First, we consider the decision version: we are given two sequences P = p1, . . . ,
pn and Q = q1, . . . , qn of n points in the plane, and we would like to decide whether
the discrete Fréchet distance between P and Q is at most 1. Katz and Sharir [53] showed
that we can compute a representation of the set of pairs ( pi , q j ) with pi − q j ≤ 1 in
O(n4/3) steps. This information suffices to complete the reachability matrix without
further comparisons. As shown by Agarwal et al. [1], one can then solve the
optimization problem at the cost of another O(log n)factor, which is absorbed into the
Onotation.
Given our results above, we prove an analogous statement for the continuous Fréchet
distance.
Theorem 9.2 There exists an algebraic decision tree for the Fréchet problem (decision
version) of depth O(n2−ε), for a fixed constant ε > 0.
Proof We reconsider the steps of our algorithm. The only phases that actually involve
the input are the second preprocessing phase and the traversal of the elementary boxes.
The reason of our choice for τ was to keep the time for the first preprocessing phase
small. This is no longer a problem. By Lemmas 4.2 and 5.1, the remaining cost is
bounded by O(nτ 5 + n2(log τ )/τ ). Choosing τ = n1/6, we get a decision tree of
depth n · n5/6 + n2−1/6 log n. This is O(n2−(1/6) log n) = O(n2−ε), for any fixed
0 < ε < 1/6.
10 Weak Fréchet Distance
The weak Fréchet distance is a variant of the Fréchet distance where we are allowed to
walk backwards along the curves [8]. More precisely, let P and Q be two polygonal
curves, each with n edges, and let be the set of all continuous functions α : [0, 1] →
[0, n] with α(0) = 0 and α(1) = n. The weak Fréchet distance between P and Q is
defined as
Compared to the regular Fréchet distance, the set now also contains nonmonotone
functions. The weak Fréchet distance was also introduced by Alt and Godau [8],
who showed how to compute it in O(n2 log n) worstcase time. We will now use our
framework to obtain an algorithm that runs in o(n2) expected time on a word RAM.
10.1 A Decision Algorithm for the Pointer Machine
As usual, we start with the decision version: given two polygonal curves P and Q,
each with n edges, decide whether dwF( P, Q) ≤ 1. This has an easy interpretation in
Fig. 12 The polygonal curves P and Q have weak Fréchet distance at most 1, but Fréchet distance larger
than 1: the point (n, n) is not in reach(P, Q), but the vertices for the cells C(0, 0) and C(4, 4) are in the
same connected component of G. The reachable region is shown dark blue, the edges if G are shown light
blue
terms of the freespace diagram. Define an undirected graph G = (V , E ) with vertex
set V = {0, 1, . . . , n − 1}2. The vertex (i, j ) ∈ V corresponds to the cell C (i, j ), and
there is an edge between two vertices (i, j ) and (i , j ) if and only if the two cells
C (i, j ) and C (i , j ) are neighboring (i.e., if i − i  +  j − j  = 1) and the door
between them is open. Then, dwF( P, Q) ≤ 1 if and only if (i)  P(0) − Q(0) ≤ 1;
(ii)  P(n) − Q(n) ≤ 1; and (iii) the vertices (0, 0) and (n − 1, n − 1) are in the same
connected component of G, see Fig. 12.
Let τ, ρ ∈ N be parameters, to be determined later. We subdivide the cells into
k = O(n/τ ) vertical strips L1, . . . , Lk , each consisting of τ consecutive columns.
Each strip Li corresponds to a subcurve Pi of P with τ edges. For each such subcurve
Pi , we define two arrangements Ai and Bi . To obtain Ai , we take for each edge e of Pi
the “stadium” ce of points with distance exactly 1 from e, and we compute the resulting
arrangement. Since two distinct curves ce, ce cross in O(1) points, the complexity of
Ai is O(τ 2), see Fig. 13. The arrangement Bi is the arrangement B described in the
proof of Lemma 4.1, i.e., the arrangement of the curves dual to the tangent lines for
the unit circles around the vertices of Pi .
Next, we subdivide each strip into = O(n/ρ) elementary boxes, each consisting of
ρ consecutive rows. We label the elementary boxes as Bi j , with 1 ≤ i ≤ k, 1 ≤ j ≤ .
The rows of an elementary box Bi j correspond to a subcurve Q j of Q with ρ edges.
The signature of Bi j consists of (i) the index i of the corresponding strip; (ii) for each
vertex of Q j the face of Ai that contains it; and (ii) for each edge e of Q j the face
of Bi that contains the point that is dual to the supporting line of e, plus two indices
a, b ∈ {1, . . . , τ } that indicate the first and the last unit circle around a vertex of Pi
that e intersects, as we walk from one endpoint to another.
Given an elementary box B, the connection graph G B of B has τρ vertices, one
for each cell in B, and an edge between two cells C , C of B if and only if C and C
share a (horizontal or vertical) edge with an open door. The connectivity list of B is a
linked list with 2τ + 2ρ − 4 entries that stores for each cell C on the boundary of B a
pointer to a record that represents the connected component of the connection graph
G B that contains C .
Lemma 10.1 There are O(nτ 8ρ+1) different signatures. The connection graph of an
elementary box Bi j depends only on its signature, and the connectivity list can be
computed in O(τρ) time on a pointer machine, given the signature.
Proof First, we count the signatures. There are O(n/τ ) strips. Once the strip index i
is fixed, a signature consists of ρ + 1 faces of Ai , ρ faces of Bi , and ρ pairs of indices
a, b ∈ {1, . . . , τ }. The arrangement Ai has O(τ 2) faces, and the arrangement Bi has
O(τ 4) faces, as explained in the proof of Lemma 4.1. Finally, there are τ 2 pairs of
indices. Thus, the number of possible signatures in one strip is τ 8ρ+2. In total, we get
O(nτ 8ρ+1) signatures.
Next, the connection graph of an elementary box Bi j is determined solely by which
doors are open and which doors are closed. We explain how to deduce this information
from the signature. A horizontal edge of Bi j corresponds to an edge e of Pi and a
vertex q of Q j . The door is open if and only if q has distance at most 1 from e.
This is determined by the face of Ai containing q. Similarly, a vertical edge of Bi j
corresponds to a vertex p of Pi and an edge e of Q j . The door is open if and only if
e intersects the unit circle with center p. As in Lemma 4.1, this can be inferred from
the face of Bi that contains the point dual to the supporting line of e, together with the
indices (a, b) of the first and last circle intersected by e.
Finally, given the signature, we can build the connection graph G Bi j in O(τρ) time,
assuming that the arrangements Ai and Bi provide suitable data structures. With G Bi j
at hand, the connection list can be found in O(τρ) steps, using breadth first search.
As usual, our strategy now is to preprocess all possible signatures and to determine
the signature of each elementary box. Using this information, we can then process
the elementary boxes quickly in our main algorithm. The next lemma describes the
preprocessing steps.
Lemma 10.2 We can determine for each elementary box Bi j a pointer to its
connectivity list in total time O(nτ 8ρ+2ρ + (n2/τ ) log τ ) on a pointer machine.
Proof First, we compute the arrangements Ai and Bi for each vertical strip. By
Lemma 4.1, this takes O(τ 6) steps per strip, for a total of O(τ 6 · n/τ ) = O(nτ 5).
Then, we enumerate all possible signatures, and we compute the connectivity list for
each of them. By Lemma 10.1, this needs O(nτ 8ρ+2ρ) time. We store the signatures
and their connectivity lists
Next, we determine the signature for each elementary box. Fix a strip Li . For each
vertex q and each edge e of Q, we determine the containing faces of Ai and Bi and the
pair (a, b) that represents the first and least intersection of e with the unit circles around
the vertices of Pi . This takes O(n log τ ) steps in total, using appropriate point location
structures for the arrangements. Now, we use the data structure from the preprocessing
to connect each elementary box to its connectivity list. A simple pointerbased structure
supports one lookup in O(ρ log τ ) time. Since there are O(n/ρ) elementary boxes in
one strip, we get a total running time of O(n log τ ) per strip. Since there are O(n/τ )
strips, the resulting running time is O((n2/τ ) log τ ).
With the information from the preprocessing phase, we can easily solve the decision
problem with a unionfind data structure.
Proof Let S be the set that contains the boundary cells of all elementary boxes. Then,
S = O n2(ττρ+ρ) . We create a UnionFind data structure for S. Then, we go through
all elementary boxes Bi j , and for each Bi j , we use the connectivity list to connect
those subsets of boundary cells that are connected inside Bi j . This can be done with
The following theorem summarizes our algorithm for the decision problem.
Theorem 10.4 Let P and Q be two polygonal curves with n edges each. We can decide
whether dwF ( p, q) ≤ 1 in O((n2/ log n)α(n) log log n) time on a pointer machine,
where α(·) denotes the inverse Ackermann function.
Proof We set τ = log n and ρ = λ log n/ log log n, for a suitable constant λ > 0.
If λ is small enough, then τ 8ρ+2 = O(n4/3), and the algorithm from Lemma 10.2
runs in time O((n2/ log n) log log n). After the preprocessing is finished, we can use
Lemma 10.3 to obtain the final result in time O((n2/ log n)α(n) log log n).
10.2 A Faster Decision Algorithm on the Word RAM
Theorem 10.4 is too weak to obtain a subquadratic algorithm for the weak Fréchet
distance. This requires the full power of the word RAM.
Theorem 10.5 Let P and Q be two polygonal curves with n edges each. We can
decide whether dwF ( p, q) ≤ 1 in O((n2/ log2 n)(log log n)5) time on a word RAM.
Proof We set τ = λ log2 n/ log log n and ρ = λ log n/ log log n. If λ is small enough,
then τ 8ρ+2 = O(n4/3), and we can perform the algorithm from Lemma 10.2 in time
O((n2/ log2 n) log log n).
We modify the algorithm from Lemma 10.2 slightly. Instead of a pointerbased
connectivity list, we compute a packed connectivity list. It consists of O(log n) words
of log n bits each. Each word stores (log n/ log log n) entries of O(log log n) bits.
An entry consists of O(1) fields, each with O(log log n) bits. As in Lemma 10.2,
the entries in the packed connectivity list represent the connected components of the
connection graph G B for the boundary cells of a given elementary box B. This is done
as follows: each entry of the connectivity list corresponds to a boundary cell C of B.
In the first field, we store a unique identifier from {1, . . . , 2τ + 2ρ − 4} that identifies
C . In the second field, we store the smallest identifier of any boundary cell of B that
lies in the same connected component as C . The remaining fields of the entry are
initialized to 0. Furthermore, we compute for each elementary box B a sequence of
O((τ / log n) log log n) words that indicates for each boundary cell of B whether the
corresponding door to the neighboring elementary box is open. This information can
be obtained in the same time by adapting the algorithm from Lemma 10.2.
Next, we group the elementary boxes into clusters. A cluster consists of log n
vertically adjacent elementary boxes from a single strip. The first set of clusters
come from bottommost log n elementary boxes, the second set of clusters from
following log n elementary boxes, etc. The boundary and the connectivity list of a
cluster are defined analogously as for an elementary box. Below, in Lemma 10.6,
we show that we can compute the (pointerbased) connectivity list of a cluster in time
O((τ + ρ)(log log n)4). Then, the lemma follows: there are O(n2/(τρ log n)) clusters,
so the total time to find connectivity lists for all clusters is O(n2/(ρ log n)(log log n)4+
n2/(τ log n)(log log n)4) = O((n2/ log2 n)(log log n)5). After that, we can solve the
decision problem in time
= O
as in Lemma 10.3.
Lemma 10.6 Given a cluster, we can compute the connectivity list for its boundary
in total time O((τ + ρ)(log log n)4) on a word RAM.
Proof Our algorithm makes extensive use of the word RAM capabilities. Data is
processed in packed form. That is, we usually deal with a sequence of packed words,
each with O(log n/ log log n) entries. Each entry has O(1) fields of O(log log n) bits.
We restrict the number of entries in each word so that the total number of bits used
is at most, say, (1/3) log n. Then, we can support any reasonable binary operation on
two packed words in O(1) time after o(n) preprocessing time. Indeed, we only need
to build a lookuptable during the preprocessing phase. First, we observe that we can
sort packed data efficiently.
Claim 10.7 Suppose we are given a sequence of μ packed words, representing a
sequence of O(μ(log n/ log log n)) entries. We can obtain a sequence of μ packed
words that represents the same entries, sorted according to any given field, in
O(μ log μ) time.
Proof We adapt usual techniques for packed integer sorting [5]. We precompute a
mergeoperation that receives two packed words, each with their entries sorted
according to a given field, and returns two packed words that represent the sorted sequence
of all entries. We also precompute an operation that receives one packed word and
returns a packed word with the same entries, sorted according to a given field. Then,
we can perform a merge sort in O(μ log μ) time, because in each level of the recursion,
the total time for merging is O(μ). Once we are down to a single word, we can sort
its entries in one step. More details can be found in the literature on packed sorting
[5,27].
Now we describe the strategy of the main algorithm. As explained above, a
cluster consists of log n vertically adjacent elementary boxes. We number the boxes
B1, . . . , Blog n, from bottom to top. For i = 1, . . . , log n, we denote by Gi the
connection graph for the boxes B1, . . . , Bi . That is, Gi is obtained by taking the union
ij=1 G B j of the individual connection graphs and adding edges for adjacent cells in
neighboring boxes that share an open door. Our algorithm proceeds in two phases.
In the first phase, we propagate the connectivity information upwards. That is, for
i = 1, . . . , log n, we compute a sequence Wi of O((τ + ρ) log log n/ log n) packed
words. Each entry in Wi corresponds to a cell C on the boundary of Bi . The first field
stores a unique identifier for the cell, and the second entry stores an identifier of the
connected component of C in Gi . In the second phase, we go in the reverse
direction. For i = log n, . . . , 1, we compute a sequence Wi of O((τ + ρ) log log n/ log n)
packed words that store for each cell C on the boundary of Bi an identifier of the
connected component of C in Glog n. Once this information is available, the
(pointerbased) connectivity list for the cluster boundary can be extracted in O(τ + ρ) time,
using appropriate precomputed operations on packed words, see Fig. 14 .
We begin with the first phase. The upper boundary of an elementary box are the
τ cells in the topmost row of the box, the lower boundary are the τ cells in the
bottommost row. From the preprocessing phase, we have a pointer to the packed
connectivity lists for all elementary boxes B1, . . . , Blog n. We make local copies of
Fig. 14 The algorithm for computing the connected components within a cluster. In the beginning, we know
the components inside each elementary box (left). In the first phase (top right), we update the connectivity
information from bottom to top, by considering the upper and lower boundaries of the vertically adjacent
elementary boxes. In the second phase (bottom right), we propagate the connected components from top to
bottom
these packed connectivity lists, and we modify them so that the identifiers of all cells
and connected components are unique. For example, we can increase all identifiers in
the connectivity list for Bi by i log3 n. Using lookuptables for these wordoperations,
this step can be carried out in O((τ + ρ) log log n) time for the whole cluster.
The sequence W1 for B1 is exactly the packed connectivity list of G1. Now, suppose
that we have computed for Bi a sequence Wi of packed words that represent the
connected component in Gi for each boundary cell of Bi . We need to compute a similar
sequence Wi+1 for Bi+1. For this, we need to determine how the doors between Bi
and Bi+1 affect the connected components in G Bi+1 .
Let H = (VH , E H ) be the graph that has one vertex for each connected component
of Gi that contains a cell at the upper boundary of Bi and one vertex for each connected
component of G Bi+1 that contains a cell at the lower boundary of Bi+1. The vertices
of H are labeled by the identifiers of their corresponding components. Suppose that
v ∈ VH represents a component Dv in Gi and w ∈ VH represents a component Dw
in G Bi+1 . There is an edge between v and w in E H if and only if Dv contains a cell
on the upper boundary of Bi and Dw contains a cell on the lower boundary of Bi+1
such that these two cells share an open door. The graph H has O(τ ) vertices and O(τ )
edges. Our goal is to obtain the connected components of H . From this, we will then
Fig. 15 The algorithm of Hirschberg, Chandra, and Sarwate. In each step, we either connect all nodes to
their largest neighbor or to their smallest neighbor. Then, we contract the resulting rooted trees and label
the contracted nodes with the index of the smallest node
derive the next sequence Wi+1. To obtain the components of H , we adapt the parallel
connectivity algorithm of Hirschberg, Chandra, and Sarwate [50].
Proof First, we obtain a list of O((τ / log n) log log n) packed words that
represent the vertices of H , as follows: we extract from Wi and from the connectivity
list of Bi+1 the entries for the upper boundary of Bi and for the lower
boundary of Bi+1. If all these entries are stored together in the respective packed
lists, this takes O((τ / log n) log log n) time to copy the relevant data words.
Then, we sort these packed sequences according to the component identifier, in
time O(τ / log n) log log2 n) (Claim 10.7). Finally, we go over the sorted packed
sequence to extract a sorted lists of the distinct component identifiers. This takes
O((τ / log n) log log n) time, using an appropriate operation on packed words.
The edges of E H can also be represented by a sequence of O((τ / log n) log log n)
words. As explained above, we have from the preprocessing phase for each elementary
box a sequence of O((τ / log n) log log n) words whose entries indicate the open doors
along its upper boundary. From this, we can obtain in O((τ / log n) log log n) time a
sequence of O((τ / log n)) log log n) words whose entries represent the edges of H ,
where each entry stores the component identifiers of the two endpoints of the edge.
We are now ready to implement the algorithm of Hirschberg, Chandra, and Sarwate.
The main steps of the algorithm are as follows, see Fig. 15.
Step 1: Find for each vertex of H the neighbor with the smallest and with the
largest identifier. This can be done in O((τ / log n) log log2 n) time by sorting the
edge lists twice, once in lexicographic order and once in reverse lexicographic
order of the identifiers of the endpoints. From these sorted lists, we can extract the
desired information in the claimed time, using appropriate word operations.
Step 2: Let VH be the vertices of H with at least one neighbor, and let VH ⊆ VH be
the vertices v ∈ VH having a neighbor with a smaller identifier than the identifier
of v. If VH  ≥ VH /2, we set the successor of each v ∈ VH to the neighbor
of v with the smallest identifier. Otherwise, at least half of the nodes in VH have
a neighbor with a larger identifier. In this case, we let VH be the set of these
nodes, and we set the successor of each v ∈ VH to the neighbor of v with the
largest identifier. The successor relation defines a directed forest F on VH such
that at least half of the vertices in VH are not a root in F . Given the information
available from Step 1 and appropriate word operations, this step can carried out in
O((τ / log n) log log n) time.
Step 3: Use pointer jumping to determine for each vertex v ∈ VH the identifier of
the root of the tree in F that contains v. For this, we set the successor of each v ∈ VH
that does not yet have a successor to v itself. Then, for log VH  = O(log log n)
rounds, we set simultaneously for each v ∈ VH the new successor of v to the old
successor of the old successor of v (pointer jumping). Each step at least halves
the distance of v to its root in F , so it takes O(log VH ) rounds until each vertex
in VH has found the root of its tree in F . Each round can be implemented in
O((τ / log n) log log2 n) time by sorting the vertices according to their successors.
Thus, this step takes O((τ / log n) log log3 n) time in total.
Step 4: Contract each tree of F into a single vertex whose identifier is the smallest
identifier in the tree.. Maintain for the original vertices of H a list that gives the
identifier of the contracted node that represents it. Again, this step can be carried
out in O(τ / log n) log log2 n) time using sorting.
After Steps 1–4, the number of nonsingleton components in H has at least
halved. Thus, by repeating the steps O(log VH ) = O(log log n) times, we can
identify the connected components of H . The total time of the algorithm is
O((τ / log n) log log4 n), as claimed.
Given the connected components of H , we can find the desired sequence Wi+1
for the boundary Bi+1. Indeed, the procedure from Claim 10.8 outputs a sequence
of O((τ / log n) log log n) packed words that gives for each vertex in H an identifier
of the component in H that contains it. We can use this list as a lookup table to
update the identifiers of the components in the connectivity list of Bi+1. This takes
O((τ + ρ)/ log n) log log2 n) time, using sorting.
In summary, since we consider log n elementary boxes, the total time for the first
phase if O((τ + ρ) log log4 n). The second phase is much easier. For i = logn, . . . , 2,
we propagate the connectivity information from Bi+1 to Bi . For this, we need to update
the identifiers of the connected components for the cells on the upper boundary of Bi
using the identifiers of the connected components on the lower boundary of Bi+1,
and then adjust the connectivity list Bi+1 with these new indices. Again, this takes
O((τ + ρ)/ log n) log log2 n) time, using sorting.
To actually compute the weak Fréchet distance, we use a simplified version of the
procedure from Sect. 7. In particular, for the weak Fréchet distance, there are only
critical values of the type vertex–vertex and vertex–edge, i.e., there are only O(n2)
critical values. However, we aim for a subquadratic running time, so we need to perform
the sampling procedure in a slightly different way.
Theorem 10.9 Suppose we can answer the decision problem for the weak Fréchet
distance in time T (n), for input curves P and Q with n edges each. Then, we can
compute the weak Fréchet distance of P and Q in expected time O(n3/2 logc n +
T (n) log n), for some fixed constant c > 0.
Proof First, we sample a set S of K = 6n1/2 critical values uniformly at random.
Then, we find a, b such that the weak Fréchet distance lies between a and b and such
that the interval [a, b] contains no other element from S. This takes O(K + T (n) log n)
time, using median finding.
Similarly to HarPeled and Raichel [49, Lem. 6.2], we see that the probability, that
the interval [a, b] has more than 2γ n3/2 ln n critical values is at most 1/nγ . Indeed,
there are at most n2 vertex–vertex and at most 2n2 vertex–edge events. Thus, the total
number of critical values is at most L ≤ 3n2. Let U + be the next γ n3/2 ln n larger
critical values after dwF( P, Q). Then, the probability that S contains no value from
U + is
L−U +
K
L
K
K −1
i=0
≤ 1 −
≤ exp
Analogously, the probability that S contains none of the next γ n3/2 ln n smaller critical
values is also at most 1/2nγ , so the claim follows. For γ > 0 large enough, the
contribution of this event to the expected running time is negligible.
Next, we find all critical values in the interval [a, b]. For this, we must determine
all vertex–vertex and vertex–edge pairs with distance in [a, b]. We report for every
vertex p of P or Q the set of vertices of the other curve that lie in the annulus with
radii a and b around p. Furthermore, we report for every edge e of P or Q the set
of vertices of the other curve that lie in the stadium with radii a and b. This can
be done efficiently with a rangesearching structure for semialgebraic sets. Agarwal,
Matoušek and Sharir [3] show that we can preprocess the vertices of P and the vertices
of Q into a data structure that can answer our desired range reporting queries in time
O(n1/2 logc n + k), where k is the output size and c > 0 is some fixed constant.
The expected preprocessing time is O(n1+ε), where ε > 0 can be made arbitrarily
small. We perform O(n) queries, and the total expected output size of our queries
is O(n3/2 ln n), so it takes expected time O(n3/2 logc n) to find the critical values in
[a, b]. Finally, we perform a binary search on theses critical values to compute the
weak Fréchet distance. This takes O(T (n) log n) time.
The following theorem summarizes our results on the weak Fréchet distance.
Theorem 10.10 The weak Fréchet distance of two polygonal curves, each with n
edges, can be computed by a randomized algorithm in time O (n2α(n) log log n) on a
pointer machine and in time O ((n2/ log n)(log log n)5) on a word RAM.
We have broken the longstanding quadratic upper bound for the decision version of
the Fréchet problem. Moreover, we have shown that this problem has an algebraic
decision tree of depth O (n2−ε), for some ε > 0 and where n is the number of vertices
of the polygonal curves. We have shown how our faster algorithm for the decision
version can be used for a faster algorithm to compute the Fréchet distance. If we allow
constanttime tablelookup, we obtain a running time in close reach of O (n2). This
leaves us with intriguing open research questions. Can we devise a quadratic or even
a slightly subquadratic algorithm for the optimization version? Can we devise such an
algorithm on the word RAM, that is, with constanttime tablelookup? What can be
said about approximation algorithms?
Acknowledgements We would like to thank Natan Rubin for pointing out to us that Fréchetrelated papers
require a witty title involving a dog, Bettina Speckmann for inspiring discussions during the early stages of
this research, Günter Rote for providing us with his Ipelet for the freespace diagram, and Otfried Cheong for
Ipe. We would also like to thank anonymous referees for suggesting simpler proofs of Lemmas 4.1 and 7.1,
for pointing out [53] for Theorem 9.1, and for suggesting a section on the weak Fréchet distance. K. Buchin
and M. Buchin in part supported by European Cooperation in Science and Technology (COST) ICT Action
IC0903 MOVE. M. Buchin supported in part by the Netherlands Organisation for Scientific Research
(NWO) under Project No. 612.001.106. W. Meulemans supported by the Netherlands Organisation for
Scientific Research (NWO) under Project No. 639.022.707. W. Mulzer supported in part by DFG Projects
MU/3501/1 and MU/3501/2.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0
International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution,
and reproduction in any medium, provided you give appropriate credit to the original author(s) and the
source, provide a link to the Creative Commons license, and indicate if changes were made.
Appendix: Computational Models
The standard machine model in computational geometry is the real RAM [56]. Here,
data is represented as a (countably) infinite sequence of storage cells. These cells
can be of two different types: they can store real numbers or integers. The model
supports standard operations on these numbers in constant time, including addition,
multiplication, and elementary functions like squareroot, sine or cosine. Furthermore,
the integers can be used as indices to memory locations. Integers can be converted to
real numbers in constant time, but we need to be careful about the reverse direction.
The floor function can be used to truncate a real number to an integer, but if we were
allowed to use it arbitrarily, the real RAM could solve PSPACEcomplete problems
The word RAM is essentially a real RAM without support for real numbers. However,
on a real RAM, the integers are usually treated as atomic, whereas the word RAM
allows for powerful bitmanipulation tricks [42]. More precisely, the word RAM
represents the data as a sequence of wbit words, where w = (log n). Data can be
accessed arbitrarily, and standard operations, such as Boolean operations (and, xor,
shl, . . .), addition, or multiplication take constant time. There are many variants of
the word RAM, depending on precisely which instructions are supported in constant
time [27, 42]. The general consensus seems to be that any function in AC0 is
acceptable.7 However, it is always preferable to rely on a set of operations as small, and as
nonexotic, as possible. Note that multiplication is not in AC0 [43], but nevertheless
is often included in the word RAM instruction set [42].
The pointer machine model disallows the use of constant time table lookup, and is
therefore a restriction of the (real) RAM model [59, 61]. The data structure is modeled
as a directed graph G with bounded outdegree. Each node in G represents a record,
with a bounded number of pointers to other records and a bounded number of (real or
integer) data items. The algorithm can access data only by following pointers from the
inputs (and a bounded number of global entry records); random access is not possible.
The data can be manipulated through the usual real RAM operations, but without
support for the floor function, for reasons mentioned above.
Algebraic Computation Tree
Algebraic computation trees (ACTs) [11] are the computational geometry analogue
of binary decision trees, and like these they are mainly used for proving lower bounds.
Let x1, . . . , xn ∈ R be the inputs. An ACT is a binary tree with two kinds of nodes:
computation and branch nodes. A computation node v has one child and is labeled
with an expression of the type yv = yu ⊕ yw, where ⊕ ∈ {+, −, ∗, /, √·} is a operator
and yu , yw is either an input variable x1, . . . , xn or corresponds to a computation node
that is an ancestor of v. A branch node has degree 2 and is labeled by yu = 0 or
yu > 0, where again yu is either an input or a variable for an ancestor. A family of
algebraic computation trees (Tn )n∈N solves a computational problem (like Delaunay
7 AC0 is the class of all functions f : {0, 1}∗ → {0, 1}∗ that can be computed by a family of circuits
(Cn )n∈N with the following properties [11]: (i) each Cn has n inputs; (ii) there exist constants a, b, such
that Cn has at most anb gates, for n ∈ N; (iii) there is a constant d such that for all n the length of the
longest path from an input to an output in Cn is at most d (i.e., the circuit family has bounded depth); (iv)
each gate has an arbitrary number of incoming edges (i.e., the fanin is unbounded).
triangulation or convex hull computation), if for each n ∈ N, the tree Tn accepts inputs
of size n, and if for any such input x1, . . . , xn the corresponding path in Tn (where the
children of the branch nodes are determined according the conditions they represent)
constitutes a computation that represents the answer in the variables yv encountered
during the path.
Algebraic decision trees are defined as follows: we allow only branch nodes.
Each branch node is labeled with a predicate of the form p(x1, . . . , xn ) = 0 or
p(x1, . . . , xn ) > 0. The leaves are labeled yes or no. Fix some r ∈ {1, . . . , n}. If
p is restricted to be of the form p(x1, . . . , xn ) = in=1 ai xi − b, with at most r
coefficients ai = 0, we call the decision tree r linear. Erickson [40] showed that
any 3linear decision tree for 3SUM has depth (n2). However, Grønlund and Pettie
showed that there is a 4linear decision tree of depth O (n3/2√log n) for the problem.
In geometric problems, linear predicates are often much too restrictive. For example,
there is no r linear decision tree for the Fréchet problem, no matter the choice of r :
with r linear decision trees, we cannot even decide whether two given points p and q
have Euclidean distance at most 1.
1. Agarwal , P.K. , Ben Avraham , R., Kaplan , H. , Sharir , M. : Computing the discrete Fréchet distance in subquadratic time . SIAM J. Comput . 43 ( 2 ), 429  449 ( 2014 )
2. Agarwal , P.K. , HarPeled , S. , Mustafa , N.H. , Wang , Y. : Nearlinear time approximation algorithms for curve simplification . Algorithmica 42 ( 34 ), 203  219 ( 2005 )
3. Agarwal , P.K. , Matoušek , J. , Sharir , M. : On range searching with semialgebraic sets . II. SIAM J. Comput . 42 ( 6 ), 2039  2062 ( 2013 )
4. Ailon , N. , Chazelle , B. : Lower bounds for linear degeneracy testing . J. ACM 52 ( 2 ), 157  171 ( 2005 )
5. Albers , S. , Hagerup, T.: Improved parallel integer sorting without concurrent writing . Inf. Comput . 136 ( 1 ), 25  51 ( 1997 )
6. Alt , H. : The computational geometry of comparing shapes . In: Albers, S. , Alt , H. , Näher , S. (eds.) Efficient Algorithms. Lecture Notes in Computer Science , vol. 5760 , pp. 235  248 . Springer, Berlin ( 2009 )
7. Alt , H. , Buchin , M. : Can we compute the similarity between surfaces? Discrete Comput . Geom. 43 ( 1 ), 78  99 ( 2010 )
8. Alt , H. , Godau , M. : Computing the Fréchet distance between two polygonal curves . Int. J. Comput. Geom. Appl . 5 ( 12 ), 75  91 ( 1995 )
9. Alt , H. , Knauer , C. , Wenk , C. : Comparison of distance measures for planar curves . Algorithmica 38 ( 1 ), 45  58 ( 2003 )
10. Aronov , B. , HarPeled , S. , Knauer , C. , Wang , Y. , Wenk , C. : Fréchet distances for curves, revisited . In: Azar, Y. , Erlebach, T. (eds.) AlgorithmsESA 2006. Lecture Notes in Computer Science , vol. 4168 , pp. 52  63 . Springer, Berlin ( 2006 )
11. Arora , S. , Barak , B. : Computational Complexity . Cambridge University Press , Cambridge ( 2009 )
12. Baran , I., Demaine , E.D. , Pa˘tras¸cu, M.: Subquadratic algorithms for 3SUM . Algorithmica 50 ( 4 ), 584  596 ( 2008 )
13. Bellman , R. , Kalaba , R.: On adaptive control processes . IRE Trans. Autom. Control 4 ( 2 ), 1  9 ( 1959 )
14. Ben Avraham , R., Filtser , O. , Kaplan , H. , Katz , M.J. , Sharir , M. : He discrete and semicontinuous Fréchet distance with shortcuts via approximate distance counting and selection . ACM Trans. Algorithms 11 ( 4 ), 29 ( 2015 )
15. Brakatsoulas , S. , Pfoser , D. , Salas , R. , Wenk , C.: On mapmatching vehicle tracking data . In: Böhm, K. , et al. (eds.) Proceedings of the 31st VLDB Conference , pp. 853  864 . ACM, New York ( 2005 )
16. Bremner , D. , Chan, T.M., Demaine , E.D. , Erickson , J. , Hurtado , F. , Iacono , J. , Langerman , S. , Pa˘tras¸cu, M. , Taslakian , P. : Necklaces, convolutions, and X + Y . Algorithmica 69 ( 2 ), 294  314 ( 2012 )
17. Bringmann , K. : Why walking the dog takes time: Fréchet distance has no strongly subquadratic algorithms unless SETH fails . In: 55th Annual IEEE Symposium on Foundations of Computer ScienceFOCS 2014 , pp. 661  670 . IEEE Computer Society, Los Alamitos ( 2014 )
18. Bringmann , K. , Mulzer , W.: Approximability of the discrete Fréchet distance . J. Comput. Geom . 7 ( 2 ), 46  76 ( 2016 )
19. Buchin , K. , Buchin , M. , Gudmundsson , J. : Constrained free space diagrams: a tool for trajectory analysis . Int. J. GIS 24 ( 7 ), 1101  1125 ( 2010 )
20. Buchin , K. , Buchin , M. , Gudmundsson , J. , Löffler , M. , Luo , J. : Detecting commuting patterns by clustering subtrajectories . Int. J. Comput. Geom. Appl . 21 ( 3 ), 253  282 ( 2011 )
21. Buchin , K. , Buchin , M. , Knauer , C. , Rote , G. , Wenk , C. : How difficult is it to walk the dog ? In: Aichholzer, O. , Hackl , T. (eds.) 23rd EuroCG/FWCG , pp. 170  173 . Technischen Universität Graz, Graz ( 2007 )
22. Buchin , K. , Buchin , M. , Meulemans , W. , Speckmann , B. : Locally correct Fréchet matchings . In: Epstein, L. , Ferragina , P. (eds.) AlgorithmsESA 2012. Lecture Notes in Computer Science , vol. 7501 , pp. 229  240 . Springer, Heidelberg ( 2012 )
23. Buchin , K. , Buchin , M. , Schulz , A. : Fréchet distance of surfaces: some simple hard cases . In: de Berg, M. , Meyer , U. (eds.) AlgorithmsESA 2010. Part II. Lecture Notes in Computer Science , vol. 6347 , pp. 63  74 . Springer, Berlin ( 2010 )
24. Buchin , K. , Buchin , M. , van Leusden , R. , Meulemans , W. , Mulzer , W.: Computing the Fréchet distance with a retractable leash . Discrete Comput. Geom . 56 ( 2 ), 315  336 ( 2016 )
25. Buchin , K. , Buchin , M. , Wang , Y. : Exact algorithms for partial curve matching via the Fréchet distance . In: Proceedings of the Twentieth Annual ACMSIAM Symposium on Discrete Algorithms , pp. 645  654 . SIAM, Philadelphia ( 2009 )
26. Buchin , K. , Buchin , M. , Wenk , C. : Computing the Fréchet distance between simple polygons . Comput. Geom . 41 ( 12 ), 2  20 ( 2008 )
27. Buchin , K. , Mulzer , W.: Delaunay triangulations in O(sort(n)) time and more . J. ACM 58 ( 2 ), 6 ( 2011 )
28. Buchin , M. : On the Computability of the Fréchet Distance Between Triangulated Surfaces . PhD thesis , Free University Berlin, Berlin ( 2007 ). http://www.diss.fuberlin. de/diss/receive/FUDISS_thesis_ 000000002618
29. Chambers , E.W. , Colin de Verdière , É. , Erickson , J. , Lazard , S. , Lazarus , F. , Thite , S. : Homotopic Fréchet distance between curves or, walking your dog in the woods in polynomial time . Comput. Geom . 43 ( 3 ), 295  311 ( 2010 )
30. Chan, T.M.: Allpairs shortest paths with real weights in O(n3/ log n) time . Algorithmica 50 ( 2 ), 236  243 ( 2008 )
31. Chan, T.M.: More algorithms for allpairs shortest paths in weighted graphs . SIAM J. Comput . 39 ( 5 ), 2075  2089 ( 2010 )
32. Chazelle , B.M. , Lee , D.T.: On a circle placement problem . Computing 36 ( 12 ), 1  16 ( 1986 )
33. Cook IV , A.F., Driemel , A. , HarPeled , S. , Sherette , J. , Wenk , C. : Computing the Fréchet distance between folded polygons . In: Dehne, F. , Iacono , J. , Sack , J.R . (eds.) Algorithms and Data Structures . Lecture Notes in Computer Science , vol. 6844 , pp. 267  278 . Springer, Heidelberg ( 2011 )
34. Cook IV , A.F., Wenk , C. : Geodesic Fréchet distance inside a simple polygon . ACM Trans. Algorithms 7 ( 1 ), 9 ( 2010 )
35. de Berg, M. , Cook IV , A.F., Gudmundsson , J. : Fast Fréchet queries . Comput. Geom . 46 ( 6 ), 747  755 ( 2013 )
36. Driemel , A. , HarPeled , S. : Jaywalking your dog: computing the Fréchet distance with shortcuts . SIAM J. Comput . 42 ( 5 ), 1830  1866 ( 2013 )
37. Driemel , A. , HarPeled , S. , Wenk , C. : Approximating the Fréchet distance for realistic curves in near linear time . Discrete Comput. Geom . 48 ( 1 ), 94  127 ( 2012 )
38. Efrat , A. , Guibas , L.J. , HarPeled , S. , Mitchell, J.S.B., Murali, T.M.: New similarity measures between polylines with applications to morphing and polygon sweeping . Discrete Comput. Geom . 28 ( 4 ), 535  569 ( 2002 )
39. Eiter , T. , Mannila , H. : Computing Discrete Fréchet Distance . Technical report CDTR 94 /65, Christian Doppler Laboratory ( 1994 )
40. Erickson , J. : Bounds for linear satisfiability problems . Chic. J. Theor. Comput. Sci . 1999 , 8 ( 1999 )
41. Fredman , M.L. : How good is the information theory bound in sorting? Theor . Comput. Sci. 1 ( 4 ), 355  361 ( 1975 /76)
42. Fredman , M.L. , Willard , D.E. : Surpassing the informationtheoretic bound with fusion trees . J. Comput. Syst. Sci . 47 ( 3 ), 424  436 ( 1993 )
43. Fredman , M.L. , Saxe , J.B. , Sipser , M. : Parity, circuits, and the polynomialtime hierarchy . Math. Syst. Theory 17 ( 1 ), 13  27 ( 1984 )
44. Gajentaan , A. , Overmars , M.H.: On a class of O(n2) problems in computational geometry . Comput. Geom . 5 ( 3 ), 165  185 ( 1995 )
45. Godau , M. : A natural metric for curvescomputing the distance for polygonal chains and approximation algorithms . In: Choffrut, C. , Jantzen , M. (eds.) STACS 91. Lecture Notes in Computer Science , vol. 480 , pp. 127  136 . Springer, Berlin ( 1991 )
46. Godau , M. : On the Complexity of Measuring the Similarity Between Geometric Objects in Higher Dimensions . PhD thesis , Free University Berlin, Berlin ( 1998 )
47. Gudmundsson , J. , Wolle , T.: Towards automated football analysis: algorithms and data structures . In: Proceedings of the 10th Australasian Conference on Mathematics and Computers in Sport 2010 . ANZIAM, Darwin ( 2010 )
48. HarPeled , S. , Nayyeri , A. , Salavatipour , M. , Sidiropoulos , A. : How to walk your dog in the mountains with no magic leash . In: Dey, T. , Whitesides , S. (eds.) Computational Geometry (SCG'12) , pp. 121  130 . ACM, New York ( 2012 )
49. HarPeled , S. , Raichel , B. : The Fréchet distance revisited and extended . ACM Trans. Algorithms 10 ( 1 ), 3 ( 2014 )
50. Hirschberg , D.S. , Chandra , A.K. , Sarwate , D.V. : Computing connected components on parallel computers . Commun. ACM 22 ( 8 ), 461  464 ( 1979 )
51. Indyk , P. : Approximate nearest neighbor algorithms for Frechet distance via product metrics . In: Computational Geometry (SCG'02) , pp. 102  106 . ACM, New York ( 2002 )
52. Grønlund , A. , Pettie , S. : Threesomes, degenerates, and love triangles . In: 55th Annual IEEE Symposium on Foundations of Computer ScienceFOCS 2014 , pp. 621  630 . IEEE Computer Society, Los Alamitos ( 2014 )
53. Katz , M.J. , Sharir , M. : An expanderbased approach to geometric optimization . SIAM J. Comput . 26 ( 5 ), 1384  1408 ( 1997 )
54. Maheshwari , A. , Sack , J.R. , Shahbaz , K. , ZarrabiZadeh , H. : Fréchet distance with speed limits . Comput. Geom . 44 ( 2 ), 110  120 ( 2011 )
55. Maheshwari , A. , Sack , J.R. , Shahbaz , K. , ZarrabiZadeh , H. : Improved algorithms for partial curve matching . In: Demetrescu, C. , Halldórsson , M.M. (eds.) AlgorithmsESA 2011. Lecture Notes in Computer Science , vol. 6942 , pp. 518  529 . Springer, Heidelberg ( 2011 )
56. Preparata , F.P. , Shamos , M.I.: Computational Geometry . Texts and Monographs in Computer Science. Springer, New York ( 1985 )
57. Pa˘tras¸cu, M.: Towards polynomial lower bounds for dynamic problems . In: STOC'10 , pp. 603  610 . ACM, New York ( 2010 )
58. Schönhage , A. : On the power of random access machines . In: Maurer, H.A. (ed.) Automata, Languages and Programming. Lecture Notes in Computer Science , vol. 71 , pp. 520  529 . Springer, Berlin ( 1979 )
59. Schönhage , A. : Storage modification machines . SIAM J. Comput . 9 ( 3 ), 490  508 ( 1980 )
60. Sharir , M. , Agarwal , P.K. : DavenportSchinzel Sequences and Their Geometric Applications . Cambridge University Press, Cambridge ( 1995 )
61. Tarjan , R.E.: Efficiency of a good but not linear set union algorithm . J. Assoc. Comput. Mach . 22 ( 2 ), 215  225 ( 1975 )
62. Thorup , M. : Randomized sorting in O(n log log n) time and linear space using addition, shift, and bitwise Boolean operations . J. Algorithms 42 ( 2 ), 205  230 ( 2002 )
63. Wenk , C. , Salas , R. , Pfoser , D. : Addressing the need for mapmatching speed: localizing global curvematching algorithms . In: Proceedings of the 18th International Conference on Scientific and Statistical Database Management , pp. 379  388 . IEEE Computer Society, Los Alamitos ( 2006 )