# Efficient algorithm to generate two diffuse, deranged permutations of a multiset at random

### Background


Let $\msS$ be the multiset $\small\{\overbrace{\po,…,\po}^{n_1},\;\overbrace{\pt,…,\pt}^{n_2},\;…,\;\overbrace{\vphantom 1\pc,…,\pc}^{n_c}\}$ representing one batch. In frequency representation, $\msS$ can also be written as $(\po^{n_1} \;\pt^{n_2}\; … \;\pc^{n_c})$.

The number of distinct permutations of $\msS$ is given by the multinomial: $$\left|\mfS_{\msS}\right|=\binom{n}{n_1,n_2,\dots,n_c}=\frac{n!}{n_1!\,n_2!\cdots n_c!}=n! \prod_{i=1}^c \frac1{n_i!}.$$

### Question

Is there an efficient algorithm to generate two diffuse, deranged permutations $P$ and $Q$ of $\msS$ at random? (The distribution should be uniform.)

• A permutation $P$ is diffuse if for every distinct element $i$ of $P$, the instances of $i$ are spaced out roughly evenly in $P$.

For example, suppose $\msS=(\po^4\;\pt^4)=\{\po,\po,\po,\po,\pt,\pt,\pt,\pt\}$.

• $\{\po, \po, \po, \pt, \pt, \pt, \pt, \po\}$ is not diffuse
• $\{\po, \pt, \po, \pt, \po, \pt, \po, \pt\}$ is diffuse

More rigorously:

• If $n_i=1$, there is only one instance of $i$ to “space out” in $P$, so let $\Delta(i)=0$.
• Otherwise, let $d(i,j)$ be the distance between instance $j$ and instance $j+1$ of $i$ in $P$. Subtract from it the expected distance between instances of $i$, defining the following: $$\delta(i,j)=d(i,j)-\frac n{n_i}\qquad\qquad\Delta(i)=\sum_{j=1}^{n_i-1} \delta(i,j)^2$$ If $i$ is evenly spaced in $P$, then $\Delta(i)$ should be zero, or very close to zero if $n_i\nmid n$.

Now define the statistic $s(P)=\sum_{i=1}^c\Delta(i)$ to measure how much every $i$ is evenly spaced in $P$. We call $P$ diffuse if $s(P)$ is close to zero, or roughly $s(P)\ll n^2$. (One can choose a threshold $k\ll1$ specific to $\msS$ so that $P$ is diffuse if $s(P)<kn^2$.)

This constraint recalls a stricter real-time scheduling problem called the pinwheel problem with multiset $\ms A=n/\msS$ (so that $a_i=n/n_i$) and density $\rho=\sum_{i=1}^c n_i/n=1$. The objective is to schedule a cyclic infinite sequence $P$ such that any subsequence of length $a_i$ contains at least one instance of $i$. In other words, a feasible schedule requires all $d(i,j)≤a_i$; if $\ms A$ is dense ($\rho= 1$), then $d(i,j)=a_i$ and $s(P)=0$. The pinwheel problem appears to be NP-complete.

• Two permutations $P$ and $Q$ are deranged if $P$ is a derangement of $Q$; that is, $P_i ≠ Q_i$ for every index $i\in[n]$.

For example, suppose $\msS=(\po^2\;\pt^2)=\{\po,\po,\pt,\pt\}$.

• $\{\po, \pt, \po, \pt\}$ and $\{\po, \po, \pt, \pt\}$ are not deranged
• $\{\po, \pt, \po, \pt\}$ and $\{\pt, \po, \pt, \po\}$ are deranged

### Exploratory analysis

I am interested in the family of multisets with $n=20$ and $n_i=4$ for $i\lesssim4$. In particular, let $\msD=(\gr1^4\,\gr2^4\,\gr3^4\,\gr4^3\,\gr5^2\,\gr6^1\,\gr7^1\,\gr8^1)$.

• The probability that two random permutations $P$ and $Q$ of $\msD$ are deranged is about 3%.

This can be calculated as follows, where $L_k$ is the $k$th Laguerre polynomial: \begin{align*} \left|{\mathfrak D}_{\msD}\right| &=\int_0^\infty \!\!dt\; e^{-t}\, \prod_{i=1}^c L_{n_i}(t) =\int_0^\infty \!\!dt\; e^{-t}\, \bigl(L_4(t)\bigr)^3\bigl(L_3(t)\bigr)\bigl(L_2(t)\bigr)\bigl(L_1(t)\bigr)^3\\ &=4.5\times10^{11}\\ \left|\mfS_{\msD}\right| &=n!\prod_{i=1}^c \frac1{n_i!} =\frac{20!}{(4!)^3\,(3!)\,(2!)\,(1!)^3} =1.5\times10^{13}\\ p&=\left|{\mathfrak D}_{\msD}\right|/ \left|\mfS_{\msD}\right|\approx0.03\end{align*} See here for an explanation.

