Permutation ranking / unranking

permuation_rank_unrank.py

from functools import lru_cache, reduce
from collections import Counter

f_memo = {0: 1}


@lru_cache(maxsize=None)
def F(n):
  '''Compute Factorial'''
  if n in f_memo:
    return f_memo[n]
  r = 1
  for i in range(2, n+1):
    r *= i
    f_memo[i] = r  # store intermediate result
  return r


@lru_cache(maxsize=None)
def choose(m, a):
  '''Compute multinomial coefficient
    => m! / (a_1! x a_2! x ... x a_k!)
  '''
  num = F(m)
  denominator = reduce(lambda x, y: x * y, map(F, a), 1)
  return num // denominator


@lru_cache(maxsize=None)
def count_unique_permutations(partition):
  '''Count the number of unique permutations of a given partition
  '''
  freq = Counter(partition)
  return choose(len(partition), freq.values())


def nth_permutation(sequence, index):
  '''Returns the nth permutation of a given sequence, ignoring duplicates.
  Sorts the passed sequence, assuming it is the first permuation
  '''
  unique = sorted(set(sequence))
  N = len(sequence)
  freq = Counter(sequence)

  if index < 1 or choose(N, freq.values()) < index:
    return ''

  P = 0
  permutation = []
  while P != index:
    P = 0
    for el in unique:
      if not freq[el]:
        continue
      el_remaining = N - len(permutation)  # elements remaining

      # remove this character
      freq[el] -= 1

      # compute total possibilities now
      total_p = choose(el_remaining - 1, freq.values())
      P += total_p

      if P >= index:
        permutation += [el]
        index -= P - total_p
        break

      if P < index:
        freq[el] += 1

  # get reamining elements in reversed order
  remaining = [el for el in unique for _ in range(freq[el]) if freq[el]]
  permutation += reversed(remaining)
  return permutation


def permuation_index(sequence):
  '''Return the lexicographic index of a given sequence in the list of possible permuations.
  Assumes the sorted sequence is the first in the list
  '''
  index = 1
  length = len(sequence)
  freq = Counter(sequence)
  sorted_uniq_seq = sorted(freq.keys())

  for n in range(length):
    el = sequence[n]

    fsum = 0
    for i in range(length):
      ch = sorted_uniq_seq[i]
      if ch == el:
        break
      fsum += freq[ch]

    fprod = reduce(lambda x, y: y*x, [F(v) for v in freq.values()])
    freq[el] -= 1

    index += (fsum * F(length - n - 1)) // fprod
  return index