Generating permutations lazily

Question

I'm looking for an algorithm to generate permutations of a set in such a way that I could make a lazy list of them in Clojure. i.e. I'd like to iterate over a list of permutations where each permutation is not calculated until I request it, and all of the permutations don't have to be stored in memory at once.

Alternatively I'm looking for an algorithm where given a certain set, it will return the "next" permutation of that set, in such a way that repeatedly calling the function on its own output will cycle through all permutations of the original set, in some order (what the order is doesn't matter).

Is there such an algorithm? Most of the permutation-generating algorithms I've seen tend to generate them all at once (usually recursively), which doesn't scale to very large sets. An implementation in Clojure (or another functional language) would be helpful but I can figure it out from pseudocode.

Solution 1

Yes, there is a "next permutation" algorithm, and it's quite simple too. The C++ standard template library (STL) even has a function called next_permutation.

The algorithm actually finds the next permutation -- the lexicographically next one. The idea is this: suppose you are given a sequence, say "32541". What is the next permutation?

If you think about it, you'll see that it is "34125". And your thoughts were probably something this: In "32541",

• there is no way to keep the "32" fixed and find a later permutation in the "541" part, because that permutation is already the last one for 5,4, and 1 -- it is sorted in decreasing order.
• So you'll have to change the "2" to something bigger -- in fact, to the smallest number bigger than it in the "541" part, namely 4.
• Now, once you've decided that the permutation will start as "34", the rest of the numbers should be in increasing order, so the answer is "34125".

The algorithm is to implement precisely that line of reasoning:

1. Find the longest "tail" that is ordered in decreasing order. (The "541" part.)
2. Change the number just before the tail (the "2") to the smallest number bigger than it in the tail (the 4).
3. Sort the tail in increasing order.

You can do (1.) efficiently by starting at the end and going backwards as long as the previous element is not smaller than the current element. You can do (2.) by just swapping the "4" with the '2", so you'll have "34521". Once you do this, you can avoid using a sorting algorithm for (3.), because the tail was, and is still (think about this), sorted in decreasing order, so it only needs to be reversed.

