I've been looking though my books, and I can't find any mention of the
shuffler which I came up with. Could someone tell me if my idea is the
same as some already existing algorithm, or can I name it "the goldberg
shuffler"? :)
Here it is in psuedo-code:
array shuffle(unordered_set x) {
unordered_set left, right;
if( x.size <= 1 ) return x.to_array();
for i in x:
if( rand_bit ) left.insert(i)
else right.insert(i)
return shuffle(left).append( shuffle(right) );
}
Here it is in [untested] C code:
void swap(int *a, int *b) { int t = *a; *a = *b; *b = t; }
void shuffle(int *start, int *end) {
while( end-start > 2 ) {
int *front = start, *back = end;
while( front < back )
rand_bit() ? swap( front, --back ) : ++front;
if( front-start < end-back ) {
shuffle( start, front );
start = front;
} else {
shuffle( back, end );
end = back;
}
}
if( end-start == 2 && rand_bit() )
swap(start, start+1);
}
/*
int x[5] = {1,2,3,4,5};
shuffle(&x[0], &x[5]);
*/
It's recursive, but by changing the tail call into a loop, and arranging
to only recurse on the smaller half of the array, the depth of recursion
is limited to up to a maximum of log2(N).
I'm trying to figure out how many calls to rand_bit are made on average
for a size N array of ints, in big-O notation, but errno=ENOCAFFEINE.
PS: Hey Tom, you took down your 'Remote Win32 C Compiler'! How come?
--
print reverse( ",rekcah", " lreP", " rehtona", " tsuJ" )."\n";
It is obvious that any sorting algorithm must be able to realize all
possible permutations on its input. You can turn an abitrary
comparison-based sorting algorithm into a shuffling algorithm by
randomizing the results of comparisons. The problem is that the
resulting permutations are generally not uniformly distributed--some
permutations occur more often than others.
--
Nicol So
Disclaimer: Views expressed here are casual comments and should
not be relied upon as the basis for decisions of consequence.
Benjamin Goldberg wrote:
>
> Ok, anti sorting is just another name for shuffling :)
[snip]
I haven't closely studied your code. But you have to
demonstrate two things in my view: (1) it has no bias,
(2) it is more efficient than known methods.
M. K. Shen
Several weeks ago, I suggested a very similar shuffling function to
Dabiker. Since then, I have implemented the algorithm in C++, and I
have trimmed that version to C. (The C++ code has been tested; the C
code has not.) I'm not that smart, so I doubt that it's an original
idea with me. It seems to have its Genesis in a sort-merge algorithm
that I invented more than 20 years ago.
I like its elegence, but am unable to prove that the algorithm is
capable of producing every possible permutation with equal probability
(given a uniformly distributed PRNG of large period). Subjectively, it
is very fast even with very large tables and it appears to produce very
good shuffling in one pass.
In the C code below, "p_aulTable" is a simple array of unsigned long.
Each element of the array is initialized to contain the value of the
subscript that points to it. The array is shuffled "in place", meaning
that no temporary storage is required. "ClockPRNG()" is a Boolean
function that clocks a PRNG and returns a single bit--true or
false--depending on the result.
void MakeRandomPermutation(
unsigned long p_aulTable[],
int nBase, int nLen)
{
if (nLen >= 2)
{
int nSplit, ii;
nSplit = nLen / 2;
MakeRandomPermutation(p_aulTable,
nBase, nSplit);
MakeRandomPermutation(p_aulTable,
nBase + nSplit, nLen - nSplit);
for (ii = 0; ii < nLen; ii++)
{
int jj;
jj = (nLen - 1) - ii;
if (ClockPRNG())
{
unsigned long nTemp;
nTemp = p_aulTable[ii];
p_aulTable[ii] = p_aulTable[jj];
p_aulTable[jj] = nTemp;
}
}
}
}
Did you mean to use 'nBase' (instead of 0) as the initial value of 'ii'?
> {
> int jj;
>
> jj = (nLen - 1) - ii;
Did you mean 'jj = (nBase + nLen - 1) - ii'?
> if (ClockPRNG())
> {
> unsigned long nTemp;
>
> nTemp = p_aulTable[ii];
> p_aulTable[ii] = p_aulTable[jj];
> p_aulTable[jj] = nTemp;
> }
> }
--
> The problem is that the
> resulting permutations are generally not uniformly distributed--some
> permutations occur more often than others.
This one works.
Proof, by induction:
We can observe that for one element, all permutations are equally
likely. Now assuming the Goldberg shuffler works for all integers
less than an arbitrary but fixed n > 1, we show it also works for n.
We are considering the case that the algorithm does terminate, so we
can assume we apply the first partitioning step until we get at
least one element on each side. Thus there are 2^n - 2 candidates
for the first partitioning.
Consider any permutation of the n elements, call it P. That
permutation can be induced by a number of different combination of
the first partitioning into left and right sides, followed by
recursive sorting each side. In particular, for each i, 1 < i < n,
permutation P is induced by:
Partitioning the first i elements of P on the left, then
shuffling the left i elements to match the first i elements of
P, and the right n-i elements to match the last n-i elements
of P.
We can see from the pseudo-code that all 2^n - 2 partitionings are
equally likely (the two that leave all the elements on the same side
are impossible since we repeat until we get a non-trivial
partitioning). Thus each partitioning has probability 1/(2^n - 2).
By our induction hypothesis, the algorithm works for all partitions
of size less than n, thus a recursive call on a partition of size i
< n has a 1/i! chance of inducing each particular permutations of
the i elements. Thus given one of the n-2 workable partitionings,
the probability of subsequent sorting inducing P is 1/(i! * (n-i)!)
Thus the probability of inducing permutation P is
P = sum for i in [1..n-1]: 1/(2^n - 2) * 1/(i! * (n-i)!).
That is, we're summing over the n-2 partitionings (for which P is
possible). The terms we sum are the probability of that
partitioning times the probability that given that partitioning, the
two halves will shuffle in just the right way to induce P.
Taking the constant factor outside the sum:
P = 1/(2^n - 2) * (sum for i in [1..n-1]: 1/(i! * (n-i)!)
We adjusting the index range ,which creates two new terms, and then
subtract these terms back out to get:
P = 1/(2^n - 2) * ((sum for i in [0..n]: 1/(i! * (n-i)!) - 2/n!)
Now we need,
Lemma:
sum for k in [0..n]: 1/(k! * (n-k)!) = 2^n / n!.
Proof:
From the binomial theorem:
sum for k in [0..n]: n!/(k! * (n-k)!) = 2^n
Dividing both sides by n! gives us the Lemma.
Applying the Lemma gives us:
P = 1/(2^n - 2) * ((2^n / n!) - 2/n!)),
= 1/(2^n - 2) * (2^n - 2)/n!,
= 1/n!.
Thus for an arbitrary n-element permutation P, the probability that
the Goldberg shuffler induces P is 1/n!, which completes the proof.
And it's just that easy. There may be a much simpler way to prove the
algorithm is correct, a fact that is almost obvious by symmetry. But I
wrote a complicated proof because I did not have time to write a simple
one.
--Bryan
shuffle gets called 2n times and has to produce 2n unordered sets.
The random number generator is called O(n log n) times. It does
appear to have the desired properties (all permutations are possible,
all permutations are equally likely), but at a very high runtime
cost of both memory and operations.
For comparison, Fisher-Yates is O(n) and only calls the random number
generator O(n) times. It can be coded in a only a few lines in C
or Perl. It is not recursive either.
> print reverse( ",rekcah", " lreP", " rehtona", " tsuJ" )."\n";
Then you might be interested in this example:
perldoc -q shuffle
Trammell
--
-----|----- hud...@swcp.com H 240-476-1373
*>=====[]L\ Trammel...@celera.com W 240-453-3317
' -'-`- http://www.swcp.com/~hudson/ KC5RNF
Bryan Olson has shown a proof indicating that it has no bias.
However, I cannot show that it is more efficient than known methods,
because you haven't said what kind of efficiency you want.
Number of swaps? Number of random bits used? Amount of computation?
The most efficient algorithm in terms of number of random bits is to
calculate a single number in the range [0,N!), and then use that number
in place of a call to a prng in a manner similar to a fisher-yates
shuffle.
Nothing can possibly beat that in terms of random bits consumed or in
terms of number of swaps. The cost of this is that you have to do large
integer arithmetic.
A fisher-yates uses the same number of swaps, but since the prng must
select random numbers in small ranges, it will be wasteful of random
bits if your rand_mod function is truly unbiased.
I don't know whether it's more efficient in terms of bits consumed than
fisher-yates. To be honest, I don't know how to formulate how efficient
f-y is, let alone my own algorithm.
The reasons are as follows:
Consider that f-y uses a rand_mod function. When f-y makes a call of
rand_mod(4), it consumes exactly two bits to produce an unbiased number
in the range [0,4). When f-y makes a call rand_mod(3), it will consume
an unknown number of bits to produce an unbiased number in the range
[0,3). When shuffling an N element array, f-y makes calls like
rand_mod(N), rand_mod(N-1) ... rand_mod(2).
I don't know how many bits get consumed by f-y on average because the
variability of how many bits get consumed by rand_mod().
And my own algorithm is rather more complicated to analyse than
rand_mod.
In practice, the input 'unordered_set' would actually be the same data
type as the output 'array', and the elements would be shuffled in-place,
as I showed in the C code. So there's really only a trivial cost in
terms of memory.
> For comparison, Fisher-Yates is O(n) and only calls the random number
> generator O(n) times. It can be coded in a only a few lines in C
> or Perl. It is not recursive either.
F-Y needs to get unbiased random numbers in varying ranges... that is,
you need one number from each of the ranges: [0,2), [0,3), [0,4) ... up
to [0,N).
So the N calls to the random-number-in-range generator in F-Y are doing
much more work than the calls to the random-bit-generator that mine is
doing.
In fact, it's possible that the total number of random bits needed by
F-Y is also O(N log N), if the random-number-in-range generator produces
unbiased numbers.
> > print reverse( ",rekcah", " lreP", " rehtona", " tsuJ" )."\n";
>
> Then you might be interested in this example:
>
> perldoc -q shuffle
I know of it, but IMHO, mine has a higher coolness factor :)
--
Benjamin Goldberg wrote:
One thing to remember is that to prevent a sort, there must
be a relationship present to prevent a particular sort, or
the random ordering has caused the sort.
A shuffle is applied to cause the shuffling relation
to the ordering. It is a good relation in that the
outcome has no relation to the initial state of
the ordering. This is the cryptographic
application of the shuffle, simple complexity.
And safe also!
The shuffle in a compact form is very
efficient, except yours is a simulation of the
shuffle code. The real shuffle causes the cpu
to shuffle directly on the memory states.
A thing done using a particular property
of C arrays called a implied definition.
The real shuffling then, is a means of causing
the array's runtime existence to cause the
relations outcome.
And so the "real" shuffle is then performed on the
definition of the array and not the contents
of the array. A simulation is simply one abstraction
shy of the real thing. Although, it is a good
demonstration in the fashion you listed.
Douglas Eagleson
Gaithersburg, MD USA
There is one way in which your algorithm is not as good as Fisher-Yates
shuffle: worst-case time-complexity. In the worst case your algorithm
takes O(n^2) time to execute.
> There is one way in which your algorithm is not as good as Fisher-Yates
> shuffle: worst-case time-complexity. In the worst case your algorithm
> takes O(n^2) time to execute.
No, the worst case doesn't terminate.
I don't think that's a significant criticism. Any method for generating
permutations (of more than two elements) from uniform independent bits
must either have a bias in the permutations it generates, or be a non-
terminator.
Furthermore, the chance of a really bad case is negligible, and
independent of the data set.
--Bryan
I'm not quite sure what I meant. I think you're right but I'll check it
out more carefully sometime today. Thanks for the flag.
Nice catch. I forgot the two degenerate cases with one the partitions
being empty, which do not occur in quicksort.
> I don't think that's a significant criticism. ...
That's not intended to be.
If you intended to write what I thought you did, your shuffling
algorithm has a flaw--it cannot implement all possible permutations.
Your shuffling algorithm bears some superficial resemblance to
mergesort. It starts by dividing the input array into two (almost) equal
halves. It then recursively applies itself to permute the two subarrays,
and somehow combines the results into one. As I noted in another
message, any comparison-based sorting algorithm can be converted into a
shuffling algorithm (with no guarantee that it will permute its input in
all possible ways with equal probability). Your algorithm is not exactly
the result of a converted mergesort.
There is a problem with the way you combine the permuted subarrays.
Let's consider a top-level call to your MakeRandomPermutation procedure.
Notice that after the upper and lower halves of the array are permuted,
each position in the lower half can only be swapped with one
*corresponding* position (which I'll call its peer) in the upper half.
Now consider the case in which two values, say X and Y, both occur in
the same half of the original array. X and Y can never end up in peer
positions in the final array.
Sorry about the nBase bug.
You may want to look at this again (in the corrected version below). If I
understood your previous paragraph, you are saying that, given the following
sequence {1,2,3,4}, 1 and 4 are in "peer positions" and 2 can never be
placed in the position occupied by 4. During the first shuffle, {1,2} could
become {2,1} and following the outer shuffle, 2 could be swapped with 4
leaving the sequence {4,1,3,2}. I am pretty sure that the algorithm can
generate all permutations, but I can't say that it will generate each with
equal probability, particularly when the Table length is not an integral
power of 2.
Corrected version follows. nSplit is rounded up for odd length arrays.
Shuffle only processes half (approximately) of array, allowing one shuffle
between peer positions instead of two.
void MakeRandomPermutation(unsigned long p_aulTable[],
int nBase, int nLen)
{
if (nLen >= 2)
{
int nSplit, ii;
nSplit = (nLen + 1) / 2;
MakeRandomPermutation(p_aulTable,
nBase, nSplit);
MakeRandomPermutation(p_aulTable,
nBase + nSplit, nLen - nSplit);
for (ii = 0; ii < nSplit; ii++)
{
int jj, kk;
kk = nBase + ii;
jj = nBase + (nLen - 1) - ii;
if (ClockPRNG())
{
unsigned long nTemp;
nTemp = p_aulTable[kk];
p_aulTable[kk] = p_aulTable[jj];
p_aulTable[jj] = nTemp;
}
}
}
}
That's not exactly the problem I described.
Let a[4] = {1, 2, 3, 4}. '3' and '4' are elements in the same half. Now
a[0] and a[3] are peer positions in the array. No matter what random
choices are made, '3' and '4' (as a pair) will not occupy those two
positions. In other words,
{3, ?, ?, 4}, and
{4, ?, ?, 3}
will never be the final value of the array. (This is because two
elements in peer positions must come from different halves of the
original array, *and* permuting the sub-arrays do not move elements from
one half to the other).
If your idea is to turn mergesort into a shuffling algorithm, why not do
it the straightforward way i.e. just randomize the results of
comparisons in a standard mergesort?
>
>I've been looking though my books, and I can't find any mention of the
>shuffler which I came up with. Could someone tell me if my idea is the
>same as some already existing algorithm, or can I name it "the goldberg
>shuffler"? :)
>
>Here it is in psuedo-code:
>
>array shuffle(unordered_set x) {
> unordered_set left, right;
> if( x.size <= 1 ) return x.to_array();
> for i in x:
> if( rand_bit ) left.insert(i)
> else right.insert(i)
> return shuffle(left).append( shuffle(right) );
>}
>
[snip happens]
>It's recursive, but by changing the tail call into a loop, and arranging
>to only recurse on the smaller half of the array, the depth of recursion
>is limited to up to a maximum of log2(N).
>
>I'm trying to figure out how many calls to rand_bit are made on average
>for a size N array of ints, in big-O notation, but errno=ENOCAFFEINE.
>
I've never seen this method before.
Might I suggest the following as an improvement;
>array shuffle(unordered_set x) {
> unordered_set left, right;
> if( x.size <= 1 ) return x.to_array();
if ( x.size == 2 ) {
if( rand_bit ) return x(1) append x(2)
return x(2) append x(1)
}
> for i in x:
> if( rand_bit ) left.insert(i)
> else right.insert(i)
> return shuffle(left).append( shuffle(right) );
>}
That would reduce the number of random bits needed.
Scott Nelson <sc...@helsbreth.org>
[...]
> In fact, it's possible that the total number of random bits needed by
> F-Y is also O(N log N), if the random-number-in-range generator produces
> unbiased numbers.
The average-case time is definitely Theta(N log N). Analysis of
randomized-quicksort applies. The number of random bits used
must therefore also be O(N log N), since they're generated one
at a time.
--Bryan
You're absolutely right. The algorithm only generates 16 out of a possible
24 permutations. Thanks for pointing out the flaws.
Yes, that is a good idea. In fact, I do that in the C version.
> > for i in x:
> > if( rand_bit ) left.insert(i)
> > else right.insert(i)
> > return shuffle(left).append( shuffle(right) );
> >}
>
> That would reduce the number of random bits needed.
If my C code were identical to my psuedo-code, I'd have had:
void shuffle(int *start, int *end) {
while( end-start > 1 ) {
....
}
}
But instead I have:
void shuffle(int *start, int *end) {
while( end-start > 2 ) {
....
}
if( end-start == 2 && rand_bit() )
swap(start, start+1);
}
, which is essentially the same as your suggestion.
How would one go about proving that F-Y will use O(N log N) bits?
> How would one go about proving that F-Y will use O(N log N) bits?
Ah, I see I was answering the wrong question.
When did F&Y publish the algorithm? Knuth credits it to Moses and
Oakford, and also to Durstenfeld.
Using a uniform random integer generator that is optimal in the number
of bits used, we can show that generating a uniform selection in [1..n]
takes Theta(lg n) bits in the average case. Such generators have
appeared here several times. My code for it is here:
http://groups.google.com/groups?selm=32F15F82.4F09%40ix.netcom.com
Note that on each test of:
if( rand < target_range )
we have a better than 1/2 chance of terminating. Before the first such
test we generate ceil(lg n) bits, and we never generate more that
ceil(lg n) bits between tests. Thus the expected number of bits used,
call it 'b' is bounded by:
ceil(lg n) < b < ceil(lg n) + ceil(lg n)/2 + ceil(lg n)/4 + ...
ceil(lg n) < b < ceil(lg n) * (1 + 1/2 + 1/4 + ...)
ceil(lg n) < b < 2 ceil(lg n)
b = Theta(lg n)
In the entire algorithm the average number of random bits generated,
call it B is:
B = sum for i in [1..N]: Theta(lg i)
random bits. Theta implies a constant factor bound above and below,
so we can take it outside the summation.
B = Theta(sum for i in [1..N]: lg i)
There are many known ways to prove the summation is Theta(N lg N), such
as showing that just the last half of the terms sum to more than N/2 *
lg(N/2) which is already Theta(N lg N), and we can bound the entire
series from above by N times it's largest term.
B = Theta(Theta(N lg N))
B = Theta(N lg N)
--Bryan
That's actually quite simple. Based on what you wrote in another
message, it seems that you're not sure how many random bits are consumed
by the operation you named rand_mod(). There are smart ways to implement
rand_mod() that wastes very few bits, but they are not necessary for
achieving the O(N log N) (average-case) bound for the overall algorithm.
Let's consider the following straightforward method for implementing
rand_mod(m):
1. Generate ceiling(log_2(m)) random bits.
2. If the random bits encode an integer in the desired range,
output it.
3. Otherwise, repeat steps 1 to 2 until success.
We note that the probability of success in step 1 is at least 0.5. You
can upper bound the expect # of tries before success by an infinite
series, which converges (quickly) to a constant. (The details are easy
and omitted here). You can then upper bound the number of random bits
consumed per rand_mod() call by (C * log N), where C is a positive
constant.
The O(N log N) bound for the expected # of random bits for the entire
algorithm then follows almost immediately.
Thank you, very cool. I suppose I'll have to figure out on my own which
one has a smaller constant factors. I *suspect* that mine will use
fewer bits on average, but I need to do some calculations to be sure.
> I suppose I'll have to figure out on my own which
> one has a smaller constant factors. I *suspect* that mine will use
> fewer bits on average, but I need to do some calculations to be sure.
Hmmm, I'd bet the other way, but I don't really know. The F-Y method is
also easily optimized for bit usage by generating more than one
iteration's bits at a time. In the extreme, generate a uniform choice
in the integers [1..N!], and that's optimal in bit usage.
If you want to beat F-Y at something, look at the time to randomly
permute an array of integers that's several times larger than the
physical memory available.
--Bryan
Ok, that's easy, using a variant of my method.
Here's a perl program which reads in newline seperated records from
stdin, shuffles them, and prints the results to stdout.
#!/usr/bin/perl
use IO::File;
sub permute {
my $filehandle = shift;
while( 1 ) {
my $left = IO::File->new_tempfile;
my $right = IO::File->new_tempfile;
my $items = 0;
while( my $item = <$filehandle> ) {
print { rand()>0.5 ? $left : $right } $item;
++$items;
}
close $filehandle; # frees it from filesystem if it was a temp file.
seek $_, 0, 0 for $left, $right;
if( $items <= 1 ) { # there was only one item
for my $fh ( $left, $right ) {
print STDOUT while <$fh>;
}
return;
}
if( -s $left < -s $right ) {
permute($left);
$filehandle = $right;
} else {
permute($right);
$filehandle = $left;
}
}
}
permute(\*STDIN);
__END__
This could be changed to C fairly easily...
To change it to work on 4 byte records, rather than newline delimited
records, add $/ = \4 at the top, and if you're on dos/win binmode on all
filehandles (binmode isn't needed on *nix).
There will be at most ceil(log2(N))+1 filehandles open... This could be
changed by using normal files (rather than ones which get unlinked
immediately after opening), which get closed after they're written to,
and reopened just before being read from: This way, there'd be at most 3
filehandles open (with the added cost of closing/opening files).
Whoops, that shoulda been new_tmpfile, sorry.
Otherwise it works ok.
> Bryan Olson wrote:
>
>>If you want to beat F-Y at something, look at the time to randomly
>>permute an array of integers that's several times larger than the
>>physical memory available.
>
> Ok, that's easy, using a variant of my method.
I just meant in ordinary virtual memory.
--Bryan
I did and it works like a champ. Preliminary analysis with a very small
number of random permutations suggests that certain permutations may show up
more often than others. I am in the process of writing a test harness that
will generate statistics on thousands of random permutations (for small
Tables).
By the way, it (C version) generates a randomized permutation of 32768
elements in about 0.25 second on a 300 MHz Pentium.
Hmm? When you said, "the physical memory available," I assumed you
meant bigger than main memory plus swap space.
F-Y needs everything in one big array; you could easily do it by using a
file is a virtual array, by fseek()ing forwards and backwards... it
won't be especially efficient, especially if the records are variable
sized, but there's nothing at all preventing you from doing it.
The variant of my method that I showed demonstrates that one can
partition the data using *just* temp files, keeping almost nothing in
memory but a single record, and one only needs to move forward through
those files, only seeking backwards *once* to get back to the beginning.
I will admit, it's not especially efficient for small files, *but*, one
can always use a hybrid approach; specifically, one can do partitioning
using files with my approach as long as those files are too big for
memory, and then switch to F-Y once they're small enough to do
efficiently in-memory.
Herewith, the modification based on your suggestion (and check out the test
results below):
#define HIGH_BIT (RAND_MAX & ~(RAND_MAX >> 1))
void random_permutation(int nA[], int nLength, int nWork[])
{
if (nLength >= 2)
{
int nSplit, ii, jj, kk;
nSplit = (nLength + 1) / 2;
random_permutation(&nA[0], nSplit, nWork);
random_permutation(&nA[nSplit], nLength - nSplit, nWork);
ii = 0;
jj = nSplit;
kk = 0;
while (kk < nLength)
{
if (ii < nSplit)
{
if (jj < nLength)
{
if (rand() & HIGH_BIT)
{
nWork[kk++] = nA[ii++];
}
else
{
nWork[kk++] = nA[jj++];
}
}
else
{
nWork[kk++] = nA[ii++];
}
}
else
{
nWork[kk++] = nA[jj++];
}
}
for (ii = 0; ii < nLength; ii++) nA[ii] = nWork[ii];
}
}
And this table shows the result of calling "random_permutation" for
2,400,000 consecutive random permutations. Since there are 24 (4x3x2x1)
permutations, you would expect each permutation to be formed about 100,000
times. The randomizer is the high-order bit of the rand() function from the
MS C-RTL that ships with Visual Studio 6.
( 0, 1, 2, 3) 100317
( 0, 1, 3, 2) 100245
( 0, 2, 1, 3) 100044
( 0, 2, 3, 1) 100045
( 0, 3, 1, 2) 99824
( 0, 3, 2, 1) 99935
( 1, 0, 2, 3) 99498
( 1, 0, 3, 2) 100347
( 1, 2, 0, 3) 99819
( 1, 2, 3, 0) 100149
( 1, 3, 0, 2) 99433
( 1, 3, 2, 0) 99479
( 2, 0, 1, 3) 100027
( 2, 0, 3, 1) 100122
( 2, 1, 0, 3) 100528
( 2, 1, 3, 0) 100330
( 2, 3, 0, 1) 100047
( 2, 3, 1, 0) 99994
( 3, 0, 1, 2) 99869
( 3, 0, 2, 1) 100700
( 3, 1, 0, 2) 99454
( 3, 1, 2, 0) 99829
( 3, 2, 0, 1) 100022
( 3, 2, 1, 0) 99943
> Herewith, the modification based on your suggestion
[code snipped]
> And this table shows the result of calling "random_permutation" for
> 2,400,000 consecutive random permutations.
[... amazingly uniform counts snipped]
Hmmm, I find that hard to believe. My theory tells me that not all
merges are equally likely. Merging single-element lists is perfect, but
when merging two-element lists, there are six possible merges, two with
probability 1/4, and four with probability 1/8.
Can you show the code that ran the tests?
--Bryan
Hmm. Any algorithm which produce all permutations without bias will
take a non-determinsitic amount of time (in an extreme case, for certain
outputs of the rng, it will loop 'forever'). In the F-Y shuffle, this
behavior is hidden inside of rand_mod. In my shuffle, this occurs if
rand_bit consistantly partitions all elements to the left side, or all
to the right side.
What kind of rng output causes your algorithm to loop forever?
That is non-obvious and if it can be proven, the proof is likely to be
well beyond my ability to comprehend. I simply don't believe that's a
true statement.
> What kind of rng output causes your algorithm to loop forever?
There is no sequence from the rng that will cause my algorithm to run
forever. You can easily satisfy yourself that this is true. The
looping condition depends on filling the "work" array. When the "work"
array is full, the contents are transferred back to the input array and
the function terminates. The rng only determines which half of the
input array will supply a word for filling the "work" array. When one
half of the input array has been exhausted, the remaining elements are
selected from the other half.
Besides the distribution test, whose results I posted for a 4-element
permutation, there is another test I run. I compute the mean sum of
moments for each permutation, which I define as
Sum(J*A[J]) / Sum(A[J])
averaged over a large number of permutations.
In each case up to N=32768, the result comes out to be about (N-1)/2, as
expected. I have not tested larger permutation tables.
For every permutation of N elements, one could define an integer in the
range [0,N!). To shuffle N elements, without bias, is equivilant to
selecting a number in [0,N!), without bias. I think this part should be
obvious so far.
Consider N = 3. Selecting one of the 6 permutations is exactly
eqivilant to selecting an integer in [0,6). There is no way to do this
without bias in a deterministic amount of time if all you have is an rng
which produces 0s and 1s.
> > What kind of rng output causes your algorithm to loop forever?
>
> There is no sequence from the rng that will cause my algorithm to run
> forever.
Indeed. Therefor, I conclude that your algorithm produces biased
output.
It doesn't quite work as other shuffles i've seen, the algorithm below take
an argument of
0-255 pseude random integers with no reocurrences and start shuffle, (the
deck). Not the bits then the deck is used to tell how to exchange the bits.
So there is no bytes or bits shuffled, however the output from the shuffle
go into a swap.
I know that the shuffle is weak for the last perm[255] but it would be easy
doing
perm[255]<->perm[0].
*until out of bits*
for (i=0;i<255;i++){
slot=perm[i+1];
slot2=perm[slot];
perm[slot]=perm[i];
perm[i]=slot2;
save=i;
}
then do
n=0;
for (i=0;i<255;i++){
swapval=+perm[i];
outbits[i]=bit[swapval];
n++;
}
}
Any comment appreciated.
"Benjamin Goldberg" <gol...@earthlink.net> skrev i meddelandet
news:3CDE55AF...@earthlink.net...
> Ok, anti sorting is just another name for shuffling :) I put that
> subject in to attempt to catch people's attentions.
>
> I've been looking though my books, and I can't find any mention of the
> shuffler which I came up with. Could someone tell me if my idea is the
> same as some already existing algorithm, or can I name it "the goldberg
> shuffler"? :)
>
> Here it is in psuedo-code:
>
> array shuffle(unordered_set x) {
> unordered_set left, right;
> if( x.size <= 1 ) return x.to_array();
> for i in x:
> if( rand_bit ) left.insert(i)
> else right.insert(i)
> return shuffle(left).append( shuffle(right) );
> }
>
> Here it is in [untested] C code:
>
> void swap(int *a, int *b) { int t = *a; *a = *b; *b = t; }
>
> void shuffle(int *start, int *end) {
> while( end-start > 2 ) {
> int *front = start, *back = end;
> while( front < back )
> rand_bit() ? swap( front, --back ) : ++front;
> if( front-start < end-back ) {
> shuffle( start, front );
> start = front;
> } else {
> shuffle( back, end );
> end = back;
> }
> }
> if( end-start == 2 && rand_bit() )
> swap(start, start+1);
> }
> /*
> int x[5] = {1,2,3,4,5};
> shuffle(&x[0], &x[5]);
> */
>
> It's recursive, but by changing the tail call into a loop, and arranging
> to only recurse on the smaller half of the array, the depth of recursion
> is limited to up to a maximum of log2(N).
>
> I'm trying to figure out how many calls to rand_bit are made on average
> for a size N array of ints, in big-O notation, but errno=ENOCAFFEINE.
>
> PS: Hey Tom, you took down your 'Remote Win32 C Compiler'! How come?
> > Hmm. Any algorithm which produce all permutations without bias will
> > take a non-determinsitic amount of time (in an extreme case, for
> certain
> > outputs of the rng, it will loop 'forever').
> >
>
> That is non-obvious and if it can be proven, the proof is likely to be
> well beyond my ability to comprehend. I simply don't believe that's a
> true statement.
>
He's assuming a random source of bits (as you use).
Given that, then, to produce a random unbiased permutation on, say, three
elements, you need to produce each permutation with probability 1/6.
If there is no sequence that will make it loop forever, then that means it
takes at most N random bits for some N. With N bits, you can only get a
multiple of 2**(-N), which cannot equal 1/6; this is because 1/6 in binary
is 0.0010101010101... a repeating decimal (base-2).
One possible solution for 3 elements is:
random_permute_3: roll_die, then pick one of the six permutations.
roll_die {
x = get_bit()
loop {
y = get_bit(), z = get_bit()
if (y is 0 or z is 0) then
return 4*y + 2*z + x + 1 and stop
else
keep going
}
}
The sequence X111111111 of random bits makes it loop forever.
> > What kind of rng output causes your algorithm to loop forever?
>
> There is no sequence from the rng that will cause my algorithm to run
> forever. You can easily satisfy yourself that this is true.
Then it's not quite equally distributed... but perhaps this effect will
only become noticeable for larger permutations; what is the maximum number
of bits it takes?
--
To email me, remove my head.
> Benjamin Goldberg wrote:
>>Any algorithm which produce all permutations without bias will
>>take a non-determinsitic amount of time (in an extreme case, for
>>certain outputs of the rng, it will loop 'forever').
>
> That is non-obvious and if it can be proven, the proof is likely to be
> well beyond my ability to comprehend.
No, no. Math is your friend. None of us were born spitting-up
theorems, and if you can program merge-sort you can understand this one.
If the program always terminates, then for each particular N (the number
of elements to permute) there must be some largest possible number of
bits the random generator could be called upon to produce, call that
number B. We could obviously generate B bits at the start of the
algorithm, and simply take them from the list as needed. If the
algorithm terminates before reading all B bits, that's fine. Each
possible list of B bits has probability 1/2^B.
For any of any particular output, the probability of that output is the
number of B-bit lists that induce that output divided by 2^B. Thus the
probability of each output must be of the form M/2^B for some integer M.
For all permutations to be equally likely, each must have probability
1/N!. For N > 2, there's no way to express 1/N! as an integer divided
by 2^B. The denominator of M/2^B will always be a power of two.
Thus if each permutation has probability 1/N!, there cannot be a B that
is the maximum number of random bits ever needed.
[...]
> Besides the distribution test, whose results I posted for a 4-element
> permutation
The output of the distribution test sure looked uniform. But Goldberg
is right, and I suspect a mistake in test.
--Bryan
If you've got the numbers 0..255 with no number occuring twice, and
they're in random order, then they are already a random permutation.
Once you *have* a random permutation, there's no gain in shuffling it
further.
The goal of shuffling is to start with a *non*random array of N distinct
elements, and a source of random bits, and have your algorithm process
those together to produce any of the N! possible orderings, with equal
likelyhood of each ordering.
Mine does that. Fisher-Yates does that. Yours does not do that.
Your code seems to be correct but, as others have pointed out, something
is not right with your experimental data.
For the test case you used, your algorithm will use up to 5 random bits.
Assuming that your RNG is (almost) perfect, the probability of all
outcomes are multiples of 1/2^5 = 0.03125. Partitioning 32 equiprobable
5-bit sequences to 4! = 24 events cannot result in a close-to-uniform
probability distribution. One would expect very noticeable variations
among the probability of different permutations.
I'd suggest a second look at the test code.
Well that would mean to apply a static permutation that is repeated, if not
the permutation work along all bits in the message?
Since my knowledge about common use of things is rather limited my
approaches may sometimes seem a bit naive but once in a while innovative. My
reason to keep on sorting is to create a shuffle stream that is reversible
once you know the *perm key* 256!. That way it is possible to shuffle bits
in O(n), by just keep on shuffle the deck and then apply the swap of bits
after each shuffle round (the block approach) or create the whole shuffle
and then apply swaps (the stream approach). Please correct me if i am wrong.
Of course you have to produce the initial shuffle state (0-255 random int)
by using a key and a generator that produce the *perm key*, i have created
such to. I can post it if anyone keen on have a look at it.
> The goal of shuffling is to start with a *non*random array of N distinct
> elements, and a source of random bits, and have your algorithm process
> those together to produce any of the N! possible orderings, with equal
> likelyhood of each ordering.
Well that would certainly be a much stronger approach then that of mine, my
just shuffle 256 bits among themselves at a time but it is a dynamic shuffle
between blocks though, while yours shuffle them among N elements is that
possible without having the whole message in memory?
> Mine does that. Fisher-Yates does that. Yours does not do that.
Ok i see you have a stronger approach without a limited deck of course i
could build out my deck with N elements to do what you do but i would still
be forced to create the *perm key* with a generator.
But i see the difference your approach is much stronger, but is it necessary
to shuffle all over file/message when bits, i can see the neccisity with
bytes though. But isn't a shuffle among 256 bits close to impossible to
break if not repeated?
Good luck with your ALG.
> Consider N = 3. Selecting one of the 6 permutations is exactly
> eqivilant to selecting an integer in [0,6). There is no way to do this
> without bias in a deterministic amount of time if all you have is an rng
> which produces 0s and 1s.
Exactly true, if all the bits are unbiased.
for N=3 I need exactly ln(3!)/ln(2) bits of information =
2.5849625007211561814537389439478 (approximately)
That's two unbiased bits, plus one bit that has probability p of being true
such that:
p*lg(p)+(1-p)*lg(1-p) = -.5849625007211561814537389439478 (approximately)
(Of course, getting your random number generator to produce a random bit
whose bias is exactly equal to some irrational number near 85.8% is
probably moving the problem rather than solving it.)
Note that if we are willing to accept a permutation bias so small that it
could barely be detected with, say, 10^100 samples, then a deterministic,
finite procedure becomes possible, and for most purposes practical. It
would be very difficult to prove that any random bitstream was so unbiased
that it would not be the primary source of error.
Anyway, the chances than an optimally functioning permutation generator
would require more than two hundred unbiased bits to generate an unbiased
permutation of ten distinct items are far exceeded by your chances of being
struck by lightning while reading this sentence.
>
> Your code seems to be correct but, as others have pointed out, something
> is not right with your experimental data.
>
> For the test case you used, your algorithm will use up to 5 random bits.
> Assuming that your RNG is (almost) perfect, the probability of all
> outcomes are multiples of 1/2^5 = 0.03125. Partitioning 32 equiprobable
> 5-bit sequences to 4! = 24 events cannot result in a close-to-uniform
> probability distribution. One would expect very noticeable variations
> among the probability of different permutations.
>
> I'd suggest a second look at the test code.
>
5 bits is far too small. First thing I'd check is: after calling
random_permutation() do you reset the array to something fixed, or do you
re-permute the permuted array?
Every time the algorithm executes, it makes up to 5 two-way random
decisions. Are you disputing that?
If the algorithm consumes at most n bits, and 2**n is a feasible
number of trials, Ursus, try running one trial with each of 2**n cooked
random number generators. Each cooked random number generator returns
the bits from a different sequence 00000, 00001, 00010, ... 11111. That
should get you definitive results a lot faster than running enough truly
random trials that you can distinguish a probability 0.03 outcome from a
probability 0.04 outcome.
--Mike Amling
I must apologize up front for the length of this post. It includes the
original "random_permutation" function as well as the test harness.
Moreover, I have modified the test harness to collect distribution
statistics for N=3, 4, 5, and 6. I have done this using a very
hammer-and-tongs approach. I will leave it to interested parties to
reproduce my distributions. Needless to say, I have found the same uniform
distributions of permutations for all values of N between 3 and 6.
The manifest constant HOWMANY specifies the number of elements in the array
to be permuted. It can be larger than 6, but if it is, only the mean moment
will be computed.
#include "stdlib.h"
#include "stdio.h"
#define HOWMANY 4
#define HIGH_BIT (RAND_MAX & ~(RAND_MAX >> 1))
void random_permutation(int nA[], int nLength, int nWork[])
{
if (nLength >= 2)
{
int nSplit, ii, jj, kk;
nSplit = (nLength + 1) / 2;
random_permutation(&nA[0], nSplit, nWork);
random_permutation(&nA[nSplit], nLength - nSplit, nWork);
ii = 0;
jj = nSplit;
kk = 0;
while (kk < nLength)
{
if (ii < nSplit)
{
if (jj < nLength)
{
if (rand() & HIGH_BIT)
{
nWork[kk++] = nA[ii++];
}
else
{
nWork[kk++] = nA[jj++];
}
}
else
{
nWork[kk++] = nA[ii++];
}
}
else
{
nWork[kk++] = nA[jj++];
}
}
for (ii = 0; ii < nLength; ii++) nA[ii] = nWork[ii];
}
}
int main(int argc, char* argv[])
{
int A[HOWMANY], WORK[HOWMANY];
double dX, dY, dMom;
int nP = 1;
#if HOWMANY == 3
int Stats[HOWMANY][HOWMANY][HOWMANY];
int i0, i1, i2;
#endif
#if HOWMANY == 4
int Stats[HOWMANY][HOWMANY][HOWMANY][HOWMANY];
int i0, i1, i2, i3;
#endif
#if HOWMANY == 5
int Stats[HOWMANY][HOWMANY][HOWMANY][HOWMANY][HOWMANY];
int i0, i1, i2, i3, i4;
#endif
#if HOWMANY == 6
int Stats[HOWMANY][HOWMANY][HOWMANY][HOWMANY][HOWMANY][HOWMANY];
int i0, i1, i2, i3, i4, i5;
#endif
for (int ii = 0; ii < HOWMANY; ii++) A[ii] = ii;
#if HOWMANY == 3
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++)
Stats[i0][i1][i2] = 0;
for (ii = 0; ii < 600000; ii++)
{
random_permutation(A, HOWMANY, WORK);
Stats[A[0]][A[1]][A[2]] += 1;
}
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++)
if (Stats[i0][i1][i2])
printf("%5d (%3d,%3d,%3d) %5d\n", nP++, i0, i1, i2, Stats[i0][i1][i2]);
printf("\n");
#endif
#if HOWMANY == 4
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++) for (i3 = 0; i3 < HOWMANY; i3++)
Stats[i0][i1][i2][i3] = 0;
for (ii = 0; ii < 2400000; ii++)
{
random_permutation(A, HOWMANY, WORK);
Stats[A[0]][A[1]][A[2]][A[3]] += 1;
}
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++) for (i3 = 0; i3 < HOWMANY; i3++)
if (Stats[i0][i1][i2][i3])
printf("%5d (%3d,%3d,%3d,%3d) %5d\n", nP++, i0, i1, i2, i3,
Stats[i0][i1][i2][i3]);
printf("\n");
#endif
#if HOWMANY == 5
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++) for (i3 = 0; i3 < HOWMANY; i3++)
for (i4 = 0; i4 < HOWMANY; i4++)
Stats[i0][i1][i2][i3][i4]= 0;
for (ii = 0; ii < 12000000; ii++)
{
random_permutation(A, HOWMANY, WORK);
Stats[A[0]][A[1]][A[2]][A[3]][A[4]] += 1;
}
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++) for (i3 = 0; i3 < HOWMANY; i3++)
for (i4 = 0; i4 < HOWMANY; i4++)
if (Stats[i0][i1][i2][i3][i4])
printf("%5d (%3d,%3d,%3d,%3d,%3d) %5d\n", nP++, i0, i1, i2, i3, i4,
Stats[i0][i1][i2][i3][i4]);
printf("\n");
#endif
#if HOWMANY == 6
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++) for (i3 = 0; i3 < HOWMANY; i3++)
for (i4 = 0; i4 < HOWMANY; i4++) for (i5 = 0; i5 < HOWMANY; i5++)
Stats[i0][i1][i2][i3][i4][i5] = 0;
for (ii = 0; ii < 72000000; ii++)
{
random_permutation(A, HOWMANY, WORK);
Stats[A[0]][A[1]][A[2]][A[3]][A[4]][A[5]] += 1;
}
for (i0 = 0; i0 < HOWMANY; i0++) for (i1 = 0; i1 < HOWMANY; i1++)
for (i2 = 0; i2 < HOWMANY; i2++) for (i3 = 0; i3 < HOWMANY; i3++)
for (i4 = 0; i4 < HOWMANY; i4++) for (i5 = 0; i5 < HOWMANY; i5++)
if (Stats[i0][i1][i2][i3][i4][i5])
printf("%5d (%3d,%3d,%3d,%3d,%3d,%3d) %5d\n", nP++, i0, i1, i2, i3, i4,
i5,
Stats[i0][i1][i2][i3][i4][i5]);
printf("\n");
#endif
dMom = 0;
for (int jj = 0; jj < 2000; jj++)
{
dX = 0.;
dY = 0.;
random_permutation(A, HOWMANY, WORK);
for (ii = 0; ii < HOWMANY; ii++)
{
dX += (double) ii * (double) A[ii];
dY += (double) A[ii];
}
dMom += (dX / dY);
}
printf("\nAverage Moment: %.1f\n", dMom / (double) jj);
return 0;
}
That question makes no sense. I would suggest you think again about
what you want to ask, try translating a bit differently, and double
check to be certain that the english sentence translates sensibly back
into your question.
> Since my knowledge about common use of things is rather limited my
> approaches may sometimes seem a bit naive but once in a while
> innovative. My reason to keep on sorting is to create a shuffle
> stream that is reversible once you know the *perm key* 256!.
To create such a key for N elements takes O(N log N) work. Shuffling is
equivlant to creating such keys. Honestly, I don't care what one does
once one has shuffled / what one does once one has such a key. It's the
creation of the shuffle in the first place that's difficult.
[snip]
> > The goal of shuffling is to start with a *non*random array of N
> > distinct elements, and a source of random bits, and have your
> > algorithm process those together to produce any of the N! possible
> > orderings, with equal likelyhood of each ordering.
>
[snip]
> while yours shuffle them among N elements is that
> possible without having the whole message in memory?
Yes. In fact, mine only needs to have one
element in memory, and can have the rest on disk.
> for (ii = 0; ii < 2400000; ii++)
> {
> random_permutation(A, HOWMANY, WORK);
>
> Stats[A[0]][A[1]][A[2]][A[3]] += 1;
> }
>
Bingo. Make that:
for (ii = 0; ii < 2400000; ii++)
{
for (int jj = 0; jj < HOWMANY; jj++) A[jj] = jj;
random_permutation(A, HOWMANY, WORK);
Stats[A[0]][A[1]][A[2]][A[3]] += 1;
}
> Daniel Mehkeri wrote:
> >
> > On Tue, 14 May 2002, Nicol So wrote:
> >
> > > For the test case you used, your algorithm will use up to 5 random bits.
> >
> > 5 bits is far too small.
>
> Every time the algorithm executes, it makes up to 5 two-way random
> decisions. Are you disputing that?
>
No. I'm saying it's far too small to make a distribution close to 1/24.
No, that's not how I use it. Why should I want to reset the value to a
fixed starting point when all I care about is getting a new random
permutation?
Because without Daniel's suggestion, you are shuffling an already
shuffled array.
Imagine you were doing stastical tests on an RNG. Which would you do:
for( i = 0; i < 2400000; ++i ) {
int j = MY_RNG_NEXT();
process(j);
}
or:
int j;
for( i = 0; i < 2400000; ++i ) {
j += MY_RNG_NEXT();
process(j);
}
Your code is like the second. Daniel is suggesting you change it to be
like the first.
> Ursus Horibilis wrote:
> >
> > "Daniel Mehkeri" <fo...@kvackmyhead.org> wrote in message
> > news:Pine.GSO.4.05.10205151310230.25533-100000@columbia...
> > >
> > >
> > > > for (ii = 0; ii < 2400000; ii++)
> > > > {
> > > > random_permutation(A, HOWMANY, WORK);
> > > >
> > > > Stats[A[0]][A[1]][A[2]][A[3]] += 1;
> > > > }
> > > >
> > >
> > >
> > > Bingo. Make that:
> > >
> > >
> > > for (ii = 0; ii < 2400000; ii++)
> > > {
> > > for (int jj = 0; jj < HOWMANY; jj++) A[jj] = jj;
> > >
> > > random_permutation(A, HOWMANY, WORK);
> > >
> > > Stats[A[0]][A[1]][A[2]][A[3]] += 1;
> > > }
> > >
> > >
> >
> > No, that's not how I use it. Why should I want to reset the value to
> > a fixed starting point when all I care about is getting a new random
> > permutation?
>
> Because without Daniel's suggestion, you are shuffling an already
> shuffled array.
>
If you want to see this point made drastically, replace
random_permutation() with:
broken_random_permutation:
if (rand_bit())
swap A[0] and A[1]
else
swap A[1], A[2], and A[3]
If you don't reset the array, then this should pass your test; all
permutations on four elements will be produced with equal frequency, but
there are definite correlations between one permutation and the next!
And your point is...? Look, when you set the array back to an ordered state
between shuffles, you are actually forcing the implicit creation of a
sequence of shuffles in which every other point in the sequence is a fixed,
known value. It is clear that we are looking for two different things here.
What I am testing is not the model you envision and the model you envision
is not particularly attractive to me, mainly because it has really lousy
uniformity of distribution.
If uniformity of distribution is all you want, you can just produce each
possible permutation in sequence. You probably also want independence
among the generated permutations, though.
My point is, by not resetting the array, you are not really testing your
shuffle algorithm properly.
Suppose one were to replace your shuffler with a 'next_permutation'
function... which produces the next lexigraphically greater ordering of
the array if the array is not in order, and reverses the array
otherwise.
This is distinctly nonrandom, but if you call it numerous times without
resetting the array between calls, it will eventually result in every
possible permutation, and will produce each permutation an exactly equal
number of times.
But if you reset the array between calls, the non-randomness of this is
quite visible.
Let's put this another way:
Suppose I've got an RNG which is supposed to produce numbers in the
range 0..255 with uniform frequency. Let's suppose I test it as
follows:
unsigned char test = 0;
for( i = 0; i < 60000; ++i ) {
test += my_rng();
samples[test] = samples[test] + 1;
}
This is essentially the same as not resetting the array between
shuffles. For almost any function you stick into my_rng, your
statistical tests will indicate that all outputs are equally likely,
even if my_rng always produces a '1'.
A much better test would be:
for( i = 0; i < 60000; ++i ) {
unsigned char test = 0;
test += my_rng();
samples[test] = samples[test] + 1;
}
This will perform actual testing. See the difference?
[snip]
> It is clear that we are looking for two different things here.
> What I am testing is not the model you envision and the model you
> envision is not particularly attractive to me, mainly because it has
> really lousy uniformity of distribution.
If you have a working shuffler, then it will produce a good distribution
even if the array is reset to an ordered state before every shuffle.
If your shuffler does not produce a good distribution if the array is
reset before each shuffler, then it is broken.
> If you've got the numbers 0..255 with no number occuring twice, and
> they're in random order, then they are already a random permutation.
A shuffled deck is a shuffled deck let us stay with that. A shuffle of an
element in the deck, could be on a subset or using the whole deck. If given
a file using 50 000 000 bytes, you tell me that each bit could be swapped
with any of the 50 000 000 * 8 bits i doubt that.
And if you actaully do then it would be madness, so i guess you actually
permutate a subset in a round, that would mean a repeated shuffle. But then
one may ask is your shuffle dynamic or static between the blocks.
> Once you *have* a random permutation, there's no gain in shuffling it
> further.
That pure bull given a shuffle of n < text, if you do not like to have your
blocks static shuffled. Of course you could go for a new shuffle for each
round ->*but that look real slow to me* and really doesn't make any sense.
; )...
> The goal of shuffling is to start with a *non*random array of N distinct
> elements, and a source of random bits, and have your algorithm process
> those together to produce any of the N! possible orderings, with equal
> likelyhood of each ordering.
> Mine does that. Fisher-Yates does that. Yours does not do that.
Oooh so you do move each bit all over the file using a sort algorithm, but
why settle for such a slow one? Why don't you use a O(n).
And actually my ALG could easely do what your does in O(n), using a key
working on a precomputed random array *but of course it would still be a DB
shuffle*.
>To create such a key for N elements takes O(N log N) work. Shuffling is
>equivlant to creating such keys. Honestly, I don't care what one does
>once one has shuffled / what one does once one has such a key. It's the
>creation of the shuffle in the first place that's difficult.
I had no idea that single loop that is not nested could work in O(N log N)
that would certainly be news to me, i must tell all my friends about it.
Ooooh i see you mean that the actual swap-> perm[slot]=perm[i]<- is a O(N
log N) operation mmmmm, you not only create better shuffles then me, i see
you are also a more experienced within how to calculate the work a
permutation demands *congratulations*.
>Here it is in psuedo-code:
>array shuffle(unordered_set x) {
> unordered_set left, right;
> if( x.size <= 1 ) return x.to_array();
> for i in x:
> if( rand_bit ) left.insert(i)
> else right.insert(i)
> return shuffle(left).append( shuffle(right) );
>}
Well looks like bad news to me that shuffle isn't keyed and wouldn't be
reversible, you really consider to use it?
And what about using an unordered set as "*non*random" how does the set
become unordered?
Indeed. My algorithm takes an unshuffled deck and produces a shuffled
one.
Once one already has a shuffled deck, one is finished. Perhaps that
shuffled deck might be used for something else, but that is beyond the
scope of the discussion.
[snip]
> > Once you *have* a random permutation, there's no gain in shuffling
> > it further.
>
> That pure bull [snip]
Think about this. Suppose that I have a random, unbiased, integer in
the range [0,X). Then I add another (independent) number to it, modulo
X. Is the resulting number more random than the original number? No.
Shuffling is like taking a fixed [nonrandom] integer in [0,N!) and
adding to it a random integer in the range [0,N!). The resulting number
is a random number in the range [0,N!). Doing more shuffling would be
pointless -- fully random is fully random
> > The goal of shuffling is to start with a *non*random array of N
> > distinct elements, and a source of random bits, and have your
> > algorithm process those together to produce any of the N! possible
> > orderings, with equal likelyhood of each ordering.
> > Mine does that. Fisher-Yates does that. Yours does not do that.
>
> Oooh so you do move each bit all over the file using a sort algorithm,
I move each record, whatever kind of record it might be, bits, bytes,
CRLF-delimited-lines, anything.
> but why settle for such a slow one? Why don't you use a O(n).
Because it is impossible to do an unbiased shuffle in O(N). You can
only do it in O(N log N).
> And actually my ALG could easely do what your does in O(n), using a
> key working on a precomputed random array *but of course it would
> still be a DB shuffle*.
Sure, once you have such a key. But it takes O(N log N) to make the
key.
> >To create such a key for N elements takes O(N log N) work. Shuffling
> >is equivlant to creating such keys. Honestly, I don't care what one
> >does once one has shuffled / what one does once one has such a key.
> >It's the creation of the shuffle in the first place that's difficult.
>
> I had no idea that single loop that is not nested could work in O(N
> log N) that would certainly be news to me, i must tell all my friends
> about it.
Show me the code for creating your key, and I will show you that either
it is biased, or else it takes O(N log N) time.
> Ooooh i see you mean that the actual swap-> perm[slot]=perm[i]<- is a
> O(N log N) operation mmmmm, you not only create better shuffles then
> me, i see you are also a more experienced within how to calculate the
> work a permutation demands *congratulations*.
No, moron. If you have a pre-shuffled array, then performing a fixed
permutation of that array takes O(1) time. That's all your algorithm
does. But the preshuffling takes O(N log N) time.
> >Here it is in psuedo-code:
> >array shuffle(unordered_set x) {
> > unordered_set left, right;
> > if( x.size <= 1 ) return x.to_array();
> > for i in x:
> > if( rand_bit ) left.insert(i)
> > else right.insert(i)
> > return shuffle(left).append( shuffle(right) );
> >}
>
> Well looks like bad news to me that shuffle isn't keyed and wouldn't
> be reversible, you really consider to use it?
Actually, it is reversible. It's just not obvious to you.
function unshuffle( array x ) {
tmp = rng.copy;
left = 0, right = 0;
for( i = 0; i < x.length ; ++x )
rng.rand_bit() ? ++left : ++right;
left_arr = unshuffle( x[0 .. left-1] );
right_arr = unshuffle( x[left .. x.length] );
right = left; left = 0;
for( i = 0; i < x.length; ++i )
tmp.rand_bit() ?
results.append(x[left++]) :
results.append(x[right++]) ;
return results;
}
> And what about using an unordered set as "*non*random" how does the
> set become unordered?
Unordered merely means not having an *explicit* order.
For example, if one has a hash table, then outputs all the keys, they
will have an order [based on the hash function], but it's not an order
which is relevant to anything.
Benjamin Goldberg wrote:
> [[ something ]]
No offense meant, and lord knows I'm not the newsgroup mom, but having read
your series of posts I really think you guys are arguing about two different
things.
Bob H
This is the part i find interesting, how many bits do you actually need to
make the array of elements random in O(n log n) could you make an array of
let us say 8 000 000 bits perfect shuffled just using a source of 2 random
bits. Or must your random bit array(that would be your key wouldn't it?) be
of equal length as the element array, that would make quite a difference.
It takes O(N log N) bits on average to select, without equal
probability, any of the N-factorial possible permutations of N distinct
items.
> could you make an array of let us say 8 000 000 bits perfect shuffled
Generally, when discussing shuffling algorithms, one assumes that all
elements are distinct. So one would not use a generic shuffle algorithm
for 8e6 BITS, since bits can only have values of 0 and 1. One might
shuffle 8e6 distinct integers, or distinct strings, etc., but if there
are duplicate elements, there are more efficient ways to shuffle
(assuming that one does not need to reverse the order).
For example, to shuffle an array of bits, which can only be 0 or 1,
consider this:
function shuffle(array x) {
int ones = 0, zeros = 0;
for(bit b in x) { b ? ++ones : ++zeros; }
while( ones && zeros ) {
if( rand_mod(ones+zeros) < ones ) {
result.append(1); --ones;
} else {
result.append(0); --zeros;
}
}
while( ones ) { result.append(1); --ones; }
while( zeros ) { result.append(0); --zeros; }
return result;
}
Note, however, that although this will be much more efficient
computationally, as it does no swaps, it will still require O(N log N)
bits on average.
> just using a source of 2 random bits.
To shuffle 8000000 elements, it will take 8000000 * 22.9 bits on
average. Far greater than 2.
After all, with 2 random bits, there are only 4 possible results.
With 8000000 distinct elements, there are 8000000! different
permutations. If your elements are bits, rather than distinct elements,
then if there are A ones and B zeros, then there are (A+B)!/(A!*B!)
distinct permutations. Even with A=7999999 and B=1 (one zero and the
rest ones), then there are still 8000000 possible different arrangement
of those bits; far more than 4.
> Or must your random bit array(that would be your key wouldn't it?) be
> of equal length as the element array, that would make quite a
> difference.
If your random bit source is a true random number generator, then it's
quite possible that the algorithm will never terminate, consuming an
infinite number of bits. As a simplified example of why this is so,
consider when one needs to get a random integer in [0,3) :
int rand_mod3() {
int x, y;
do {
x = rand_bit(), y = rand_bit();
} while( x && y );
return x*2 + y;
}
Although it's highly likely that this will terminate quickly, returning
a number evenly distributed in [0,3), it's *possible* that it will
continue forever.
That said, the number of bits needs to be quite a bit greater than the
number of elements, since if you have N elements, and N random bits,
then the most you could do with those N random bits is to randomly place
your N elements into two groups (a 0 bit indicating that the
corresponding element goes into the first group, and a 1 bit indicating
that the corresponding element goes into the second group).
This will divide your input array roughly in half. If you randomly
shuffle these two halves, and concatenate them together, the result will
be perfectly shuffled. How to shuffle each group? Well, if there's
only 0 or 1 elements in the group, it's shuffled. If there's two, pick
one random bit and either swap them or don't swap them. If there's more
than two elements, one can use recursion.
Since the partitions are generally approximiately halves, there will be
an average of log2(n) recursive calls. And the total number of random
bits used on each level of recursion is N. So the total number of bits
used of all levels of recursion put together is the number of bits used
per level times the number of levels. Or, N times log2 N.
The shuffler shown below may still have bias but, even following your rules,
it's hard to see in the distributions for N=3, 4, 5, and 6. Unfortunately,
it is completely unrelated to the shuffling algorithm you originally posted
and so is basically irrelevant to the original thread.
void random_permutation(int nA[], int nLength)
{
for (int jj = 0; jj < nLength; jj++)
{
int nFrom = jj + (int) (((double) rand() /
((double) RAND_MAX + 1.0))
* (double) (nLength - jj));
int nTemp = nA[jj];
nA[jj] = nA[nFrom];
nA[nFrom] = nTemp;
}
}
I took my own advice and ran it with cooked random functions that
illustrate the exact distribution of permutations the
random_permutation() function produces. There are six possible
permutations of 3 elements, and random_permutations produces them in
proportion to 1,1,1,2,1,2.
Distribution for random_permutation(3 elements) over 8 "cooked" random functions:
1,1,1,2,1,2,
Distribution for random_permutation(4 elements) over 32 "cooked" random functions:
2,2,1,1,1,1,2,2,1,1,1,1,1,1,1,1,2,2,1,1,1,1,2,2,
More results are at http://www.isaurian.com/UHRP.text, with source at http://www.isaurian.com/UHRP.java.text.
--Mike Amling
Those results are approximately what I found when I ran it with the
modification suggested by Mehkeri and Goldberg: resetting the
permutation table to a fixed state between calls to
"random_permutation". Nevertheless, given that you say you "cooked" the
PRNGs, you have tested a different algorithm from the one I posted (and
possibly under different operating conditions--it is not clear from your
description whether you modified it to reset your table to a fixed state
between calls).
B. Goldberg's (and others') analysis is flawed and my algorithm
demonstrates this. Setting the permutation table to a fixed state
between shuffles is neither necessary nor relevant. Denying this is
akin to the apocryphal story about the aeronautical engineer who claims
he can prove that the Bumble Bee can't fly; it's just that nobody
explained this to the Bumble Bee. Here's why.
Suppose that I wrap the "random_permutation()" function in a loop that
simply calls it M times. I redefine the "random_permutation()" function
to be this new wrapped function consisting of an iterated call to the
original "random_permutation" (which is now called something else). I
choose M to be large enough so that you when you look at the
distribution, you can't tell the difference between bias and random
fluctuations caused by the PRNG.
Now let's review. We have a deterministic, terminating, O(n log n)
algorithm that produces random permutations with a distribution that you
can't distinguish from uniform. What's wrong with this picture? I
don't know, but I suspect that it involves a misapplication of
probability theory that is taken from the continuous domain and applied
to the discrete domain.
This tends to be confirmed if you look further down in the thread at my
post that has the subject changed to "Shuffling by Brute Force". In it,
I give an example of another algorithm that, by actual experiment,
produces highly uniform distributions of permutations, this time using
"Goldberg's rules". Again, it is deterministic, terminating and O(n log
n). (It "looks" like O(n), but you have to consider that it requires
approximately log(n) random bits per decision.)
No, he hasn't. A real RNG will give every sequence of random numbers
with equal probability. So, testing your algorithm *once* with each
possible bitstring that an RNG might produce will have the same results
as testing X million times with a real rng.
> (and possibly under different operating conditions--it is not clear
> from your description whether you modified it to reset your table to a
> fixed state> between calls).
I think that it is clear that he did modify the tester so that your
shuffler always began with an ordered set.
> B. Goldberg's (and others') analysis is flawed and my algorithm
> demonstrates this.
Eh? An analysis can show an algorithm to be flawed, but an algorithm
cannot show an analysis to be flawed.
> Setting the permutation table to a fixed state between shuffles is
> neither necessary nor relevant.
Oh, yes it is necessary and relevant for testing.
> Denying this is akin to the apocryphal story about the aeronautical
> engineer who claims he can prove that the Bumble Bee can't fly; it's
> just that nobody explained this to the Bumble Bee.
Heh. Yes, you are denying the need to reset the table to a fixed state.
So, you are indeed like the aeronautical engineer.
> Here's why.
>
> Suppose that I wrap the "random_permutation()" function in a loop that
> simply calls it M times. I redefine the "random_permutation()"
> function to be this new wrapped function consisting of an iterated
> call to the original "random_permutation" (which is now called
> something else).
Hmm. Let's suppose I have a biased RNG, called X(). Now suppose I
define a new RNG Y() which is merely (X() xor X()...M times... xor X()).
> I choose M to be large enough so that you when you look at the
> distribution, you can't tell the difference between bias and random
> fluctuations caused by the PRNG.
I choose M large enough that when .... by the PRNG.
> Now let's review. We have a deterministic, terminating, O(n log n)
> algorithm that produces random permutations with a distribution that
> you can't distinguish from uniform. What's wrong with this picture?
> I don't know, but I suspect that it involves a misapplication of
> probability theory that is taken from the continuous domain and
> applied to the discrete domain.
There are two problems:
One: A flawed algorithm whos bias is hidden is still a flawed algorithm.
Two: The constant factor hidden by the big-oh in *my* O(n log n) is
close to one, and the constant factor hidden by the big-oh in yours is
some ridiculously large M.
Although my algorithm can theoretically take more time to complete, in
practice, mine will finish first.
> This tends to be confirmed if you look further down in the thread at
> my post that has the subject changed to "Shuffling by Brute Force".
> In it, I give an example of another algorithm that, by actual
> experiment, produces highly uniform distributions of permutations,
> this time using "Goldberg's rules". Again, it is deterministic,
> terminating and O(n log n). (It "looks" like O(n), but you have to
> consider that it requires approximately log(n) random bits per
> decision.)
What does "tends to be confirmed" mean?
> Suppose that I wrap the "random_permutation()" function in a loop that
> simply calls it M times. I redefine the "random_permutation()" function
> to be this new wrapped function consisting of an iterated call to the
> original "random_permutation" (which is now called something else). I
> choose M to be large enough so that you when you look at the
> distribution, you can't tell the difference between bias and random
> fluctuations caused by the PRNG.
>
> Now let's review. We have a deterministic, terminating, O(n log n)
> algorithm that produces random permutations with a distribution that you
> can't distinguish from uniform. What's wrong with this picture? I
> don't know, but I suspect that it involves a misapplication of
> probability theory that is taken from the continuous domain and applied
> to the discrete domain.
>
It's because by repeating the random_permutation a certain number of
times, you actually reduce the bias. If you do it enough, the bias is not
noticeable, so although not theoretically perfect, it can be made as close
to perfect as desired, which is all that you need. A while after I posted,
it did actually occur to me to do that, though it might turn out to need
too many iterations to be efficient.
If you don't iterate the procedure and just re-use the shuffled array, you
will, as you've discovered, get a uniformly distributed sequence of
permutations, but they will not be independant. So for instance, the
chance of any given permutation being X is about 1/24, but the chance of
seeing Y immediately after an X is either 1/16 or 1/32, depending on the X
and the Y. This is probably not a good thing for whatever your application
is. Then again, it might not matter.
Ursus Horibilis wrote:
>
> "Michael Amling" <nos...@nospam.com> wrote in message
> news:3CE482F1...@nospam.com...
> >
> > I took my own advice and ran it with cooked random functions that
> > illustrate the exact distribution of permutations the
> > random_permutation() function produces. There are six possible
> > permutations of 3 elements, and random_permutations produces them in
> > proportion to 1,1,1,2,1,2.
> >
> > Distribution for random_permutation(3 elements) over 8 "cooked" random
> functions:
> > 1,1,1,2,1,2,
> > Distribution for random_permutation(4 elements) over 32 "cooked"
> random functions:
> > 2,2,1,1,1,1,2,2,1,1,1,1,1,1,1,1,2,2,1,1,1,1,2,2,
> >
>
> Those results are approximately what I found when I ran it with the
> modification suggested by Mehkeri and Goldberg: resetting the
> permutation table to a fixed state between calls to
> "random_permutation". Nevertheless, given that you say you "cooked" the
> PRNGs, you have tested a different algorithm from the one I posted (and
> possibly under different operating conditions--it is not clear from your
> description whether you modified it to reset your table to a fixed state
> between calls).
I did transliterate the code from C to Java. I certainly intended to
keep the algorithm the same, although coding errors are always possible.
Each invocation permuted an array in which initially a[x]==x for all x.
>
> B. Goldberg's (and others') analysis is flawed and my algorithm
> demonstrates this. Setting the permutation table to a fixed state
> between shuffles is neither necessary nor relevant.
I think they're right. The presumption is that when used in production
it would just be called once per desired permutation.
> ...
> Suppose that I wrap the "random_permutation()" function in a loop that
> simply calls it M times. I redefine the "random_permutation()" function
> to be this new wrapped function consisting of an iterated call to the
> original "random_permutation" (which is now called something else). I
> choose M to be large enough so that you when you look at the
> distribution, you can't tell the difference between bias and random
> fluctuations caused by the PRNG.
Actually, that's an interesting idea. For fixed M, the O(...) would
still be the same. As with Rijndael, 10 rounds could a good algorithm
even if 1 round by itself is not.
>
> Now let's review. We have a deterministic, terminating, O(n log n)
> algorithm that produces random permutations with a distribution that you
> can't distinguish from uniform. What's wrong with this picture?
Well, for one thing, I'd want to know how M should vary, if at all,
with the number of elements being shuffled. Also, it uses the auxiliary
array, which F-Y doesn't.
Is your random_permutation based on a particular sorting algorithm?
Can shuffling algorithms be derived from other sorting algorithms? Wait,
wait, don't tell me....
I coded up a QuickSort with "comparison" results supplied by a cooked
RNG. M/L is the probability of the most probable output permutation
divided by the probability of the least probable output permutation.
Algorithm #Elems M/L Max(bits used)
UH's r_p 3 2 3
" 4 2 5
" 5 8 8
" 6 16 11
" 7 32 14
QuickShuf 3 2 3
" 4 2 5
" 5 4 8
" 6 11 1/3 16
" 7 25.99 25
Looking at the results, I don't think I'm going to replace the F-Y
shuffle, for which M/L is 1, with either of these algorithms.
--Mike Amling
A flawed algorithm is one in which a flaw is demonstrable. In other
words, it demonstrably does not conform to a specification. If the flaw
is truly hidden, then not only is it not a flaw, it's a feature of the
algorithm. It requires a peculiar intellectual dishonesty to claim that
an algorithm is flawed even though no one can demonstrate it.
>
> Two: The constant factor hidden by the big-oh in *my* O(n log n) is
> close to one, and the constant factor hidden by the big-oh in yours is
> some ridiculously large M.
>
> Although my algorithm can theoretically take more time to complete, in
> practice, mine will finish first.
>
If your algorithm can't be guaranteed to terminate, and that is what you
claimed, then big-oh has no meaning. Your algorithm doesn't
"theoretically take more time to complete"; under some conditions, it
*doesn't* complete--you said so. What happens when you average one
infinity of operations in with your O(n log n)?
Let's suppose that 128 bit AES could be broken with 2**120 work.
If this were so, we would not just called it flawed, we'd call it
fatally broken.
> > Two: The constant factor hidden by the big-oh in *my* O(n log n) is
> > close to one, and the constant factor hidden by the big-oh in yours
> > is some ridiculously large M.
> >
> > Although my algorithm can theoretically take more time to complete,
> > in practice, mine will finish first.
>
> If your algorithm can't be guaranteed to terminate, and that is what
> you claimed, then big-oh has no meaning.
Actually, it does have meaning. It describes the *average* time to
complete.
Often, algorithms (especially ones whose behavior depends upon a random
number source) will have three different big-ohs: The minimum time to
complete, an average time to complete, and a maximum time to complete.
In mine, the min and average times are both O(N log N), and the maximum
is unbounded (not infinite, but there is no deterministic number which
the algorithm will take less time than to complete).
> Your algorithm doesn't "theoretically take more time to complete";
> under some conditions, it *doesn't* complete--you said so.
Hmm, actually, it *does* eventually complete, but for any finite length
of time you site before running the algorithm, it's possible for it to
take more than that finite length of time. That's what nondeterministic
means.
Put it this way, suppose I have a true random bit source (every bit
produced has 50/50 odds of being zero/one), and someone asks, is it
possible for the first million bits to be all zeros? Yes, though it's
highly unlikely. Is it possible for the first billion bits to be all
zeros? Still yes, though even more unlikely. Could it produce zeros
literally forever and no ones? No... It will eventually produce ones,
but for any finite number of bits you ask about, it *could* produce all
zeros, or all ones, or any fixed bitstring you might suggest.
> What happens when you average one infinity
> of operations in with your O(n log n)?
The algorithm *always* completes in a finite amount of time. But that
time might be nondeterministically large. However, this is *already
taken into account* in that O(n log n).
I give up! You win! There's nothing more to squeeze out of this
thread. And believe it or not, you have actually been helpful in some
respects.