]> mj.ucw.cz Git - saga.git/blob - ram.tex
Call subgraphs clusters.
[saga.git] / ram.tex
1 \ifx\endpart\undefined
2 \input macros.tex
3 \fi
4
5 \chapter{Fine Details of Computation}
6 \id{ramchap}
7
8 \section{Models and machines}
9
10 Traditionally, computer scientists use a~variety of computational models
11 as a~formalism in which their algorithms are stated. If we were studying
12 NP-completeness, we could safely assume that all the models are equivalent,
13 possibly up to polynomial slowdown which is negligible. In our case, the
14 differences between good and not-so-good algorithms are on a~much smaller
15 scale. In this chapter, we will replace the usual ``tape measure'' by a~micrometer,
16 state our computation models carefully and develop a repertoire of basic
17 data structures taking advantage of the fine details of the models.
18
19 We would like to keep the formalism close enough to the reality of the contemporary
20 computers. This rules out Turing machines and similar sequentially addressed
21 models, but even the remaining models are subtly different from each other. For example, some of them
22 allow indexing of arrays in constant time, while on the others,
23 arrays have to be emulated with pointer structures, requiring $\Omega(\log n)$
24 time to access a~single element of an~$n$-element array. It is hard to say which
25 way is superior --- while most ``real'' computers have instructions for constant-time
26 indexing, it seems to be physically impossible to fulfil this promise regardless of
27 the size of memory. Indeed, at the level of logical gates inside the computer,
28 the depth of the actual indexing circuits is logarithmic.
29
30 In recent decades, most researchers in the area of combinatorial algorithms
31 have been considering two computational models: the Random Access Machine and the Pointer
32 Machine. The former is closer to the programmer's view of a~real computer,
33 the latter is slightly more restricted and ``asymptotically safe.''
34 We will follow this practice and study our algorithms in both models.
35
36 \para
37 The \df{Random Access Machine (RAM)} is not a~single coherent model, but rather a~family
38 of closely related machines, sharing the following properties.
39 (See Cook and Reckhow \cite{cook:ram} for one of the usual formal definitions
40 and Hagerup \cite{hagerup:wordram} for a~thorough description of the differences
41 between the RAM variants.)
42
43 The \df{memory} of the machine is represented by an~array of \df{memory cells}
44 addressed by non-negative integers, each of them containing a~single non-negative integer.
45 The \df{program} is a~finite sequence of \df{instructions} of two basic kinds: calculation
46 instructions and control instructions.
47
48 \df{Calculation instructions} have two source arguments and one destination
49 argument, each \df{argument} being either an~immediate constant (not available
50 as destination), a~directly addressed memory cell (specified by its number)
51 or an~indirectly addressed memory cell (its address is stored in a~directly
52 addressed memory cell).
53
54 \df{Control instructions} include branches (to a~specific instruction in
55 the program), conditional branches (e.g., jump if two arguments specified as
56 in the calculation instructions are equal) and an~instruction to halt the program.
57
58 At the beginning of the computation, the memory contains the input data
59 in specified cells and arbitrary values in all other cells.
60 Then the program is executed one instruction at a~time. When it halts,
61 specified memory cells are interpreted as the program's output.
62
63 \para\id{wordsize}%
64 In the description of the RAM family, we have omitted several details
65 on~purpose, because different members of the family define them differently.
66 These are: the size of the available integers, the time complexity of a~single
67 instruction, the space complexity assigned to a~single memory cell and the set
68 of operations available in calculation instructions.
69
70 If we impose no limits on the magnitude of the numbers and we assume that
71 arithmetic and logical operations work on them in constant time, we get
72 a~very powerful parallel computer --- we can emulate an~exponential number
73 of parallel processors using arithmetics and suddenly almost everything can be
74 computed in constant time, modulo encoding and decoding of input and output.
75 Such models are unrealistic and there are two basic possibilities how to
76 avoid this behavior:
77
78 \numlist\ndotted
79 \:Keep unbounded numbers, but increase costs of instructions: each instruction
80   consumes time proportional to the number of bits of the numbers it processes,
81   including memory addresses. Similarly, space usage is measured in bits,
82   counting not only the values, but also the addresses of the respective memory
83   cells.
84 \:Place a~limit on the size of the numbers ---define the \df{word size~$W$,}
85   the number of bits available in the memory cells--- and keep the cost of
86   instructions and memory cells constant. The word size must not be constant,
87   since we can address only~$2^W$ cells of memory. If the input of the algorithm
88   is stored in~$N$ cells, we need~$W\ge\log N$ just to be able to read the input.
89   On the other hand, we are interested in polynomial-time algorithms only, so $\Theta(\log N)$-bit
90   numbers should be sufficient. In practice, we pick~$W$ to be the larger of
91   $\Theta(\log N)$ and the size of integers used in the algorithm's input and output.
92   We will call an integer that fits in a~single memory cell a~\df{machine word.}
93 \endlist
94
95 Both restrictions easily avoid the problems of unbounded parallelism. The first
96 choice is theoretically cleaner and Cook et al.~show nice correspondences to the
97 standard complexity classes, but the calculations of time and space complexity tend
98 to be somewhat tedious. What more, when compared with the RAM with restricted
99 word size, the complexities are usually exactly $\Theta(W)$ times higher.
100 This does not hold in general (consider a~program that uses many small numbers
101 and $\O(1)$ large ones), but it is true for the algorithms we are interested in.
102 Therefore we will always assume that the operations have unit cost and we make
103 sure that all numbers are limited by the available word size.
104
105 \para
106 As for the choice of RAM operations, the following three instruction sets are often used:
107
108 \itemize\ibull
109 \:\df{Word-RAM} --- allows the ``C-language operators'', i.e., addition,
110   subtraction, multiplication, division, remainder, bitwise $\band$, $\bor$, exclusive
111   $\bor$ ($\bxor$) and negation ($\bnot$), and bitwise shifts ($\shl$ and~$\shr$).
112 \:\df{${\rm AC}^0$-RAM} --- allows all operations from the class ${\rm AC}^0$, i.e.,
113   those computable by constant-depth polynomial-size boolean circuits with unlimited
114   fan-in and fan-out. This includes all operations of the Word-RAM except for multiplication,
115   division and remainders, and also many other operations like computing the Hamming
116   weight (number of bits set in a~given number).
117 \:Both restrictions at once.
118 \endlist
119
120 Thorup discusses the usual techniques employed by RAM algorithms in~\cite{thorup:aczero}
121 and he shows that they work on both Word-RAM and ${\rm AC}^0$-RAM, but the combination
122 of the two restrictions is too weak. On the other hand, the intersection of~${\rm AC}^0$
123 with the instruction set of modern processors is already strong enough (e.g., when we
124 add some floating-point operations and multimedia instructions available on the Intel's
125 Pentium~4~\cite{intel:pentium}).
126
127 We will therefore use the Word-RAM instruction set, mentioning differences from the
128 ${\rm AC}^0$-RAM where necessary.
129
130 \nota
131 When speaking of the \df{RAM,} we implicitly mean the version with numbers limited
132 by a~specified word size of $W$~bits, unit cost of operations and memory cells and the instruction
133 set of the Word-RAM. This corresponds to the usage in recent algorithmic literature,
134 although the authors rarely mention the details. In some cases, a~non-uniform variant
135 of the Word-RAM is considered as well (e.g., in~\cite{hagerup:dd}):
136
137 \defn\id{nonuniform}%
138 A~Word-RAM is called \df{weakly non-uniform,} if it is equipped with $\O(1)$-time
139 access to a~constant number of word-sized constants, which depend only on the word
140 size. These are called \df{native constants} and they are available in fixed memory
141 cells when the program starts. (By analogy with the high-level programming languages,
142 these constants can be thought of as computed at ``compile time.'')
143
144 \para
145 The \df{Pointer Machine (PM)} also does not have any well established definition. The
146 various kinds of pointer machines are mapped by Ben-Amram in~\cite{benamram:pm},
147 but unlike the RAM's they turn out to be equivalent up to constant slowdown.
148 Our definition will be closely related to the \em{linking automaton} proposed
149 by Knuth in~\cite{knuth:fundalg}, we will only adapt it to use RAM-like
150 instructions instead of an~opaque control unit.
151
152 The PM works with two different types of data: \df{symbols} from a~finite alphabet
153 and \df{pointers}. The memory of the machine consists of a~fixed amount of \df{registers}
154 (some of them capable of storing a~single symbol, each of the others holds a~single pointer)
155 and an~arbitrary amount of \df{cells}. The structure of all cells is the same: each of them
156 again contains a~fixed number of fields for symbols and pointers. Registers can be addressed
157 directly, the cells only via pointers --- by using a~pointer stored either in a~register,
158 or in a~cell pointed to by a~register (longer chains of pointers cannot be followed in
159 constant time).
160
161 We can therefore view the whole memory as a~directed graph, whose vertices
162 correspond to the cells (the registers are stored in a~single special cell).
163 The outgoing edges of each vertex correspond to pointer fields of the cells and they are
164 labeled with distinct labels drawn from a~finite set. In addition to that,
165 each vertex contains a~fixed amount of symbols. The program can directly access
166 vertices within distance~2 from the register vertex.
167
168 The program is a~finite sequence of instructions of the following kinds:
169
170 \itemize\ibull
171 \:\df{symbol instructions,} which read a~pair of symbols, apply an~arbitrary
172   function on them and write the result to a~symbol register or field;
173 \:\df{pointer instructions} for assignment of pointers to pointer registers/fields
174   and for creation of new memory cells (a~pointer to the new cell is assigned
175   immediately);
176 \:\df{control instructions} --- similarly to the RAM; conditional jumps can decide
177   on~arbitrary unary relations on symbols and compare pointers for equality.
178 \endlist
179
180 Time and space complexity are defined in the straightforward way: all instructions
181 have unit cost and so do all memory cells.
182
183 Both input and output of the machine are passed in the form of a~linked structure
184 pointed to by a~designated register. For example, we can pass graphs back and forth
185 without having to encode them as strings of numbers or symbols. This is important,
186 because with the finite alphabet of the~PM, symbolic representations of graphs
187 generally require super-linear space and therefore also time.\foot{%
188 The usual representation of edges as pairs of vertex labels uses $\Theta(m\log n)$ bits
189 and as a~simple counting argument shows, this is asymptotically optimal for general
190 sparse graphs. On the other hand, specific families of sparse graphs can be stored
191 more efficiently, e.g., by a~remarkable result of Tur\'an~\cite{turan:succinct},
192 planar graphs can be encoded in~$\O(n)$ bits. Encoding of dense graphs is of
193 course trivial as the adjacency matrix has only~$\Theta(n^2)$ bits.}
194
195 \para
196 Compared to the RAM, the PM lacks two important capabilities: indexing of arrays
197 and arithmetic instructions. We can emulate both with poly-logarithmic slowdown,
198 but it will turn out that they are rarely needed in graph algorithms. We are
199 also going to prove that the RAM is strictly stronger, so we will prefer to
200 formulate our algorithms in the PM model and use RAM only when necessary.
201
202 \thm
203 Every program for the Word-RAM with word size~$W$ can be translated to a~PM program
204 computing the same with $\O(W^2)$ slowdown (given a~suitable encoding of inputs and
205 outputs, of course). If the RAM program does not use multiplication, division
206 and remainder operations, $\O(W)$~slowdown is sufficient.
207
208 \proofsketch
209 Represent the memory of the RAM by a~balanced binary search tree or by a~radix
210 trie of depth~$\O(W)$. Values are encoded as~linked lists of symbols pointed
211 to by the nodes of the tree. Both direct and indirect accesses to the memory
212 can therefore be done in~$\O(W)$ time. Use standard algorithms for arithmetic
213 on big numbers: $\O(W)$ per operation except for multiplication, division and
214 remainders which take $\O(W^2)$.\foot{We could use more efficient arithmetic
215 algorithms, but the quadratic bound is good enough for our purposes.}
216 \qed
217
218 \FIXME{Add references, especially to the unbounded parallelism remark.}
219
220 \thm
221 Every program for the PM running in polynomial time can be translated to a~program
222 computing the same on the Word-RAM with only $\O(1)$ slowdown.
223
224 \proofsketch
225 Encode each cell of the PM's memory to $\O(1)$ integers. Store the encoded cells to
226 the memory of the RAM sequentially and use memory addresses as pointers. As the symbols
227 are finite and there is only a~polynomial number of cells allocated during execution
228 of the program, $\O(\log N)$-bit integers suffice ($N$~is the size of the program's input).
229 \qed
230
231 \para
232 There are also \df{randomized} versions of both machines. These are equipped
233 with an~additional instruction for generating a~single random bit. The standard
234 techniques of design and analysis of randomized algorithms apply (see for
235 example Motwani and Raghavan~\cite{motwani:randalg}).
236
237 \FIXME{Consult sources. Does it make more sense to generate random words at once on the RAM?}
238
239 \rem
240 There is one more interesting machine: the \df{Immutable Pointer Machine} (see
241 the description of LISP machines in \cite{benamram:pm}). It differs from the
242 ordinary PM by the inability to modify existing memory cells. Only the contents
243 of the registers are allowed to change. All cell modifications thus have to
244 be performed by creating a~copy of the particular cell with some fields changed.
245 This in turn requires the pointers to the cell to be updated, possibly triggering
246 a~cascade of further cell copies. For example, when a~node of a~binary search tree is
247 updated, all nodes on the path from that node to the root have to be copied.
248
249 One of the advantages of this model is that the states of the machine are
250 persistent --- it is possible to return to a~previously visited state by recalling
251 the $\O(1)$ values of the registers (everything else could not have changed
252 since that time) and ``fork'' the computations. This corresponds to the semantics
253 of pure functional languages, e.g., Haskell~\cite{jones:haskell}.
254
255 Unless we are willing to accept a~logarithmic penalty in execution time and space
256 (in fact, our emulation of the Word-RAM on the PM can be easily made immutable),
257 the design of efficient algorithms for the immutable PM requires very different
258 techniques. Therefore, we will concentrate on the imperative models instead
259 and refer the interested reader to the thorough treatment of purely functional
260 data structures in the Okasaki's monograph~\cite{okasaki:funcds}.
261
262 %--------------------------------------------------------------------------------
263
264 \section{Bucket sorting and unification}\id{bucketsort}%
265
266 The Contractive Bor\o{u}vka's algorithm (\ref{contbor}) needed to contract a~given
267 set of edges in the current graph and flatten it afterwards, all this in time $\O(m)$.
268 We have spared the technical details for this section, in which we are going to
269 explain several techniques based on bucket sorting. These will be useful in further
270 algorithms, too.
271
272 As already suggested in the proof of Lemma \ref{contbor}, contractions can be performed
273 in linear time by building an~auxiliary graph and finding its connected components.
274 We will thus take care only of the subsequent flattening.
275
276 \paran{Flattening on RAM}%
277 On the RAM, we can view the edges as ordered pairs of vertex identifiers with the
278 smaller of the identifiers placed first and sort them lexicographically. This brings
279 parallel edges together, so that a~simple linear scan suffices to find each bunch
280 of parallel edges and remove all but the lightest one.
281 Lexicographic sorting of pairs can be accomplished in linear time by a~two-pass
282 bucket sort with $n$~buckets corresponding to the vertex identifiers.
283
284 However, there is a~catch in this. Suppose that we use the standard representation
285 of graphs by adjacency lists whose heads are stored in an array indexed by vertex
286 identifiers. When we contract and flatten the graph, the number of vertices decreases,
287 but if we inherit the original vertex identifiers, the arrays will still have the
288 same size. We could then waste a~super-linear amount of time by scanning the increasingly
289 sparse arrays, most of the time skipping unused entries.
290
291 To avoid this problem, we have to renumber the vertices after each contraction to component
292 identifiers from the auxiliary graph and create a~new vertex array. This helps to
293 keep the size of the representation of the graph linear with respect to its current
294 size.
295
296 \paran{Flattening on PM}%
297 The pointer representation of graphs does not suffer from sparsity since the vertices
298 are always identified by pointers to per-vertex structures. Each such structure
299 then contains all attributes associated with the vertex, including the head of its
300 adjacency list. However, we have to find a~way how to perform bucket sorting
301 without indexing of arrays.
302
303 We will keep a~list of the per-vertex structures that defines the order of~vertices.
304 Each such structure will be endowed with a~pointer to the head of the list of items in
305 the corresponding bucket. Inserting an~edge to a~bucket can be then done in constant time
306 and scanning all~$n$ buckets takes $\O(n+m)$ time.
307
308 \paran{Tree isomorphism}%
309 Another nice example of pointer-based radix sorting is a~pointer algorithm for
310 deciding whether two rooted trees are isomorphic. Let us assume for a~moment that
311 the outdegree of each vertex is at most a~fixed constant~$k$. We sort the subtrees
312 of both trees by their depth by running the depth-first search to calculate the
313 depths and bucket-sorting them with $n$~buckets afterwards.
314
315 Then we proceed from depth~1 to the maximum depth and for each of them we identify
316 the isomorphism equivalence classes of subtrees of that particular depth. We will assign some
317 sort of identifiers to the classes; at most~$n+1$ of them are needed as there are
318 $n+1$~subtrees in the tree (including the empty subtree). As the PM does not
319 have numbers as a~first-class type, we just create a~``\df{yardstick}'' ---
320 a~list of $n+1$~distinct items and we use pointers to these items as identifiers.
321 Isomorphism of the whole trees can be finally decided by comparing the
322 identifiers assigned to their roots.
323
324 Suppose that classes of depths $1,\ldots,d-1$ are already computed and we want
325 to identify those of depth~$d$. We will denote their number of~$n_d$. We take
326 a~root of every such tree and label it with an~ordered $k$-tuple of identifiers
327 of its subtrees; when it has less than $k$ sons, we pad the tuple with empty
328 subtrees. Tuples corresponding to isomorphic subtrees are identical up to
329 reordering of elements. We therefore sort the codes inside each tuple and then
330 sort the tuples, which brings the equivalent tuples together.
331
332 The first sort (inside the tuples) would be easy on the RAM, but on the PM we
333 have no means of comparing two identifiers for anything else than equality.
334 To work around this, we sort the set $\{ (x,i,j) \mid \hbox{$x$ is the $i$-th
335 element of the $j$-th tuple} \}$ on~$x$, reset all tuples and insert the elements
336 back in the increasing order of~$x$, ignoring the original positions. The second
337 sort is a~straightforward $k$-pass bucket sort.
338
339 If we are not careful, a~single sorting pass takes $\O(n_d + n)$ time, because
340 while we have only $n_d$~items to sort, we have to scan all $n$~buckets. This can
341 be easily avoided if we realize that the order of buckets does not need to be
342 fixed --- in every pass, we can use a~completely different order and it still
343 does bring the equivalent tuples together. Thus we can keep a~list of buckets
344 which were used in the current pass and look only inside these buckets. This way,
345 the pass takes $\O(n_d)$ time only and the whole algorithm is $\O(\sum_d n_d) = \O(n)$.
346
347 Our algorithm can be easily modified for trees with unrestricted degrees.
348 We replace the fixed $d$-tuples by general sequences of identifiers. The first
349 sort does not need any changes. In the second sort, we proceed from the first
350 position to the last one and after each bucket sort pass we put aside the sequences
351 that have just ended. (They are obviously not equivalent to any other sequences.)
352 The second sort is linear in the sum of the lengths of the sequences, which is
353 $n_{d+1}$ for depth~$d$. We can therefore decide isomorphism of the whole trees
354 in time $\O(\sum_d (n_d + n_{d+1})) = \O(n)$.
355
356 The unification of sequences by bucket sorting will be useful in many
357 other situations, so we will state it as a~separate lemma:
358
359 \lemman{Unification of sequences}\id{uniflemma}%
360 Partitioning of a~collection of sequences $S_1,\ldots,S_n$, whose elements are
361 arbitrary pointers and symbols from a~finite alphabet, to equality classes can
362 be performed on the Pointer machine in time $\O(n + \sum_i \vert S_i \vert)$.
363
364 \rem
365 The first linear-time algorithm that partitions all subtrees to isomorphism equivalence
366 classes is probably due to Zemlayachenko \cite{zemlay:treeiso}, but it lacks many
367 details. Dinitz et al.~\cite{dinitz:treeiso} have recast this algorithm in modern
368 terminology and filled the gaps. Our algorithm is easier to formulate than those,
369 because it replaces the need for auxiliary data structures by more elaborate bucket
370 sorting.
371
372 \paran{Topological graph computations}%
373 Many graph algorithms are based on the idea of so called \df{micro/macro decomposition:}
374 We decompose a~graph to subgraphs on roughly~$k$ vertices and solve the problem
375 separately inside these ``micrographs'' and in the ``macrograph'' obtained by
376 contraction of the micrographs. If $k$~is small enough, many of the micrographs
377 are isomorphic, so we can compute the result only once for each isomorphism class
378 and recycle it for all micrographs in that class. On the other hand, the macrograph
379 is roughly $k$~times smaller than the original graph, so we can use a~less efficient
380 algorithm and it will still run in linear time with respect to the size of the original
381 graph.
382
383 This type of decomposition was traditionally used for trees, especially in the
384 algorithms for the Lowest Common Ancestor problem (cf.~Section \ref{verifysect}
385 and the survey paper \cite{alstrup:nca}) and for online maintenance of marked ancestors
386 (cf.~Alstrup et al.~\cite{alstrup:marked}). Let us take a~glimpse at what happens when
387 we set~$k$ to $1/4\cdot\log n$. There are at most $2^{2k} = \sqrt n$ non-isomorphic subtrees,
388 because each isomorphism class is uniquely determined by the sequence of $2k$~up/down steps
389 performed by depth-first search. Suppose that we are able to decompose the input and identify
390 the equivalence classes of microtrees in linear time, then solve the problem in time $\O(\poly(k))$ for
391 each microtree and finally in $\O(n'\log n')$ for the macrotree. When we put these pieces
392 together, we get an~algorithm with time complexity $\O(n + \sqrt{n}\cdot\poly(\log n) + n/\log n\cdot\log(n/\log n)) = \O(n)$
393 for the whole problem.
394
395 Decompositions are usually implemented on the RAM, because subgraphs can be easily
396 encoded in numbers, which can be then used to index arrays containing precomputed
397 results. As the previous algorithm for subtree isomorphism shows, indexing is not
398 required for identifying and equivalent microtrees and it can be replaced by bucket
399 sorting on the Pointer machine. Buchsbaum et al.~\cite{buchsbaum:verify} have extended
400 this technique to general graphs in form of topological graph computations.
401
402 \defn
403 A~\df{graph computation} takes a~\df{labeled undirected graph} as its input. The labels of both
404 vertices and edges can be arbitrary symbols drawn from a~finite alphabet. The output
405 of the computation is another labeling of the same graph. This time, the vertices and
406 edges can be labeled with not only symbols of the alphabet, but also with alphabet pointers to the vertices
407 and edges of the graph given as the input, and possibly also with pointers to outside objects.
408 A~graph computation is called \df{topological} iff it produces isomorphic
409 outputs for isomorphic inputs. The isomorphism has of course preserve not only
410 the structure of the graph, but also the labels in the obvious way.
411
412 \obs
413 The topological graph computations cover a~great variety of graph problems, ranging
414 from searching for spanning trees or Eulerian tours to the Traveling Salesman Problem.
415 The MST problem itself however does not lie in this class, because we do not have any means
416 of representing the edge weights, unless of course there is only a~fixed amount
417 of possible weights.
418
419 As in the case of tree decompositions, we would like to identify the equivalent subgraphs
420 and process only a~single instance from each equivalence class. The obstacle is that
421 graph isomorphism is known to be computationally hard (it is one of the few
422 problems that are neither known to lie in~$\rm P$ nor to be $\rm NP$-complete,
423 see Arvind and Kurur \cite{arvind:isomorph} for recent results on its complexity).
424 We will therefore manage with a~weaker form of equivalence, based on some sort
425 of graph encodings:
426
427 \defn
428 A~\df{canonical encoding} of a~given labeled graph represented by adjancency lists
429 can be obtained by running the depth-first search on the graph. When we enter
430 a~vertex, we assign an~identifier to it (again using a~yardstick to represent numbers)
431 and we append the label of this vertex to the encoding. Then we scan all back edges
432 going from this vertex and append the identifiers of their destinations, accompanied
433 by the edges' labels. Finally we append a~special terminator to mark the boundary
434 between the code of this vertex and its successor.
435
436 \obs
437 The canonical encoding is well defined in the sense that two non-isomorphic graphs
438 receive different encodings. Obviously, encodings of isomorphic graphs can differ,
439 depending on the order of vertices and also of the adjacency lists. A~graph
440 on~$n$ vertices with $m$~edges is assigned an~encoding of length at most $2n+2m$
441 (for each vertex, we record its label and a~single terminator; edges contribute
442 by identifiers and labels). This encoding can be constructed in linear time.
443 Let us use the encodings for unification of graphs:
444
445 \lemman{Unification of graphs}\id{uniflemma}%
446 A~collection~$\C$ of labeled graphs can be partitioned into classes which
447 share the same canonical encoding in time $\O(\Vert\C\Vert)$, where $\Vert\C\Vert$
448 is the total size of the collection, i.e., $\sum_{G\in\C} n(G) + m(G)$.
449
450 \para
451 When we have to perform a~topological computation on a~collection of graphs
452 on $k$~vertices, we first precompute its result for all possible canonical
453 encodings (the \df{generic graphs}) and then we use unification to match
454 the actual graphs to the generic graphs. This gives us the following theorem:
455
456 \thmn{Batched topological computations, Buchsbaum et al.~\cite{buchsbaum:verify}}\id{topothm}%
457 Suppose we have a~topological graph computation~$\cal T$ that can be performed in time
458 $T(k)$ for graphs on $k$~vertices. Then we can run~$\cal T$ on a~collection~$\C$
459 of graphs on~$k$ vertices in time $\O(\Vert\C\Vert + (k+s)^{k(2k+1)}\cdot (T(k)+k^2))$,
460 where~$s$ is the number of symbols used as vertex/edge labels.
461
462 \proof
463 A~graph on~$k$ vertices has less than~$k^2/2$ edges, so the canonical encodings of
464 all such graphs are shorter than $2k + 2k^2/2 = k(2k+1)$. Each element of the encoding
465 is either a~vertex identifier or a~symbol, so it can attain at most $k+s$ possible values.
466 We can therefore enumerate all possible encodings, convert them to the collection $\cal G$
467 of all generic graphs and run the computation on them in time $\O(\vert{\cal G}\vert \cdot T(k))
468 = \O((k+s)^{k(2k+1)}\cdot T(k))$.
469
470 We then use the Unification lemma (\ref{uniflemma}) on the union of the collections
471 $\C$ and~$\cal G$ to match the generic graphs with the equivalent graphs in~$\C$
472 in time $\O(\Vert\C\Vert + \Vert{\cal G}\Vert) = \O(\Vert\C\Vert + \vert{\cal G}\vert \cdot k^2)$.
473 Finally we create a~copy of the generic result for each of the actual graphs.
474 If the computation uses pointers to input vertices in its output, this involves
475 redirecting them to the actual input vertices, but we can do that by associating
476 the output vertices that refer to an~input vertex with the corresponding places
477 in the encoding of the input graph. The whole output can be generated in time
478 $\O(\Vert\C\Vert + \Vert{\cal G}\Vert)$.
479 \qed
480
481 \rem
482 The topological computations and the Unification lemma will play important
483 roles in Sections \ref{verifysect} and \ref{optalgsect}.
484
485 %--------------------------------------------------------------------------------
486
487 \section{Data structures on the RAM}
488 \id{ramdssect}
489
490 There is a~lot of data structures designed specifically for the RAM, taking
491 advantage of both indexing and arithmetics. In many cases, they surpass the known
492 lower bounds for the same problem on the~PM and they often achieve constant time
493 per operation, at least when either the magnitude of the values or the size of
494 the data structure are suitably bounded.
495
496 A~classical result of this type are the trees of van Emde Boas~\cite{boas:vebt},
497 which represent a~subset of the integers $\{0,\ldots,U-1\}$, allowing insertion,
498 deletion and order operations (minimum, maximum, successor etc.) in time $\O(\log\log U)$,
499 regardless of the size of the subset. If we replace the heap used in the Jarn\'\i{}k's
500 algorithm (\ref{jarnik}) by this structure, we immediately get an~algorithm
501 for finding the MST in integer-weighted graphs in time $\O(m\log\log w_{max})$,
502 where $w_{max}$ is the maximum weight. We will show later that it is even
503 possible to achieve linear time complexity for arbitrary integer weights.
504
505 A~real breakthrough has been made by Fredman and Willard, who introduced
506 the Fusion trees~\cite{fw:fusion} which again perform membership and predecessor
507 operation on a~set of $n$~integers, but this time with complexity $\O(\log_W n)$
508 per operation on a~Word-RAM with $W$-bit words. This of course assumes that
509 each element of the set fits in a~single word. As $W$ must at least~$\log n$,
510 the operations take $\O(\log n/\log\log n)$ and we are able to sort $n$~integers
511 in time~$o(n\log n)$. This was a~beginning of a~long sequence of faster and
512 faster sorting algorithms, culminating with the work by Thorup and Han.
513 They have improved the time complexity of integer sorting to $\O(n\log\log n)$ deterministically~\cite{han:detsort}
514 and expected $\O(n\sqrt{\log\log n})$ for randomized algorithms~\cite{hanthor:randsort},
515 both in linear space.
516
517 Despite the recent progress, the corner-stone of most RAM data structures
518 is still the representation of data structures by integers introduced by Fredman
519 and Willard. It will also form a~basis for the rest of this chapter.
520
521 \FIXME{Add more history.}
522
523 %--------------------------------------------------------------------------------
524
525 \section{Bits and vectors}\id{bitsect}
526
527 In this rather technical section, we will show how the RAM can be used as a~vector
528 computer to operate in parallel on multiple elements, as long as these elements
529 fit in a~single machine word. At the first sight this might seem useless, because we
530 cannot require large word sizes, but surprisingly often the elements are small
531 enough relative to the size of the algorithm's input and thus also relative to
532 the minimum possible word size. Also, as the following lemma shows, we can
533 easily emulate slightly longer words:
534
535 \lemman{Multiple-precision calculations}
536 Given a~RAM with $W$-bit words, we can emulate all calculation and control
537 instructions of a~RAM with word size $kW$ in time depending only on the~$k$.
538 (This is usually called \df{multiple-precision arithmetics.})
539
540 \proof
541 We split each word of the ``big'' machine to $W'$-bit blocks, where $W'=W/2$, and store these
542 blocks in $2k$ consecutive memory cells. Addition, subtraction, comparison and
543 bitwise logical operations can be performed block-by-block. Shifts by a~multiple
544 of~$W'$ are trivial, otherwise we can combine each block of the result from
545 shifted versions of two original blocks.
546 To multiply two numbers, we can use the elementary school algorithm using the $W'$-bit
547 blocks as digits in base $2^{W'}$ --- the product of any two blocks fits
548 in a~single word.
549
550 Division is harder, but Newton-Raphson iteration (see~\cite{ito:newrap})
551 converges to the quotient in a~constant number of iterations, each of them
552 involving $\O(1)$ multiple-precision additions and multiplications. A~good
553 starting approximation can be obtained by dividing the two most-significant
554 (non-zero) blocks of both numbers.
555
556 Another approach to division is using the improved elementary school algorithm as described
557 by Knuth in~\cite{knuth:seminalg}. It uses $\O(k^2)$ steps, but the steps involve
558 calculation of the most significant bit set in a~word. We will show below that it
559 can be done in constant time, but we have to be careful to avoid division instructions.
560 \qed
561
562 \notan{Bit strings}\id{bitnota}%
563 We will work with binary representations of natural numbers by strings over the
564 alphabet $\{\0,\1\}$: we will use $\(x)$ for the number~$x$ written in binary,
565 $\(x)_b$ for the same padded to exactly $b$ bits by adding leading zeroes,
566 and $x[k]$ for the value of the $k$-th bit of~$x$ (with a~numbering of bits such that $2^k[k]=1$).
567 The usual conventions for operations on strings will be utilized: When $s$
568 and~$t$ are strings, we write $st$ for their concatenation and
569 $s^k$ for the string~$s$ repeated $k$~times.
570 When the meaning is clear from the context,
571 we will use $x$ and $\(x)$ interchangeably to avoid outbreak of symbols.
572
573 \defn
574 The \df{bitwise encoding} of a~vector ${\bf x}=(x_0,\ldots,x_{d-1})$ of~$b$-bit numbers
575 is an~integer~$x$ such that $\(x)=\(x_{d-1})_b\0\(x_{d-2})_b\0\ldots\0\(x_0)_b$, i.e.,
576 $x = \sum_i 2^{(b+1)i}\cdot x_i$. (We have interspersed the elements with \df{separator bits.})
577
578 \notan{Vectors}\id{vecnota}%
579 We will use boldface letters for vectors and the same letters in normal type
580 for their encodings. The elements of a~vector~${\bf x}$ will be written as
581 $x_0,\ldots,x_{d-1}$.
582
583 \para
584 If we want to fit the whole vector in a~single word, the parameters $b$ and~$d$ must satisty
585 the condition $(b+1)d\le W$.
586 By using multiple-precision arithmetics, we can encode all vectors satisfying $bd=\O(W)$.
587 We will now describe how to translate simple vector manipulations to sequences of $\O(1)$ RAM operations
588 on their codes. As we are interested in asymptotic complexity only, we prefer clarity
589 of the algorithms over saving instructions. Among other things, we freely use calculations
590 on words of size $\O(bd)$, assuming that the Multiple-precision lemma comes to save us
591 when necessary.
592
593 \para
594 First of all, let us observe that we can use $\band$ and $\bor$ with suitable constants
595 to write zeroes or ones to an~arbitrary set of bit positions at once. These operations
596 are usually called \df{bit masking}. Also, any element of a~vector can be extracted or
597 replaced by a~different value in $\O(1)$ time by masking and shifts.
598
599 \newdimen\slotwd
600 \def\setslot#1{\setbox0=#1\slotwd=\wd0}
601 \def\slot#1{\hbox to \slotwd{\hfil #1\hfil}}
602 \def\[#1]{\slot{$#1$}}
603 \def\9{\rack{\0}{\hss$\cdot$\hss}}
604
605 \def\alik#1{%
606 \medskip
607 \halign{\hskip 0.15\hsize\hfil $ ##$&\hbox to 0.6\hsize{${}##$ \hss}\cr
608 #1}
609 \medskip
610 }
611
612 \algn{Operations on vectors with $d$~elements of $b$~bits each}\id{vecops}
613
614 \itemize\ibull
615
616 \:$\<Replicate>(\alpha)$ --- creates a~vector $(\alpha,\ldots,\alpha)$:
617
618 \alik{\<Replicate>(\alpha)=\alpha\cdot(\0^b\1)^d. \cr}
619
620 \:$\<Sum>(x)$ --- calculates the sum of the elements of~${\bf x}$, assuming that
621 the result fits in $b$~bits:
622
623 \alik{\<Sum>(x) = x \bmod \1^{b+1}. \cr}
624
625 This is correct because when we calculate modulo~$\1^{b+1}$, the number $2^{b+1}=\1\0^{b+1}$
626 is congruent to~1 and thus $x = \sum_i 2^{(b+1)i}\cdot x_i \equiv \sum_i 1^i\cdot x_i \equiv \sum_i x_i$.
627 As the result should fit in $b$~bits, the modulo makes no difference.
628
629 If we want to avoid division, we can use double-precision multiplication instead:
630
631 \setslot{\hbox{~$\0x_{d-1}$}}
632 \def\.{\[]}
633 \def\z{\[\0^b\1]}
634 \def\dd{\slot{$\cdots$}}
635 \def\vd{\slot{$\vdots$}}
636 \def\rule{\noalign{\medskip\nointerlineskip}$\hrulefill$\cr\noalign{\nointerlineskip\medskip}}
637
638 \alik{
639 \[\0x_{d-1}] \dd \[\0x_2] \[\0x_1] \[\0x_0] \cr
640 *~~ \z \dd \z\z\z \cr
641 \rule
642 \[x_{d-1}] \dd \[x_2] \[x_1] \[x_0] \cr
643 \[x_{d-1}] \[x_{d-2}] \dd \[x_1] \[x_0] \. \cr
644 \[x_{d-1}] \[x_{d-2}] \[x_{d-3}] \dd \[x_0] \. \. \cr
645 \vd\vd\vd\vd\.\.\.\cr
646 \[x_{d-1}] \dd \[x_2]\[x_1]\[x_0] \. \. \. \. \cr
647 \rule
648 \[r_{d-1}] \dd \[r_2] \[r_1] \[s_d] \dd \[s_3] \[s_2] \[s_1] \cr
649 }
650
651 This way, we even get the vector of all partial sums:
652 $s_k=\sum_{i=0}^{k-1}x_i$, $r_k=\sum_{i=k}^{d-1}x_i$.
653
654 \:$\<Cmp>(x,y)$ --- element-wise comparison of~vectors ${\bf x}$ and~${\bf y}$,
655 i.e., a~vector ${\bf z}$ such that $z_i=1$ if $x_i<y_i$ and $z_i=0$ otherwise.
656
657 We replace the separator zeroes in~$x$ by ones and subtract~$y$. These ones
658 change back to zeroes exactly at the positions where $x_i<y_i$ and they stop
659 carries from propagating, so the fields do not interact with each other:
660
661 \setslot{\vbox{\hbox{~$x_{d-1}$}\hbox{~$y_{d-1}$}}}
662 \def\9{\rack{\0}{\hss ?\hss}}
663 \alik{
664    \1 \[x_{d-1}] \1 \[x_{d-2}] \[\cdots] \1 \[x_1] \1 \[x_0]  \cr
665 -~ \0 \[y_{d-1}] \0 \[y_{d-2}] \[\cdots] \0 \[y_1] \0 \[y_0]  \cr
666 \rule
667    \9 \[\ldots]  \9 \[\ldots]  \[\cdots] \9 \[\ldots] \9 \[\ldots] \cr
668 }
669
670 It only remains to shift the separator bits to the right positions, negate them
671 and mask out all other bits.
672
673 \:$\<Rank>(x,\alpha)$ --- returns the number of elements of~${\bf x}$ which are less than~$\alpha$,
674 assuming that the result fits in~$b$ bits:
675
676 \alik{
677 \<Rank>(x,\alpha) = \<Sum>(\<Cmp>(x,\<Replicate>(\alpha))). \cr
678 }
679
680 \:$\<Insert>(x,\alpha)$ --- inserts~$\alpha$ into a~sorted vector $\bf x$:
681
682 We calculate the rank of~$\alpha$ in~$x$ first, then we insert~$\alpha$ as the $k$-th
683 field of~$\bf x$ using masking operations and shifts.
684
685 \algo
686 \:$k\=\<Rank>(x,\alpha)$.
687 \:$\ell\=x \band \1^{(b+1)(n-k-1)}\0^{(b+1)(k+1)}$. \cmt{``left'' part of the vector}
688 \:$r=x \band \1^{(b+1)k}$. \cmt{``right'' part}
689 \:Return $(\ell\shl (b+1)) \bor (\alpha\shl ((b+1)k)) \bor r$.
690 \endalgo
691
692 \:$\<Unpack>(\alpha)$ --- creates a~vector whose elements are the bits of~$\(\alpha)_d$.
693 In other words, inserts blocks~$\0^b$ between the bits of~$\alpha$. Assuming that $b\ge d$,
694 we can do it as follows:
695
696 \algo
697 \:$x\=\<Replicate>(\alpha)$.
698 \:$y\=(2^{b-1},2^{b-2},\ldots,2^0)$. \cmt{bitwise encoding of this vector}
699 \:$z\=x \band y$.
700 \:Return $\<Cmp>(z,y)$.
701 \endalgo
702
703 Let us observe that $z_i$ is either zero or equal to~$y_i$ depending on the value
704 of the $i$-th bit of the number~$\alpha$. Comparing it with~$y_i$ normalizes it
705 to either zero or one.
706
707 \:$\<Unpack>_\pi(\alpha)$ --- like \<Unpack>, but changes the order of the
708 bits according to a~fixed permutation~$\pi$: The $i$-th element of the
709 resulting vector is equal to~$\alpha[\pi(i)]$.
710
711 Implemented as above, but with a~mask $y=(2^{\pi(b-1)},\ldots,2^{\pi(0)})$.
712
713 \:$\<Pack>(x)$ --- the inverse of \<Unpack>: given a~vector of zeroes and ones,
714 it produces a~number whose bits are the elements of the vector (in other words,
715 it crosses out the $\0^b$ blocks).
716
717 We interpret the~$x$ as an~encoding of a~vector with elements one bit shorter
718 and we sum these elements. For example, when $n=4$ and~$b=4$:
719
720 \setslot{\hbox{$x_3$}}
721 \def\z{\[\0]}
722 \def\|{\hskip1pt\vrule height 10pt depth 4pt\hskip1pt}
723 \def\.{\hphantom{\|}}
724
725 \alik{
726 \|\z\.\z\.\z\.\z\.\[x_3]\|\z\.\z\.\z\.\z\.\[x_2]\|\z\.\z\.\z\.\z\[x_1]\|\z\.\z\.\z\.\z\.\[x_0]\|\cr
727 \|\z\.\z\.\z\.\z\|\[x_3]\.\z\.\z\.\z\|\z\.\[x_2]\.\z\.\z\|\z\.\z\[x_1]\.\z\|\z\.\z\.\z\.\[x_0]\|\cr
728 }
729
730 However, this ``reformatting'' does not produce a~correct encoding of a~vector,
731 because the separator zeroes are missing. For this reason, the implementation
732 of~\<Sum> using modulo does not work correctly (it produces $\0^b$ instead of $\1^b$).
733 We therefore use the technique based on multiplication instead, which does not need
734 the separators. (Alternatively, we can observe that $\1^b$ is the only case
735 affected, so we can handle it separately.)
736
737 \endlist
738
739 \para
740 We can use the aforementioned tricks to perform interesting operations on individual
741 numbers in constant time, too. Let us assume for a~while that we are
742 operating on $b$-bit numbers and the word size is at least~$b^2$.
743 This enables us to make use of intermediate vectors with $b$~elements
744 of $b$~bits each.
745
746 \algn{Integer operations in quadratic workspace}\id{lsbmsb}
747
748 \itemize\ibull
749
750 \:$\<Weight>(\alpha)$ --- computes the Hamming weight of~$\alpha$, i.e., the number of ones in~$\(\alpha)$.
751
752 We perform \<Unpack> and then \<Sum>.
753
754 \:$\<Permute>_\pi(\alpha)$ --- shuffles the bits of~$\alpha$ according
755 to a~fixed permutation~$\pi$.
756
757 We perform $\<Unpack>_\pi$ and \<Pack> back.
758
759 \:$\<LSB>(\alpha)$ --- finds the least significant bit of~$\alpha$,
760 i.e., the smallest~$i$ such that $\alpha[i]=1$.
761
762 By a~combination of subtraction with $\bxor$, we create a~number
763 that contains ones exactly at the position of $\<LSB>(\alpha)$ and below:
764
765 \alik{
766 \alpha&=                        \9\9\9\9\9\1\0\0\0\0\cr
767 \alpha-1&=                      \9\9\9\9\9\0\1\1\1\1\cr
768 \alpha\bxor(\alpha-1)&=         \0\9\9\9\0\1\1\1\1\1\cr
769 }
770
771 Then we calculate the \<Weight> of the result and subtract~1.
772
773 \:$\<MSB>(\alpha)$ --- finds the most significant bit of~$\alpha$ (the position
774 of the highest bit set).
775
776 Reverse the bits of the number~$\alpha$ first by calling \<Permute>, then apply \<LSB>
777 and subtract the result from~$b-1$.
778
779 \endlist
780
781 \rem
782 As noted by Brodnik~\cite{brodnik:lsb} and others, the space requirements of
783 the \<LSB> operation can be reduced to linear. We split the input to $\sqrt{b}$
784 blocks of $\sqrt{b}$ bits each. Then we determine which blocks are non-zero and
785 identify the lowest such block (this is a~\<LSB> of a~number whose bits
786 correspond to the blocks). Finally we calculate the \<LSB> of this block. In
787 both calls to \<LSB,> we have a $\sqrt{b}$-bit number in a~$b$-bit word, so we
788 can use the previous algorithm. The same trick of course works for finding the
789 \<MSB> as well.
790
791 The following algorithm shows the details.
792
793 \algn{LSB in linear workspace}
794
795 \algo\id{lsb}
796 \algin A~$w$-bit number~$\alpha$.
797 \:$b\=\lceil\sqrt{w}\,\rceil$. \cmt{size of a~block}
798 \:$\ell\=b$. \cmt{the number of blocks is the same}
799 \:$x\=(\alpha \band (\0\1^b)^\ell) \bor (\alpha \band (\1\0^b)^\ell)$.
800 \hfill\break
801 \cmt{encoding of a~vector~${\bf x}$ such that $x_i\ne 0$ iff the $i$-th block is non-zero}%
802 \foot{Why is this so complicated? It is tempting to take $\alpha$ itself as a~code of this vector,
803 but we unfortunately need the separator bits between elements, so we create them and
804 relocate the bits we have overwritten.}
805 \:$y\=\<Cmp>(0,x)$. \cmt{$y_i=1$ if the $i$-th block is non-zero, otherwise $y_0=0$}
806 \:$\beta\=\<Pack>(y)$. \cmt{each block compressed to a~single bit}
807 \:$p\=\<LSB>(\beta)$. \cmt{the index of the lowest non-zero block}
808 \:$\gamma\=(\alpha \shr bp) \band \1^b$. \cmt{the contents of that block}
809 \:$q\=\<LSB>(\gamma)$. \cmt{the lowest bit set there}
810 \algout $\<LSB>(\alpha) = bp+q$.
811 \endalgo
812
813 \rem
814 We have used a~plenty of constants that depend on the format of the vectors.
815 Either we can write non-uniform programs (see \ref{nonuniform}) and use native constants,
816 or we can observe that all such constants can be easily manufactured. For example,
817 $(\0^b\1)^d = \1^{(b+1)d} / \1^{b+1} = (2^{(b+1)d}-1)/(2^{b+1}-1)$. The only exceptions
818 are the~$w$ and~$b$ in the LSB algorithm \ref{lsb}, which we are unable to produce
819 in constant time. In practice we use the ``bit tricks'' as frequently called subroutines
820 in an~encompassing algorithm, so we usually can spend a~lot of time on the precalculation
821 of constants performed once during algorithm startup.
822
823 %--------------------------------------------------------------------------------
824
825 \section{Q-Heaps}\id{qheaps}%
826
827 We have shown how to perform non-trivial operations on a~set of values
828 in constant time, but so far only under the assumption that the number of these
829 values is small enough and that the values themselves are also small enough
830 (so that the whole set fits in $\O(1)$ machine words). Now we will show how to
831 lift the restriction on the magnitude of the values and still keep constant time
832 complexity. We will describe a~slightly simplified version of the Q-heaps developed by
833 Fredman and Willard in~\cite{fw:transdich}.
834
835 The Q-heap represents a~set of at most~$k$ word-sized integers, where $k\le W^{1/4}$
836 and $W$ is the word size of the machine. It will support insertion, deletion, finding
837 of minimum and maximum, and other operations described below, in constant time, provided that
838 we are willing to spend~$\O(2^{k^4})$ time on preprocessing.
839
840 The exponential-time preprocessing may sound alarming, but a~typical application uses
841 Q-heaps of size $k=\log^{1/4} N$, where $N$ is the size of the algorithm's input.
842 This guarantees that $k\le W^{1/4}$ and $\O(2^{k^4}) = \O(N)$. Let us however
843 remark that the whole construction is primarily of theoretical importance
844 and that the huge constants involved everywhere make these heaps useless
845 in practical algorithms. Many of the tricks used however prove themselves
846 useful even in real-life implementations.
847
848 Spending the time on reprocessing makes it possible to precompute tables for
849 almost arbitrary functions and then assume that they can be evaluated in
850 constant time:
851
852 \lemma\id{qhprecomp}%
853 When~$f$ is a~function computable in polynomial time, $\O(2^{k^4})$ time is enough
854 to precompute a~table of the values of~$f$ for all arguments whose size is $\O(k^3)$ bits.
855
856 \proof
857 There are $2^{\O(k^3)}$ possible combinations of arguments of the given size and for each of
858 them we spend $\poly(k)$ time on calculating the function. It remains
859 to observe that $2^{\O(k^3)}\cdot \poly(k) = \O(2^{k^4})$.
860 \qed
861
862 \para
863 We will first show an~auxiliary construction based on tries and then derive
864 the real definition of the Q-heap from it.
865
866 \nota
867 Let us introduce some notation first:
868 \itemize\ibull
869 \:$W$ --- the word size of the RAM,
870 \:$k = \O(W^{1/4})$ --- the limit on the size of the heap,
871 \:$n\le k$ --- the number of elements stored in the heap,
872 \:$X=\{x_1, \ldots, x_n\}$ --- the elements themselves: distinct $W$-bit numbers
873 indexed in a~way that $x_1 < \ldots < x_n$,
874 \:$g_i = \<MSB>(x_i \bxor x_{i+1})$ --- the position of the most significant bit in which $x_i$ and~$x_{i+1}$ differ,
875 \:$R_X(x)$ --- the rank of~$x$ in~$X$, that is the number of elements of~$X$ which are less than~$x$
876 (where $x$~itself need not be an~element of~$X$).\foot{We will dedicate the whole chapter \ref{rankchap} to the
877 study of various ranks.}
878 \endlist
879
880 \defn
881 A~\df{trie} for a~set of strings~$S$ over a~finite alphabet~$\Sigma$ is
882 a~rooted tree whose vertices are the prefixes of the strings in~$S$ and there
883 is an~edge going from a~prefix~$\alpha$ to a~prefix~$\beta$ iff $\beta$ can be
884 obtained from~$\alpha$ by appending a~single symbol of the alphabet. The edge
885 will be labeled with the particular symbol. We will also define the~\df{letter depth}
886 of a~vertex as the length of the corresponding prefix. We mark the vertices
887 which match a~string of~$S$.
888
889 A~\df{compressed trie} is obtained from the trie by removing the vertices of outdegree~1
890 except for the root and marked vertices.
891 Whereever is a~directed path whose internal vertices have outdegree~1 and they carry
892 no mark, we replace this path by a~single edge labeled with the contatenation
893 of the original edge's labels.
894
895 In both kinds of tries, we order the outgoing edges of every vertex by their labels
896 lexicographically.
897
898 \obs
899 In both tries, the root of the tree is the empty word and for every other vertex, the
900 corresponding prefix is equal to the concatenation of edge labels on the path
901 leading from the root to that vertex. The letter depth of the vertex is equal to
902 the total size of these labels. All leaves correspond to strings in~$S$, but so can
903 some internal vertices if there are two strings in~$S$ such that one is a~prefix
904 of the other.
905
906 Furthermore, the labels of all edges leaving a~common vertex are always
907 distinct and when we compress the trie, no two such labels have share their initial
908 symbols. This allows us to search in the trie efficiently: when looking for
909 a~string~$x$, we follow the path from the root and whenever we visit
910 an~internal vertex of letter depth~$d$, we test the $d$-th character of~$x$,
911 we follow the edge whose label starts with this character, and we check that the
912 rest of the label matches.
913
914 The compressed trie is also efficient in terms of space consumption --- it has
915 $\O(\vert S\vert)$ vertices (this can be easily shown by induction on~$\vert S\vert$)
916 and all edge labels can be represented in space linear in the sum of the
917 lengths of the strings in~$S$.
918
919 \defn
920 For our set~$X$, we define~$T$ as a~compressed trie for the set of binary
921 encodings of the numbers~$x_i$, padded to exactly $W$~bits, i.e., for $S = \{ \(x)_W ; x\in X \}$.
922
923 \obs
924 The trie~$T$ has several interesting properties. Since all words in~$S$ have the same
925 length, the leaves of the trie correspond to these exact words, that is to the numbers~$x_i$.
926 The inorder traversal of the trie enumerates the words of~$S$ in lexicographic order
927 and therefore also the~$x_i$'s in the order of their values. Between each
928 pair of leaves $x_i$ and~$x_{i+1}$ it visits an~internal vertex whose letter depth
929 is exactly~$W-1-g_i$.
930
931 \para
932 Let us now modify the algorithm for searching in the trie and make it compare
933 only the first symbols of the edges. In other words, we will test only the bits~$g_i$
934 which will be called \df{guides} (as they guide us through the tree). For $x\in
935 X$, the modified algorithm will still return the correct leaf. For all~$x$ outside~$X$
936 it will no longer fail and instead it will land on some leaf~$x_i$. At the
937 first sight the number~$x_i$ may seem unrelated, but we will show that it can be
938 used to determine the rank of~$x$ in~$X$, which will later form a~basis for all
939 Q-heap operations:
940
941 \lemma\id{qhdeterm}%
942 The rank $R_X(x)$ is uniquely determined by a~combination of:
943 \itemize\ibull
944 \:the trie~$T$,
945 \:the index~$i$ of the leaf found when searching for~$x$ in~$T$,
946 \:the relation ($<$, $=$, $>$) between $x$ and $x_i$,
947 \:the bit position $b=\<MSB>(x\bxor x_i)$ of the first disagreement between~$x$ and~$x_i$.
948 \endlist
949
950 \proof
951 If $x\in X$, we detect that from $x_i=x$ and the rank is obviously~$i-1$.
952 Let us assume that $x\not\in X$ and imagine that we follow the same path as when
953 searching for~$x$,
954 but this time we check the full edge labels. The position~$b$ is the first position
955 where~$\(x)$ disagrees with a~label. Before this point, all edges not taken by
956 the search were leading either to subtrees containing elements all smaller than~$x$
957 or all larger than~$x$ and the only values not known yet are those in the subtree
958 below the edge that we currently consider. Now if $x[b]=0$ (and therefore $x<x_i$),
959 all values in that subtree have $x_j[b]=1$ and thus they are larger than~$x$. In the other
960 case, $x[b]=1$ and $x_j[b]=0$, so they are smaller.
961 \qed
962
963 \para
964 The preceding lemma shows that the rank can be computed in polynomial time, but
965 unfortunately the variables on which it depends are too large for a~table to
966 be efficiently precomputed. We will carefully choose an~equivalent representation
967 of the trie which is compact enough.
968
969 \lemma\id{citree}%
970 The trie is uniquely determined by the order of the guides~$g_1,\ldots,g_{n-1}$.
971
972 \proof
973 We already know that the letter depths of the trie vertices are exactly
974 the numbers~$W-1-g_i$. The root of the trie must have the smallest of these
975 letter depths, i.e., it must correspond to the highest numbered bit. Let
976 us call this bit~$g_i$. This implies that the values $x_1,\ldots,x_i$
977 must lie in the left subtree of the root and $x_{i+1},\ldots,x_n$ in its
978 right subtree. Both subtrees can be then constructed recursively.\foot{This
979 construction is also known as the \df{cartesian tree} for the sequence
980 $g_1,\ldots,g_n$ and it is useful in many other algorithms as it can be
981 built in $\O(n)$ time. A~nice application on the Lowest Common Ancestor
982 and Range Minimum problems has been described by Bender et al.~in \cite{bender:lca}.}
983 \qed
984
985 \para
986 Unfortunately, the vector of the $g_i$'s is also too long (is has $k\log W$ bits
987 and we have no upper bound on~$W$ in terms of~$k$), so we will compress it even
988 further:
989
990 \nota\id{qhnota}%
991 \itemize\ibull
992 \:$B = \{g_1,\ldots,g_n\}$ --- the set of bit positions of all the guides, stored as a~sorted array,
993 \:$G : \{1,\ldots,n\} \rightarrow \{1,\ldots,n\}$ --- a~function mapping
994 the guides to their bit positions in~$B$: $g_i = B[G(i)]$,
995 \:$x[B]$ --- a~bit string containing the bits of~$x$ originally located
996 at the positions given by~$B$, i.e., the concatenation of bits $x[B[1]],
997 x[B[2]],\ldots, x[B[n]]$.
998 \endlist
999
1000 \obs\id{qhsetb}%
1001 The set~$B$ has $\O(k\log W)=\O(W)$ bits, so it can be stored in a~constant number
1002 of machine words in form of a~sorted vector. The function~$G$ can be also stored as a~vector
1003 of $\O(k\log k)$ bits. We can change a~single~$g_i$ in constant time using
1004 vector operations: First we delete the original value of~$g_i$ from~$B$ if it
1005 is not used anywhere else. Then we add the new value to~$B$ if it was not
1006 there yet and we write its position in~$B$ to~$G(i)$. Whenever we insert
1007 or delete a~value in~$B$, the values at the higher positions shift one position
1008 up or down and we have to update the pointers in~$G$. This can be fortunately
1009 accomplished by adding or subtracting a~result of vector comparison.
1010
1011 In this representation, we can reformulate our lemma on ranks as follows:
1012
1013 \lemma\id{qhrank}%
1014 The rank $R_X(x)$ can be computed in constant time from:
1015 \itemize\ibull
1016 \:the function~$G$,
1017 \:the values $x_1,\ldots,x_n$,
1018 \:the bit string~$x[B]$,
1019 \:$x$ itself.
1020 \endlist
1021
1022 \proof
1023 Let us prove that all ingredients of Lemma~\ref{qhdeterm} are either small
1024 enough or computable in constant time.
1025
1026 We know that the shape of the trie~$T$ is uniquely determined by the order of the $g_i$'s
1027 and therefore by the function~$G$ since the array~$B$ is sorted. The shape of
1028 the trie together with the bits in $x[B]$ determine the leaf~$x_i$ found when searching
1029 for~$x$ using only the guides. This can be computed in polynomial time and it
1030 depends on $\O(k\log k)$ bits of input, so according to Lemma~\ref{qhprecomp}
1031 we can look it up in a~precomputed table.
1032
1033 The relation between $x$ and~$x_i$ can be obtained directly as we know the~$x_i$.
1034 The bit position of the first disagreement can be calculated in constant time
1035 using the LSB/MSB algorithm (\ref{lsb}).
1036
1037 All these ingredients can be stored in $\O(k\log k)$ bits, so we may assume
1038 that the rank can be looked up in constant time as well.
1039 \qed
1040
1041 \para
1042 In the Q-heap we would like to store the set~$X$ as a~sorted array together
1043 with the corresponding trie, which will allow us to determine the position
1044 for a~newly inserted element in constant time. However, the set is too large
1045 to fit in a~vector and we cannot perform insertion on an~ordinary array in
1046 constant time. This can be worked around by keeping the set in an~unsorted
1047 array together with a~vector containing the permutation that sorts the array.
1048 We can then insert a~new element at an~arbitrary place in the array and just
1049 update the permutation to reflect the correct order.
1050
1051 We are now ready for the real definition of the Q-heap and for the description
1052 of the basic operations on it.
1053
1054 \defn
1055 A~\df{Q-heap} consists of:
1056 \itemize\ibull
1057 \:$k$, $n$ --- the capacity of the heap and the current number of elements (word-sized integers),
1058 \:$X$ --- the set of word-sized elements stored in the heap (an~array of words in an~arbitrary order),
1059 \:$\varrho$ --- a~permutation on~$\{1,\ldots,n\}$ such that $X[\varrho(1)] < \ldots < X[\varrho(n)]$
1060 (a~vector of $\O(n\log k)$ bits; we will write $x_i$ for $X[\varrho(i)]$),
1061 \:$B$ --- a~set of ``interesting'' bit positions
1062 (a~sorted vector of~$\O(n\log W)$ bits),
1063 \:$G$ --- the function that maps the guides to the bit positions in~$B$
1064 (a~vector of~$\O(n\log k)$ bits),
1065 \:precomputed tables of various functions.
1066 \endlist
1067
1068 \algn{Search in the Q-heap}\id{qhfirst}%
1069 \algo
1070 \algin A~Q-heap and an~integer~$x$ to search for.
1071 \:$i\=R_X(x)+1$, using Lemma~\ref{qhrank} to calculate the rank.
1072 \:If $i\le n$ return $x_i$, otherwise return {\sc undefined.}
1073 \algout The smallest element of the heap which is greater or equal to~$x$.
1074 \endalgo
1075
1076 \algn{Insertion to the Q-heap}
1077 \algo
1078 \algin A~Q-heap and an~integer~$x$ to insert.
1079 \:$i\=R_X(x)+1$, using Lemma~\ref{qhrank} to calculate the rank.
1080 \:If $x=x_i$, return immediately (the value is already present).
1081 \:Insert the new value to~$X$:
1082 \::$n\=n+1$.
1083 \::$X[n]\=x$.
1084 \::Insert~$n$ at the $i$-th position in the permutation~$\varrho$.
1085 \:Update the $g_j$'s:
1086 \::Move all~$g_j$ for $j\ge i$ one position up. \hfil\break
1087    This translates to insertion in the vector representing~$G$.
1088 \::Recalculate $g_{i-1}$ and~$g_i$ according to the definition.
1089    \hfil\break Update~$B$ and~$G$ as described in~\ref{qhsetb}.
1090 \algout The updated Q-heap.
1091 \endalgo
1092
1093 \algn{Deletion from the Q-heap}
1094 \algo
1095 \algin A~Q-heap and an~integer~$x$ to be deleted from it.
1096 \:$i\=R_X(x)+1$, using Lemma~\ref{qhrank} to calculate the rank.
1097 \:If $i>n$ or $x_i\ne x$, return immediately (the value is not in the heap).
1098 \:Delete the value from~$X$:
1099 \::$X[\varrho(i)]\=X[n]$.
1100 \::Find $j$ such that~$\varrho(j)=n$ and set $\varrho(j)\=\varrho(i)$.
1101 \::$n\=n-1$.
1102 \:Update the $g_j$'s like in the previous algorithm.
1103 \algout The updated Q-heap.
1104 \endalgo
1105
1106 \algn{Finding the $i$-th smallest element in the Q-heap}\id{qhlast}%
1107 \algo
1108 \algin A~Q-heap and an~index~$i$.
1109 \:If $i<1$ or $i>n$, return {\sc undefined.}
1110 \:Return~$x_i$.
1111 \algout The $i$-th smallest element in the heap.
1112 \endalgo
1113
1114 \para
1115 The heap algorithms we have just described have been built from primitives
1116 operating in constant time, with one notable exception: the extraction
1117 $x[B]$ of all bits of~$x$ at positions specified by the set~$B$. This cannot be done
1118 in~$\O(1)$ time on the Word-RAM, but we can implement it with ${\rm AC}^0$
1119 instructions as suggested by Andersson in \cite{andersson:fusion} or even
1120 with those ${\rm AC}^0$ instructions present on real processors (see Thorup
1121 \cite{thorup:aczero}). On the Word-RAM, we need to make use of the fact
1122 that the set~$B$ is not changing too much --- there are $\O(1)$ changes
1123 per Q-heap operation. As Fredman and Willard have shown, it is possible
1124 to maintain a~``decoder'', whose state is stored in $\O(1)$ machine words,
1125 and which helps us to extract $x[B]$ in a~constant number of operations:
1126
1127 \lemman{Extraction of bits}\id{qhxtract}%
1128 Under the assumptions on~$k$, $W$ and the preprocessing time as in the Q-heaps,\foot{%
1129 Actually, this is the only place where we need~$k$ to be as low as $W^{1/4}$.
1130 In the ${\rm AC}^0$ implementation, it is enough to ensure $k\log k\le W$.
1131 On the other hand, we need not care about the exponent because it can
1132 be arbitrarily increased using the Q-heap trees described below.}
1133 it is possible to maintain a~data structure for a~set~$B$ of bit positions,
1134 which allows~$x[B]$ to be extracted in $\O(1)$ time for an~arbitrary~$x$.
1135 When a~single element is inserted to~$B$ or deleted from~$B$, the structure
1136 can be updated in constant time, as long as $\vert B\vert \le k$.
1137
1138 \proof
1139 See Fredman and Willard \cite{fw:transdich}.
1140 \qed
1141
1142 \para
1143 This was the last missing bit of the mechanics of the Q-heaps. We are
1144 therefore ready to conclude this section by the following theorem
1145 and its consequences:
1146
1147 \thmn{Q-heaps, Fredman and Willard \cite{fw:transdich}}\id{qh}%
1148 Let $W$ and~$k$ be positive integers such that $k=\O(W^{1/4})$. Let~$Q$
1149 be a~Q-heap of at most $k$-elements of $W$~bits each. Then the Q-heap
1150 operations \ref{qhfirst} to \ref{qhlast} on~$Q$ (insertion, deletion,
1151 search for a~given value and search for the $i$-th smallest element)
1152 run in constant time on a~Word-RAM with word size~$W$, after spending
1153 time $\O(2^{k^4})$ on the same RAM on precomputing of tables.
1154
1155 \proof
1156 Every operation on the Q-heap can be performed in a~constant number of
1157 vector operations and calculations of ranks. The ranks are computed
1158 in $\O(1)$ steps involving again $\O(1)$ vector operations, binary
1159 logarithms and bit extraction. All these can be calculated in constant
1160 time using the results of section \ref{bitsect} and Lemma \ref{qhxtract}.
1161 \qed
1162
1163 \rem
1164 We can also use the Q-heaps as building blocks of more complex structures
1165 like Atomic heaps and AF-heaps (see once again \cite{fw:transdich}). We will
1166 show a~simpler, but useful construction, sometimes called the \df{Q-heap tree.}
1167 Suppose we have a~Q-heap of capacity~$k$ and a~parameter $d\in{\bb N}^+$. We
1168 can build a~balanced $k$-ary tree of depth~$d$ such that its leaves contain
1169 a~given set and every internal vertex keeps the minimum value in the subtree
1170 rooted in it, together with a~Q-heap containing the values in all its sons.
1171 This allows minimum to be extracted in constant time (it is placed in the root)
1172 and when any element is changed, it is sufficient to recalculate the values
1173 from the path from this element to the root, which takes $\O(d)$ Q-heap
1174 operations.
1175
1176 \corn{Q-heap trees}\id{qhtree}%
1177 For every positive integer~$r$ and $\delta>0$ there exists a~data structure
1178 capable of maintaining the minimum of a~set of at most~$r$ word-sized numbers
1179 under insertions and deletions. Each operation takes $\O(1)$ time on a~Word-RAM
1180 with word size $W=\Omega(r^{\delta})$, after spending time
1181 $\O(2^{r^\delta})$ on precomputing of tables.
1182
1183 \proof
1184 Choose $\delta' \le \delta$ such that $r^{\delta'} = \O(W^{1/4})$. Build
1185 a~Q-heap tree of depth $d=\lceil \delta/\delta'\rceil$ containing Q-heaps of
1186 size $k=r^{\delta'}$. \qed
1187
1188 \rem\id{qhtreerem}%
1189 When we have an~algorithm with input of size~$N$, the word size is at least~$\log N$
1190 and we can spend time $\O(N)$ on preprocessing, so we can choose $r=\log N$ and
1191 $\delta=1$ in the above corollary and get a~heap of size $\log N$ working in
1192 constant time per operation.
1193
1194 \endpart