Generating Permutations

Background: Radix Method for Generating Permutations

In Volume 4A of his classic work, The Art of Computer Programming, Donald Knuth turns to the question of generating permutations for a given combinatoric system.

The book opens with him jumping immediately into the problem of generating permutations, as a distinct problem from counting permutations. The seemingly simple problem of actually iterating through every possible combination has many subtleties and interesting features.

Algorith M is the simplest algorithm, the one that Knuth opens with. It is, in essence, equivalent to labeling every possible combination with a number, from 0 to $N-1$, and generating each combination in a loop from 0 to $N-1$.

Beginning with the example of generating a sequence of $n$ binary digits, Knuth shows that generating every possible outcome is equivalent to counting from 0 to $N-1$ in binary. If we were generating a sequence of $n$ ternary digits (0, 1, or 2), this would be equivalent to counting from 0 to $N-1$ in base 3.

The method Knuth uses, and the algorithm he provides, works for the case where each digit position may have a different number of items to choose from. (In the examples above, all $n$ positions had either 2 or 3 things to choose from, but the more general approach allows us to generate, say, "all 5-digit numbers where the first three digits are odd and the last 2 digits are even, and the last digit is divisible by 4").

To do this, consider the system of $n$ independent variables as completely speciied by the $n$-tuple

$$ (a_1, a_2, \dots, a_n) $$

where $a_i$ takes on any possible vaules of independent variable from 0 to the maximum $0 \leq a_i \leq m_i$. Now this can be considered as a number in a variable radix system:

The number is $\left[ a_1, a_2, \dots, a_n \right]$

The radix is $\left[ m_1, m_2, \dots, m_n \right]$

By repeatedly adding 1 to the number $\left[ a_1, a_2, \dots, a_n \right]$, we iterate through every possible tuple $(a_1, a_2, \dots, a_n)$ and therefore through every possible combination of the independent variables.


Algorithm M (Mixed-radix generation). This algorithm visits all $n$-tuples that satisfy the number/radix expressions above, by repeatedly adding 1 to the mixed-radix number until overflow occurs. (Aux. variables $a_0$ and $m_0$ introduced for convenience only.)

M1. [Initialize.] Set $a_j \rightarrow 0$ for $0 \leq j \leq n$, set $m_0 \rightarrow 2$.

M2. [Visit.] Visit the $n$-tuple $(a_1, \dots, a_n)$. (The program that wants to examine all $n$-tulpes now does its thing.)

M3. [Prepare to add one.] Set $j \rightarrow n$.

M4. [Carry if necessary.] If $a_j = m_j - 1$, set $a_j \rightarrow 0, j \rightarrow j-1$ and repeat this step.

M5. [Increase, unless done.] If $j=0$, terminate algorithm. Otherwise set $a_j = a_j + 1$ and go back to step $M2$.


Implementing Algorithm M

Unfortunately, this pseudocode takes some work to translate. Fortunately, that's already done in the method below.

The method below implements Algorithm M in Python to generate random sequences on a 4x4 Rubik's Cube (called the Rubik's Revenge, or RR for short). The RR cube has six faces that can each be rotated clockwise or counterclockwise by a single layer, denoted by the uppercase letters U, D, B, F, L, R (up, down, back, front, left, right, respectively) for clockwise rotations and U', D', B', F', L', R' for counterclockwise rotations, for 12 total moves.

