Java (n=8)
import java.util.*;
import java.util.concurrent.*;
public class HankelCombinatorics {
public static final int NUM_THREADS = 8;
private static final int[] FACT = new int[13];
static {
FACT[0] = 1;
for (int i = 1; i < FACT.length; i++) FACT[i] = i * FACT[i-1];
}
public static void main(String[] args) {
long prevElapsed = 0, start = System.nanoTime();
for (int i = 1; i < 12; i++) {
long count = count(i), elapsed = System.nanoTime() - start;
System.out.format("%d in %dms, total elapsed %dms\n", count, (elapsed - prevElapsed) / 1000000, elapsed / 1000000);
prevElapsed = elapsed;
}
}
@SuppressWarnings("unchecked")
private static long count(int n) {
int[][] perms = new int[FACT[n]][];
genPermsInner(0, 0, new int[n], perms, 0);
// We partition by canonical representation of the row sum multiset, discarding any with a density > 50%.
Map<CanonicalMatrix, Map<CanonicalMatrix, Integer>> part = new HashMap<CanonicalMatrix, Map<CanonicalMatrix, Integer>>();
for (int m = 0; m < 1 << (2*n-1); m++) {
int density = 0;
int[] key = new int[n];
for (int i = 0; i < n; i++) {
key[i] = Integer.bitCount((m >> i) & ((1 << n) - 1));
density += key[i];
}
if (2 * density <= n * n) {
CanonicalMatrix _key = new CanonicalMatrix(key);
Map<CanonicalMatrix, Integer> map = part.get(_key);
if (map == null) part.put(_key, map = new HashMap<CanonicalMatrix, Integer>());
map.put(new CanonicalMatrix(m, perms[0]), m);
}
}
List<Job> jobs = new ArrayList<Job>();
ExecutorService pool = Executors.newFixedThreadPool(NUM_THREADS);
for (Map.Entry<CanonicalMatrix, Map<CanonicalMatrix, Integer>> e : part.entrySet()) {
Job job = new Job(n, perms, e.getKey().sum() << 1 == n * n ? 0 : 1, e.getValue());
jobs.add(job);
pool.execute(job);
}
pool.shutdown();
try {
pool.awaitTermination(1, TimeUnit.DAYS); // i.e. until it's finished - inaccurate results are useless
}
catch (InterruptedException ie) {
throw new IllegalStateException(ie);
}
long total = 0;
for (Job job : jobs) total += job.subtotal;
return total;
}
private static int genPermsInner(int idx, int usedMask, int[] a, int[][] perms, int off) {
if (idx == a.length) perms[off++] = a.clone();
else for (int i = 0; i < a.length; i++) {
int m = 1 << (a[idx] = i);
if ((usedMask & m) == 0) off = genPermsInner(idx+1, usedMask | m, a, perms, off);
}
return off;
}
static class Job implements Runnable {
private volatile long subtotal = 0;
private final int n;
private final int[][] perms;
private final int shift;
private final Map<CanonicalMatrix, Integer> unseen;
public Job(int n, int[][] perms, int shift, Map<CanonicalMatrix, Integer> unseen) {
this.n = n;
this.perms = perms;
this.shift = shift;
this.unseen = unseen;
}
public void run() {
long result = 0;
int[][] perms = this.perms;
Map<CanonicalMatrix, Integer> unseen = this.unseen;
while (!unseen.isEmpty()) {
int m = unseen.values().iterator().next();
Set<CanonicalMatrix> equiv = new HashSet<CanonicalMatrix>();
for (int[] perm : perms) {
CanonicalMatrix canonical = new CanonicalMatrix(m, perm);
if (equiv.add(canonical)) {
result += canonical.weight() << shift;
unseen.remove(canonical);
}
}
}
subtotal = result;
}
}
static class CanonicalMatrix {
private final int[] a;
private final int hash;
public CanonicalMatrix(int m, int[] r) {
this(permuteRows(m, r));
}
public CanonicalMatrix(int[] a) {
this.a = a;
Arrays.sort(a);
int h = 0;
for (int i : a) h = h * 37 + i;
hash = h;
}
private static int[] permuteRows(int m, int[] perm) {
int[] cols = new int[perm.length];
for (int i = 0; i < perm.length; i++) {
for (int j = 0; j < cols.length; j++) cols[j] |= ((m >> (perm[i] + j)) & 1L) << i;
}
return cols;
}
public int sum() {
int sum = 0;
for (int i : a) sum += i;
return sum;
}
public int weight() {
int prev = -1, count = 0, weight = FACT[a.length];
for (int col : a) {
if (col == prev) weight /= ++count;
else {
prev = col;
count = 1;
}
}
return weight;
}
@Override public boolean equals(Object obj) {
// Deliberately unsuitable for general-purpose use, but helps catch bugs faster.
CanonicalMatrix that = (CanonicalMatrix)obj;
for (int i = 0; i < a.length; i++) {
if (a[i] != that.a[i]) return false;
}
return true;
}
@Override public int hashCode() {
return hash;
}
}
}
Save as HankelCombinatorics.java, compile as javac HankelCombinatorics.java, run as java -Xmx2G HankelCombinatorics.
With NUM_THREADS = 4 on my quad-core machine it gets 20420819767436 for n=8 in 50 to 55 seconds elapsed, with a fair amount of variability between runs; I expect that it should easily manage the same on your octa-core machine but will take an hour or more to get n=9.
How it works
Given n, there are 2^(2n-1) binary nxn Hankel matrices. The rows can be permuted in n! ways, and the columns in n! ways. All we need to do is to avoid double-counting...
If you calculate the sum of each row, then neither permuting the rows nor permuting the columns changes the multiset of sums. E.g.
0 1 1 0 1
1 1 0 1 0
1 0 1 0 0
0 1 0 0 1
1 0 0 1 0
has row sum multiset {3, 3, 2, 2, 2}, and so do all Hankelable matrices derived from it. This means that we can group the Hankel matrices by these row sum multisets and then handle each group independently, exploiting multiple processor cores.
There's also an exploitable symmetry: the matrices with more zeroes than ones are in bijection with the matrices with more ones than zeroes.
Double-counting occurs when Hankel matrix M_1 with row permutation r_1 and column permutation c_1 matches Hankel matrix M_2 with row permutation r_2 and column permutation c_2 (with up to two but not all three of M_1 = M_2, r_1 = r_2, c_1 = c_2). The row and column permutations are independent, so if we apply row permutation r_1 to M_1 and row permutation r_2 to M_2, the columns as multisets must be equal. So for each group, I calculate all of the column multisets obtained by applying a row permutation to a matrix in the group. The easy way to get a canonical representation of the multisets is to sort the columns, and that is also useful in the next step.
Having obtained the distinct column multisets, we need to find how many of the n! permutations of each are unique. At this point, double-counting can only occur if a given column multiset has duplicate columns: what we need to do is to count the number of occurrences of each distinct column in the multiset and then calculate the corresponding multinomial coefficient. Since the columns are sorted, it's easy to do the count.
Finally we add them all up.
The asymptotic complexity isn't trivial to calculate to full precision, because we need to make some assumptions about the sets. We evaluate on the order of 2^(2n-2) n! column multisets, taking n^2 ln n time for each (including the sorting); if grouping doesn't take more than a ln n factor, we have time complexity Theta(4^n n! n^2 ln n). But since the exponential factors completely dominate the polynomial ones, it's Theta(4^n n!) = Theta((4n/e)^n).
n=6the total is260357434. I think memory pressure is a bigger issue than CPU time. \$\endgroup\$