Practical Algorithms for Selection on Coarse-Grained Parallel Computers.
AbstractIn this paper, we consider the problem of selection on coarse-grained distributed memory parallel computers. We discuss several deterministic and randomized algorithms for parallel selection. We also consider several algorithms for load balancing needed to keep a balanced distribution of data across processors during the execution of the selection algorithms. We have carried out detailed implementations of all the algorithms discussed on the CM-5 and report on the experimental results. The results clearly demonstrate the role of randomization in reducing communication overhead.
better in practice than its deterministic counterpart due to the low constant associated with the
algorithm.
Parallel selection algorithms are useful in such practical applications as dynamic distribution
of multidimensional data sets, parallel graph partitioning and parallel construction of multidimensional
binary search trees. Many parallel algorithms for selection have been designed for the PRAM
model [2, 3, 4, 9, 14] and for various network models including trees, meshes, hypercubes and re-configurable
architectures [6, 7, 13, 16, 22]. More recently, Bader et.al. [5] implement a parallel
deterministic selection algorithm on several distributed memory machines including CM-5, IBM
SP-2 and INTEL Paragon. In this paper, we consider and evaluate parallel selection algorithms
for coarse-grained distributed memory parallel computers. A coarse-grained parallel computer consists
of several relatively powerful processors connected by an interconnection network. Most of
the commercially available parallel computers belong to this category. Examples of such machines
include CM-5, IBM SP-1 and SP-2, nCUBE 2, INTEL Paragon and Cray T3D.
The rest of the paper is organized as follows: In Section 2, we describe our model of parallel
computation and outline some primitives used by the algorithms. In Section 3, we present two deterministic
and two randomized algorithms for parallel selection. Selection algorithms are iterative
and work by reducing the number of elements to consider from iteration to iteration. Since we can
not guarantee that the same number of elements are removed on every processor, this leads to load
imbalance. In Section 4, we present several algorithms to perform such a load balancing. Each of
the load balancing algorithms can be used by any selection algorithm that requires load balancing.
In Section 5, we report and analyze the results we have obtained on the CM-5 by detailed implementation
of the selection and load balancing algorithms presented. In section 6, we analyze the
selection algorithms for meshes and hypercubes. Section 7 discusses parallel weighted selection.
We conclude the paper in Section 8.
Preliminaries
2.1 Model of Parallel Computation
We model a coarse-grained parallel machine as follows: A coarse-grained machine consists of several
relatively powerful processors connected by an interconnection network. Rather than making specific
assumptions about the underlying network, we assume a two-level model of computation. The
two-level model assumes a fixed cost for an off-processor access independent of the distance between
the communicating processors. Communication between processors has a start-up overhead of - ,
while the data transfer rate is 1
- . For our complexity analysis we assume that - and - are constant
and independent of the link congestion and distance between two processors. With new techniques,
such as wormhole routing and randomized routing, the distance between communicating processors
seems to be less of a determining factor on the amount of time needed to complete the communica-
tion. Furthermore, the effect of link contention is eased due to the presence of virtual channels and
the fact that link bandwidth is much higher than the bandwidth of node interface. This permits us
to use the two-level model and view the underlying interconnection network as a virtual crossbar
network connecting the processors. These assumptions closely model the behavior of the CM-5 on
which our experimental results are presented. A discussion on other architectures is presented in
Section 6.
2.2 Parallel Primitives
In the following, we describe some important parallel primitives that are repeatedly used in our
algorithms and implementations. We state the running time required for each of these primitives
under our model of parallel computation. The analysis of the run times for the primitives described
is fairly simple and is omitted in the interest of brevity. The interested reader is referred to [15].
In what follows, p refers to the number of processors.
1. Broadcast
In a Broadcast operation, one processor has an element of data to be broadcast to all other
processors. This operation can be performed in O((-) log p) time.
2. Combine
Given an element of data on each processor and a binary associative and commutative op-
eration, the Combine operation computes the result of combining the elements stored on all
the processors using the operation and stores the result on every processor. This operation
can also be performed in O((-) log p) time.
3. Parallel Prefix
Suppose that x are p data elements with processor P i containing x i .
Let\Omega be a
binary associative operation. The Parallel Prefix operation stores the value of x
on processor P i . This operation can be be performed in O((-) log p) time.
4. Gather
Given an element of data on each processor, the Gather operation collects all the data and
stores it in one of the processors. This can be accomplished in O(- log
5. Global Concatenate
This is same as Gather except that the collected data should be stored on all the processors.
This operation can also be performed in O(- log
6. Transportation Primitive
The transportation primitive performs many-to-many personalized communication with possibly
high variance in message size. If the total length of the messages being sent out or
received at any processor is bounded by t, the time taken for the communication is 2-t (+
lower order terms) when t - O(p p-). If the outgoing and incoming traffic bounds are
r and c instead, the communication takes time c) (+ lower order terms) when either
3 Parallel Algorithms for Selection
Parallel algorithms for selection are also iterative and work by reducing the number of elements
to be considered from iteration to iteration. The elements are distributed across processors and
each iteration is performed in parallel by all the processors. Let n be the number of elements
and p be the number of processors. To begin with, each processor is given d n
Otherwise, this can be easily achieved by using one of the load balancing techniques to be described
in Section 4. Let n (j)
i be the number of elements in processor P i at the beginning of iteration j.
Algorithm 1 Median of Medians selection algorithm
Total number of elements
Total number of processors labeled from 0 to
List of elements on processor P i , where jL
desired rank among the total elements
On each processor P i
Step 1. Use sequential selection to find median m i of list L i [l; r]
2.
Step 3. On P0
Find median of M , say MoM , and broadcast it to all processors.
Step 4. Partition L i into - MoM and ? MoM to give index i , the split index
Step 5. count = Combine(index i , add) calculates the number of elements
Step 6. If (rank - count )
else
Step 7. LoadBalance(L
Step 8.
Step 9. On P0
Perform sequential selection to find element q of rank in L

Figure

1: Median on Medians selection Algorithm
. Let k (j) be the rank of the element we need to identify among these n (j)
elements. We use this notation to describe all the selection algorithms presented in this paper.
3.1 Median of Medians Algorithm
The median of medians algorithm is a straightforward parallelization of the deterministic sequential
algorithm [8] and has recently been suggested and implemented by Bader et. al. [5]. This algorithm

Figure

load balancing at the beginning of each iteration.
At the beginning of iteration j, each processor finds the median of its n (j)
elements using the sequential deterministic algorithm. All such medians are gathered on one pro-
cessor, which then finds the median of these medians. The median of medians is then estimated
to be the median of all the n (j) elements. The estimated median is broadcast to all the processors.
Each processor scans through its set of points and splits them into two subsets containing elements
less than or equal to and greater than the estimated median, respectively. A Combine operation
and a comparison with k (j) determines which of these two subsets to be discarded and the value of
k (j+1) needed for the next iteration.
Selecting the median of medians as the estimated median ensures that the estimated median
will have at least a guaranteed fraction of the number of elements below it and at least a guaranteed
fraction of the elements above it, just as in the sequential algorithm. This ensures that the worst
case number of iterations required by the algorithm is O(log n). Let n (j)
. Thus,
finding the local median and splitting the set of points into two subsets based on the estimated
median each requires O(n (j)
in the j th iteration. The remaining work is one Gather, one
Broadcast and one Combine operation. Therefore, the worst-case running time of this algorithm is
log
p ), the running time is O( n
log n+ - log p log n+
-p log n).
This algorithm requires the use of load balancing between iterations. With load balancing,
. Assuming load balancing and ignoring the cost of load balancing itself, the running
time of the algorithm reduces to P log
3.2 Bucket-Based Algorithm
The bucket-based algorithm [17] attempts to reduce the worst-case running time of the above algorithm
without requiring load balance. This algorithm is shown in Figure 2. First, in order to keep
the algorithm deterministic without a balanced number of elements on each processor, the median
of medians is replaced by the weighted median of medians. As before, local medians are computed
on each processor. However, the estimated median is taken to be the weighted median of the local
medians, with each median weighted by the number of elements on the corresponding processor.
This will again guarantee that a fixed fraction of the elements is dropped from consideration every
iteration. The number of iterations of the algorithm remains O(log n).
The dominant computational work in the median of medians algorithm is the computation of
the local median and scanning through the local elements to split them into two sets based on the
estimated median. In order to reduce this work which is repeated every iteration, the bucket-based
approach preprocesses the local data into O(log p) buckets such that for any 0 -
every element in bucket i is smaller than any element in bucket j. This can be accomplished by
finding the median of the local elements, splitting them into two buckets based on this median
and recursively splitting each of these buckets into log pbuckets using the same procedure. Thus,
preprocessing the local data into O(log p) buckets requires O( n
log log p) time.
Bucketing the data simplifies the task of finding the local median and the task of splitting the
local data into two sets based on the estimated median. To find the local median, identify the
bucket containing the median and find the rank of the median in the bucket containing the median
Algorithm 2 Bucket-based selection algorithm
Total number of elements
Total number of processors labeled from 0 to
List of elements on processor P i , where jL
desired rank among the total elements
On each processor P i
Step 0. Partition L i on P i into log p buckets of equal size such that if r 2 bucket j , and s 2 bucketk , then r ! s if
whilen ? C (a constant)
Step 1. Find the bucket bktk containing the median element using a binary search on the remaining
buckets. This is followed by finding the appropriate rank in bktk to find the median m i . Let N i be the
number of remaining keys on P i .
2.
Step 3. On P0
Find the weighted median of M , say WM and broadcast it.
Step 4. Partition L i into - WM and ? WM using the buckets to give index i ; the split index
Step 5. count = Combine(index i , add) calculates the number of elements less than WM
Step 6. If (rank - count )
else
Step 7.
Step 8. On P0
Perform sequential selection to find element q of rank in L

Figure

2: Bucket-based selection algorithm
in O(log log p) time using binary search. The local median can be located in the bucket by the
sequential selection algorithm in O( n
time. The cost of finding the local median reduces from
O( n
log p ). To split the local data into two sets based on the estimated median,
first identify the bucket that should contain the estimated median. Only the elements in this bucket
need to be split. Thus, this operation also requires only O(log log
log p
time.
After preprocessing, the worst-case run time for selection is O(log log p log
log p log
log p log n+ -p log n) = O( n
log p
log log log p. Therefore, the
worst-case run time of the bucket-based approach is O( n
log
without any load balancing.
Algorithm 3 Randomized selection algorithm
Total number of elements
Total number of processors labeled from 0 to
List of elements on processor P i , where jL
desired rank among the total elements
On each processor P i
whilen ? C (a constant)
Step
Step 1.
Step 2. Generate a random number nr (same on all processors) between 0 and
Step 3. On Pk where (nr
Step 4. Partition L i into - mguess and ? mguess to give index i , the split index
Step 5. count = Combine(index i , add) calculates the number of elements less than mguess
Step 6. If (rank - count )
else
Step 7.
Step 8. On P0
Perform sequential selection to find element q of rank in L

Figure

3: Randomized selection algorithm
3.3 Randomized Selection Algorithm
The randomized median finding algorithm (Figure 3) is a straightforward parallelization of the
randomized sequential algorithm described in [12]. All processors use the same random number
generator with the same seed so that they can produce identical random numbers. Consider the
behavior of the algorithm in iteration j. First, a parallel prefix operation is performed on the
's. All processors generate a random number between 1 and n (j) to pick an element at random,
which is taken to be the estimate median. From the parallel prefix operation, each processor can
determine if it has the estimated median and if so broadcasts it. Each processor scans through
its set of points and splits them into two subsets containing elements less than or equal to and
greater than the estimated median, respectively. A Combine operation and a comparison with k (j)
determines which of these two subsets to be discarded and the value of k (j+1) needed for the next
iteration.
Since in each iteration approximately half the remaining points are discarded, the expected
number of iterations is O(log n) [12]. Let n (j)
. Thus, splitting the set of points into
two subsets based on the median requires O(n (j)
in the j th iteration. The remaining work
is one Parallel Prefix, one Broadcast and one Combine operation. Therefore, the total expected
running time of the algorithm is P log
p ), the expected running time is O( n
log n). In practice,
one can expect that n (j)
max reduces from iteration to iteration, perhaps by half. This is especially
true if the data is randomly distributed to the processors, eliminating any order present in the
input. In fact, by a load balancing operation at the end of every iteration, we can ensure that for
every iteration j, n (j)
. With load balancing and ignoring the cost of it, the running
time of the algorithm reduces to P log
log n). Even
without this load balancing, assuming that the initial data is randomly distributed, the running
time is expected to be O( n
log n).
3.4 Fast Randomized Selection Algorithm
The expected number of iterations required for the randomized median finding algorithm is O(log n).
In this section we discuss an approach due to Rajasekharan et. al. [17] that requires only
O(log log n) iterations for convergence with high probability (Figure 4).
Suppose we want to find the k th smallest element among a given set of n elements. Sample a set
S of o(n) keys at random and sort S. The element with rank
e in S will have an expected
rank of k in the set of all points. Identify two keys l 1 and l 2 in S with ranks
ffi is a small integer such that with high probability the rank of l 1 is ! k and the rank of l 2 is ? k in
the given set of points. With this, all the elements that are either ! l 1 or ? l 2 can be eliminated.
Recursively find the element with rank in the remaining elements. If the number
of elements is sufficiently small, they can be directly sorted to find the required element.
If the ranks of l 1 and l 2 are both ! k or both ? k, the iteration is repeated with a different
sample set. We make the following modification that may help improve the running time of the
algorithm in practice. Suppose that the ranks of l 1 and l 2 are both ! k. Instead of repeating the
iteration to find element of rank k among the n elements, we discard all the elements less than l 2
and find the element of rank in the remaining elements. If the ranks of
l 1 and l 2 are both ? k, elements greater than l 1 can be discarded.
Rajasekharan et. al. show that the expected number of iterations of this median finding
algorithm is O(log log n) and that the expected number of points decreases geometrically after each
iteration. If n (j) is the number of points at the start of the j th iteration, only a sample of o(n (j) )
keys is sorted. Thus, the cost of sorting, o(n (j) log n (j) ) is dominated by the O(n (j) ) work involved
in scanning the points.
Algorithm 4 Fast randomized selection algorithm
Total number of elements
Total number of processors labeled from 0 to
List of elements on processor P i , where jL
desired rank among the total elements
On each processor P i
whilen ? C (a constant)
Step
Step 1. Collect a sample S i from L i [l; r] by picking n i n ffl
n elements at random on P i between l and r.
Step 2.
On P0
Step 3. Pick k1 , k2 from S with ranks d ijSj
jSjlogne and d ijSj
jSjlogne
Step 4. Broadcast k1 and k2.The rank to be found will be in [k1 , k2 ] with high probability.
Step 5. Partition L i between l and r into ! k1 , [k1 , k2 ] and ? k2 to give counts less, middle
and high and splitters s1 and s2 .
Step 6.
Step 7. cless = Combine(less, add)
Step 8. If (rank 2 (cless ; cmid ])
else
else
Step 9.
Step 10. On P0
Perform sequential selection to find element q of rank in L

Figure

4: Fast Randomized selection Algorithm
In iteration j, Processor P (j)
randomly selects n (j)
n (j) of its n (j)
elements. The selected elements
are sorted using a parallel sorting algorithm. Once sorted, the processors containing the elements
l (j)and l (j)broadcast them. Each processor finds the number of elements less than l (j)and greater
than l (j)
contained by it. Using Combine operations, the ranks of l (j)
1 and l (j)
are computed and
the appropriate action of discarding elements is undertaken by each processor. A large value of
ffl increases the overhead due to sorting. A small value of ffl increases the probability that both
the selected elements (l (j)
1 and l (j)
lie on one side of the element with rank k (j) , thus causing an
unsuccessful iteration. By experimentation, we found a value of 0:6 to be appropriate.
As in the randomized median finding algorithm, one iteration of the median finding algorithm
takes O(n (j)
log log n iterations are required, median
finding requires O( n
log log n + (-) log p log log n) time.
As before, we can do load balancing to ensure that n (j)
reduces by half in every iteration.
Assuming this and ignoring the cost of load balancing, the running time of median finding reduces
to
log log
log log n). Even without this load balancing,
the running time is expected to be O( n
log log n).
4 Algorithms for load balancing
In order to ensure that the computational load on each processor is approximately the same during
every iteration of a selection algorithm, we need to dynamically redistribute the data such that
every processor has nearly equal number of elements. We present three algorithms for performing
such a load balancing. The algorithms can also be used in other problems that require dynamic
redistribution of data and where there is no restriction on the assignment of data to processors.
We use the following notation to describe the algorithms for load balancing: Initially, processor
is the total number of elements on all the processors, i.e.
4.1 Order Maintaining Load Balance
Suppose that each processor has its set of elements stored in an array. We can view the n elements
as if they were globally sorted based on processor and array indices. For any i ! j, any element
in processor P i appears earlier in this sorted order than any element in processor P j . The order
maintaining load balance algorithm is a parallel prefix based algorithm that preserves this global
order of data after load balancing.
The algorithm first performs a Parallel Prefix operation to find the position of the elements it
contains in the global order. The objective is to redistribute data such that processor P i contains
Algorithm 5 Modified order maintaining load balance
- Number of total elements
Total number of processors labeled from 0 to
List of elements on processor P i of size n i
On each processor P i
Step
increment navg
Step 1.
Step 2. diff
Step 3. If diff [j] is positive P j is labeled as a source. If diff [j] is negative P j is labeled as
a sink.
Step 4. If P i is a source calculate the prefix sum of the positive diff [ ] in an array p src, else calculate
the prefix sums for sinks using negative diff [ ] in p snk.
Step 5. l
Step 7. Calculate the range of destination processors [P l ; Pr ] using a binary search on p snk.
Step 8. while(l - r)
elements to P l and increment l
Step 5. l
Step 7. Calculate the range of source processors [P l ; Pr ] using a binary search on
src.
Step 8. while( l - r)
Receive elements from P l and increment l

Figure

5: Modified order maintaining load balance
the elements with positions n avg in the global order. Using the parallel prefix
operation, each processor can figure out the processors to which it should send data and the amount
of data to send to each processor. Similarly, each processor can figure out the amount of data it
should receive, if any, from each processor. Communication is generated according to this and the
data is redistributed.
In our model of computation, the running time of this algorithm only depends on the maximum
communication generated/received by a processor. The maximum number of messages sent out by
a processor is d nmax
navg e+1 and the maximum number of elements sent is n max . The maximum number
of elements received by a processor is n avg . Therefore, the running time is O(- nmax
The order maintaining load balance algorithm may generate much more communication than
necessary. For example, consider the case where all processors have n avg elements except that P 0
has one element less and P p\Gamma1 has one element more than n avg . The optimal strategy is to transfer
the one extra element from P p to P 0 . However, this algorithm transfers one element from P i to
messages.
Since preserving the order of data is not important for the selection algorithm, the following
modification is done to the algorithm: Every processor retains minfn of its original elements.
the processor has (n excess and is labeled a source. Otherwise,
the processor needs (n avg \Gamma n i ) elements and is labeled a sink. The excessive elements in the source
processors and the number of elements needed by the sink processors are ranked separately using
two Parallel Prefix operations. The data is transferred from sources to sinks using a strategy similar
to the order maintaining load balance algorithm. This algorithm (Figure 5), which we call modified
order maintaining load balance algorithm (modified OMLB), is implemented in [5].
The maximum number of messages sent out by a processor in modified OMLB is O(p) and the
maximum number of elements sent is (n The maximum number of elements received
by a processor is n avg . The worst-case running time is O(-p
4.2 Dimension Exchange Method
The dimension exchange method (Figure 6) is a load balancing technique originally proposed for
hypercubes [11][21]. In each iteration of this method, processors are paired to balance the load
locally among themselves which eventually leads to global load balance. The algorithm runs in
log p iterations. In iteration i (0 processors that differ in the i th least significant
bit position of their id's exchange and balance the load. After iteration i, for any 0 -
processors have the same number of elements.
In each iteration, ppairs of processors communicate in parallel. No processor communicates
more than nmaxelements in an iteration. Therefore, the running time is O(- log
However, since 2 j processors hold the maximum number of elements in iteration j, it is likely that
either n max is small or far fewer elements than nmaxare communicated. Therefore, the running
time in practice is expected to be much better than what is predicated by the worst-case.
4.3 Global Exchange
This algorithm is similar to the modified order maintaining load balance algorithm except that
processors with large amounts of data are directly paired with processor with small amounts of
data to minimize the number of messages (Figure 7).
As in the modified order maintaining load balance algorithm, every processor retains minfn
of its original elements. If the processor has (n excess and is la-
Algorithm 6 Dimension exchange method
- Number of total elements
Total number of processors labeled from 0 to
List of elements on processor P i of size n i
On each processor P i
Step 1. P
Step 2. Exchange the count of elements between P
Step 3.
Step 4. Send elements from L i [navg ] to processor P l
Step 5. n
else
Step 4. Receive n l \Gamma navg elements from processor P l at
Step 5. Increment n i by n l \Gamma navg

Figure

exchange method for load balancing
beled a source. Otherwise, the processor needs (n avg \Gamma n i ) elements and is labeled a sink. All
the source processors are sorted in non-increasing order of the number of excess elements each
processor holds. Similarly, all the sink processors are sorted in non-increasing order of the number
of elements each processor needs. The information on the number of excessive elements in each
source processor is collected using a Global Concatenate operation. Each processor locally ranks
the excessive elements using a prefix operation according to the order of the processors obtained by
the sorting. Another Global Concatenate operation collects the number of elements needed by each
sink processor. These elements are then ranked locally by each processor using a prefix operation
performed using the ordering of the sink processors obtained by sorting.
Using the results of the prefix operation, each source processor can find the sink processors to
which its excessive elements should be sent and the number of element that should be sent to each
such processor. The sink processors can similarly compute information on the number of elements
to be received from each source processor. The data is transferred from sources to sinks. Since the
sources containing large number of excessive elements send data to sinks requiring large number of
elements, this may reduce the total number of messages sent.
In the worst-case, there may be only one processor containing all the excessive elements and thus
the total number of messages sent out by the algorithm is O(p). No processor will send more than
data and the maximum number of elements received by any processor is
n avg . The worst-case run time is O(-p
Algorithm 7 Global Exchange load balance
- Number of total elements
Total number of processors labeled from 0 to
List of elements on processor P i of size n i
On each processor P i
Step
increment navg
Step 1.
for j /0 to
Step 2. diff
Step 3. If diff [j] is positive P j is labeled as a source. If diff [j] is negative P j is labeled as
a sink.
Step 4. For k 2 [0; sources in descending order maintaining appropriate
processor indices. Also sort diff [k] for sinks in ascending order.
Step 5. If P i is a source calculate the prefix sum of the positive diff [ ] in an array p src, else calculate
the prefix sums for sinks using negative diff [ ] in p snk.
Step 6. If P i is a source calculate the prefix sum of the positive diff [ ] in an array p src, else calculate
the prefix sums for sinks using negative diff [ ] in p snk.
Step 7. l
8. r
Step 9. Calculate the range of destination processors [P l ; Pr ] using a binary search on p snk.
Step 10. while(l - r)
elements to P l and increment l
Step 7. l
8. r
Step 9. Calculate the range of source processors [P l ; Pr ] using a binary search on
src.
Step 10. while( l - r)
Receive elements from P l and increment l

Figure

7: Global exchange method for load balancing
Selection Algorithm Run-time
Median of Medians O( n
Randomized O( n
log n)
Fast randomized O( n
log log n)

Table

1: The running times of various selection algorithm assuming but not including the cost of
load balancing
Selection Algorithm Run-time
Median of Medians O( n
Bucket-based O( n
log
Randomized O( n
log log n)
Fast randomized O( n
log log n + (-) log p log log n)

Table

2: The worst-case running times of various selection algorithms
5 Implementation Results
The estimated running times of various selection algorithms are summarized in Table 1 and Table 2.

Table

1 shows the estimated running times assuming that each processor contains approximately
the same number of elements at the end of each iteration of the selection algorithm. This can be
expected to hold for random data even without performing any load balancing and we also observe
this experimentally. Table 2 shows the worst-case running times in the absence of load balancing.
We have implemented all the selection algorithms and the load balancing techniques on the CM-
5. To experimentally evaluate the algorithms, we have chosen the problem of finding the median of a
given set of numbers. We ran each selection algorithm without any load balancing and with each of
the load balancing algorithms described (except for the bucket-based approach which does not use
load balancing). We have run all the resulting algorithms on 32k, 64k, 128k, 256k, 512k, 1024k and
2048k numbers using 2, 4, 8, 16, 32, 64 and 128 processors. The algorithms are run until the total
number of elements falls below p 2 , at which point the elements are gathered on one processor and
the problem is solved by sequential selection. We found this to be appropriate by experimentation,
to avoid the overhead of communication when each processor contains only a small number of
elements. For each value of the total number of elements, we have run each of the algorithms on
two types of inputs - random and sorted. In the random case, n
p elements are randomly generated on
each processor. To eliminate peculiar cases while using the random data, we ran each experiment on
five different random sets of data and used the average running time. Random data sets constitute
close to the best case input for the selection algorithms. In the sorted case, the n numbers are
chosen to be the numbers containing the numbers i n
The sorted input is a close to the worst-case input for the selection algorithms. For example, after
the first iteration of a selection algorithm using this input, approximately half of the processors lose
all their data while the other half retains all of their data. Without load balancing, the number of
active processors is cut down by about half every iteration. The same is true even if modified order
maintaining load balance and global exchange load balancing algorithms are used. After every
iteration, about half the processors contain zero elements leading to severe load imbalance for the
load balancing algorithm to rectify. Only some of the data we have collected is illustrated in order
to save space.
The execution times of the four different selection algorithms without using load balancing for
random data (except for median of medians algorithm requiring load balancing for which global
exchange is used) with 128k, 512k and 2048k numbers is shown in Figure 8. The graphs clearly
demonstrate that all four selection algorithms scale well with the number of processors. An immediate
observation is that the randomized algorithms are superior to the deterministic algorithms by
an order of magnitude. For example, with the median of medians algorithm ran
at least 16 times slower and the bucket-based selection algorithm ran at least 9 times slower than
either of the randomized algorithms. Such an order of magnitude difference is uniformly observed
even using any of the load balancing techniques and also in the case of sorted data. This is not
surprising since the constants involved in the deterministic algorithms are higher due to recursively
finding the estimated median. Among the deterministic algorithms, the bucket-based approach
consistently performed better than the median of medians approach by about a factor of two for
random data. For sorted data, the bucket-based approach which does not use any load balancing
ran only about 25% slower than median of medians approach with load balancing.
In each iteration of the parallel selection algorithm, each processor also performs a local selection
algorithm. Thus the algorithm can be split into a parallel part where the processors combine the
results of their local selections and a sequential part involving executing the sequential selection
locally on each processor. In order to convince ourselves that randomized algorithms are superior
in either part, we ran the following hybrid experiment. We ran both the deterministic parallel
selection algorithms replacing the sequential selection parts by randomized sequential selection. The
running time of the hybrid algorithms was in between the deterministic and randomized parallel
selection algorithms. We made the following observation: The factor of improvement in randomized
parallel selection algorithms over deterministic parallel selection is due to improvements in both the
sequential and parallel parts. For large n, much of the improvement is due to the sequential part.
For large p, the improvement is due to the parallel part. We conclude that randomized algorithms
are faster in practice and drop the deterministic algorithms from further consideration.
Time
(in
seconds)
Number of Processors
Median of Medians
Bucket Based
Randomized
Fast Randomized0.0150.0250.0350.0450.0550.065
Time
(in
seconds)
Number of Processors
Randomized
Fast Randomized0.51.52.53.5
Time
(in
seconds)
Number of Processors
Median of Medians
Bucket Based
Randomized
Fast Randomized0.040.080.120.160.20.24
Time
(in
seconds)
Number of Processors
Randomized
Fast Randomized261014
Time
(in
seconds)
Number of Processors
Median of Medians
Bucket Based
Randomized
Fast Randomized0.10.30.50.7
Time
(in
seconds)
Number of Processors
Randomized
Fast Randomized

Figure

8: Performance of different selection algorithms without load balancing (except for median
of medians selection algorithm for which global exchange is used) on random data sets.
Time
(in
seconds)
Number of Processors
Random data, n=512k
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global Exchange0.10.30.50.7
Time
(in
seconds)
Number of Processors
Random data, n=2M
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global Exchange0.050.150.250.350.45
Time
(in
seconds)
Number of Processors
Sorted data, n=512k
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global Exchange0.20.611.41.8
Time
(in
seconds)
Number of Processors
Sorted data, n=2m
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global Exchange

Figure

9: Performance of randomized selection algorithm with different load balancing strategies
on random and sorted data sets.
To facilitate an easier comparison of the two randomized algorithms, we show their performance
separately in Figure 8. Fast randomized selection is asymptotically superior to randomized selection
for worst-case data. For random data, the expected running times of randomized and fast
randomized algorithms are O( n
log n) and O( n
log log n), respectively.
Consider the effect of increasing n for a fixed p. Initially, the difference in log n and log log n is not
significant enough to offset the overhead due to sorting in fast randomized selection and randomized
selection performs better. As n is increased, fast randomized selection begins to outperform
randomized selection. For large n, both the algorithms converge to the same execution time since
the O( n
dominates. Reversing this point view, we find that for any fixed n, as
we increase p, randomized selection will eventually perform better and this can be readily observed
in the graphs.
The effect of the various load balancing techniques on the randomized algorithms for random
data is shown in Figure 9 and Figure 10. The execution times are consistently better without
using any load balancing than using any of the three load balancing techniques. Load balancing
Time
(in
seconds)
Number of Processors
Random data, n=512k
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global Exchange0.10.30.50.7
Time
(in
seconds)
Number of Processors
Random data, n=2m
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global
Time
(in
seconds)
Number of Processors
Sorted data, n=512k
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global
Time
(in
seconds)
Number of Processors
Sorted data, n=2m
Balance
Mod. order Maintaining Load Balance
Dimension Exchange
Global Exchange

Figure

10: Performance of fast randomized selection algorithm with different load balancing strategies
on random and sorted data sets.
for random data almost always had a negative effect on the total execution time and this effect
is more pronounced in randomized selection than in fast randomized selection. This is explained
by the fact that fast randomized selection has fewer iterations (O(log log n) vs. O(log n)) and less
data in each iteration.
The observation that load balancing has a negative effect on the running time for random
data can be easily explained: In load balancing, a processor with more elements sends some of
its elements to another processor. The time taken to send the data is justified only if the time
taken to process this data in future iterations is more than the time for sending it. Suppose that
a processor sends m elements to another processor. The processing of this data involves scanning
it in each iteration based on an estimated median and discarding part of the data. For random
data, it is expected that half the data is discarded in every iteration. Thus, the estimated total
time to process this data is O(m). The time for sending the data is (-m), which is also O(m).
By observation, the constants involved are such that load balancing is taking more time than the
reduction in running time caused by it.
Time
(in
seconds)
Number of Processors
Comparing two randomized selection algorithms using sorted data for n=512k
Randomized
Fast Randomized0.20.611.41.8
Time
(in
seconds)
Number of Processors
Comparing two randomized selection algorithm using sorted datas for n=2M
Randomized
Fast Randomized

Figure

11: Performance of the two randomized selection algorithms on sorted data sets using the
best load balancing strategies for each algorithm \Gamma no load balancing for randomized selection and
modified order maintaining load balancing for fast randomized selection.
Consider the effect of the various load balancing techniques on the randomized algorithms for
sorted data (see Figure 9 and Figure 10). Even in this case, the cost of load balancing more than
offset the benefit of it for randomized selection. However, load balancing significantly improved the
performance of fast randomized selection.
In

Figure

11, we see a comparison of the two randomized algorithms for sorted data with the
best load balancing strategies for each algorithm \Gamma no load balancing for randomized selection and
modified order maintaining load balancing for fast randomized algorithm (which performed slightly
better than other strategies). We see that, for large n, fast randomized selection is superior. We
also observe (see Figure 11 and Figure 8) that the fast randomized selection has better comparative
advantage over randomized selection for sorted data.
Finally, we consider the time spent in load balancing itself for the randomized algorithms on
both random and sorted data (see Figure 12 and Figure 13). For both types of data inputs, fast
randomized selection spends much less time than randomized selection in balancing the load. This
is reflective of the number of times the load balancing algorithms are utilized (O(log log n) vs.
O(log n)). Clearly, the cost of load balancing increases with the amount of imbalance and the
number of processors. For random data, the overhead due to load balancing is quite tolerable
for the range of n and p used in our experiments. For sorted data, a significant fraction of the
execution time of randomized selection is spent in load balancing. Load balancing never improved
the running time of randomized selection. Fast randomized selection benefited from load balancing
for sorted data. The choice of the load balancing algorithm did not make a significant difference in
the running time.
Consider the variance in the running times between random and sorted data for both the
Number of Processors0.10.30.50.7Time
(in
seconds)
Randomized selection , random data
load balancing time
O
Number of Processors0.20.61.01.41.8
Time
(in
seconds)
Randomized selection , sorted data
load balancing time
O

Figure

12: Performance of randomized selection algorithm with different load balancing strategies
balancing (N), Order maintaining load balancing (O), Dimension exchange method (D)
and Global exchange (G).
Number of Processors0.10.30.5Time
(in
seconds)
Fast Randomized selection , random data
load balancing time
O
Number of Processors0.20.61.0Time
(in
seconds)
Fast Randomized selection , sorted data
load balancing time
OD
G

Figure

13: Performance of fast randomized selection algorithm with different load balancing strategies
balancing (N), Order maintaining load balancing (O), Dimension exchange method
(D) and Global exchange (G).
Primitive Two-level model Hypercube Mesh
Broadcast O((-) log p) O((-) log p) O(- log
Combine O((-) log p) O((-) log p) O(- log
Parallel Prefix O((-) log p) O((-) log p) O(- log
Gather O(- log
Global Concatenate O(- log
Transportation O(-p

Table

3: Running time for basic communication primitives on meshes and hypercubes using cut-through
routing. For the transportation primitive, t refers to the maximum of the total size of
messages sent out or received by any processor.
randomized algorithms. The randomized selection algorithm ran 2 to 2.5 times faster for random
data than for sorted data (see Figure 12). Using any of the load balancing strategies, there is very
little variance in the running time of fast randomized selection (Figure 13). The algorithm performs
equally well on both best and worst-case data. For the case of 128 processors the stopping criterion
results in execution of one iteration in most runs. Thus, load balancing has a detrimental
effect on the overall cost. We had decided to choose the same stopping criterion to provide a fair
comparison between the different algorithms. However, an appropriate fine tuning of this stopping
criterion and corresponding increase in the number of iterations should provide time improvements
with load balancing for 2M data size on 128 processors.
6 Selection on Meshes and Hypercubes
Consider the analysis of the algorithms presented for cut-through routed hypercubes and square
meshes with p processors. The running time of the various algorithms on meshes and hypercubes is
easily obtained by substituting the corresponding running times for the basic parallel communication
primitives used by the algorithms. Table 3 shows the time required for each parallel primitive
on the two-level model of computation, a hypercube of p processors and a p p \Theta
mesh. The analysis is omitted to save space and similar analysis can be found in [15].
Load balancing can be achieved by using the communication pattern of the transportation
primitive [24] which involves two all-to-all personalized communications. Each processor has O( n
elements to be sent out. The worst-case time of order maintaining load balance is O(-p
for the hypercube and mesh respectively when n ? O(-p
exchange load balancing algorithm on the hypercube has worst-case run time of O(- log p+- n
log p)
and on the mesh it is O(- log
The global exchange load balancing algorithmm has the
same time complexities as the modified order maintaining load balancing algorithm, both on the
hypercube and the mesh. These costs must be added to the selection algorithms if analysis of the
algorithms with load balancing is desired.
From the table, the running times of all the primitives remain the same on a hypercube and
hence the analysis and the experimental results obtained for the two-level model will be valid for
hypercubes. Thus, the time complexity of all the selection algorithms is the same on the hypercube
as the two-level model discussed in this paper. If the ratio of unit computation cost to the unit
communication cost is large, i.e the processor is much faster than the underlying communication
network, cost of load balancing will offset its advantages and fast randomized algorithm without
load balancing will have superior performance for practical scenarios.
Load balancing on the mesh results in asymptotically worse time requirements. We would
expect load balancing to be useful for small number of processors. For large number of processors,
even one step of load balancing would dominate the overall time and hence would not be effective.
In the following we present results for the performance for best case and worst case data on a mesh.
1. Deterministic Algorithms: The communication primitives used in the deterministic selection
algorithms are Gather, Broadcast and Combine. Even though Broadcast and Combine require
more time than the two-level model, their cost is absorbed by the time required for the Gather
operation which is identical on the mesh and the two-level model. Hence, the complexity of
the deterministic algorithms on the mesh remains the same as the two-level model. The total
time requirements for the median of medians algorithm are O( n
log n) for
the best case and O( n
log n) for the worst case. The bucket-based
deterministic algorithm runs in to O( n
log time in the
worst case without load balancing.
2. Randomized Algorithms: The communication for the the randomized algorithm includes one
PrefixSum, one Broadcast and one Combine. The communication time on a mesh for one
iteration of the randomized algorithm is O(- log p+- p
p), making its overall time complexity
O( n
log n) for the best case and O( n
log log n) for
the worst case data.
The fast randomized algorithm involves a parallel sort of the sample for which we use bitonic
sort. A sample of n ffl (0 chosen from the n elements and sorted. On the
mesh, sorting a sample of n ffl elements using bitonic sort takes O(
log
should be acceptably small to keep the sorting phase from dominating every iteration.
The run-time of the fast randomized selection on the mesh is O( n
for the best case data. For the worst case data the time requirement would be O( n
log log p+
7 Weighted selection
having a corresponding weight w i attached to
it. The problem of weighted selection is to find an element x i such that
any x l
As an example, the weighted median
will be the element that divides the data set S with sum of weights W , into two sets S 1 and S 2
with approximately equal sum of weights. Simple modifications can be made to the deterministic
algorithms to adapt them for weighted selection. In iteration j of the selection algorithms, a set
S (j) of elements is split into two subsets S (j)
1 and S (j)
2 and a count of elements is used to choose
the subset in which the desired element can be found. Weighted selection is performed as follows:
First, the elements of S (j) are divided into two subsets S (j)
1 and S (j)
2 as in the selection algorithm.
The sum of weights of all the elements in subset S (j)
1 is computed. Let k j be the weight metric
in iteration j. If k j is greater than the sum of weights of S (j)
1 , the problem reduces to performing
weighted selection with k
(j). Otherwise, we need to
perform weighted selection with k
1 . This method retains the property
that a guaranteed fraction of elements will be discarded at each iteration keeping the worst case
number of iterations to be O(log n). Therefore, both the median of medians selection algorithm
and the bucket-based selection algorithm can be used for weighted selection without any change in
their run time complexities.
The randomized selection algorithm can also be modified in the same way. However, the same
modification to the fast randomized selection will not work. This algorithm works by sorting a
sample of the data set and picking up two elements that with high probability lie on either side
of the element with rank k in sorted order. In weighted selection, the weights determine the
position of the desired element in the sorted order. Thus, one may be tempted to select a sample of
weights. However, this does not work since the weights of the elements should be considered in the
order of the sorted data and a list of the elements sorted according to the weights does not make
sense. Hence, randomized selection without load balancing is the best choice for parallel weighted
selection.
Conclusions
In this paper, we have tried to identify the selection algorithms that are most suited for fast execution
on coarse-grained distributed memory parallel computers. After surveying various algorithms,
we have identified four algorithms and have described and analyzed them in detail. We also considered
three load balancing strategies that can be used for balancing data during the execution of
the selection algorithms.
Based on the analysis and experimental results, we conclude that randomized algorithms are
faster by an order of magnitude. If determinism is desired, the bucket-based approach is superior
to the median of medians algorithm. Of the two randomized algorithms, fast randomized selection
with load balancing delivers good performance for all types of input distributions with very little
variation in the running time. The overhead of using load balancing with well-behaved data is
insignificant. Any of the load balancing techniques described can be used without significant
variation in the running time. Randomized selection performs well for well-behaved data. There is
a large variation in the running time between best and worst-case data. Load balancing does not
improve the performance of randomized selection irrespective of the input data distribution.
9

Acknowledgements

We are grateful to Northeast Parallel Architectures Center and Minnesota Supercomputing Center
for allowing us to use their CM-5. We would like to thank David Bader for providing us a copy of
his paper and the corresponding code.



--R

Deterministic selection in O(log log N) parallel time
The design and analysis of parallel algorithms
Parallel selection in O(log log n) time using O(n
An optimal algorithm for parallel selection
Practical parallel algorithms for dynamic data redistribution

Technical Report CMU-CS-90-190
Time bounds for selection
A parallel median algorithm
Introduction to algorithms.
Dynamic load balancing for distributed memory multiprocessors
Expected time bounds for selection
Selection on the Reconfigurable Mesh
An introduction to parallel algorithms
Introduction to Parallel Computing: Design and Analysis of Algorithms
Efficient computation on sparse interconnection networks
Unifying themes for parallel selection
Derivation of randomized sorting and selection algorithms
Randomized parallel selection

Programming a Hypercube Multicomputer
Efficient parallel algorithms for selection and searching on sorted matrices
Finding the median
Random Data Accesses on a Coarse-grained Parallel Machine II
Load balancing on a hypercube
--TR

--CTR
Ibraheem Al-Furaih , Srinivas Aluru , Sanjay Goil , Sanjay Ranka, Parallel Construction of Multidimensional Binary Search Trees, IEEE Transactions on Parallel and Distributed Systems, v.11 n.2, p.136-148, February 2000
David A. Bader, An improved, randomized algorithm for parallel selection with an experimental study, Journal of Parallel and Distributed Computing, v.64 n.9, p.1051-1059, September 2004
Marc Daumas , Paraskevas Evripidou, Parallel Implementations of the Selection Problem: A Case Study, International Journal of Parallel Programming, v.28 n.1, p.103-131, February 2000
parallel algorithms;selection;parallel computers;coarse-grained;median finding;randomized algorithms;load balancing;meshes;hypercubes