- Email: [email protected]
Unranking and ranking spanning trees of a graph
JOURNAL
OF ALGONTW~S
Unranking
10,271-286 (1989)
and Ranking Spanning Trees of a Graph CHARLES J. COLBOURN
Department
of Combinatorics
and Optimization, University Ontario, Canada N2L 3Gl
of Waterloo,
Waterloo,
AND ROBERT P. J. DAY AND LOUIS D. NEL Department
of Computer
Science, University of Waterloo, Cattaah N2L 3Gl
Waterloo,
Ontario,
Received February 17,1988
The set S of spanning trees of an n-vertex graph G can be placed in one-to-one correspondence with the integers in the interval [l, s], where s = 1S I. We develop O(n3) unranking and ranking functions for the spanning trees of an arbitrary graph. The uua&ing function maps any interval [l, S] to the corresponding tree, while the ranking function maps a spaMing tree to the appropriate index in the interval. The unranking function provides an O(n3) method for generating a random spanning tree of a graph with uniform distribution. 0 1989 AC&& PWS.S, Inc.
1. INTRODUCTION One main interest in the field of combinatorics lies in counting the number of structures with a particular property and in the complexity of the enumeration of such structures. From the realm of graph theory, we find some counting problems which are efficiently solvable and others which are apparently intractable. An example of the first type is the number of spanning trees of a general graph, while examples of the second type include the Hamiltonian circuits and the spanning connected subgraphs of a graph. If, as in the first case, this enumeration can be done efficiently (that is, in time which is polynomial in the size of the graph), an associated problem can be considered-if we impose some initial total ordering on the set S of these structures, is there some way to . given one of the structures, determine its position in the predefined ordering (can we rank the element)? 271 01%-6774/89
$3.00
Copyright Q 1989 by Academic Press, Inc. All rights of repmduction in any foml mcrvcd.
212
COLBOURN,
DAY,
AND
NEL
. given an number in the interval [ 1, s], where s = ]S ], can we determine the corresponding element of S, again based on the given ordering (can we unsunk the number to the element of S)? As a trivial example, again from graph theory, given some initial order on the edges of a graph, one can quickly perform the rank and unrank operations between edges and numbers. Some structures with more complicated properties that still have efficient ranking and unranking operations include spanning trees of complete, undirected graphs [lo], binary trees (and, due to a simple l-l correspondence, ordered forests) on a given number of nodes [6], spanning rooted aborescences of directed, acyclic graphs [8], and k-subsets of a subset [4,7]. A major motivation for these kinds of functions is as a basis for random selection of the structures in question. Given some ordering of the set S for which an unranking function can be defined, the problem of random selection from S reduces quite naturally to random selection of an integer from the interval [l, s].
2. SPANNING
TREES OF UNDIRECTED
GRAPHS
The motivation for developing an unranking function for spanning trees is primarily for use in a random tree selection algorithm (for applications for such a procedure, refer to [l, 91). GuCnoche [3] reports an 0(n5) algorithm for this problem, and Kulkarni [8] improves this to O(n3m). This problem is not as elementary as it might seem at first glance. Given the graph shown in Fig. 1, if one were to iteratively construct a spanning tree by choosing edges at random and discard any that introduced cycles (essentially Kruskal’s algorithm for minimum-weight spanning trees with equally-weighted edges), the eight distinct spanning trees of the graph would be chosen with the probabilities shown in Fig. 2.
El 1
2
3
4
FIGURE
1
UNRANKING
SPANNING
FIGURE
2
273
TREES
z
P’E
8
Each of the spanning trees should, of course, be chosen with equal probability. The flaw in the above algorithm is that each edge should be chosen, not equally probably, but based on the fraction of spanning trees of the graph that include that particular edge. We first define some notation. Given some graph G = (V, E), with e E E, we can form the graph G - e = (V, E - {e}), read “G delete e.” Also, given the same graph G containing edge e = (x, JJ), we can simplify G by deleting (x, y) and identifying the two vertices. The resulting graph is denoted G * e, read “G contract e.” The second operation may introduce multiple edges; these are retained. If we define the function st(G) as the number of spanning trees of graph G, it is easy to see that st(G) = st(G - e) + st(G . e), an observation due originally to Feussner [2]. This implies that, for a graph G containing edge e, the fraction of spanning trees of G that contain edge e is st(G . e)/st(G). Thus, for any edge e, we choose with the above probability to add the edge to the spanning tree being constructed and to discard the edge otherwise. Depending on whether the edge was kept or discarded, we would continue this process of edge inspection with the graph G . e or G - e, respectively. It is not hard to prove that this technique constructs trees with a truly uniform distribution [9]. The final component of this algorithm is a procedure capable of calculating the fraction of spanning trees involving a particular edge. This can be done with the aid of an elegant method for counting spanning trees developed by Kirchoff [5]. For any undirected graph G = (V, E), with I VI = n and ]E( = rn, we define an n X n matrix D = di, j called the degree matrix, with entries as follows: d,,i = deg(i), di,j = -k,
the degree of node i k being the multiplicity
of edge (i, j) .
214
COLBOURN,
DAY,
AND
NEL
Shown below is the degree matrix D(G) for the graph in Fig. 1:
Deleting any i th row and i th column from the degree matrix leaves us with a reduced matrix D’(G), whose determinant is exactly the number of spanning trees of G. Thus, the fraction of spanning trees of G that do not contain a given edge e is simply det[D’(G
- e)]
det[ D’(G)]
.
(1)
Given that solving for the determinant of an n X n matrix using Gaussian elimination can be done in O(n3) time, an O(n3m) algorithm for the selection of a single random spanning tree is now straightforward. For each edge e in the graph G, calculate the ratio in (1) and either delete the edge from G if the edge is not chosen for the spanning tree or contract the edge (identify the endpoints of e) if the edge is selected-repeat until a spanning tree has been constructed. One property of the above algorithm is that the order in which the edges are processed is irrelevant. If we first label the edges of the graph 1 to m, we can define a natural order for the spanning trees based on the ascending lexicographic order of the set of edges in the tree (in effect, ranking the trees within the interval [l, s]). This particular ranking allows us to partition the entire sequence of trees into two classes-those trees that contain the first edge in the order and those that do not. If we are given some number K to unrank into the corresponding tree, then using the above view of the “number line” along which the trees lie, the inclusion status of this first edge should be based on which of the two sub-intervals must contain the spanning tree with rank K. Once this has been done, the trees in the resulting sub-interval can be partitioned into two further classes-those that contain the second edge and those that do not. Continuing to process the edges in the pre-defined order is equivalent to refining the interval containing the spanning tree of given rank until we have an interval of size one. The resulting O(n3m) algorithm is somewhat impractical, and hence a more efficient method is of interest, The basic algorithm employs O(m) determinant computations, each costing O(n3) time. To improve upon this, we could attempt to reduce the number of determinant computations, to reduce the cost of each determinant computation individually or to reduce the total cost of all determinant computation. Colbourn, Debroni, and
UNRANKING
SPANNING
TREES
275
Myrvold [l] develop a scheme which reduces the number of determinant computations when many random spanning trees are required, but it affords no improvement when a single random spanning tree is required. It appears that an inclusion/exclusion method may be required to examine essentially all of the edges to form a single spanning tree. To accelerate the computation of each determinant, one could in principle use improved algorithms for matrix multiplication; however, for the applications of interest, the traditional O(n3) algorithm is more practical. Finally, what remains is the option of reducing the overall cost of computing the O(m) determinants, by avoiding unnecessary recalculation; we pursue this strategy in the following sections.
3.
ORDERING
THE SPANNING
TREES
Part of the basis for the O(n3) unranking algorithm for spanning trees is a novel technique for ordering the edges of a graph that imposes the same natural ordering on the set of trees as described in the previous section. This ordering scheme is later used to decide the inclusion status, not of single edges, but of sets of edges at a time. The edge-ordering scheme works as follows. For simplicity, we assume throughout that the number of nodes is a power of 2. Partition the n nodes in the graph G = (V, E) into two equal-sized classes, PO and Pr, as shown in Fig. 3 (the basis for the partition is explained shortly). All edges with both endpoints in class PO are designated as “type 1” edges, while edges with both endpoints in Pr are “type 2” edges. Any edge straddling these two partitions is unaccounted for-such edges are classi-
n
Ro PO
Pl
FIGURES
276
COLBOURN,
DAY,
FIGURE
AND
NEL
4
fied by further subdividing both partitions into equal-sized subclasses Poe, PoI, Pro, and PII and classifying these edges as type 3, 4, 5 or 6 (see Fig. 4). Every edge in the graph has now been labelled with an edge type and, more importantly, each set of edges of any one type has all of its endpoints in a set of nodes of size n/2 (for instance, all type 1 edges have their endpoints in PO, type 3 edges have their endpoints in subclasses PO0 and PI,,, and so on). Our order on the edges is now defined by stating that edges of type 1 precede edges of type 2, which precede edges of type 3, and so on. Since this type classification does not, in general, impose a total ordering on the edges, the edge labelling must be further refined. This is done by selecting all edges of a particular type, collecting the n/2 nodes that are the endpoints of all edges of that type, and performing the partitioning operation on this subset of nodes (as an example, for edges of type 4, we choose the nodes from subclasses PO0 and PII to be partitioned further). When done for each edge type, this operation gives rise to edges of type 11,12,13 , . . . ,21,22, and so on. This partitioning and refinement of labelling is done recursively until every edge in the graph has a unique label and the edges are now totally ordered lexicographically by their label. The remaining issue is how the actual partitioning is done. Consider the binary representation of each node number, and partition based on the bit value. At the first step, nodes whose binary number representations have a high-order bit of zero are placed in PO, while the remaining nodes make up PI; partitioning these classes is done based on the second bit. To complete the partitioning recursively, a simple rule suffices: for any set of nodes, the subclasses Pm, Pal, Plo, and PII are defined by partitioning on the first
IJNRANKING
277
SPANNING TREES
two high-order bits that do not have an identical value across all nodes in the set to be partitioned. There are two different categories of edge type at the top level: types 1 and 2, whose endpoints all lie in some set of size n/2 and may have, as endpoints, any pair of nodes in the set, and types 3, 4, 5, and 6, whose endpoints also lie in some set of size n/2 but must have, as endpoints, one node in one half of the set and one node in the other half. Because of this, edges of type 1 or 2 may be further refined into all of the subtypes 1 through 6, but edges of type 3 through 6 can only be further refined into subtypes 3 through 6 (that is, there are no edges of type 31,32,41,. . . ). To unrank a number into a spanning tree with this edge ordering, the simple strategy given earlier suffices. However, the edge ordering introduces an important feature: in processing type i edges (i = 1, . . . ,6), the contractions and deletions performed involve only n/2 nodes. In fact, for any type t of length k, edges of type t involve at most n/2k nodes; all contractions and deletions performed on edges of this type affect only this (small) set of nodes. Hence our edge ordering isolates changes in the degree matrix to a small submatrix during the processing of edges of any one type. 4. SHORTCUTS FOR DETERMINANT CALCULATION REDUCED MATRIX
OF THE
The second component of the spanning tree selection algorithm is based on the observation that one does not always need to recalculate the determinant of the new reduced matrix from scratch during the computation. To motivate this, consider the deletion of an edge from the graph. Recall the reduced matrix (with row n and column n deleted) for the graph in Fig. 1 and consider its upper-triangular form after Gaussian elimination, as shown above:
Deleting edge (3,4) from the original graph would result in the new reduced matrix and its eliminated form below:
The change in both matrices involves subtracting one from the element in position (3,3), which suggests (at least for this very restricted case) the
278
COLBOURN,
DAY,
AND
NEL
\
: FIG. 5.
O[ 1
Partially
D*
eliminated
I
D matrix.
maxim, “Deletion followed by elimination is equivalent to elimination followed by deletion.” We formalize this observation and establish similar results for contraction and row/column permutation. Define matrix operations to be performed on the degree matrix D that correspond to graph operations on G; these matrix operations, and their graphical equivalents, are: (a) Bij( D): decrement Di, j and Dj, j by -D;, j, and set both Di, j and Dj, i to zero, corresponding to deleting all edges between vertices i and j in G. (b) Vij( D): add row i to row j, then column i to column j, and delete row and column i, corresponding to contracting edge (i, j) in G and labelling the new vertex j. (c) Pij( D): interchange rows i and j and columns i and j; that is, swap the labels of vertices i and j in G. (d) c?~(D): perform Gaussian elimination of columns 1,. . . , k of D, leaving a lower right remaining submatrix D* as shown in Fig. 5. This operation has no graphical equivalent. The above matrix operations are defined in terms of the degree matrix D of a graph G-they work equally well for a reduced matrix D’, where any operations to be performed on a row or column that has already been removed are ignored. LEMMA 4.1. Let D be any n - 1 X n - 1 reduced matrix for G (with vertex set { 1, . . . , n}). For vertices i, j of G and k -C i, j,
(9 ‘k;c(gij(D)) = 44.j(dk,(D))* (‘) gkk(gijtD)) = sij(gkk(D))* (iii) ‘,c;c(‘ip,j(D)) = ‘ip,j(dk(D))* Proof: The lemma is proved for a general fore holds for the restricted case of a reduced lemma holds for the case k = 1, it holds for remainder of the proof therefore deals strictly
square matrix D and therematrix for any graph. If the general k by induction; the with the case k = 1.
UNRANKING
SPANNING
279
TREES
For some n X n matrix D, it is easy to see that the three commutativity conditions hold for elements in the first row and first column of D. What remains is to prove equality under all three operations for the remaining n-lXn1 submatrix of D. (i) 8t(Sij(
D)). Easily derived.
D)) = Pij(4’l(
(ii) &l(gij(D)) = Bij(gll(D))* Since the result of the operation gij( D) is to effect the assignment dij = dij + aij for some 2-by-2 submatrix of D and some values a, it suffices to observe that (dij + aij)
- $-dlj
= dij - >dlj
11
+ aij.
11
(in) eO,(%Tij(D)) = cipj(8r( D)). We need only consider the effect on the elements in row and column j here. If we extend the domain of our operators to individual entries dij of D, then for some element dkj, k # j, we have &",('iP,j(d,j))
=
qj(
=
'ltdkj))
(d/ci
+
dkj)
-
and equality holds. A similar technique can be used to prove equality for elements djk, k + j, while proving equality for djj is only slightly more complicated. This lemma is fundamental in the timing analysis; it guarantees that in accomplishing a deletion, contraction, or permutation, we can avoid the inspection of D and consider only c?~(D), thereby avoiding repeating a large portion of the Gaussian elimination.
5.
A
GENERAL
DESCRIPTION AND
A SIMPLE
OF
THE
ALGORITHM
EXAMPLE
We now have the building blocks for an O(n3) spanning tree unranking algorithm. The edge ordering technique from Section 3 allows us to order the set of all spanning trees in the interval [l, s], while the actual unranking procedure resembles the partition refinement algorithm from Section 2-given some interval that contains the tree with the given rank, we use some as-yet-unspecified criteria to partition the interval into two sub-intervals and continue inspecting the sub-interval that contains the tree with that rank. Eventually, we refine to an interval of size one.
280
COLBOURN,
DAY,
AND
NEL
The partitioning described above is accomplished by, firstly, deciding on the inclusion/exclusion status of all edges of type 1, followed by all edges of type 2, and so on. To see how this corresponds to a partitioning refinement of the interval, observe that the trees in the interval, because of the way they have been ordered, can be partitioned based on whether or not they contain any type 1 edges, as shown below: 1 ~-----------------------I------------------------~ trees with type 1 edges
trees with no type 1 edges
It is a simple matter to remove all type 1 edges from the graph, calculate the determinant of the resulting reduced matrix, and decide which partition contains the tree of given rank. If there are no type 1 edges, the right sub-interval can be regarded as 1 ~------------------------,------------------------~ trees with type 2 edges
trees with no type 2 edges
and the original process repeated for type 2 edges. If there are type 1 edges, then we can regard the left sub-interval as follows: 1 ~----------------+----------~~~~ trees with type 11 edges
no type 11 edges
In either case, we have reduced our search space (the size of the interval). As an example, consider again the graph shown in Fig. 1, shown below with an edge labelling and all eight of its spanning trees, all of which also have the appropriate edge labelling and have been placed along the number line accordingly. (Since there is only one edge of each type, we refer to the edge of type 1 as edge 1, and so on.)
3
Ii3 1
2
1
3
6
flVlLl~K7l7Z 123
4
1
125
126
136
156
235
236
356
UNRANKING
SPANNING
TREES
281
Assume that we are searching for the spanning tree with rank seven, that is, the tree with edges 2, 3, and 6. Initially, we decide on the status of all type 1 edges (all one of them) and determine that, since there are five trees that involve edge 1, we are interested in the other sub-interval, that is, interval [6,8]. Edge 1 is deleted and our graph is now
In deciding the status of type-2 edges (again, all one of them), forced to keep edge 2 since this corresponds to the interval [6,7] number line, which contains the rank we are interested in. This contracted and the graph corresponding to the current situation
we are on our edge is is now
So far, we have decided to exclude edge 1 and include edge 2, and we are examining the interval [6,7] on the number line. We note that, in the current graph, there are two edges that have identical endpoints but different edge types. This occurs because we require an edge to retain its original edge 1abeIling throughout the life of the algorithm, regardless of any contractions that may occur.
282
COLBOURN,
DAY,
AND
NEL
The status of edge 3 is obvious-it must be included to keep the graph connected, so the edge is contracted to produce the graph
There are no edges of type 4 and, based on the resulting sub-interval, we decide to delete edge 5 (corresponding to the interval [6,6]), which leaves us with the graph
corresponding to the interval [7,7]. This forces us to keep the last edge, edge 6, and our record of contractions tells us that the spanning tree constructed consists of edges 2, 3, and 6, this being the tree with rank seven. In the (likely) event that the graph has more than one edge of type 1, we are forced to determine the inclusion/exclusion status of the entire set of type-l edges before moving on to type-2 edges. If we decide that at least one type-l edge must be included, our next question is whether we must include any edges of type 11, followed by edges of type 12, and so on. We would like to do this without recalculating the determinant of the new reduced matrix from scratch. Using Lemma 4.1(i), we can begin by permuting the rows and columns of the degree matrix to move the endpoints of all type-l edges (at most n/2 of them) to the lower right-hand submatrix, in effect moving all type-l edges into this submatrix. We are now free to perform Gaussian elimination down as far as this submatrix, knowing that any decisions about type-l edges do not affect this elimination, see Fig. 6. Once we reach this stage, determining
UNRANKING
SPANNING TREES
283
FIG. 6. Partially eliminated reduced matrix
the status of the type-l edges resembles the initial problem-six consider, but each a problem of size n/2-by-n/2 only.
subcases to
6. THE RECURSIVEUNRANKING ALGORITHM All of the basic ingredients are now present to describe the O(n3) algorithm. The method operates recursively, and the graph changes during the course of the computation. Each recursive invocation is to determine, for some type t of length k, which of the edges of type t appear in the tree (and are to be contracted) and which do not appear in the tree (and are to be deleted). At entry to the procedure, endnodes of edges of type t all appear in the lower right n k X nk submatrix (nk I n/2k). The degree matrix has had partial Gaussian elimination performed on it to eliminate subdiagonal entries in the first n - nk columns; hence at this level, the degree matrix is specified by a product p of the first n - nk diagonal entries, together with the remaining nk x nk matrix D. The procedure also takes as input the start s and end f of the interval of spanning trees under consideration, along with the rank r of the desired spanning tree. Gn return, the procedure must report which edges of type t are to be deleted and which are to be contracted. The algorithm ensures that whenever the procedure is invoked, at least one edge of type t must be contracted. The procedure itself operates as follows. For each type t’ = tl, t2,. . . , t6 (of length k + 1) in turn, first locate all edges of type t’. Form a permuted copy C of D so as to place endnodes of edges of type t ’ in the lower right submatrix of C (nk+l I n/2 k+l). Perform partial Gaussian nk+lX nk+l elimination on C to eliminate the subdiagonal in the first nk - nk+ 1 columns; retain p’, the product of the first nk - n k+l diagonal elements, and D’, the remaining nk+l X nk+l matrix. Having produced D’, form D” by deleting all edges of type t’ in D’ (that is, for each edge e = (i, j) of type t’, apply gij). Compute d = det( D”).
284
COLBOURN,
DAY,
AND
NEL
Now p . p’ . d is the number of spanning trees not containing any type t’ edges. Hence if r, the rank desired, is at least f - p . p’ - d, no edge of type t’ appears in the spanning tree. In this case, we delete all edges of type t’ in D and advance to the next edge type. This has the effect of changing the start s of the interval to f - p - p’ - d. Otherwise at least one edge of type t’ must be included, and hence the end f of the interval of interest is f - p * p’ . d. If there is only one edge of type t’, we contract it in D and proceed to edges of the next type (this has the effect of reducing nk by one as well). The most complicated case arises when some edge(s) of type t’ are to be contracted, and there is more than one type t’ edge. In this case, we recursively invoke the procedure with the current start, end, and rank, and the matrix specified by diagonal product p . p’ and remaining matrix D’. A list of deletions and contractions for edges of type t’ is returned. We apply each of the deletions and contractions in D, updating s and f accordingly. We then proceed to edges of the next type. When all types tl, t2,. . . , t6 have been processed, we return the list of all deletions and contractions for all type t edges (we must ensure that the deletions and contractions are modified to undo the effect of any permutations applied within this invocation). At the top level, the procedure is simply processing edges of types 1,2,..., 6 in turn; upon completion, the list of edges to be contracted which it returns is precisely the edge set of the spanning tree with rank r. THEOREM 6.1. The above algorithm to select a random spanning tree of an undirected graph of order n runs in O(n3) time.
Proof If the initial partial Gaussian elimination requires time O(n3), and we reconstruct the degree matrix entirely after each edge type decision, the recurrence relation for the amount of work for an n-by-n matrix is given by T(n)
5 6n3 + 6T(n/2),
and hence T(n) is O(n3). The relation is actually rather conservative since, as mentioned before, not all edge types give rise to all six subcases. In addition, it may be the case that edges of two different types actually have identical endpoints due to a series of contractions, so that no reorganization of the submatrix is necessary. There remains the problem of bookkeeping that allows us to perform the matrix permutations described above. These operations employ two operations we must be able to perform efficiently; they are (1) Given some edge type, generate the list of edges of this type. (2) Given the edges in a particular type, determine the set of nodes that are the endnodes of the edges.
UNRANKING
SPANNING TREES
285
These two operations do not necessarily constitute the same problem due to the contractions that occur as the algorithm runs. Since the contraction operation involves identifying two nodes in the graph, we should be able to, given any node label, determine its current label. The first of these operations can be implemented efficiently by a number of common data structures, either hierarchical or array-based, while the relabelling of vertices due to contractions can be handled by a reasonable implementation of a UNION-FIND abstract data type. Since by Lemma 4.1, the operations performed on the partially eliminated submatrices yield the same result as in the original matrix, we can also state THEOREM 6.2. The algorithm the corresponding spanning tree.
given above correctly unranks
7. CONCLUDING
any value to
REMARKS
An O(n3) ranking algorithm for producing the rank of a given spanning tree is now straightforward; rather than employing information about the rank to determine the inclusion/exclusion status of edges, we simply use the known status of edges of each type to refine the interval progressively until the desired rank is found. A number of practical issues arise in implementing this algorithm. Some of the determinant computations now performed can be avoided by checking for bridges (cut-edges); if a bridge is formed during the computation, it can be immediately contracted since it must appear in the spanning tree. Similarly, if a loop is formed it can be immediately deleted. These do not affect the O(n3) running time, but afford some practical improvement. Another improvement in practice can be realized by observing that our tree ordering depends on the edge ordering, which in turn depends on the node ordering, which is arbitrary. While any node ordering leads to an O(n3) worst case, the variation in actual computation time can be of practical significance. A useful heuristic is to order the nodes so that at the top level(s) in the partitioning, as few different edge types as possible are realized. This has the effect of minimizing the number of subproblems of size n/2, n/4 and so on. Finally, a practical improvement can be obtained by “batching”; if many random spanning trees are required, for example, we can sort the ranks r,, . . . , r, of the desired trees. In unranking ri and ri+l, much of the computation is duplicated typically; avoiding recomputation of the common part of the two computations yields a dramatic savings in practice.
286
COLBOURN,
DAY,
AND
NEL
It is perhaps remarkable that selecting a random spanning tree via many determinant computations is only a constant factor more expensive than simply performing one determinant computation. The sixfold divide-andconquer used here gives a powerful technique for partitioning graph problems. We expect that this technique will prove useful in other problems involving random selection of combinatorial objects. ACKNOWLEDGMENTS
Thanks to Wendy Myrvold and Bill Pulleyblank for helpful discussions about this research. The research of the first author is supported by NSERC Canada under Grant AO579.
REFERENCES
1. C. J. COLBOURN, B. M. DEBRONI, AND W. J. MYRVOLD, Estimating the coefficients of the reliability polynomial, in “Proceedings, Seventeenth Manitoba Conference on Numerical Mathematics and Computing, 1987, to appear. 2. W. FEUSSNER, iiber Stromberzeigung in nekformigen Leitem, Ann. Physik 9 (1902), 1304-1329; 15 (1904), 385-394. 3. A. GUENOCEIB, Random spanning tree, J. Algorithms 4 (1983), 214-220. 4. G. KATONA, A theorem of finite sets, “Theory of Graphs,” (P. Erd& and G. Katona, Eds.), pp. 187-207, Akademia Kiado, Budapest, 1966. 5. G. KIRCHOFF, iiber die Aufl&mg der Gleichungen auf welche man sei der Untersuchung der Linearen Verteihmg Galvanischer Strome Gefuhrt wind, Poppendorg’s Ann. Phys. Chem. 72 (1847), 497-508 (On the solution of equations obtained from the investigation of the linear distribution of galvanic currents, IRE Trans. Circuit Theory CT-5 (1958), 4-7. 6. G. D. KNOIT, A numbering system for binary trees, Comm. ACM 28 (1977), 113-115. 7. J. B. KRUSKAL, The number of simplices in a complex, in “Mathematical Optimization Techniques” (R. Bellman, Pd.), pp. 251-278, Univ. of California Press, Berkeley, 1963. 8. V. G. KULKARNI, “Generating Random Spanning Trees in a Graph,” Tech. Report UNC/ORSA/TR-85/l, University of North Carolina, 1985. 9. L. D. NEL AND C. J. COLBOURN, Combining Monte Carlo estimates and bounds for network reliability, submitted for publication, 1988. 10. H. hijFER, Neuer Reweis e-inesSatzes Uber Permutationen, Arch. Muth. Phys. 27 (1918), 742-744.
OF ALGONTW~S
Unranking
10,271-286 (1989)
and Ranking Spanning Trees of a Graph CHARLES J. COLBOURN
Department
of Combinatorics
and Optimization, University Ontario, Canada N2L 3Gl
of Waterloo,
Waterloo,
AND ROBERT P. J. DAY AND LOUIS D. NEL Department
of Computer
Science, University of Waterloo, Cattaah N2L 3Gl
Waterloo,
Ontario,
Received February 17,1988
The set S of spanning trees of an n-vertex graph G can be placed in one-to-one correspondence with the integers in the interval [l, s], where s = 1S I. We develop O(n3) unranking and ranking functions for the spanning trees of an arbitrary graph. The uua&ing function maps any interval [l, S] to the corresponding tree, while the ranking function maps a spaMing tree to the appropriate index in the interval. The unranking function provides an O(n3) method for generating a random spanning tree of a graph with uniform distribution. 0 1989 AC&& PWS.S, Inc.
1. INTRODUCTION One main interest in the field of combinatorics lies in counting the number of structures with a particular property and in the complexity of the enumeration of such structures. From the realm of graph theory, we find some counting problems which are efficiently solvable and others which are apparently intractable. An example of the first type is the number of spanning trees of a general graph, while examples of the second type include the Hamiltonian circuits and the spanning connected subgraphs of a graph. If, as in the first case, this enumeration can be done efficiently (that is, in time which is polynomial in the size of the graph), an associated problem can be considered-if we impose some initial total ordering on the set S of these structures, is there some way to . given one of the structures, determine its position in the predefined ordering (can we rank the element)? 271 01%-6774/89
$3.00
Copyright Q 1989 by Academic Press, Inc. All rights of repmduction in any foml mcrvcd.
212
COLBOURN,
DAY,
AND
NEL
. given an number in the interval [ 1, s], where s = ]S ], can we determine the corresponding element of S, again based on the given ordering (can we unsunk the number to the element of S)? As a trivial example, again from graph theory, given some initial order on the edges of a graph, one can quickly perform the rank and unrank operations between edges and numbers. Some structures with more complicated properties that still have efficient ranking and unranking operations include spanning trees of complete, undirected graphs [lo], binary trees (and, due to a simple l-l correspondence, ordered forests) on a given number of nodes [6], spanning rooted aborescences of directed, acyclic graphs [8], and k-subsets of a subset [4,7]. A major motivation for these kinds of functions is as a basis for random selection of the structures in question. Given some ordering of the set S for which an unranking function can be defined, the problem of random selection from S reduces quite naturally to random selection of an integer from the interval [l, s].
2. SPANNING
TREES OF UNDIRECTED
GRAPHS
The motivation for developing an unranking function for spanning trees is primarily for use in a random tree selection algorithm (for applications for such a procedure, refer to [l, 91). GuCnoche [3] reports an 0(n5) algorithm for this problem, and Kulkarni [8] improves this to O(n3m). This problem is not as elementary as it might seem at first glance. Given the graph shown in Fig. 1, if one were to iteratively construct a spanning tree by choosing edges at random and discard any that introduced cycles (essentially Kruskal’s algorithm for minimum-weight spanning trees with equally-weighted edges), the eight distinct spanning trees of the graph would be chosen with the probabilities shown in Fig. 2.
El 1
2
3
4
FIGURE
1
UNRANKING
SPANNING
FIGURE
2
273
TREES
z
P’E
8
Each of the spanning trees should, of course, be chosen with equal probability. The flaw in the above algorithm is that each edge should be chosen, not equally probably, but based on the fraction of spanning trees of the graph that include that particular edge. We first define some notation. Given some graph G = (V, E), with e E E, we can form the graph G - e = (V, E - {e}), read “G delete e.” Also, given the same graph G containing edge e = (x, JJ), we can simplify G by deleting (x, y) and identifying the two vertices. The resulting graph is denoted G * e, read “G contract e.” The second operation may introduce multiple edges; these are retained. If we define the function st(G) as the number of spanning trees of graph G, it is easy to see that st(G) = st(G - e) + st(G . e), an observation due originally to Feussner [2]. This implies that, for a graph G containing edge e, the fraction of spanning trees of G that contain edge e is st(G . e)/st(G). Thus, for any edge e, we choose with the above probability to add the edge to the spanning tree being constructed and to discard the edge otherwise. Depending on whether the edge was kept or discarded, we would continue this process of edge inspection with the graph G . e or G - e, respectively. It is not hard to prove that this technique constructs trees with a truly uniform distribution [9]. The final component of this algorithm is a procedure capable of calculating the fraction of spanning trees involving a particular edge. This can be done with the aid of an elegant method for counting spanning trees developed by Kirchoff [5]. For any undirected graph G = (V, E), with I VI = n and ]E( = rn, we define an n X n matrix D = di, j called the degree matrix, with entries as follows: d,,i = deg(i), di,j = -k,
the degree of node i k being the multiplicity
of edge (i, j) .
214
COLBOURN,
DAY,
AND
NEL
Shown below is the degree matrix D(G) for the graph in Fig. 1:
Deleting any i th row and i th column from the degree matrix leaves us with a reduced matrix D’(G), whose determinant is exactly the number of spanning trees of G. Thus, the fraction of spanning trees of G that do not contain a given edge e is simply det[D’(G
- e)]
det[ D’(G)]
.
(1)
Given that solving for the determinant of an n X n matrix using Gaussian elimination can be done in O(n3) time, an O(n3m) algorithm for the selection of a single random spanning tree is now straightforward. For each edge e in the graph G, calculate the ratio in (1) and either delete the edge from G if the edge is not chosen for the spanning tree or contract the edge (identify the endpoints of e) if the edge is selected-repeat until a spanning tree has been constructed. One property of the above algorithm is that the order in which the edges are processed is irrelevant. If we first label the edges of the graph 1 to m, we can define a natural order for the spanning trees based on the ascending lexicographic order of the set of edges in the tree (in effect, ranking the trees within the interval [l, s]). This particular ranking allows us to partition the entire sequence of trees into two classes-those trees that contain the first edge in the order and those that do not. If we are given some number K to unrank into the corresponding tree, then using the above view of the “number line” along which the trees lie, the inclusion status of this first edge should be based on which of the two sub-intervals must contain the spanning tree with rank K. Once this has been done, the trees in the resulting sub-interval can be partitioned into two further classes-those that contain the second edge and those that do not. Continuing to process the edges in the pre-defined order is equivalent to refining the interval containing the spanning tree of given rank until we have an interval of size one. The resulting O(n3m) algorithm is somewhat impractical, and hence a more efficient method is of interest, The basic algorithm employs O(m) determinant computations, each costing O(n3) time. To improve upon this, we could attempt to reduce the number of determinant computations, to reduce the cost of each determinant computation individually or to reduce the total cost of all determinant computation. Colbourn, Debroni, and
UNRANKING
SPANNING
TREES
275
Myrvold [l] develop a scheme which reduces the number of determinant computations when many random spanning trees are required, but it affords no improvement when a single random spanning tree is required. It appears that an inclusion/exclusion method may be required to examine essentially all of the edges to form a single spanning tree. To accelerate the computation of each determinant, one could in principle use improved algorithms for matrix multiplication; however, for the applications of interest, the traditional O(n3) algorithm is more practical. Finally, what remains is the option of reducing the overall cost of computing the O(m) determinants, by avoiding unnecessary recalculation; we pursue this strategy in the following sections.
3.
ORDERING
THE SPANNING
TREES
Part of the basis for the O(n3) unranking algorithm for spanning trees is a novel technique for ordering the edges of a graph that imposes the same natural ordering on the set of trees as described in the previous section. This ordering scheme is later used to decide the inclusion status, not of single edges, but of sets of edges at a time. The edge-ordering scheme works as follows. For simplicity, we assume throughout that the number of nodes is a power of 2. Partition the n nodes in the graph G = (V, E) into two equal-sized classes, PO and Pr, as shown in Fig. 3 (the basis for the partition is explained shortly). All edges with both endpoints in class PO are designated as “type 1” edges, while edges with both endpoints in Pr are “type 2” edges. Any edge straddling these two partitions is unaccounted for-such edges are classi-
n
Ro PO
Pl
FIGURES
276
COLBOURN,
DAY,
FIGURE
AND
NEL
4
fied by further subdividing both partitions into equal-sized subclasses Poe, PoI, Pro, and PII and classifying these edges as type 3, 4, 5 or 6 (see Fig. 4). Every edge in the graph has now been labelled with an edge type and, more importantly, each set of edges of any one type has all of its endpoints in a set of nodes of size n/2 (for instance, all type 1 edges have their endpoints in PO, type 3 edges have their endpoints in subclasses PO0 and PI,,, and so on). Our order on the edges is now defined by stating that edges of type 1 precede edges of type 2, which precede edges of type 3, and so on. Since this type classification does not, in general, impose a total ordering on the edges, the edge labelling must be further refined. This is done by selecting all edges of a particular type, collecting the n/2 nodes that are the endpoints of all edges of that type, and performing the partitioning operation on this subset of nodes (as an example, for edges of type 4, we choose the nodes from subclasses PO0 and PII to be partitioned further). When done for each edge type, this operation gives rise to edges of type 11,12,13 , . . . ,21,22, and so on. This partitioning and refinement of labelling is done recursively until every edge in the graph has a unique label and the edges are now totally ordered lexicographically by their label. The remaining issue is how the actual partitioning is done. Consider the binary representation of each node number, and partition based on the bit value. At the first step, nodes whose binary number representations have a high-order bit of zero are placed in PO, while the remaining nodes make up PI; partitioning these classes is done based on the second bit. To complete the partitioning recursively, a simple rule suffices: for any set of nodes, the subclasses Pm, Pal, Plo, and PII are defined by partitioning on the first
IJNRANKING
277
SPANNING TREES
two high-order bits that do not have an identical value across all nodes in the set to be partitioned. There are two different categories of edge type at the top level: types 1 and 2, whose endpoints all lie in some set of size n/2 and may have, as endpoints, any pair of nodes in the set, and types 3, 4, 5, and 6, whose endpoints also lie in some set of size n/2 but must have, as endpoints, one node in one half of the set and one node in the other half. Because of this, edges of type 1 or 2 may be further refined into all of the subtypes 1 through 6, but edges of type 3 through 6 can only be further refined into subtypes 3 through 6 (that is, there are no edges of type 31,32,41,. . . ). To unrank a number into a spanning tree with this edge ordering, the simple strategy given earlier suffices. However, the edge ordering introduces an important feature: in processing type i edges (i = 1, . . . ,6), the contractions and deletions performed involve only n/2 nodes. In fact, for any type t of length k, edges of type t involve at most n/2k nodes; all contractions and deletions performed on edges of this type affect only this (small) set of nodes. Hence our edge ordering isolates changes in the degree matrix to a small submatrix during the processing of edges of any one type. 4. SHORTCUTS FOR DETERMINANT CALCULATION REDUCED MATRIX
OF THE
The second component of the spanning tree selection algorithm is based on the observation that one does not always need to recalculate the determinant of the new reduced matrix from scratch during the computation. To motivate this, consider the deletion of an edge from the graph. Recall the reduced matrix (with row n and column n deleted) for the graph in Fig. 1 and consider its upper-triangular form after Gaussian elimination, as shown above:
Deleting edge (3,4) from the original graph would result in the new reduced matrix and its eliminated form below:
The change in both matrices involves subtracting one from the element in position (3,3), which suggests (at least for this very restricted case) the
278
COLBOURN,
DAY,
AND
NEL
\
: FIG. 5.
O[ 1
Partially
D*
eliminated
I
D matrix.
maxim, “Deletion followed by elimination is equivalent to elimination followed by deletion.” We formalize this observation and establish similar results for contraction and row/column permutation. Define matrix operations to be performed on the degree matrix D that correspond to graph operations on G; these matrix operations, and their graphical equivalents, are: (a) Bij( D): decrement Di, j and Dj, j by -D;, j, and set both Di, j and Dj, i to zero, corresponding to deleting all edges between vertices i and j in G. (b) Vij( D): add row i to row j, then column i to column j, and delete row and column i, corresponding to contracting edge (i, j) in G and labelling the new vertex j. (c) Pij( D): interchange rows i and j and columns i and j; that is, swap the labels of vertices i and j in G. (d) c?~(D): perform Gaussian elimination of columns 1,. . . , k of D, leaving a lower right remaining submatrix D* as shown in Fig. 5. This operation has no graphical equivalent. The above matrix operations are defined in terms of the degree matrix D of a graph G-they work equally well for a reduced matrix D’, where any operations to be performed on a row or column that has already been removed are ignored. LEMMA 4.1. Let D be any n - 1 X n - 1 reduced matrix for G (with vertex set { 1, . . . , n}). For vertices i, j of G and k -C i, j,
(9 ‘k;c(gij(D)) = 44.j(dk,(D))* (‘) gkk(gijtD)) = sij(gkk(D))* (iii) ‘,c;c(‘ip,j(D)) = ‘ip,j(dk(D))* Proof: The lemma is proved for a general fore holds for the restricted case of a reduced lemma holds for the case k = 1, it holds for remainder of the proof therefore deals strictly
square matrix D and therematrix for any graph. If the general k by induction; the with the case k = 1.
UNRANKING
SPANNING
279
TREES
For some n X n matrix D, it is easy to see that the three commutativity conditions hold for elements in the first row and first column of D. What remains is to prove equality under all three operations for the remaining n-lXn1 submatrix of D. (i) 8t(Sij(
D)). Easily derived.
D)) = Pij(4’l(
(ii) &l(gij(D)) = Bij(gll(D))* Since the result of the operation gij( D) is to effect the assignment dij = dij + aij for some 2-by-2 submatrix of D and some values a, it suffices to observe that (dij + aij)
- $-dlj
= dij - >dlj
11
+ aij.
11
(in) eO,(%Tij(D)) = cipj(8r( D)). We need only consider the effect on the elements in row and column j here. If we extend the domain of our operators to individual entries dij of D, then for some element dkj, k # j, we have &",('iP,j(d,j))
=
qj(
=
'ltdkj))
(d/ci
+
dkj)
-
and equality holds. A similar technique can be used to prove equality for elements djk, k + j, while proving equality for djj is only slightly more complicated. This lemma is fundamental in the timing analysis; it guarantees that in accomplishing a deletion, contraction, or permutation, we can avoid the inspection of D and consider only c?~(D), thereby avoiding repeating a large portion of the Gaussian elimination.
5.
A
GENERAL
DESCRIPTION AND
A SIMPLE
OF
THE
ALGORITHM
EXAMPLE
We now have the building blocks for an O(n3) spanning tree unranking algorithm. The edge ordering technique from Section 3 allows us to order the set of all spanning trees in the interval [l, s], while the actual unranking procedure resembles the partition refinement algorithm from Section 2-given some interval that contains the tree with the given rank, we use some as-yet-unspecified criteria to partition the interval into two sub-intervals and continue inspecting the sub-interval that contains the tree with that rank. Eventually, we refine to an interval of size one.
280
COLBOURN,
DAY,
AND
NEL
The partitioning described above is accomplished by, firstly, deciding on the inclusion/exclusion status of all edges of type 1, followed by all edges of type 2, and so on. To see how this corresponds to a partitioning refinement of the interval, observe that the trees in the interval, because of the way they have been ordered, can be partitioned based on whether or not they contain any type 1 edges, as shown below: 1 ~-----------------------I------------------------~ trees with type 1 edges
trees with no type 1 edges
It is a simple matter to remove all type 1 edges from the graph, calculate the determinant of the resulting reduced matrix, and decide which partition contains the tree of given rank. If there are no type 1 edges, the right sub-interval can be regarded as 1 ~------------------------,------------------------~ trees with type 2 edges
trees with no type 2 edges
and the original process repeated for type 2 edges. If there are type 1 edges, then we can regard the left sub-interval as follows: 1 ~----------------+----------~~~~ trees with type 11 edges
no type 11 edges
In either case, we have reduced our search space (the size of the interval). As an example, consider again the graph shown in Fig. 1, shown below with an edge labelling and all eight of its spanning trees, all of which also have the appropriate edge labelling and have been placed along the number line accordingly. (Since there is only one edge of each type, we refer to the edge of type 1 as edge 1, and so on.)
3
Ii3 1
2
1
3
6
flVlLl~K7l7Z 123
4
1
125
126
136
156
235
236
356
UNRANKING
SPANNING
TREES
281
Assume that we are searching for the spanning tree with rank seven, that is, the tree with edges 2, 3, and 6. Initially, we decide on the status of all type 1 edges (all one of them) and determine that, since there are five trees that involve edge 1, we are interested in the other sub-interval, that is, interval [6,8]. Edge 1 is deleted and our graph is now
In deciding the status of type-2 edges (again, all one of them), forced to keep edge 2 since this corresponds to the interval [6,7] number line, which contains the rank we are interested in. This contracted and the graph corresponding to the current situation
we are on our edge is is now
So far, we have decided to exclude edge 1 and include edge 2, and we are examining the interval [6,7] on the number line. We note that, in the current graph, there are two edges that have identical endpoints but different edge types. This occurs because we require an edge to retain its original edge 1abeIling throughout the life of the algorithm, regardless of any contractions that may occur.
282
COLBOURN,
DAY,
AND
NEL
The status of edge 3 is obvious-it must be included to keep the graph connected, so the edge is contracted to produce the graph
There are no edges of type 4 and, based on the resulting sub-interval, we decide to delete edge 5 (corresponding to the interval [6,6]), which leaves us with the graph
corresponding to the interval [7,7]. This forces us to keep the last edge, edge 6, and our record of contractions tells us that the spanning tree constructed consists of edges 2, 3, and 6, this being the tree with rank seven. In the (likely) event that the graph has more than one edge of type 1, we are forced to determine the inclusion/exclusion status of the entire set of type-l edges before moving on to type-2 edges. If we decide that at least one type-l edge must be included, our next question is whether we must include any edges of type 11, followed by edges of type 12, and so on. We would like to do this without recalculating the determinant of the new reduced matrix from scratch. Using Lemma 4.1(i), we can begin by permuting the rows and columns of the degree matrix to move the endpoints of all type-l edges (at most n/2 of them) to the lower right-hand submatrix, in effect moving all type-l edges into this submatrix. We are now free to perform Gaussian elimination down as far as this submatrix, knowing that any decisions about type-l edges do not affect this elimination, see Fig. 6. Once we reach this stage, determining
UNRANKING
SPANNING TREES
283
FIG. 6. Partially eliminated reduced matrix
the status of the type-l edges resembles the initial problem-six consider, but each a problem of size n/2-by-n/2 only.
subcases to
6. THE RECURSIVEUNRANKING ALGORITHM All of the basic ingredients are now present to describe the O(n3) algorithm. The method operates recursively, and the graph changes during the course of the computation. Each recursive invocation is to determine, for some type t of length k, which of the edges of type t appear in the tree (and are to be contracted) and which do not appear in the tree (and are to be deleted). At entry to the procedure, endnodes of edges of type t all appear in the lower right n k X nk submatrix (nk I n/2k). The degree matrix has had partial Gaussian elimination performed on it to eliminate subdiagonal entries in the first n - nk columns; hence at this level, the degree matrix is specified by a product p of the first n - nk diagonal entries, together with the remaining nk x nk matrix D. The procedure also takes as input the start s and end f of the interval of spanning trees under consideration, along with the rank r of the desired spanning tree. Gn return, the procedure must report which edges of type t are to be deleted and which are to be contracted. The algorithm ensures that whenever the procedure is invoked, at least one edge of type t must be contracted. The procedure itself operates as follows. For each type t’ = tl, t2,. . . , t6 (of length k + 1) in turn, first locate all edges of type t’. Form a permuted copy C of D so as to place endnodes of edges of type t ’ in the lower right submatrix of C (nk+l I n/2 k+l). Perform partial Gaussian nk+lX nk+l elimination on C to eliminate the subdiagonal in the first nk - nk+ 1 columns; retain p’, the product of the first nk - n k+l diagonal elements, and D’, the remaining nk+l X nk+l matrix. Having produced D’, form D” by deleting all edges of type t’ in D’ (that is, for each edge e = (i, j) of type t’, apply gij). Compute d = det( D”).
284
COLBOURN,
DAY,
AND
NEL
Now p . p’ . d is the number of spanning trees not containing any type t’ edges. Hence if r, the rank desired, is at least f - p . p’ - d, no edge of type t’ appears in the spanning tree. In this case, we delete all edges of type t’ in D and advance to the next edge type. This has the effect of changing the start s of the interval to f - p - p’ - d. Otherwise at least one edge of type t’ must be included, and hence the end f of the interval of interest is f - p * p’ . d. If there is only one edge of type t’, we contract it in D and proceed to edges of the next type (this has the effect of reducing nk by one as well). The most complicated case arises when some edge(s) of type t’ are to be contracted, and there is more than one type t’ edge. In this case, we recursively invoke the procedure with the current start, end, and rank, and the matrix specified by diagonal product p . p’ and remaining matrix D’. A list of deletions and contractions for edges of type t’ is returned. We apply each of the deletions and contractions in D, updating s and f accordingly. We then proceed to edges of the next type. When all types tl, t2,. . . , t6 have been processed, we return the list of all deletions and contractions for all type t edges (we must ensure that the deletions and contractions are modified to undo the effect of any permutations applied within this invocation). At the top level, the procedure is simply processing edges of types 1,2,..., 6 in turn; upon completion, the list of edges to be contracted which it returns is precisely the edge set of the spanning tree with rank r. THEOREM 6.1. The above algorithm to select a random spanning tree of an undirected graph of order n runs in O(n3) time.
Proof If the initial partial Gaussian elimination requires time O(n3), and we reconstruct the degree matrix entirely after each edge type decision, the recurrence relation for the amount of work for an n-by-n matrix is given by T(n)
5 6n3 + 6T(n/2),
and hence T(n) is O(n3). The relation is actually rather conservative since, as mentioned before, not all edge types give rise to all six subcases. In addition, it may be the case that edges of two different types actually have identical endpoints due to a series of contractions, so that no reorganization of the submatrix is necessary. There remains the problem of bookkeeping that allows us to perform the matrix permutations described above. These operations employ two operations we must be able to perform efficiently; they are (1) Given some edge type, generate the list of edges of this type. (2) Given the edges in a particular type, determine the set of nodes that are the endnodes of the edges.
UNRANKING
SPANNING TREES
285
These two operations do not necessarily constitute the same problem due to the contractions that occur as the algorithm runs. Since the contraction operation involves identifying two nodes in the graph, we should be able to, given any node label, determine its current label. The first of these operations can be implemented efficiently by a number of common data structures, either hierarchical or array-based, while the relabelling of vertices due to contractions can be handled by a reasonable implementation of a UNION-FIND abstract data type. Since by Lemma 4.1, the operations performed on the partially eliminated submatrices yield the same result as in the original matrix, we can also state THEOREM 6.2. The algorithm the corresponding spanning tree.
given above correctly unranks
7. CONCLUDING
any value to
REMARKS
An O(n3) ranking algorithm for producing the rank of a given spanning tree is now straightforward; rather than employing information about the rank to determine the inclusion/exclusion status of edges, we simply use the known status of edges of each type to refine the interval progressively until the desired rank is found. A number of practical issues arise in implementing this algorithm. Some of the determinant computations now performed can be avoided by checking for bridges (cut-edges); if a bridge is formed during the computation, it can be immediately contracted since it must appear in the spanning tree. Similarly, if a loop is formed it can be immediately deleted. These do not affect the O(n3) running time, but afford some practical improvement. Another improvement in practice can be realized by observing that our tree ordering depends on the edge ordering, which in turn depends on the node ordering, which is arbitrary. While any node ordering leads to an O(n3) worst case, the variation in actual computation time can be of practical significance. A useful heuristic is to order the nodes so that at the top level(s) in the partitioning, as few different edge types as possible are realized. This has the effect of minimizing the number of subproblems of size n/2, n/4 and so on. Finally, a practical improvement can be obtained by “batching”; if many random spanning trees are required, for example, we can sort the ranks r,, . . . , r, of the desired trees. In unranking ri and ri+l, much of the computation is duplicated typically; avoiding recomputation of the common part of the two computations yields a dramatic savings in practice.
286
COLBOURN,
DAY,
AND
NEL
It is perhaps remarkable that selecting a random spanning tree via many determinant computations is only a constant factor more expensive than simply performing one determinant computation. The sixfold divide-andconquer used here gives a powerful technique for partitioning graph problems. We expect that this technique will prove useful in other problems involving random selection of combinatorial objects. ACKNOWLEDGMENTS
Thanks to Wendy Myrvold and Bill Pulleyblank for helpful discussions about this research. The research of the first author is supported by NSERC Canada under Grant AO579.
REFERENCES
1. C. J. COLBOURN, B. M. DEBRONI, AND W. J. MYRVOLD, Estimating the coefficients of the reliability polynomial, in “Proceedings, Seventeenth Manitoba Conference on Numerical Mathematics and Computing, 1987, to appear. 2. W. FEUSSNER, iiber Stromberzeigung in nekformigen Leitem, Ann. Physik 9 (1902), 1304-1329; 15 (1904), 385-394. 3. A. GUENOCEIB, Random spanning tree, J. Algorithms 4 (1983), 214-220. 4. G. KATONA, A theorem of finite sets, “Theory of Graphs,” (P. Erd& and G. Katona, Eds.), pp. 187-207, Akademia Kiado, Budapest, 1966. 5. G. KIRCHOFF, iiber die Aufl&mg der Gleichungen auf welche man sei der Untersuchung der Linearen Verteihmg Galvanischer Strome Gefuhrt wind, Poppendorg’s Ann. Phys. Chem. 72 (1847), 497-508 (On the solution of equations obtained from the investigation of the linear distribution of galvanic currents, IRE Trans. Circuit Theory CT-5 (1958), 4-7. 6. G. D. KNOIT, A numbering system for binary trees, Comm. ACM 28 (1977), 113-115. 7. J. B. KRUSKAL, The number of simplices in a complex, in “Mathematical Optimization Techniques” (R. Bellman, Pd.), pp. 251-278, Univ. of California Press, Berkeley, 1963. 8. V. G. KULKARNI, “Generating Random Spanning Trees in a Graph,” Tech. Report UNC/ORSA/TR-85/l, University of North Carolina, 1985. 9. L. D. NEL AND C. J. COLBOURN, Combining Monte Carlo estimates and bounds for network reliability, submitted for publication, 1988. 10. H. hijFER, Neuer Reweis e-inesSatzes Uber Permutationen, Arch. Muth. Phys. 27 (1918), 742-744.