However, on a 4x4 cube, we can also rotate two layers at a time. (That's the limit; moving three layers at a time is equivalent to a reversed rotation of the remaining single layer.) This type of rotation is denoted with a "w".

Thus, rotating the top two layers of each of the six faces clockwise is denoted Uw, Dw, Bw, Fw, Lw, Rw, and counterclockwise rotations are denoted Uw', Dw', Bw', Fw', Lw', Rw', for 12 additional moves.

We have one more type of move, which is where the second layer only is removed. This type of move is denoted with 2, and the face whose second layer is being rotated, for six moves: 2U, 2D, 2B, 2F, 2L, 2R. The prime notation denotes again a counterclockwise rotation, for an additional six moves: 2U', 2D', 2B', 2F', 2L', 2R'. This yields another 12 moves.

There are 36 total moves that can be executed on the 4x4 Rubik's Revenge cube.

In [1]:
def algorithm_m(n):
    """
    Knuth's Algorithm M for permutation generation,
    via AOCP Volume 4 Fascile 2.
    This is a generator that returns permtuations 
    generated using the variable-radix method. 

    This generates ALL permutations.
    Many of these are rotations of one another,
    so use the get_rotations() function
    (defined later) to eliminate redundant sequences.
    """
    moves = ['U',    'D',    'B',    'F',    'L',    'R',
             'U\'',  'D\'',  'B\'',  'F\'',  'L\'',  'R\'',
             'Uw',   'Dw',   'Bw',   'Fw',   'Lw',   'Rw',
             'Uw\'', 'Dw\'', 'Bw\'', 'Fw\'', 'Lw\'', 'Rw\'',
             '2U',   '2D',   '2B',   '2F',   '2L',   '2R', 
             '2U\'', '2D\'', '2B\'', '2F\'', '2L\'', '2R\'']

    # M1 - Initialize
    a = [0,]*n
    m = [len(moves),]*n

    j = n-1

    nvisits = 0
    while True:

        # M2 - visit
        move_sequence = " ".join([ moves[int(aj)] for aj in a])
        yield move_sequence 

        nvisits += 1

        # M3 - prepare to +1
        j = n-1

        # M4 - carry
        while( a[j] == m[j]-1):
            a[j] = 0
            j = j-1
        
        # M5 - increase unless done
        if(j<0):
            break
        else:
            a[j] = a[j] + 1

Test Drive

Let's take a look at how Algorithm M looks when it is applied. We gave it 36 possible moves, so we expect to see $36^n$ permutations:

In [2]:
def count_permutations(n, p=12):
    permutations = list(algorithm_m(n))
    print("Generated %d permutations (should be 36^%d = %d)"%(len(permutations), n, 36**n))
    print("First %d permutations: %s"%(p, ", ".join(permutations[:p])))
    print("")

count_permutations(2)
count_permutations(3)
Generated 1296 permutations (should be 36^2 = 1296)
First 12 permutations: U U, U D, U B, U F, U L, U R, U U', U D', U B', U F', U L', U R'

Generated 46656 permutations (should be 36^3 = 46656)
First 12 permutations: U U U, U U D, U U B, U U F, U U L, U U R, U U U', U U D', U U B', U U F', U U L', U U R'

In [3]:
from rubikscubennnsolver.RubiksCube444 import RubiksCube444, solved_4x4x4
from pprint import pprint

def get_cube():
    """
    Get a 4x4 Rubiks Cube.
    """
    order = 'URFDLB'
    cube = RubiksCube444(solved_4x4x4, order)
    return cube
In [4]:
permutations = list(algorithm_m(3))[500:505]
for sequence in permutations:
    r = get_cube()
    for move in sequence.split(' '):
        r.rotate(move)
    print("-"*40)
    print("After sequence %s:"%(sequence))
    r.print_cube()
----------------------------------------
After sequence U Dw 2B':
         B B L F 
         B B L F 
         U U U U 
         U U U U 

D D F F  R R R R  B B U U  L B R R 
D D L L  F F F F  R R U U  L B R R 
D D B B  L L L L  F F U U  L B R R 
D D B B  L L L L  F F U U  L B R R 

         D D D D 
         D D D D 
         F F R B 
         F F R B 


----------------------------------------
After sequence U Dw 2F':
         U U U U 
         U U U U 
         B R F F 
         B R F F 

F F U U  R F L L  D D B B  L L L L 
L L U U  R F L L  D D R R  B B B B 
B B U U  R F L L  D D F F  R R R R 
B B U U  R F L L  D D F F  R R R R 

         F L B B 
         F L B B 
         D D D D 
         D D D D 


----------------------------------------
After sequence U Dw 2L':
         R R U U 
         F F U U 
         L L U U 
         L L U U 

F L B B  D D R R  B B B B  L L U U 
F L B B  D D F F  R R R R  B B U U 
F L B B  D D L L  F F F F  R R U U 
F L B B  D D L L  F F F F  R R U U 

         R R D D 
         R R D D 
         B B D D 
         L L D D 


----------------------------------------
After sequence U Dw 2R':
         U U R R 
         U U R R 
         U U B B 
         U U L L 

F F F F  R R U U  B R F F  D D L L 
L L L L  F F U U  B R F F  D D B B 
B B B B  L L U U  B R F F  D D R R 
B B B B  L L U U  B R F F  D D R R 

         D D R R 
         D D F F 
         D D L L 
         D D L L 


----------------------------------------
After sequence U Bw U:
         U U B B 
         U U R R 
         U U R R 
         U U R R 

R R R R  B B D D  B B B L  U U F F 
U U L L  F F F F  R R D D  B B B L 
U U L L  F F F F  R R D D  B B B L 
U U L L  F F F F  R R D D  B B B L 

         D D D D 
         D D D D 
         F L L L 
         F L L L