The C++ code does precisely this (look at the source in /usr/include/c++/4.0.0/bits/stl_algo.h on your system, or see this article); it should be simple to translate it to your language: [Read "BidirectionalIterator" as "pointer", if you're unfamiliar with C++ iterators. The code returns false if there is no next permutation, i.e. we are already in decreasing order.]

template <class BidirectionalIterator>
bool next_permutation(BidirectionalIterator first,
                      BidirectionalIterator last) {
    if (first == last) return false;
    BidirectionalIterator i = first;
    ++i;
    if (i == last) return false;
    i = last;
    --i;
    for(;;) {
        BidirectionalIterator ii = i--;
        if (*i <*ii) {
            BidirectionalIterator j = last;
            while (!(*i <*--j));
            iter_swap(i, j);
            reverse(ii, last);
            return true;
        }
        if (i == first) {
            reverse(first, last);
            return false;
        }
    }
}

It might seem that it can take O(n) time per permutation, but if you think about it more carefully, you can prove that it takes O(n!) time for all permutations in total, so only O(1) -- constant time -- per permutation.

The good thing is that the algorithm works even when you have a sequence with repeated elements: with, say, "232254421", it would find the tail as "54421", swap the "2" and "4" (so "232454221"), reverse the rest, giving "232412245", which is the next permutation.

Solution 2

Assuming that we're talking about lexicographic order over the values being permuted, there are two general approaches that you can use:

1. transform one permutation of the elements to the next permutation (as ShreevatsaR posted), or
2. directly compute the nth permutation, while counting n from 0 upward.

For those (like me ;-) who don't speak c++ as natives, approach 1 can be implemented from the following pseudo-code, assuming zero-based indexing of an array with index zero on the "left" (substituting some other structure, such as a list, is "left as an exercise" ;-):

1. scan the array from right-to-left (indices descending from N-1 to 0)
1.1. if the current element is less than its right-hand neighbor,
     call the current element the pivot,
     and stop scanning
1.2. if the left end is reached without finding a pivot,
     reverse the array and return
     (the permutation was the lexicographically last, so its time to start over)
2. scan the array from right-to-left again,
   to find the rightmost element larger than the pivot
   (call that one the successor)
3. swap the pivot and the successor
4. reverse the portion of the array to the right of where the pivot was found
5. return

Here's an example starting with a current permutation of CADB:

1. scanning from the right finds A as the pivot in position 1
2. scanning again finds B as the successor in position 3
3. swapping pivot and successor gives CBDA
4. reversing everything following position 1 (i.e. positions 2..3) gives CBAD
5. CBAD is the next permutation after CADB

For the second approach (direct computation of the nth permutation), remember that there are N! permutations of N elements. Therefore, if you are permuting N elements, the first (N-1)! permutations must begin with the smallest element, the next (N-1)! permutations must begin with the second smallest, and so on. This leads to the following recursive approach (again in pseudo-code, numbering the permutations and positions from 0):

To find permutation x of array A, where A has N elements:
0. if A has one element, return it
1. set p to ( x / (N-1)! ) mod N
2. the desired permutation will be A[p] followed by
   permutation ( x mod (N-1)! )
   of the elements remaining in A after position p is removed

So, for example, the 13th permutation of ABCD is found as follows:

perm 13 of ABCD: {p = (13 / 3!) mod 4 = (13 / 6) mod 4 = 2; ABCD[2] = C}
C followed by perm 1 of ABD {because 13 mod 3! = 13 mod 6 = 1}
  perm 1 of ABD: {p = (1 / 2!) mod 3 = (1 / 2) mod 2 = 0; ABD[0] = A}
  A followed by perm 1 of BD {because 1 mod 2! = 1 mod 2 = 1}
    perm 1 of BD: {p = (1 / 1!) mod 2 = (1 / 1) mod 2 = 1; BD[1] = D}
    D followed by perm 0 of B {because 1 mod 1! = 1 mod 1 = 0}
      B (because there's only one element)
    DB
  ADB
CADB

Incidentally, the "removal" of elements can be represented by a parallel array of booleans which indicates which elements are still available, so it is not necessary to create a new array on each recursive call.

So, to iterate across the permutations of ABCD, just count from 0 to 23 (4!-1) and directly compute the corresponding permutation.

Solution 3

You should check the Permutations article on wikipeda. Also, there is the concept of Factoradic numbers.

Anyway, the mathematical problem is quite hard.

In C# you can use an iterator, and stop the permutation algorithm using yield. The problem with this is that you cannot go back and forth, or use an index.

Solution 4

More examples of permutation algorithms to generate them.

Source: http://www.ddj.com/architect/201200326

1. Uses the Fike's Algorithm, that is the one of fastest known.
2. Uses the Algo to the Lexographic order.
3. Uses the nonlexographic, but runs faster than item 2.

1.

PROGRAM TestFikePerm;
CONST marksize = 5;
VAR
    marks : ARRAY [1..marksize] OF INTEGER;
    ii : INTEGER;
    permcount : INTEGER;

PROCEDURE WriteArray;
VAR i : INTEGER;
BEGIN
FOR i := 1 TO marksize
DO Write ;
WriteLn;
permcount := permcount + 1;
END;

PROCEDURE FikePerm ;
{Outputs permutations in nonlexicographic order.  This is Fike.s algorithm}
{ with tuning by J.S. Rohl.  The array marks[1..marksizn] is global.  The   }
{ procedure WriteArray is global and displays the results.  This must be}
{ evoked with FikePerm(2) in the calling procedure.}
VAR
    dn, dk, temp : INTEGER;
BEGIN
IF 
THEN BEGIN { swap the pair }
    WriteArray;
    temp :=marks[marksize];
    FOR dn :=  DOWNTO 1
    DO BEGIN
        marks[marksize] := marks[dn];
        marks [dn] := temp;
        WriteArray;
        marks[dn] := marks[marksize]
        END;
    marks[marksize] := temp;
    END {of bottom level sequence }
ELSE BEGIN
    FikePerm;
    temp := marks[k];
    FOR dk :=  DOWNTO 1
    DO BEGIN
        marks[k] := marks[dk];
        marks[dk][ := temp;
        FikePerm;
        marks[dk] := marks[k];
        END; { of loop on dk }
    marks[k] := temp;l
    END { of sequence for other levels }
END; { of FikePerm procedure }

BEGIN { Main }
FOR ii := 1 TO marksize
DO marks[ii] := ii;
permcount := 0;
WriteLn ;
WrieLn;
FikePerm ; { It always starts with 2 }
WriteLn ;
ReadLn;
END.

2.

PROGRAM TestLexPerms;
CONST marksize = 5;
VAR
    marks : ARRAY [1..marksize] OF INTEGER;
    ii : INTEGER;
    permcount : INTEGER;

PROCEDURE WriteArray;
VAR i : INTEGER;
BEGIN
FOR i := 1 TO marksize
DO Write ;
permcount := permcount + 1;
WriteLn;
END;

PROCEDURE LexPerm ;
{ Outputs permutations in lexicographic order.  The array marks is global   }
{ and has n or fewer marks.  The procedure WriteArray () is global and  }
{ displays the results.                             }
VAR 
    work : INTEGER:
    mp, hlen, i : INTEGER;
BEGIN
IF 
THEN BEGIN { Swap the pair }
    work := marks[1];
    marks[1] := marks[2];
    marks[2] := work;
    WriteArray ;
    END
ELSE BEGIN
    FOR mp :=  DOWNTO 1
    DO BEGIN
        LexPerm<>;
        hlen :=  DIV 2;
        FOR i := 1 TO hlen
        DO BEGIN { Another swap }
            work := marks[i];
            marks[i] := marks[n - i];
            marks[n - i] := work
            END;
        work := marks[n]; { More swapping }
        marks[n[ := marks[mp];
        marks[mp] := work;
        WriteArray;
        END;
    LexPerm<>
    END;
END;

BEGIN { Main }
FOR ii := 1 TO marksize
DO marks[ii] := ii;
permcount := 1; { The starting position is permutation }
WriteLn <   Starting position:  >;
WriteLn
LexPerm ;
WriteLn <   PermCount is   , permcount>;
ReadLn;
END.

3.

PROGRAM TestAllPerms;
CONST marksize = 5;
VAR
    marks : ARRAY [1..marksize] of INTEGER;
    ii : INTEGER;
    permcount : INTEGER;

PROCEDURE WriteArray;
VAR i : INTEGER;
BEGIN
FOR i := 1 TO marksize
DO Write ;
WriteLn;
permcount := permcount + 1;
END;

PROCEDURE AllPerm (n : INTEGER);
{ Outputs permutations in nonlexicographic order. The array marks is    }
{ global and has n or few marks.  The procedure WriteArray is global and    }
{ displays the results.                             }
VAR
    work : INTEGER;
    mp, swaptemp : INTEGER;
BEGIN
IF 
THEN BEGIN { Swap the pair }
    work := marks[1];
    marks[1] := marks[2];
    marks[2] := work;
    WriteArray;
    END
ELSE BEGIN
    FOR mp :=  DOWNTO 1
    DO BEGIN
        ALLPerm<< n - 1>>;
        IF > 
        THEN swaptemp := 1
        ELSE swaptemp := mp;
        work := marks[n];
        marks[n] := marks[swaptemp};
        marks[swaptemp} := work;
        WriteArray;
    AllPerm< n-1 >;
    END;
END;

BEGIN { Main }
FOR ii := 1 TO marksize
DO marks[ii] := ii
permcount :=1;
WriteLn < Starting position;   >;
WriteLn;
Allperm < marksize>;
WriteLn < Perm count is  , permcount>;
ReadLn;
END.

View more solutions

Comments

[sidgeon smythe]

"Anyway, the mathematical problem is quite hard." No it's not :-)

[Chris Conway]

This will work, assuming you have a total order on the elements.

[sidgeon smythe]

If you start with a set, you can arbitrarily define a total order on the elements; map the elements to distinct numbers. :-)

[Bogdan Maxim]

Well, it is.. if you don't know about Factoradic numbers there is no way you could come up with a proper algorithm in an acceptable time. It's like trying to solve a 4th degree equation without knowing the method.

[sidgeon smythe]

You're welcome; I enjoyed revisiting it and thinking about the details too :-)

[sidgeon smythe]

Oh sorry, I thought you were talking about the original problem. I still don't see why you need "Factoradic numbers" anyway... it's pretty simple to assign a number to each of the n! permutations of a given set, and to construct a permutation from a number. [Just some dynamic programming/counting..]

[Die in Sente]

++ Your answer is underappreciated. Not to take away from the accepted answer, but the second approach is more powerful because it can be generalized to combinations as well. A complete discussion would show the reverse function from sequence to index.

[Jason S]

whee! neat. +1 -- I'm not an algorithms expert, but I've seen most of the basic ones and it's been a while since I've seen an algorithm that was both new to me, and easy to understand.

[Daniel C. Sobral]

This answer just doesn't get enough upvotes, but I can only upvote it once... :-)

[Stephen Swensen]

Very nice, I've implemented it in F#: projecteulerfun.blogspot.com/2010/05/…

[sidgeon smythe]

Indeed. I agree with the previous comment — even though my answer does slightly fewer operations for the specific question asked, this approach is more general, since it works for e.g. finding the permutation that is K steps away from a given one.

[Stephen Swensen]

And here is a further developed version which is generic, takes a function for defining a total ordering on the elements, and is wrapped by a immutable functional sequence expression: stackoverflow.com/questions/286427/…

[Masse]

Probably a stupid question, but if at (1) you look for the longest tail of decreasing order, and in (3) you make the tail increasing order, isn't the longest tail afterwards always of length 1? (Probably a stupid question, but let's put it on account of sleep deprivation)

[sidgeon smythe]

@Masse: Not exactly... roughly, you can go from 1 to a larger number. Using the example: Start with 32541. The tail is 541. After doing the necessary steps, the next permutation is 34125. Now the tail is just 5. Incrementing 3412 using the 5 and swapping, the next permutation is 34152. Now the tail is 52, of length 2. Then it becomes 34215 (tail length 1), 34251 (tail length 2), 34512 (length 1), 34521 (length 3), 35124 (length 1), etc. You are right that the tail is small most of the time, which is why the algorithm has good performance over multiple calls.

[Drew Noakes]

In idiomatic C#, an iterator is more correctly referred to as an enumerator.

[Harsh Pathak]

@ShreevatsaR: How would you do that short of generating all permutations? E.g. if you needed to generate the n!th permutation.

[Sam Stoelinga]

Wouldn't it take O(n!) for all permutation if we can generate the next one in O(1) time?

[sidgeon smythe]

@SamStoelinga: You are right, actually. O(n log n) is O(log n!). I should have said O(n!).

[sidgeon smythe]

@Jacob: Sorry, seeing this only now. You can generate the kth permutation (0≤k<n!) without generating all; this is the answer of joel.neely above. The kth permutation of n elements a[0], a[1], ...a[n-1] is a[k/(n-1)!] (here '/' is integer division, i.e. floor) followed by the k%((n-1)!)'th permutation of the other n-1 elements. (This incidentally turns out to be related to the factorial number system, but it is straightforward to derive without explicitly going that route — it's a general approach that works for generating the kth combinatorial structure whenever counting is easy.)

[Harsh Pathak]

@ShreevatsaR: Thanks! It makes sense since I just used the factoradic system for my problem.

[Nikunj Banka]

@ShreevatsaR for generating the next permutation, we will have to look at all the numbers in the worst case. So doesn't the algorithm take O(N) time for generating the next permutation.

[sidgeon smythe]

@NikunjBanka: Yes there's an error in the answer; see previous comment above. Generating the next permutation can take O(N) time in the worst case for some particular permutations, but when amortized over all n! successive calls to the function, IIRC it takes O(1) time on average. That's what I was trying to get at.

[Nikunj Banka]

Can you please present an analysis as to how it takes O(1) time on average.

[Nikunj Banka]

@ShreevatsaR I got the answer to the question. Thanks.

[Apriori]

Is there any reason to not use a binary search when finding the smallest number bigger than the number preceding the tail in the tail? I know that finding the tail was possibly O(n) (but O(1) average), so pairing it with another linear step can be swept under the rug. But a binary search would still be less operations for large strings, so why not use it?

[sidgeon smythe]

I didn't know it at the time I posted the answer, but this algorithm is the one given by Knuth in The Art of Computer Programming Volume 4A, section 7.2.1.2 "Generating all permutations". (Earlier published in Vol 4 Fasc 2, and before that in pre-fascicle 2B.) It is "Algorithm L" on page 319, the first of 37 pages of sections 7.2.1.2 on this and related fun questions (not counting the 21 pages of answers to exercises). In the exercises and solutions he analyzes this algorithm to a good level of detail.

[sidgeon smythe]

@Apriori Sorry it took me more than 2 years to reply to your question. The answer is that just as finding the tail is O(1) on average, finding the smallest bigger number in the tail is also O(1) on average. Knuth (see my comment just above) in Exercise 5(b) of 7.2.1.2 analyzes that this step (step L3 in his book) performs exactly k comparisons with probability sum_(j=k+1 to n) (1/j!). For large n this gives approximately e/2 ≈ 1.359 comparisons on average, with variance of ≈ 0.418. So surprisingly enough, doing binary search, which takes O(log n) comparisons, would actually make it slower!