Zählen Sie die Anzahl der Hankelable-Matrizen

12

Hintergrund

Eine binäre Hankel-Matrix ist eine Matrix mit konstanten Schrägdiagonalen (positiv abfallenden Diagonalen), die nur 0s und 1s enthält. ZB sieht eine 5x5 binäre Hankelmatrix so aus

a b c d e
b c d e f
c d e f g
d e f g h
e f g h i

wo a, b, c, d, e, f, g, h, isind entweder 0oder 1.

Definieren wir eine Matrix M als Hankelbar, wenn es eine Permutation der Reihenfolge der Zeilen und Spalten von M gibt, so dass M eine Hankel-Matrix ist. Dies bedeutet, dass man eine Permutation auf die Reihenfolge der Zeilen und eine möglicherweise andere auf die Spalten anwenden kann.

Die Herausforderung

Die Herausforderung besteht darin , zu zählen , wie viele Hankelable n von nMatrizen für alle da sindn bis zu so großen Wert wie möglich.

Ausgabe

Geben Sie für jede ganze Zahl n ab 1 die Anzahl der Hankelablen von nMatrizen mit Einträgen aus, die 0oder sind 1.

Dafür sollten n = 1,2,3,4,5die Antworten sein2,12,230,12076,1446672 . (Danke an orlp für den Code, der diese erstellt hat.)

Zeitlimit

Ich werde Ihren Code auf meinem Computer ausführen und ihn nach 1 Minute stoppen. Der Code, der die richtigen Antworten bis zum größten Wert von n ausgibt, gewinnt. Die Frist ist für alles von n = 1bis zum größten Wert vonn auf den Sie eine Antwort geben.

Der Gewinner wird die beste Antwort bis Ende Samstag, den 18. April sein.

Kabelbinder

Bei einem Unentschieden nwerde ich mal sehen, wie lange es dauert, bis die Ausgänge erreicht sind n+1und der Schnellste gewinnt. In dem Fall, dass sie in der gleichen Zeit innerhalb einer Sekunde bis zu laufenn+1 , gewinnt die erste Einreichung.

Sprachen und Bibliotheken

Sie können jede Sprache verwenden, die über einen frei verfügbaren Compiler / Interpreter / etc. Verfügt. für Linux und alle Bibliotheken, die auch für Linux frei verfügbar sind.

Meine Maschine

Die Timings werden auf meinem Computer ausgeführt. Dies ist eine Ubuntu-Standardinstallation auf einem AMD FX-8350-Prozessor mit acht Kernen auf einem Asus M5A78L-M / USB3-Motherboard (Sockel AM3 +, 8 GB DDR3). Dies bedeutet auch, dass ich in der Lage sein muss, Ihren Code auszuführen. Verwenden Sie daher nur leicht verfügbare kostenlose Software und fügen Sie vollständige Anweisungen zum Kompilieren und Ausführen Ihres Codes bei.

Anmerkungen

Ich empfehle, nicht alle n-mal-n-Matrizen zu iterieren und zu ermitteln, ob jede die von mir beschriebene Eigenschaft hat. Erstens gibt es zu viele und zweitens scheint es keine schnelle Möglichkeit zu geben, diese Erkennung durchzuführen .

Bisher führende Einträge

  • n = 8 von Peter Taylor. Java
  • n = 5 von orlp. Python
Gemeinschaft
quelle
4
Ich mag das Wort "Hankelable" wirklich.
Alex A.
3
Für n=6die Summe ist 260357434. Ich denke, der Speicherdruck ist ein größeres Problem als die CPU-Zeit.
Peter Taylor
Das ist eine wunderbare Frage. Ich wurde gründlich von Nerds geschnappt.
Alexander-Brett

Antworten:

7

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;
        }
    }
}

Speichern als HankelCombinatorics.java, Kompilieren als javac HankelCombinatorics.java, Ausführen als java -Xmx2G HankelCombinatorics.

Mit NUM_THREADS = 4auf meiner Quad-Core - Maschine wird es 20420819767436für n=8in 50 bis 55 Sekunden verstrichen, mit einem fairen Betrag von Variabilität zwischen den Läufen; Ich gehe davon aus, dass es auf Ihrem Octa-Core-Computer problemlos funktioniert, aber es dauert mindestens eine Stunde, bis es verfügbar istn=9 .

Wie es funktioniert

Vorausgesetzt n, es gibt 2^(2n-1)binäre nx nHankel-Matrizen. Die Zeilen können auf verschiedene n!Arten und die Spalten auf verschiedene n!Arten permutiert werden. Alles was wir tun müssen, ist Doppelzählungen zu vermeiden ...

Wenn Sie die Summe jeder Zeile berechnen, ändert weder das Permutieren der Zeilen noch das Permutieren der Spalten die Summenmehrfachmenge. Z.B

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

hat Zeilensummen-Multiset {3, 3, 2, 2, 2} und alle daraus abgeleiteten Hankelable-Matrizen. Dies bedeutet, dass wir die Hankel-Matrizen nach diesen Zeilensummen-Multisätzen gruppieren und dann jede Gruppe unabhängig behandeln können, wobei mehrere Prozessorkerne ausgenutzt werden.

Es gibt auch eine ausnutzbare Symmetrie: Die Matrizen mit mehr Nullen als Einsen stehen im Widerspruch zu den Matrizen mit mehr Einsen als Nullen.

