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