• The probability that a random permutation $P$ of $\msD$ is diffuse is about 0.01%, setting the arbitrary threshold at roughly $s(P)<25$.

Below is an empirical probability plot of 100,000 samples of $s(P)$ where $P$ is a random permutation of $\msD$.

At medium sample sizes, $s(P)\sim \text{Gamma}(\alpha\approx8,\beta\approx18)$.

\begin{array}{ccl}\renewcommand\mfm[1]{\textbf{#1}} \hline P & s(P) & \text{cdf}(s(P)) \\ \hline \{\po, \ss8, \pt, \pth, \pf, \po, \pv, \pt, \pth, \ss6, \po, \pf, \pt, \pth, \ss7, \po, \pv, \pt, \pf, \pth\} & \frac{11}9\approx1\, & <10^{-5} \\ \{\ss8, \pt, \pth, \pf, \po, \ss6, \pv, \pt, \pth, \pf, \po, \ss7, \po, \pt, \pth, \pv, \pf, \po, \pt, \pth\} & \frac{140}9\approx16 & <10^{-4} \\ \{\pth, \ss6, \pv, \po, \pth, \pf, \pt, \po, \pt, \ss7, \ss8, \pv, \pt, \pf, \po, \pth, \pth, \pt, \po, \pf\} & \frac{650}9\approx72 & \phantom{<1}0.05 \\ \{\pth, \po, \pth, \pf, \ss8, \pt, \pt, \po, \po, \pv, \pth, \pth, \pt, \ss6, \pf, \pf, \pt, \po, \ss7, \pv\} & \frac{1223}9\approx136 & \phantom{<1}0.45 \\ \{\pf, \po, \po, \pf, \pv, \pv, \po, \pth, \pth, \ss7, \po, \pt, \pt, \pf, \pth, \pth, \ss8, \pt, \pt, \ss6\} & \frac{1697}9\approx189 & \phantom{<1}0.80 \\ \hline \end{array}

The probability that two random permutations are valid (both diffuse and deranged) is around $v\approx(0.03)(0.0001)^2\approx10^{-10}$.

### Inefficient algorithms

A common “fast” algorithm to generate a random derangement of a set is rejection-based:

do
P ← random_permutation(D)
until is_derangement(D, P)
return P


which takes approximately $e$ iterations, since there are roughly $n!/e$ possible derangements. However a rejection-based randomized algorithm would not be efficient for this problem, as it would take on the order of $1/v\approx10^{10}$ iterations.

In the algorithm used by Sage, a random derangement of a multiset “is formed by choosing an element at random from the list of all possible derangements.” Yet this too is inefficient, as there are $v\,|\mfS_{\msD}|^2\approx10^{16}$ valid permutations to enumerate, and besides, one would need an algorithm just to do that anyway.

### Further questions

What is the complexity of this problem? Can it be reduced to any familiar paradigm, such as network flow, graph coloring, or linear programming?

• Regarding your definition of "spaced out", don't you want $d(i,j) - n/(n_i + 1)$ for $0 \leq i \leq j \leq n+1$ with $P_0 = P_{n+1} = i$ as sentinels? That is to say, a single element should be in the middle, two should partition the permutation in thirds, and so on.
– Raphael
May 4, 2015 at 17:43
• What happens if $S = \{ 1^{n-k}, 2^k\}$ for evil $k$ (small, but large enough); do we even have diffuse permutations than? We certainly don't stand a change to find two deranged ones! Seems that no element can occur more than $n/2$ times.
– Raphael
May 4, 2015 at 17:49
• What is the ratio of all pairs of deranged permutations among all pairs of diffuse permutations? Similarly, out of all pairs of deranged permutations, how many consist of two diffuse ones? (If either ratio is "high", we can concentrate our effort on one half of the process, leaving the other at rejection.)
– Raphael
May 4, 2015 at 17:54
• @Raphael (#3a) Of 1 million random permutations of $\mathsf D$, these 561 diffuse ones had $s(P)\le 30$. $6118/\binom{561}{2}=6118/157080\approx3.9\%$ of pairs are deranged.
– hftf
May 5, 2015 at 8:10
• @Raphael (#3b) Of 10 million random pairs of permutations of $\mathsf D$, 306893 pairs were deranged. Only 29 of those pairs had both permutations with $s(P)\le50$. Here is a histogram (values).
– hftf
May 5, 2015 at 9:11

One approach: you can reduce this to the following problem: Given a boolean formula $\varphi(x)$, choose an assignment $x$ uniformly at random from among all the satisfying assignments of $\varphi(x)$. This problem is NP-hard, but there are standard algorithms for generating an $x$ that is approximately uniformly distributed, borrowing methods from #SAT algorithms. For instance, one technique is to pick a hash function $h$ whose range has a carefully chosen size (about the same size as the number of satisfying assignments of $\varphi$), choose uniformly at random a value $y$ from within the range of $h$, and then use a SAT solver to find a satisfying assignment to the formula $\varphi(x) \land (h(x)=y)$. To make it efficient, you can choose $h$ to be a sparse linear map.

This might be shooting a flea with a cannon, but if you have no other approaches that seem workable, this is one you could try.

• finding this hard to follow. $\varphi(x)$ is a boolean value and $h(x)$ is a binary string (set of binary variables)? so the final equation means...?
– vzn
May 4, 2015 at 20:18

some extended discussion/ analysis of this problem commenced in cs chat with further background which uncovered some subjectivity in the complex requirements of the problem but did not find any outright errors or oversights.1

here is some tested/ analyzed code which, compared to the other solution based on SAT is relatively "quick and dirty" but was nontrivial/ tricky to debug. its loosely conceptually based on a local pseudorandom/ greedy optimization scheme somewhat similar to eg 2-OPT for TSP. the basic idea is to start with a random solution that fits some constraint, and then perturb it locally to look for improvements, greedily searching for improvements and iterating through them, and terminating when all local improvements have been exhausted. a design criteria was that the algorithm should be as efficient/ avoid rejection as much as possible.

there is some research into derangement algorithms[4] eg used in SAGE[5] but they are not oriented around multisets.

the simple perturbation is only "swaps" of two positions in the tuple(s). the implementation is in ruby. following is some overview/ notes with refs to line numbers.

qb2.rb (gist-github)

the approach here is to start with two deranged tuples (#106) and then locally/ greedily improve the dispersion (#107), combined in a concept called derangesperse (#97), preserving the derangement. note that swapping two same positions in the tuple pair preserves derangement and can improve dispersion and that is (part of) the disperse method/ strategy.

the derange subroutine works left to right on the array (multiset) and swaps with elements later in the array where the swap is not with the same element (#10). the algorithm succeeds if, with no further swaps at the last position, the two tuples are still deranged (#16).

there are 3 different approaches to deranging the initial tuples. the 2nd tuple p2 is always shuffled. one can start with tuple 1 (p1) ordered by a. "highest powers 1st order" (#128), b. shuffled order (#127), c. and "lowest powers 1st order" ("highest powers last order") (#126).

the dispersion routine disperse is more involved but conceptually not so difficult. again it uses swaps. rather than trying to optimize dispersion in general over all elements it simply tries to iteratively alleviate the current worst case. the idea is to find the 1st least disperse elements, left to right. the perturbation is to swap either the left or right elements (x, y indexes) of the least disperse pair with other elements but never any between the pair (which would always decrease dispersion), and also skip attempting to swap with same elements (select in #71). m is the midpoint index of the pair (#65).

however the dispersion is measured/ optimized over both tuples in the pair (#40) using the "least/ leftmost" dispersion in each pair (#25, #44).

the algorithm attempts to swap "farthest" elements 1st (sort_by / reverse #71).

there are two different strategies true, false for deciding whether to swap the left or right element of the least disperse pair (#80), either the left element for swap position to the left/ right element to the right side, or the farthest left or right element in the disperse pair from the swap element.

the algorithm finishes (#91) when it can no longer improve dispersion (either moving worst disperse location right or increasing max dispersion over the entire tuple pair (#85)).

statistics are output for rejects over c=1000 derangements over the 3 approaches (#116) and c=1000 derangesperses (#97), looking at the 2 disperse algorithms for a deranged pair from rejection (#19, #106). the latter tracks total average dispersion (after guaranteed derangement). an example run is as follows

c       0.661000
b       0.824000
a       0.927000
[2.484, 2, 4]
[2.668, 2, 4]


this shows that the a-true algorithm gives best results with ~92% nonrejection and an average worst disperse distance of ~2.6, and a guaranteed minimum of 2 over 1000 trials, ie at least 1 nonequal intervening element between all same-element pairs. it found solutions as high as 3 nonequal intervening elements.

the derangement algorithm is linear time pre-rejection, and the dispersion algorithm (running on deranged input) seems to be possibly ~$O(n \log n)$.

1 the problem is to find quizbowl packet arrangements that satisfy so-called "feng shui"[1] or a "nice" random ordering where "nice" is somewhat subjective and not yet "officially" quantified; the author of the problem has analyzed/ reduced it to the derange/ dispersion criteria based on research wrt the quizbowl community and "feng shui experts".[2] there are different ideas on "feng shui rules." some "published" research has been done on algorithms but it appears in early stages.[3]

[1] Packet feng shui / QBWiki

[2] Quizbowl packets and feng shui / Lifshitz

[3] Question Placement, HSQuizbowl resource center forum

[4] Generating Random Derangements / Martinez, Panholzer, Prodinger

[5] Sage derangement algorithm (python) / McAndrew

• fyi further thought, theres a glitch in the derange routine & it doesnt always derange. the swap position may advance without swapping anything. theres a simple fix to test success correctly.
– vzn
May 13, 2015 at 20:53