Doppelzählung auftritt , wenn Hankel - Matrix M_1mit Zeilen Permutation r_1und Spaltenpermutation c_1einstimmt Hankel - Matrix M_2mit Zeilen Permutation r_2und Spaltenpermutation c_2(mit bis zu zwei , aber nicht alle drei von M_1 = M_2, r_1 = r_2, c_1 = c_2). Die Zeilen- und Spaltenpermutationen sind unabhängig. Wenn Sie also die Zeilenpermutation r_1auf M_1und die Zeilenpermutation r_2auf anwenden M_2, müssen die Spalten als Multisets gleich sein. Daher berechne ich für jede Gruppe alle Spalten-Multisets, die durch Anwenden einer Zeilenpermutation auf eine Matrix in der Gruppe erhalten werden. Der einfache Weg, eine kanonische Darstellung der Multisets zu erhalten, besteht darin, die Spalten zu sortieren. Dies ist auch im nächsten Schritt nützlich.

Nachdem wir die verschiedenen Spalten-Multisets erhalten haben, müssen wir herausfinden, wie viele der n!Permutationen von jeder einzigartig sind. Zu diesem Zeitpunkt kann eine Doppelzählung nur stattfinden, wenn ein gegebenes Spalten-Multiset doppelte Spalten enthält. Wir müssen lediglich die Anzahl der Vorkommen jeder einzelnen Spalte im Multiset zählen und dann den entsprechenden Multinomialkoeffizienten berechnen. Da die Spalten sortiert sind, ist es einfach, sie zu zählen.

Schließlich addieren wir sie alle.

Die asymptotische Komplexität ist nicht einfach mit voller Genauigkeit zu berechnen, da wir einige Annahmen über die Mengen treffen müssen. Wir bewerten die Reihenfolge der 2^(2n-2) n!Spalten-Multisets, wobei wir uns jeweils n^2 ln nZeit nehmen (einschließlich der Sortierung). Wenn die Gruppierung nicht mehr als einen ln nFaktor ausmacht, haben wir zeitliche Komplexität Theta(4^n n! n^2 ln n). Aber da die Exponentialfaktoren die Polynomfaktoren vollständig dominieren, ist dies der Fall Theta(4^n n!) = Theta((4n/e)^n).

Peter Taylor
quelle
Das ist sehr beeindruckend. Können Sie etwas über den von Ihnen verwendeten Algorithmus sagen?
3

Python2 / 3

Ziemlich naiver Ansatz in einer langsamen Sprache:

import itertools

def permute_rows(m):
    for perm in itertools.permutations(m):
        yield perm

def permute_columns(m):
    T = zip(*m)
    for perm in itertools.permutations(T):
        yield zip(*perm)

N = 1
while True:
    base_template = ["abcdefghijklmnopqrstuvwxyz"[i:i+N] for i in range(N)]

    templates = set()
    for c in permute_rows(base_template):
        for m in permute_columns(c):
            templates.add("".join("".join(row) for row in m))

    def possibs(free, templates):
        if free == 2*N - 1:
            return set(int(t, 2) for t in templates)

        s = set()
        for b in "01":
            new_templates = set(t.replace("abcdefghijklmnopqrstuvwxyz"[free], b) for t in templates)
            s |= possibs(free + 1, new_templates)

        return s

    print(len(possibs(0, templates)))
    N += 1

Führen Sie durch Eingabe python script.py.

orlp
quelle
Sie haben die Sprache als Python 2/3 aufgeführt, aber damit sie in Python 2 funktioniert, benötigen Sie from __future__ import print_function(oder so etwas) nicht?
Alex A.
2
@AlexA. Normalerweise ja, aber in diesem Fall nicht. Berücksichtigen Sie bei der Eingabe das Verhalten von Python2 return(1). Jetzt ersetzen returndurch print:)
orlp
Cool! Ich lerne jeden Tag etwas Neues. :)
Alex A.
2

Haskell

import Data.List
import Data.Hashable
import Control.Parallel.Strategies
import Control.Parallel
import qualified Data.HashSet as S

main = mapM putStrLn $ map (show.countHankellable) [1..]

a§b=[a!!i|i<-b]

hashNub :: (Hashable a, Eq a) => [a] -> [a]
hashNub l = go S.empty l
    where
      go _ []     = []
      go s (x:xs) = if x `S.member` s then go s xs
                                    else x : go (S.insert x s) xs

pmap = parMap rseq

makeMatrix :: Int->[Bool]->[[Bool]]
makeMatrix n vars = [vars§[i..i+n-1]|i<-[0..n-1]]

countHankellable :: Int -> Int
countHankellable n = let
    s = permutations [0..n-1]
    conjugates m = concat[permutations[r§q|r<-m]|q<-s]
    variableSets = sequence [[True,False]|x<-[0..2*(n-1)]]
 in
    length.hashNub.concat.pmap (conjugates.makeMatrix n ) $ variableSets

Nirgendwo so schnell wie bei Peter - das ist ein ziemlich beeindruckendes Setup, das er dort hat! Jetzt mit viel mehr Code aus dem Internet kopiert. Verwendung:

$ ghc -threaded hankell.hs
$ ./hankell
Alexander-Brett
quelle
Eine Antwort von Haskell ist immer willkommen. Vielen Dank.
@Lembik - wie geht es mir auf deiner Maschine?
Alexander-Brett