/******************************************************************************
* File: FactoradicPermutation.hh
* Author: Keith Schwarz (htiek@cs.stanford.edu)
*
* An algorithm for manipulating and generating permutations in lexicographic
* order by using the factoradic number system.
*
* The factoradic number system (also called the factorial number system) is a
* way of writing out numbers in a base system that depends on factorials,
* rather than powers of numbers. In the factoradic system, the last digit is
* always 0 and is in base 0!. The digit before that is either 0 or 1 and is
* in base 1!. The digit before that is either 0, 1, or 2 and is in base 2!.
* More generally, the nth-to-last digit in always 0, 1, 2, ..., or n and is in
* base n!. When writing a number, we usually append a subscript ! to indicate
* that it's written in factoradic base. In this (plain-text) C++ source file,
* I'll denote this by suffixing the number with _!.
*
* As an example, the factoradic number 2110_! is equal to
*
* 2 x 3! + 1 x 2! + 1 x 1! + 0 x 0!
* = 2 x 6 + 1 x 2 + 1 x 1 + 0 x 1
* = 12 + 2 + 1 + 0
* = 15
*
* Similarly, the first few factoradic numbers are
*
* 0_! 10_! 100_! 110_! 200_! 210_! 1000_! 1010_! 1100_! 1110_!
*
* One reason that the factoradic numbers are so useful is that there's a close
* connection between factoradic numbers and permutations of distinct values.
* To understand where this comes from, suppose that we're interested in
* generating all of the permutations of some collection in lexicographic
* order. One way to do this might be to think about generating the first such
* permutation, then the second, then the third, etc. In other words, we want
* to look for some way of mapping the integers onto permutations in a way that
* allows us to efficiently answer the question "what is the nth permutation of
* this set?"
*
* Let's consider permutations of four elements - a, b, c, and d - and the
* order in which they're generated. Here's a table of these permutations:
*
* 0 abcd 6 bacd 12 cabd 18 dabc
* 1 abdc 7 badc 13 cadb 19 dacb
* 2 acbd 8 bcad 14 cbad 20 dbac
* 3 acdb 9 bcda 15 cbda 21 dbca
* 4 adbc 10 bdac 16 cdab 22 dcab
* 5 adcb 11 bdca 17 cdba 23 dcba
*
* Initially, there might not seem to be any pattern relating the numbers to
* the permutations, but there are a few trends we can spot here. For example,
* every sixth permutation starts with a new letter, and 6 = 3!. The second
* letter of each permutation changes every two permutations, and 2 = 2!.
*
* If we rewrite the above table with the integers in factoradic base, then a
* more obvious trend starts to appear:
*
* 0000_! abcd 1000_! bacd 2000_! cabd 3000_! dabc
* 0010_! abdc 1010_! badc 2010_! cadb 3010_! dacb
* 0100_! acbd 1100_! bcad 2100_! cbad 3100_! dbac
* 0110_! acdb 1110_! bcda 2110_! cbda 3110_! dbca
* 0200_! adbc 1200_! bdac 2200_! cdab 3200_! dcab
* 0210_! adcb 1210_! bdca 2210_! cdba 3210_! dcba
*
* We noted above that the first letter changes every six permutations, and now
* we can see that this coincides with the times at which the first digit of
* the factoradic representation changes. Similarly, when we noted that the
* second letter of the permutation changed every two permutations, we can see
* this is occurring whenever the second digit of the factoradic representation
* of the index is changing.
*
* The connection here is due to a construction called "Lehmer codes," a way of
* labeling permutations in a way that also allows their construction. The
* idea behind Lehmer codes is as follows. Suppose that you are given a
* permutation of a totally-ordered set of elements, such as the permutation
* BEDAC of the letters A-E. We start off by looking at the first letter of
* this permutation - in this case, B - and seeing how many elements it's
* bigger than (in this case, 1). We then split this letter off of the
* permutation and record the fact that it was bigger than one element:
*
* B EDAC
* 1
*
* Now, we repeat this process on EDAC. Here, E is larger than three of the
* elements of EDAC (everything except itself), so we write out a three to get
*
* BE DAC
* 13
*
* Repeating this on D gives us
*
* BED AC
* 132
*
* And on A gives us
*
* BEDA C
* 1320 0
*
* And finally on the last C gives us
*
* BEDAC
* 13200
*
* So the Lehmer code for this permutation is (1, 3, 2, 0, 0). Note that the
* length of a Lehmer code is equal to the number of elements in the
* permutation.
*
* Given a permutation's Lehmer code, we can generate the permutation it stands
* for by running this process backwards. Suppose, for example, that we're
* given the Lehmer code (3, 1, 0, 0) and want to use it to recover some given
* permutation of ABCD. To do this, we'd start off by taking the element of
* this collection of letters that's bigger than three other elements. This is
* equivalent to asking which element is at position 3 using zero-indexing, and
* gives us D. So now we have
*
* D ABC
* 3
*
* and the remaining Lehmer code (1, 0, 0). Looking at the 1 in the first
* position of this Lehmer code, we should take the element out of the unused
* elements at position 1 (which, zero-indexed, is B), giving us
*
* DB AC
* 31
*
* The rest of our code is (0, 0), so we take the zeroth-largest element of the
* unused letters (A) and put it next:
*
* DBA C
* 310
*
* And finally, we use the last part of our code, (0), to select the last
* element of the permutation, which is C:
*
* DBAC
* 3100
*
* The clincher in all this is that there's a simple, one-to-one mapping
* between the factoradic numbers and Lehmer codes. The idea is simple - given
* a Lehmer code (d0, d1, ..., dn), we just pick the number
*
* d0 d1 d2 ... dn_!
*
* For example, given Lehmer code (3, 1, 0, 0), we'd pick the factoradic number
* 3100_!. Similarly, given factoradic number 210_!, we'd pick the Lehmer code
* (2, 1, 0).
*
* In order to convince ourselves that this is even mathematically well-defined
* we need to see that each of the digits we assign are within the range
* specified for that digit. For example, if it were possible that we could
* have a Lehmer code (1, 4), then we couldn't convert it to the factoradic
* number 14_!, since the last digit of any factoradic number must be 0. This
* is not particularly difficult to see. Remember that the digits of the
* Lehmer code are zero-based indices into the list of elements we haven't yet
* used yet. Thus the last digit must be zero, since there's only one element
* remaining; the penultimate digit can be either one or zero since there's
* two possible choices; the antepenultimate digit can be zero, one or two;
* etc.
*
* Given the bijiective correspondence between Lehmer codes and the factoradic
* numbers, we have a way of taking a natural number N and mapping it into a
* permutation of some number of elements. However, as of now we have no
* reason to believe that this will list every permutation in lexicographic
* order. After all, given some random encoding system for permutations, why
* should we expect anything special out of our ordering? Fortunately, though
* we have the following lemma:
*
* Lemma: Given two permutations p_1 and p_2 with Lehmer codes L_1 and L_2,
* p_1 lexicographically precedes p_2 iff L_1 lexicographically precedes L_2.
*
* Corollary: Every permutation has a unique Lehmer code.
*
* Before we go into the proof, I should remark that this guarantees that the
* following algorithm will list all permutations of n elements in order:
*
* For i = 0 to n! - 1:
* Represent i in factoradic base.
* Convert this representation to a Lehmer code.
* Generate the permutation based on this Lehmer code.
*
* This algorithm works correctly because every permutation has a Lehmer code,
* and by iterating over all n! different codes we'll hit every permutation.
* Moreover, our lemma guarantees that no two permutations have the same
* Lehmer code, we know that we'll get every permutation in order and exactly
* once.
*
* Proof of lemma: We first prove that if p_1 < p_2, then L_1 < L_2. Given the
* two permutations, consider the first elements at which the permutations
* disagree, call this index i. Note that i can't be the last element of the
* permutations, because that would mean that the permutations agree in every
* position except the last, and since the last element doesn't match it means
* that the permutations aren't of the same elements. Since the permutations
* agree on the first i - 1 elements, the Lehmer codes for the first i - 1
* elements must agree. This means that at the time the ith digit of the
* Lehmer code is being chosen, in both permutations the same set of leftover
* elements will remain. Since p_1 < p_2, the ith element of p_1 must be
* smaller than the ith element of p_2, and so the index of that element in
* what's left must be smaller than the index of the corresponding element from
* p_2 in what's left. Consequently, the digit chosen for L_1 will be less
* than the digit for L_2, and since the codes agree on their prefix before
* this we have that L_1 < L_2.
*
* To prove the opposite direction, suppose that L_1 < L_2 and find the first
* spot at which they disagree; say it's index i. i cannot be the last spot
* in the code, since the very last digit must be a 0, so it has to be some
* place before that. Since the prefixes of the codes match up to that point,
* at the spot where L_1 chooses a lower number than L_2, the partial
* permutations built up so far must match and the remaining elements must be
* the same. Consequently, when L_1 has a lower value than L_2, the resulting
* permutation will have a smaller value in p_1 than in p_2, and since the
* prefixes match the result will be that p_1 < p_2. QED.
*
* The functions exported by this file allow for the generation and
* manipulation of permutations using factoradic number representations. They
* run in quadratic time, though faster implementations are possible using more
* complex algorithms.
*/
#ifndef FactoradicPermutation_Included
#define FactoradicPermutation_Included
#include <iterator> // For std::iterator_traits
#include <algorithm> // For std::sort, std::rotate, std::count_if
#include <functional> // For std::less, std::bind2nd
/**
* Function: NthPermutation(RandomIterator begin, RandomIterator end,
* Integer n);
* Usage: NthPermutation(v.begin(), v.end(), n);
* ----------------------------------------------------------------------------
* Given a range of iterators [begin, end) and an integer n, rearranges the
* values in [begin, end) so that they are in the nth lexicographically
* smallest permutation. It is assumed that n is large enough to hold
* (end - begin)! without overflow.
*/
template <typename RandomIterator, typename Integer>
void NthPermutation(RandomIterator begin, RandomIterator end, Integer n);
/**
* Function: PermutationIndex(RandomIterator begin, RandomIterator end);
* Usage: size_t n = PermutationIndex<size_t>(begin, end);
* ----------------------------------------------------------------------------
* Given a range of iterators spanning a permutation of distinct values,
* returns the number of that permutation in the lexicographical ordering of
* permutations.
*/
template <typename Integer, typename RandomIterator>
Integer PermutationIndex(RandomIterator begin, RandomIterator end);
/**
* Function: NthPermutation(RandomIterator begin, RandomIterator end,
* Integer n, Comparator comp);
* Usage: NthPermutation(v.begin(), v.end(), n, std::greater<int>());
* ----------------------------------------------------------------------------
* Given a range of iterators [begin, end) and an integer n, rearranges the
* values in [begin, end) so that they are in the nth lexicographically
* smallest permutation. It is assumed that n is large enough to hold
* (end - begin)! without overflow. Comparisons are done using the provided
* comparator comp.
*/
template <typename RandomIterator, typename Integer, typename Comparator>
void NthPermutation(RandomIterator begin, RandomIterator end, Integer n,
Comparator comp);
/**
* Function: PermutationIndex(RandomIterator begin, RandomIterator end,
* Comparator comp);
* Usage: size_t n = PermutationIndex<size_t>(begin, end, std::greater<int>());
* ----------------------------------------------------------------------------
* Given a range of iterators spanning a permutation of distinct values,
* returns the number of that permutation in the lexicographical ordering of
* permutations. The specified comparator is used to determine the ordering of
* individual values in the sequence.
*/
template <typename Integer, typename RandomIterator, typename Comparator>
Integer PermutationIndex(RandomIterator begin, RandomIterator end,
Comparator comp);
/* * * * * Implementation Below This Point * * * * */
template <typename RandomIterator, typename Integer, typename Comparator>
void NthPermutation(RandomIterator begin, RandomIterator end, Integer n,
Comparator comp) {
/* Sort the elements into ascending order according to the comparator. The
* logic for Lehmer encoding requires this.
*/
std::sort(begin, end, comp);
/* We now need to convert the integer n into something in factoradic base.
* To do this, we attempt to recover the factoradic digits one at a time much
* in the way that you would recover the digits of, say, a digital or binary
* number. In particular, if the range has size k, we compute (k - 1)!, then
* divide n by that number to yield the most significant factoradic digit.
* We then divide (k - 1)! by k - 1 to get (k - 2)! and repeat this process
* one digit at a time.
*
* As a first step, compute k and (k - 1)!.
*/
const Integer k(end - begin);
/* Compute (k - 1)!. */
Integer radix = Integer(1);
for (Integer i = Integer(2); i < k; ++i)
radix *= i;
/* This loop holds the invariant that on entry, radix is equal to i!, where
* i is the factoradic digit being considered right now.
*
* We also hold the invariant that the range [out, out + k) is the unsorted
* part of the permutation, while the elements before this are the generated
* portion of the permutation.
*
* Again, to avoid integer underflow, we count upward rather than downward.
*/
for (Integer i = 0; i < k; ++i) {
/* Recover the next factoradic digit by dividing n by the current radix. */
Integer digit = n / radix;
/* We now want to pick that element as the next term of the permutation.
* After we do this, we need to shuffle all the elements after that element
* down a spot. This is equivalent to rotating all of the elements in the
* range [out, out + digit) one spot to the right. For example, given this
* setup where we want to put 3 as the first element of the permutation:
*
* 1 2 3 4 5
* ^
* we do a cyclic right shift to get
*
* 3 1 2 4 5
*/
std::rotate(begin, begin + digit, begin + digit + 1);
/* Set up for the next iteration. We need to move the output iterator
* forward a step so that it writes to the next location.
*/
++ begin;
/* We also need to remove the contribution of the current radix to n by
* modding it out.
*/
n %= radix;
/* Finally, we need to divide the radix by our index so that on the next
* iteration it has the correct value. However, if this is the final
* iteration, then this would be a divide-by-zero, and so we don't do the
* divide.
*/
if (i + 1 != k)
radix /= (k - i - 1);
}
}
/* To determine the index of a permutation, we determine the order of each
* element in the range - that is, the number of elements smaller than it -
* then scale those values up by the appropriate base in the factoradic
* representation.
*/
template <typename Integer, typename RandomIterator, typename Comparator>
Integer PermutationIndex(RandomIterator begin, RandomIterator end,
Comparator comp) {
/* As in the generating permutations case, compute k and (k - 1)!. */
const Integer k(end - begin);
Integer radix = Integer(1);
for (Integer i = Integer(2); i < k; ++i)
radix *= i;
/* Loop across the elements, computing the order of each and scaling it by
* the appropriate factorial base.
*/
Integer result = Integer(0);
for (; begin != end; ++begin) {
/* Compute the order of the element by seeing how many elements in the
* range are less than the value.
*/
result += radix * std::count_if(begin, end, std::bind2nd(comp, *begin));
/* Scale down the radix, unless this is the very last iteration. */
if (begin + 1 != end)
radix /= (end - begin) - 1;
}
return result;
}
/* Non-comparator versions implemented in terms of comparator. */
template <typename RandomIterator, typename Integer>
void NthPermutation(RandomIterator begin, RandomIterator end, Integer n) {
NthPermutation
(begin, end, n,
std::less<typename std::iterator_traits<RandomIterator>::value_type>());
}
template <typename Integer, typename RandomIterator>
Integer PermutationIndex(RandomIterator begin, RandomIterator end) {
return PermutationIndex<Integer>
(begin, end,
std::less<typename std::iterator_traits<RandomIterator>::value_type>());
}
#endif