13

What is the fastest way to write a function nthPermutation[xs_List, n_Integer] which is equivalent to Permutations[xs][[n]] but does not explicitly compute all the permutations of xs? This feature seems conspicuously absent from Mathematica.

I have a basic implementation below which almost does what I want, but does not preserve the original ordering of xs (it is instead equivalent to Permutations[Sort[xs]][[n]]), does not check that n is in bounds, and is almost certainly not the fastest method.

nthPermutation[{}, 1] = {};

nthPermutation[xs_List, n_Integer] := Block[{u = Union[xs], p, i = 1},
   p = Accumulate[numPermutations[DeleteCases[xs, #, {1}, 1]] & /@ u];
   While[n > p[[i]], i++];
   Prepend[nthPermutation[DeleteCases[xs, u[[i]], {1}, 1],
     If[i == 1, n, n - p[[i - 1]]]], u[[i]]]];

Edit: I should clarify that, if possible, nthPermutation[xs, n] should be equivalent to Permutations[xs][[n]] whenever it is defined, even with repeated and out-of-order elements in xs. For example, we should have that Permutations[{1, 3, 1, 2, 3}][[20]] === nthPermutation[{1, 3, 1, 2, 3}, 20] === {3, 2, 1, 3, 1}.

J. M.'s missing motivation
  • 124,525
  • 11
  • 401
  • 574
David Zhang
  • 2,316
  • 15
  • 25
  • If you have $n$ arbitrary permutations, you must step through each before you have the final state. You can break this up in a number of ways (tracing the sequence of positions of each individual element, for instance), but you must still step through all in the arbitrary case. If, however, there are cycles or identity permutations or other non-random permutations, you can speed up the calculation. Then of course there are ways to parallelize the computation. – David G. Stork Mar 16 '15 at 23:07
  • This may be useful: (6884708) – Mr.Wizard Mar 17 '15 at 00:49
  • @Mr.Wizard I think the StackExchange question you linked concerns the ranking of permutations on $n$ distinct symbols. Ranking permutations on $n$ symbols, possibly with duplicates, is a rather different problem. – David Zhang Mar 17 '15 at 01:10
  • Okay, how about: (14212030) – Mr.Wizard Mar 17 '15 at 01:14
  • @Mr.Wizard If I understand the accepted answer correctly, that (un)ranking function is precisely the one I've implemented in the nthPermutation code above. The problem is, that ranking function gives the permutations in lexicographic order, while the ranking function Mathematica uses internally is decidedly not lexicographic--examine the output of Permutations[{1, 3, 1, 2, 3}] and you'll see what I mean. – David Zhang Mar 17 '15 at 01:19

2 Answers2

12

A simple recursive function I've used many times:

permutation[{e_}, _] := {e}
permutation[l_, n_] := 
 With[{m = (Length[l] - 1)!}, 
  permutation[Delete[l, Ceiling[n/m]], Mod[n, m, 1]]~Prepend~
   l[[Ceiling[n/m]]]]

Note that this function assumes that n is a valid permutation number and that all elements are unique.

Update

Here is a solution that meets your criteria. It's not pretty, but it works on the examples I tested.

permutation[l_List, n_Integer] := Module[
  {i = n, remaining = l, tally = Tally[l], m = 1, len, perm = {}},
  If[Multinomial @@ Last /@ tally >= n > 0,
   While[Length[remaining] > 0,
    tally[[m, 2]]--;
    len = Multinomial @@ Last /@ tally;
    If[len >= i,
     AppendTo[perm, tally[[m, 1]]];
     remaining = DeleteCases[remaining, tally[[m, 1]], {1}, 1];
     tally = Tally[remaining];
     m = 1,
     i -= len;
     tally[[m, 2]]++;
     m++;
     ]
    ];
   perm,
   {}
   ]
  ]

Note that it's iterative, not recursive. This lets it handle very large cases, e.g. permutation[RandomInteger[100, 10000], 100!] (only takes around seven seconds on my computer). I think it's worst-case time is $\mathrm{O}(n^2)$ in the number of elements.

Mr.Wizard
  • 271,378
  • 34
  • 587
  • 1,371
2012rcampion
  • 7,851
  • 25
  • 44
  • Cute. Gets tragically slow on large problems, and should check for invalid permutation number. +1 nonetheless. – ciao Mar 17 '15 at 00:14
  • @rasher Slow, really? For me, event at the recursion limit (Length[l] == 1022) it's pretty fast, less than 30ms on my computer. I know manipulating lists with large, non-numeric elements is slow, maybe just try permuting the indices and then selecting from your list? – 2012rcampion Mar 17 '15 at 00:33
  • This is quite elegant, but doesn't seem to handle repeated elements correctly. Try permutation[{1, 3, 1, 2, 3}, #] & /@ Range[30]--you'll only get 18 of the 30 permutations. – David Zhang Mar 17 '15 at 00:38
  • @DavidZhang I didn't see your edit, sorry. This function assumes all elements are unique. You should use Range[120] to see all the permutations, including the repeats (four repeats for each permutation). – 2012rcampion Mar 17 '15 at 00:43
  • @2012rcampion Ah, that's alright. I've added a little more emphasis to ensure that part is clear. – David Zhang Mar 17 '15 at 00:46
  • 1
    @2012rcampion: Try yours on a target with say 10k elements... it rapidly degrades in performance (and eats memory like Oprah at an all-you-can-eat) once target lists get larger... still a neat and clean piece of code, no slight intended. – ciao Mar 17 '15 at 00:53
  • @rasher I greatly underestimated the sizes you were looking at! Yah, with many thousands of elements the list manipulation at each step just kills it. – 2012rcampion Mar 17 '15 at 01:04
  • 1
    Wow, I'm impressed with your new solution! It seems to do everything I hoped for, and is more than fast enough for my needs. – David Zhang Mar 17 '15 at 01:46
  • Nice solution. I would like to see if I can refine it further. Would you prefer that I include my version in this answer or post it separately? – Mr.Wizard Mar 17 '15 at 02:17
  • @Mr.Wizard Go ahead and post it separately. I'd love to see this cleaned up =) – 2012rcampion Mar 17 '15 at 02:27
  • Okay, I've give it a go later tonight. – Mr.Wizard Mar 17 '15 at 03:06
  • @Mr.Wizard: Combine mixed-radix with the ol' shuffle-product, will destroy all-comers performance-wise, though it seems OP is satisfied with current performance. – ciao Mar 17 '15 at 04:46
  • @rasher You make me feel ignorant again. I haven't looked at mixed-radix in years and all I remember now is confusion, and I had to look up shuffle product. By all means if you have the time and interest implement it and teach us something. – Mr.Wizard Mar 17 '15 at 05:14
  • 2
    I failed to make any meaningful improvement to this. I am going to edit your answer to use ++ and -- instead of += 1 and -= 1, then leave it be for now. If you dislike the edit feel free to revert it. – Mr.Wizard Mar 17 '15 at 08:05
  • @rasher If there's a more efficient algorithm I'd love to know about it. How does one employ mixed-radix representation and shuffle products for this problem? – David Zhang Mar 17 '15 at 11:27
  • @Mr.Wizard Do you see any way of rewriting this function to be Compile-able? I can't get Multinomial @@ Last /@ tally into any form Compile will accept. – David Zhang Mar 17 '15 at 11:30
  • @David Multinomial is not in the list of compilable functions. – Mr.Wizard Mar 17 '15 at 11:45
  • @Mr.Wizard Right, what I mean is, do you see any way of computing Multinomial from inside a compiled function? Using factorials directly will almost certainly cause integer overflow, so I'm wondering if there's some clever algorithm I could use. – David Zhang Mar 17 '15 at 11:49
  • @DavidZhang From what I can tell, you can: Compile[{{x, _Integer, 1}}, Module[{r = 1, i = 1}, Do[r = Quotient[r i, j]; i++, {e, x}, {j, e}]; r]]. However, it won't work when large lists are involved, since the answer can only hold up to 64 bits of result. Plus it's only faster by half. – 2012rcampion Mar 17 '15 at 22:20
5
factNum[n_] := 
  Module[{cnt = 1}, 
   Reverse@NestWhileList[
      QuotientRemainder[#[[1]], cnt++] &, {n, 1}, #[[1]] != 0 &][[2 ;;,2]]];

pickPerm[fn_] := Module[{p = Range@Length@fn, fnx = fn + 1},
   Reap[Scan[(Sow[p[[#]]]; p = Drop[p, {#}]) &, fnx]][[2, 1]]];

permN[target_, n_] := 
  If[n > Length[target]!, {}, 
   target[[pickPerm[PadLeft[factNum[n - 1], Length@target]]]]];

Gives permutations in lexicographic order, 1 being identity, empty list for invalid partition number.

v = Permutations[Range@5];
pn = permN[Range@5, #] & /@ Range[120];
pn == v

(* True *)

permN[{1, 2, 3, 3, 2, 1, 4, 4}, #] & /@ {1, 100, 1000}

(* {{1, 2, 3, 3, 2, 1, 4, 4}, {1, 2, 3, 4, 3, 1, 4, 2}, {1, 3, 2, 3, 4, 1, 4, 2}} *)
ciao
  • 25,774
  • 2
  • 58
  • 139
  • This doesn't seem to handle repeated elements correctly. Try permN[{1, 3, 1, 2, 3}, #] & /@ Range[30]--you'll only get 18 distinct permutations, of which there should be 30. – David Zhang Mar 17 '15 at 00:40
  • No, it handles them fine. Use a range of 120... – ciao Mar 17 '15 at 00:52
  • But that's exactly the point, if you have to use 120 it clearly doesn't account for repeated elements. The function I'm looking for ought to return the nth distinct permutation (in the order Mathematica does). – David Zhang Mar 17 '15 at 00:56
  • @DavidZhang: Not going to spend more time on this, it does precisely what title asked for. Probably might want to change your title to reflect your needs. – ciao Mar 17 '15 at 01:03
  • 1
    To be fair, the question body quite clearly states what I'm looking for. Still, changing the title is a good idea. – David Zhang Mar 17 '15 at 01:05