Jon Louis Bentley is an American computer scientist. He received a ...
Robert Floyd was an eminent American Computer Scientist who receive...
#### Random Number Generators
Random numbers are a very interest...
If you want to appreciate just how hard it can be to generate rando...
Here is the 1986 column Bentley is talking about: [Literate Program...
Let me attempt to simplify Floyd's algorithm by first teasing out s...
This is the modern version of the Fisher-Yates shuffling algorithm,...
programming
by
jon
Bentley
with Special Guest Oyster Bob Floyd
pearls
A SAMPLE OF BRILLIANCE
How can a computer deal a poker hand? If we assign
each card in the deck its own integer between 1 and 52,
then we can make a hand from a “random sample” of
5 integers in the range 1..52, for instance,
4 a 31 46 47
(It is important that no number appear twice; holding
more than one ace of spades can jeopardize a card-
player’s health.) Random samples also arise in appli-
cations such as simulation, program testing, and
statistics.
The first section of this column reviews several
standard algorithms for random sampling. The next
section describes an elegant new algorithm by Floyd.
The third section describes how Floyd extends his
algorithm to generate random permutations.
A Sampling of Sampling Algorithms
Before we can generate a random sample, we have to
be able to generate a single random number. That topic
is worth a column all its own, so we will simply assume
that we have a function
RandInt
( I,
J
) that returns
an integer uniformly distributed over I..J
It is easy to generate a random sequence of
M
inte-
gers in the range l..N, so long as we don’t mind dupli-
cates
for I :=
1
to M do
print
RandInt(1, N)
When I set
M
to 12 and N to 3, that code produced the
sequence
313311121131
This very sequence might come in handy for your next
tough game of rock-paper-scissors; more serious appli-
cations include testing finite state machines and sorting
programs.
Many applications, though, require a random sample
without duplicates. A statistical analysis, for instance,
might waste work by observing the same element
twice. Such samples are often referred to as “samples
without replacement” or as “combinations”; for the
0 1987 ACM OOOl-0782/0900-0754 $1.50
remainder of this column, though, the word “sample”
will denote a random sample with no duplicates.
The December 1984 column described several sam-
pling algorithms.’ The majority of the algorithms were
based on the pseudocode in Algorithm S, which stores
the sample in the set S.
en
ALGORITHM S. A Typical Sampling Algorithm
If the set S is implemented correctly and if
RandInt
produces random integers, then Algorithm S produces a
random sample. That is, each M-element subset is pro-
duced with probability I/(%). The loop invariant is that
S always contains a random sample of Size: integers in
the range l..N.
There are four operations on S: initializing it to
empty, testing an integer for membership, inserting
a new integer, and printing all the members. The
December 1984 column sketched five data structures
that could be used to implement the set S: bit vectors,
arrays (sorted and unsorted), binary search trees, and
bins. In the May
1986
column, Don Knuth presented
a “literate” sampling program that implements S as
an ordered hash table; David Gries described the
same data structure using abstract data types in the
April 1987 column.
Floyd’s Algorithm
Algorithm S has many virtues: it is correct, fairly effi-
cient, and remarkably succinct. It is so good, as a mat-
ter of fact, that I thought one couldn’t do better. In the
’ An updated version of that column appears as Column 11 in my 1986
Addison-Wesley book
Programming Pearls.
It contains several algorithms not in
the original column.
754
Communications of
the ACM
September 1987 Volume 30 Number 9
Programming
Pearls
December 1984 column I therefore charged ahead and
described it in detail.
Unfortunately, I was wrong; fortunately, Bob Floyd
caught me sleeping. When he studied Gries’s imple-
mentation of Algorithm S in the April 1987 column, he
was bothered by a flaw that is displayed clearly when
M = N = 100. When Size = 99, the set S contains all but
one of the desired integers. The algorithm closes its
eyes and blindly guesses integers until it stumbles over
the right one, which requires 100 random numbers on
the average. That analysis assumes that the random
number generator is truly random; for some non-
random sequences, the algorithm won’t even terminate.
Floyd set out to find an algorithm that uses exactly
one call of
RandInt
for each random number in S. The
structure of Floyd’s algorithm is easy to understand re-
cursively: to generate a 5-element sample from l..lO,
we first generate a 4-element sample from 1..9, and
then add the fifth element. The recursive program is
sketched as Program Fl.
We can appreciate the correctness of Floyd’s algo-
rithm anecdotally. When M = 5 and N = 10, the algo-
rithm first recursively computes in S a 4-element
random sample from 1..9. Next it assigns to T a random
integer in the range l..lO. Of the 10 values that T can
assume, exactly 5 result in inserting 10 into S: the four
values already in S, and the value 10 itself. Thus ele-
ment 10 is inserted into the set with the correct
probability of 5/10. The next section proves that
Algorithm Fl produces each M-element sample with
equal probability.
function SampleCM, N)
if M = 0 then
return the empty set
else
S
:= SampleCM-I,
N-l)
T := RandInt(1, N)
if T is not in S then
insert T in S
else
insert N in S
return S
ALGORITHM Fl. Floyd’s Algorithm-Recursive
Because Algorithm Fl uses a restricted form of recur-
sion. Floyd was able to translate it to an iterative form
by introducing a new variable,
J.
(Problem 8 addresses
recursion removal in more general terms.) The result is
Algorithm F2, which is much more efficient than Algo-
rithm S, yet is almost as succinct. Problems 2 and 3
address the data structures that might be used to
implement the set S.
For those who might scoff that Algorithm F2 is “just
pseudocode,” Program F2 implements Floyd’s algorithm
in the Awk language. The June 1985 column sketched
initialize set S to empty
for J := N - M +
1
to N do
T
:= RandInt(1, J)
if T is not in S then
insert T in S
else
insert J in S
ALGORITHM F2. Floyd’s Algorithm-Iterative
that language, with particular emphasis on its associa-
tive arrays (which are used to implement sets in Pro-
gram F2). Complete with input and output, the program
requires only nine lines of Awk.
BEGIN {
m = ARGVLII; n = ARGV121
for Cj = n-m+l; j i= n; j++) {
t =
1
+ int(j * randO)
if (t in s) s[jl =
1
else s[tl =
1
for (i in s) print i
PROGRAM F2.
Floyd’s Algorithm in Awk
Random Permutations
Some programs that use a random sample require that
the elements of the sample occur in a random order.
Such a sequence is called a random permutation with-
out replacement of M integers from l..N. In testing a
sorting program, for instance, it is important that ran-
domly generated input occur in random order (if the
input were always sorted, for instance, it might not
fully exercise the sort code).
We could use Floyd’s Algorithm F2 to generate a
random sample, then copy it to an array, and finally
shuffle the elements of the array. This code randomly
scrambles the array X[l..M]:
for I
:= M downto 2 do
Swap(X[RandInt(l, I)], X[Il)
This three-step method uses 2M calls to
RandInt.
Floyd’s random permutation generator is very similar
to Algorithm F2; it is shown in Algorithm P, next page.’
To compute an M-element permutation from l..N, it
first computes an (M - I)-element permutation from
l..N - 1; a recursive version of the algorithm ehmi-
nates the variable I. The primary data structure of Al-
gorithm P, though, is a sequence rather than a set.
‘Flovd observes that Aleorithm P can be viewed as a nondeterministic aleo-
rithm and macroexpand~d into a backtracking procedure that generates aTI
permutations of M out of N. For details, see Floyd’s “Nondeterministic Algo-
rithms” in 1.
ACM 14. 4
(Oct. 1967).
636-644.
September 1987 Volume 30 Number 9
Communications
of
the ACM
755
Programming Pearls
initialize sequence S to empty
for J
:= N - M +
1
to N do
T := RandInt(1, J)
if T is
not in
S then
prefix T to S
else
insert J in S after T
ALGORITHM P. Floyd’s Permutation Algorithm
Problem 5 shows that Algorithm P is remarkably effi-
cient in terms of the number of random bits it uses;
Problem 6 deals with efficient implementations of the
sequence S.
We can get an intuitive feeling for Algorithm P by
considering its behavior when M = N, in which case it
generates a random permutation of N elements. It iter-
ates J from
1
to N. Before execution of the loop body,
S is a random permutation of the integers in
1.J
- 1.
The loop body maintains the invariant by inserting J
into the sequence; 1 is the first element when T = J,
otherwise ] is placed after one of the J - 1 existing
elements at random.
In general, Algorithm P generates every M-element
permutation of
1
. . N with equal probability. Floyd’s
proof of correctness uses the loop invariant that after
i times through the iteration, J = N - M + i and S can
be any permutation of i distinct integers in
l..J,
and
that there is a single way that the permutation can
be generated.
Doug McIlroy found an elegant way to phrase Floyd’s
proof: there is one and only one way to produce each
permutation, because the algorithm can be run back-
ward. Suppose, for instance, that
M
= 5, N =
10,
and
the final sequence is
72915
Because
10
(the final value of J) does not occur in S, the
previous sequence must have been
2915
and
RandInt
returned
T
= 7. Because 9 (the relevant
value of 1) occurs in the J-element sequence after 2, the
previous T was 2. Problem 4 shows that one can simi-
larly recover the entire sequence of random values.
Because all random sequences are (supposedly) equally
likely, all permutations are also.
We can now prove the correctness of Algorithm F2
by its similarity to Algorithm P. At all times, the set S
in Algorithm F2 contains exactly the elements in the
sequence S in Algorithm P. Thus each final M-element
subset of l..N is generated by M! random sequences,
so all occur equiprobably.
Principles
Algorithm S is a pretty good algorithm, but not good
enough for Bob Floyd. Not content with its inefficiency,
he developed optimal algorithms for generating random
samples and random permutations. His programs are a
model of efficiency, simplicity, and elegance. The “Fur-
ther Reading” section sketches some of the methods
that Floyd uses to achieve such programs.
Problems
1.
2.
3.
4.
5.
6.
7.
8.
How do the various algorithms behave when the
RandInt
procedure is nonrandom? (Consider, for
instance, generators that always return 0, or that
cycle over a range that is much smaller than or
much greater than
M or N.)
Describe efficient representations for the set S in
Algorithm F2.
Algorithms S and F2 both use a set S. Is a data
structure that is efficient in one algorithm neces-
sarily efficient in the other?
Complete the proof of correctness of Algorithm P by
showing how to compute from a final permutation
the sequence of random values that produced it.
How many random bits does Algorithm P consume?
Show that this is close to optimal. Perform a similar
analysis for Algorithm F2. Can you find algorithms
that are more efficient?
Describe representations for the sequence S such
that Algorithm P runs in
O(M)
expected time or in
O(M
log
M)
worst-case time. Your structures should
use
O(M)
worst-case space.
Implement Floyd’s algorithms in your :favorite pro-
gramming language.
Algorithm F2 is an iterative version of the recursive
Algorithm
Fl.
There are many general methods for
transforming recursive functions to equivalent iter-
ative programs (one method is often illustrated on
a recursive factorial function). Consider a recursive
function with this form
function A(M)
if M = 0 then
return X
else
S := ACM-I)
return G(S, M)
where
M
is an integer, S and X have the same
type
T,
and function G maps a T and an integer
to a
T.
Show how the function can be transformed
to this iterative form
function B(M)
S
:= x
for J :=
1
to M do
S := G(S, J)
return S
756 Communications
of
the ACM
September 1987
Volume 30 Number 9
Programming Pearls
Further Reading
Robert W. Floyd received the ACM Turing Award in
1978. In his Turing lecture on “The Paradigms of
Programming,” Floyd writes, “In my own experience
of designing difficult algorithms, I find a certain
technique most helpful in expanding my own capa-
bilities. After solving a challenging problem, I solve
it again from scratch, retracing only the
insight
of the
earlier solution. I repeat this until the solution is as
clear and direct as I can hope for. Then I look for a
general rule for attacking similar problems, that
would have led me to approach the given problem in
the most efficient way the first time. Often, such a
rule is of permanent value.”
Floyd’s key rule in this problem was, in his own
words, to “look for a mathematical characterization
of the solution before you think about an algorithm
to obtain it.” His key insight dealt with the probabil-
ity of the algorithm generating any particular subset,
When Floyd calculated the probabilities of certain
events in Algorithm S, he noticed that the denomi-
nators were powers of N, while the denominators in
the solution were falling factorials. His algorithms
use a simple structure to generate the correct prob-
abilities. When Floyd finally conceived Algorithm P,
he recalls, “I knew it was right even before I
proved it.”
Floyd’s 1978 Turing lecture was published origi-
nally in the August
1979 Communications.
It also
appears in the
ACM Turing Award Lectures: The First
Twenty Years: 1966-1985,
which was recently pub-
lished by the newly formed ACM Press (a collabora-
tion between ACM and Addison-Wesley).
Solutions to July’s Problems
1. The problem can be rephrased as asking how many
assignments this routine makes after the array
X[l..N] has been sprinkled with random real
numbers.
Max :=
X[l]
for I : =
2
to N do
if X[I]
> Max then
Max := X[I]
A simple argument assumes that
if
statements are
executed about half the time, so the program will
make roughly N/2 assignments. I profiled the pro-
gram for ten runs with N = 1000, and the number
of assignments were, in sorted order: 4, 4, 5, 5,
6, 7, 8, 8, 8, 9. In Section 1.2.10 of
Fundamental Algo-
rithms,
Knuth shows that the algorithm makes
HN
- 1 comparisons, on the average, where
HN =
1+1/2 +1/3 + ...
+ l/N is the Nth harmonic
number. For N = 1000, this analysis gives an expec-
tation of 6.485; the average of the ten experiments
is 6.4.
This C program implements a Sieve of Eratosthenes
for computing all primes less than
n.
main( 1
1
int i, p, n;
char
x[1000023;
n = 100000;
for (i =
1;
i <= n; i++)
x[i] =
1;
x[ll = 0; x[n+ll =
1;
p = 2;
while (p <= n) {
printf ( “%d\n” , p) ;
for (i=2*p; i<=n; i+=p)
x[il = 0;
do
1
1
100000
1
1
9593
9592
9592
256808
9592
99999
99999
I
p++:
while (x[p] == 0);
1
The profile shows that there are 9592 primes less
than 100,000, and that the algorithm made about
2.57N assignments. In general, the algorithm makes
about N log log N comparisons; the analysis involves
the density of prime numbers and the harmonic
numbers mentioned in Solution 1. For faster imple-
mentations of prime sieves,
see
the
Communications
papers by Mairson (September 1977), Gries and
Misra (December 1978) and Pritchard (January
1981).
A simple statement-count profiler increments a
counter at each statement. One can decrease mem-
ory and run time by making do with fewer
counters. For instance, one might associate a coun-
ter with each basic block in the program’s
flow graph. One can further reduce counters by
“Kirchoff’s first law”; if you have a counter for an
if-then-else statement and one for the
then
branch, then you don’t need one for the else
branch.
Problem P6 is a correct program that is hard to
prove correct. The for loop in function prime
could give an infinite loop. To show that it always
terminates, one must prove that if
P
is a prime,
then there is another prime less than
P2.
For Correspondence: Jon Bentley, AT&T Bell Laboratories, Room X-317,
600 Mountain Avenue. Murray Hill, NJ 07974.
Permission to copy without fee all or part of this material is granted
provided that the copies are not made or distributed for direct commer-
cial advantage, the ACM copyright notice and the title of the publication
and its date appear, and notice is given that copying is by permission of
the Association for Computing Machinery. To copy otherwise, or to
republish, requires a fee and/or specific permission.
September 1987 Volume 30 Number 9
Communications of the ACM
757
Here is the 1986 column Bentley is talking about: [Literate Programming](http://fermatslibrary.com/s/literate-programming).
In that column Bentley defines a *"Literate"* program as a piece of code that reads like literature.
This is the modern version of the Fisher-Yates shuffling algorithm, widely used in computer science. https://en.wikipedia.org/wiki/Fisher%E2%80%93Yates_shuffle
If you want to appreciate just how hard it can be to generate randomness, you might want to try to play a [game of Rock-Paper-Scissors with this bot](http://www.nytimes.com/interactive/science/rock-paper-scissors.html?_r=0)
)](http://i.imgur.com/GelhVnF.png)
Jon Louis Bentley is an American computer scientist. He received a B.S. in Mathematical sciences from Stanford in 1974 and a Ph.D from the University of North Carolina at Chapel Hill in 1976. Bentley interned at the Xerox Palo Alto research and later worked at Bell Laboratories.
Bentley wrote the Programming Pearls column for the Communications of the ACM magazine. This column usually featured a relatively short, but interesting problem or concept. The entries of that column were later compiled into two books of the same name.
Let me attempt to simplify Floyd's algorithm by first teasing out similarities between Algorithm F and Algorithm S, then I'll show the key differences.
First, both F and S:
- Generate a random number, T
- Test for presence of T in the invariant
- Append T to invariant if not present
Seen this way they both do the same work. Now let's understand the genius of Algorithm F. For Algorithm S, when the test for presence of T succeeds, i.e. T is already present in invariant, it restarts the process: generate another number. Algorithm F doesn't, and that's the key departure. For Algorithm F, if the test succeeds, i.e. T is already present in the invariant, it inserts another number that is **known** to not be in the invariant. This number is the upper bound of the random numbers.
Robert Floyd was an eminent American Computer Scientist who received the Turing Award in 1978.
Born in New York on June 8, 1936, Floyd was recognized as a child prodigy at age 6; he skipped three grades and finished high school at 14. A scholarship allowed him to study at the University of Chicago, where he received a bachelor's degree in liberal arts in 1953 at age 17. After that he supported himself and earned another bachelor's degree in physics in 1958. He never received a Ph.D, but still managed to obtain full professorship at Stanford at the age of 32.
Floyd made numerous contributions to Computer Science. He is well know for the Floyd-Warshall Algorithm for finding shortest paths in weighted graphs.
#### Random Number Generators
Random numbers are a very interesting topic. Generating truly random numbers turns out to be a fairly hard problem. Having the ability to generate randomness is essential to areas such as computer simulation and cryptography.
There are two principal types of random numbers:
- **Pseudo-random numbers**: These are generated based solely on deterministic logic. In other words, computational algorithms are used to produce long sequences of apparently random numbers, which are in fact determined by a short initial value know as seed. These are appropriate for applications that only need a modest amount of unpredictability (e.g. Random Fact of the day).
- **"True" Random numbers**: These are produced by measuring some physical phenomena that is recognized to be random. Example sources include atmospheric noise, cosmic background radiation and radioactive decay.