Paul Michael Kielstra

One of the fun things about being specifically an applied mathematician is knowing when my colleagues are dealing with a combinatorial explosion. Recently, a friend of mine asked me if there were an easy way to determine how many of the partitions of a set of sixteen fixed elements obeyed a certain set of conditions. Doing this with brute force was a non-starter: there were over \(10^{10}\) of them, and if we were able to check sixteen thousand partitions per second the calculation would still have taken more than a week.

However, one of the conditions could be stated as the partitions having to be invariant under a permutation that factored into a product of eight independent \(2\)-cycles. If we could generate and test only the partitions which obeyed this rule, it would cut down the compute time dramatically. The main idea of this post is an algorithm to do just that, based around the idea that any such partition must be associated with a partition of a set of just eight elements.

More specifically, let \(A={\alpha_1,\dots,\alpha_n}\) and \(B={\beta_1,\dots,\beta_n}\). Let \(s\) be the permutation of \(A \cup B\) that sends each \(\alpha_i\) to the corresponding \(\beta_i\) and vice versa. Let \(P’\) be some partition of \(A\cup B\). Let \(P\) be \(P’\) restricted to being a partition of \(A\), and let \(\tilde{P}\) be \(P’\) restricted to being a partition of \(B\). Note in particular that, up to swapping the letter \(\alpha\) for the letter \(\beta\), \(P\) and \(\tilde{P}\) are the same partition.

Now, take some set \(T’\in P’\) and define \(T\) and \(\tilde T\) in the logical way: \(T=T’\cap A\) and \(\tilde T=T’\cap B\). One of these may be empty. Take \(U’=s(T)\), which must also be in \(P’\) by symmetry, and define \(U\) and \(\tilde{U}\) similarly. Then \(P’={T’,U’}\cup Q’\), where \(Q’\) is some partition with fewer sets than \(P’\). But \(T’=T\cup s(U)\) and \(U’=U\cup s(T)\). Therefore, by induction, \(P’\) splits into sets \({T_1’,U_1’}\cup\cdots\cup{T_j’,U_j’}\), where each \({T’_i,U’_i}\) is of the form \({T_i\cup s(U_i),U_i\cup s(T_i)}\) for some \(T,U\in P\cup{\emptyset}\), and all the non-empty \(T_i\)s and \(U_i\)s are pairwise distinct. (The \(T_i\)s are, in particular, distinct from the \(U_i\)s.)

This gives us the bones of our algorithm: we begin with a partition \(P\) of \(A\), take a set \(T\) that we assume WLOG is non-empty, take another set \(U\) that may be empty, build a set \({T’,U’}\), remove both \(T\) and \(U\) from \(P\), and recurse. By iterating over all possible choices we obtain all possible partitions. In implementation, there are two special cases that should be handled separately: the case when \(U=\emptyset\), when \({T’,U’}={T,s(T)}\), and the case when \(U=T\), when \({T’,U’}={T\cup s(T)}\). These make no difference to the mathematics of the method, but can create awkward edge cases if not dealt with.

To quantify the speed-up in the case \(n=8\), we start from the fact that there are \(4140\) partitions of a set of eight elements. (This is the eighth Bell number.) Our algorithm works on the sets in the partitions, rather than the elements in those sets, so the number of partitions it produces from any given input is entirely dependent on the size of the input as a partition. In particular, given a partition of size \(k\), it produces as many partitions as there are symmetric bipartite graphs on \(2k\) vertices with degree at most \(1\). By a lemma which we prove in another post, there are \(\leq k!\) of these.

Furthermore, out of all \(4140\) partitions of \(8\) elements, only one partition has size \(8\), and only \(8\) more have size \(7\). Therefore, the number of partitions generated by our algorithm is \(<8!+8\times7!+(4140-1-8)\times6!\), which is \(3,054,960\). This is only \(\approx 10^6\), so we’ve already saved four orders of magnitude. This bound is far from tight, and in fact the actual number is \(428,131\), which is another power of ten saved. This is coming within the realms of possibility: checking ten permutations per second would require only eleven hours.

Since the bottleneck in this code is largely iteration, recursion, and list manipulation, I implemented it in Julia instead of my usual Python.

using Combinatorics

n = 8

function s(partition)
    return partition .+ n
end

function sym_partition(base_partition)
    n = length(base_partition)
    if n == 0
        return [[]]
    end
    base_set = base_partition[1]
    backend_partitions = sym_partition(base_partition[2:end])
    full_partitions = [vcat([base_set, s(base_set)], b) for b in backend_partitions] # U is empty
    append!(full_partitions, [vcat([vcat(base_set, s(base_set))], b) for b in backend_partitions]) # U = T
    for j = 2:n
        backend_partitions = sym_partition(base_partition[1:end .∉ [[1,j]]])
        initial_element = [vcat(base_set, s(base_partition[j])), vcat(s(base_set), base_partition[j])]
        append!(full_partitions, [vcat(initial_element, b) for b in backend_partitions])
    end
    return full_partitions
end

t = time()

full_partitions = mapreduce(sym_partition, vcat, partitions(1:n))

println((length(full_partitions), time() - t))

(428131, 3.8334240913391113)

This implementation usually runs in a little under four seconds. The Python equivalent takes just under ten. (The order of growth is bad: Kotesovek 2022 shows that the number of partitions grows faster than \(e^{e^n}\), and each partition requires a recursion depth of up to possibly \(n\) to output.)

The number of permutations of length \(2n\) invariant under a permutation of \(n\) independent \(2\)-cycles is given in OEIS A002872.