Maintenant que tout le monde a développé son expertise (souvent étonnante) de codage de bas niveau pour Python, quelle lenteur? (Ou à quelle vitesse votre langue est-elle?) Et à quel point Python est-il vraiment lent (deuxième partie)? Il est temps de relever un défi qui améliorera également votre capacité à améliorer un algorithme.

Le code suivant calcule une liste de longueur 9. Position idans la liste compte le nombre de fois où au moins ides zéros consécutifs ont été trouvés lors du calcul des produits internes entre Fet S. Pour faire cela exactement, il parcourt toutes les listes Fde longueur possibles net les listes Sde longueur n+m-1.

#!/usr/bin/python
import itertools
import operator

n=8
m=n+1
leadingzerocounts = [0]*m
for S in itertools.product([-1,1], repeat = n+m-1):
    for F in itertools.product([-1,1], repeat = n):
        i = 0
        while (i<m and sum(map(operator.mul, F, S[i:i+n])) == 0):
            leadingzerocounts[i] +=1
            i+=1
print leadingzerocounts

La sortie est

[4587520, 1254400, 347648, 95488, 27264, 9536, 4512, 2128, 1064]

Si vous augmentez n à 10,12,14,16,18,20 en utilisant ce code, il devient très vite beaucoup trop lent.

Règles

• Le défi consiste à donner la sortie correcte pour le plus grand nombre possible. Seules les valeurs de n sont pertinentes.
• En cas d’égalité, la victoire revient au code le plus rapide de ma machine pour le plus grand n.
• Je me réserve le droit de ne pas tester le code qui prend plus de 10 minutes.
• Vous pouvez modifier l’algorithme de la manière qui vous convient dans la mesure où il donne le bon résultat. En fait, vous devrez changer l’algorithme pour progresser décemment vers la victoire.
• Le gagnant sera récompensé d'une semaine à partir de la date à laquelle la question a été posée.
• La prime sera remise à son échéance, peu après le vainqueur.

Ma machine Les timings seront exécutés sur ma machine. Il s’agit d’une installation ubuntu standard sur un processeur à huit cœurs AMD FX-8350. Cela signifie également que je dois pouvoir exécuter votre code. En conséquence, utilisez uniquement des logiciels libres facilement disponibles et veuillez inclure des instructions complètes sur la compilation et l’exécution de votre code.

Statut .

• C . n = 12 à 49 secondes par @Fors
• Java . n = 16 à 3:07 de @PeterTaylor
• C ++ . n = 16 en 2:21 par @ilmale
• Rpython . n = 22 en 3:11 par @primo
• Java . n = 22 à 6:56 de @PeterTaylor
• Nimrod . n = 24 en 9:28 secondes par @ReimerBehrends

Le gagnant était Reimer Behrends avec une inscription à Nimrod!

À titre de vérification, la sortie pour n = 22 devrait être [12410090985684467712, 2087229562810269696, 351473149499408384, 59178309967151104, 9975110458933248, 1682628717576192, 284866824372224, 48558946385920, 8416739196928, 1518499004416, 301448822784, 71620493312, 22100246528, 8676573184, 3897278464, 1860960256, 911646720, 451520512, 224785920, 112198656, 56062720, 28031360, 14015680].

La compétition est terminée mais ... Je vais offrir 200 points pour chaque soumission augmentant n de 2 (moins de 10 minutes sur mon ordinateur) jusqu'à ce que je manque de points. Cette offre est ouverte pour toujours .

Nimrod (N = 22)

import math, locks

const
  N = 20
  M = N + 1
  FSize = (1 shl N)
  FMax = FSize - 1
  SStep = 1 shl (N-1)
  numThreads = 16

type
  ZeroCounter = array[0..M-1, int]
  ComputeThread = TThread[int]

var
  leadingZeros: ZeroCounter
  lock: TLock
  innerProductTable: array[0..FMax, int8]

proc initInnerProductTable =
  for i in 0..FMax:
    innerProductTable[i] = int8(countBits32(int32(i)) - N div 2)

initInnerProductTable()

proc zeroInnerProduct(i: int): bool =
  innerProductTable[i] == 0

proc search2(lz: var ZeroCounter, s, f, i: int) =
  if zeroInnerProduct(s xor f) and i < M:
    lz[i] += 1 shl (M - i - 1)
    search2(lz, (s shr 1) + 0, f, i+1)
    search2(lz, (s shr 1) + SStep, f, i+1)

when defined(gcc):
  const
    unrollDepth = 1
else:
  const
    unrollDepth = 4

template search(lz: var ZeroCounter, s, f, i: int) =
  when i < unrollDepth:
    if zeroInnerProduct(s xor f) and i < M:
      lz[i] += 1 shl (M - i - 1)
      search(lz, (s shr 1) + 0, f, i+1)
      search(lz, (s shr 1) + SStep, f, i+1)
  else:
    search2(lz, s, f, i)

proc worker(base: int) {.thread.} =
  var lz: ZeroCounter
  for f in countup(base, FMax div 2, numThreads):
    for s in 0..FMax:
      search(lz, s, f, 0)
  acquire(lock)
  for i in 0..M-1:
    leadingZeros[i] += lz[i]*2
  release(lock)

proc main =
  var threads: array[numThreads, ComputeThread]
  for i in 0 .. numThreads-1:
    createThread(threads[i], worker, i)
  for i in 0 .. numThreads-1:
    joinThread(threads[i])

initLock(lock)
main()
echo(@leadingZeros)

Compiler avec

nimrod cc --threads:on -d:release count.nim

(Nimrod peut être téléchargé ici .)

Cela s'exécute dans le temps imparti pour n = 20 (et pour n = 18 si vous n'utilisez qu'un seul thread, ce qui prend environ 2 minutes dans ce dernier cas).

L'algorithme utilise une recherche récursive, élaguant l'arbre de recherche chaque fois qu'un produit interne différent de zéro est rencontré. Nous avons également réduit l’espace de recherche de moitié en observant que, pour toute paire de vecteurs, (F, -F)il suffit de prendre en compte l’un car l’autre produit les mêmes ensembles de produits internes (en inversant Ségalement).

L'implémentation utilise les installations de métaprogrammation de Nimrod pour dérouler / intégrer les premiers niveaux de la recherche récursive. Cela fait gagner un peu de temps lorsque vous utilisez gcc 4.8 et 4.9 comme back-end de Nimrod et une bonne somme pour clang.

On pourrait encore élaguer l’espace de recherche en observant que nous n'avons besoin que de considérer les valeurs de S qui diffèrent par un nombre pair des N premières positions de notre choix de F. Cependant, la complexité ou les besoins en mémoire de celui-ci ne sont pas proportionnels aux grandes valeurs. de N, étant donné que le corps de la boucle est complètement ignoré dans ces cas.

La tabulation où le produit intérieur est égal à zéro semble être plus rapide que d'utiliser une fonctionnalité de comptage de bits dans la boucle. Apparemment, accéder à la table a une très bonne localité.

Il semble que le problème doive être réglé par une programmation dynamique, tenant compte du fonctionnement de la recherche récursive, mais il n’existe aucun moyen apparent de le faire avec une quantité de mémoire raisonnable.

Exemple de sorties:

N = 16:

@[55276229099520, 10855179878400, 2137070108672, 420578918400, 83074121728, 16540581888, 3394347008, 739659776, 183838720, 57447424, 23398912, 10749184, 5223040, 2584896, 1291424, 645200, 322600]

N = 18:

@[3341140958904320, 619683355033600, 115151552380928, 21392898654208, 3982886961152, 744128512000, 141108051968, 27588886528, 5800263680, 1408761856, 438001664, 174358528, 78848000, 38050816, 18762752, 9346816, 4666496, 2333248, 1166624]

N = 20:

@[203141370301382656, 35792910586740736, 6316057966936064, 1114358247587840, 196906665902080, 34848574013440, 6211866460160, 1125329141760, 213330821120, 44175523840, 11014471680, 3520839680, 1431592960, 655872000, 317675520, 156820480, 78077440, 39005440, 19501440, 9750080, 4875040]

Pour comparer l'algorithme avec d'autres implémentations, N = 16 prend environ 7,9 secondes sur ma machine avec un seul thread et 2,3 secondes avec quatre cœurs.

N = 22 prend environ 15 minutes sur une machine à 64 cœurs avec 4.4.6 gcc, car le backend de Nimrod et déborde d'entiers 64 bits leadingZeros[0](éventuellement non signés, ils ne l'ont pas examiné).

Mise à jour: j'ai trouvé de la place pour quelques améliorations supplémentaires. Premièrement, pour une valeur donnée de F, nous pouvons énumérer les 16 premières entrées des Svecteurs correspondants avec précision, car elles doivent différer exactement à des N/2endroits. Nous avons donc Précalculer une liste de vecteurs de bits de taille Nqui ont des N/2bits ensemble et les utiliser pour tirer la partie initiale de Spartir F.

Deuxièmement, nous pouvons améliorer la recherche récursive en observant que nous connaissons toujours la valeur de F[N](car le MSB a la valeur zéro dans la représentation binaire). Cela nous permet de prédire avec précision dans quelle branche nous recurse à partir du produit intérieur. Cela nous permettrait en fait de transformer toute la recherche en une boucle récursive, mais il arrive en fait que la prédiction de branche soit bousillée un peu, nous avons donc gardé les niveaux supérieurs dans leur forme originale. Nous gagnons encore du temps, principalement en réduisant le nombre de branches que nous effectuons.

Pour certains nettoyages, le code utilise maintenant des entiers non signés et les corrige en 64 bits (juste au cas où quelqu'un voudrait l'exécuter sur une architecture 32 bits).

L'accélération globale est comprise entre un facteur x3 et x4. N = 22 a encore besoin de plus de huit cœurs pour fonctionner en moins de 10 minutes, mais sur une machine à 64 cœurs, le temps passe maintenant à environ quatre minutes (avec une numThreadsmontée en puissance correspondante). Je ne pense pas qu'il y ait beaucoup plus de place à l'amélioration sans un algorithme différent, cependant.

N = 22:

@[12410090985684467712, 2087229562810269696, 351473149499408384, 59178309967151104, 9975110458933248, 1682628717576192, 284866824372224, 48558946385920, 8416739196928, 1518499004416, 301448822784, 71620493312, 22100246528, 8676573184, 3897278464, 1860960256, 911646720, 451520512, 224785920, 112198656, 56062720, 28031360, 14015680]

Mis à jour à nouveau, en utilisant d'autres réductions possibles de l'espace de recherche. Fonctionne dans environ 9:49 minutes pour N = 22 sur ma machine quadcore.

Dernière mise à jour (je pense). De meilleures classes d'équivalence pour les choix de F, réduisant le temps d'exécution de N = 22 à 3:19 minutes 57 secondes (edit: j'avais accidentellement exécuté cela avec un seul thread) sur ma machine.

Ce changement tient compte du fait qu’une paire de vecteurs produit les mêmes zéros d’avant si l’un peut être transformé en un autre en le faisant pivoter. Malheureusement, une optimisation de bas niveau assez critique nécessite que le bit le plus haut de F dans la représentation de bit soit toujours identique, et tout en utilisant cette équivalence, a considérablement réduit l'espace de recherche et réduit le temps d'exécution d'environ un quart par rapport à l'utilisation d'un espace d'état différent. réduction sur F, les frais généraux liés à l'élimination de l'optimisation de bas niveau l'ont plus que compensée. Cependant, il s'avère que ce problème peut être éliminé en considérant également le fait que F inverses les uns des autres sont également équivalents. Bien que cela ait ajouté un peu à la complexité du calcul des classes d'équivalence, cela m'a également permis de conserver l'optimisation de bas niveau susmentionnée, conduisant à une accélération d'environ x3.

Une mise à jour supplémentaire pour prendre en charge les entiers 128 bits pour les données accumulées. Pour compiler avec des entiers 128 bits, vous aurez besoin longint.nimd' ici et avec -d:use128bit. N = 24 prend toujours plus de 10 minutes, mais j'ai inclus le résultat ci-dessous pour les personnes intéressées.

N = 24:

@[761152247121980686336, 122682715414070296576, 19793870419291799552, 3193295704340561920, 515628872377565184, 83289931274780672, 13484616786640896, 2191103969198080, 359662314586112, 60521536552960, 10893677035520, 2293940617216, 631498735616, 230983794688, 102068682752, 48748969984, 23993655296, 11932487680, 5955725312, 2975736832, 1487591936, 743737600, 371864192, 185931328, 92965664]

import math, locks, unsigned

when defined(use128bit):
  import longint
else:
  type int128 = uint64 # Fallback on unsupported architectures
  template toInt128(x: expr): expr = uint64(x)

const
  N = 22
  M = N + 1
  FSize = (1 shl N)
  FMax = FSize - 1
  SStep = 1 shl (N-1)
  numThreads = 16

type
  ZeroCounter = array[0..M-1, uint64]
  ZeroCounterLong = array[0..M-1, int128]
  ComputeThread = TThread[int]
  Pair = tuple[value, weight: int32]

var
  leadingZeros: ZeroCounterLong
  lock: TLock
  innerProductTable: array[0..FMax, int8]
  zeroInnerProductList = newSeq[int32]()
  equiv: array[0..FMax, int32]
  fTable = newSeq[Pair]()

proc initInnerProductTables =
  for i in 0..FMax:
    innerProductTable[i] = int8(countBits32(int32(i)) - N div 2)
    if innerProductTable[i] == 0:
      if (i and 1) == 0:
        add(zeroInnerProductList, int32(i))

initInnerProductTables()

proc ror1(x: int): int {.inline.} =
  ((x shr 1) or (x shl (N-1))) and FMax

proc initEquivClasses =
  add(fTable, (0'i32, 1'i32))
  for i in 1..FMax:
    var r = i
    var found = false
    block loop:
      for j in 0..N-1:
        for m in [0, FMax]:
          if equiv[r xor m] != 0:
            fTable[equiv[r xor m]-1].weight += 1
            found = true
            break loop
        r = ror1(r)
    if not found:
      equiv[i] = int32(len(fTable)+1)
      add(fTable, (int32(i), 1'i32))

initEquivClasses()

when defined(gcc):
  const unrollDepth = 4
else:
  const unrollDepth = 4

proc search2(lz: var ZeroCounter, s0, f, w: int) =
  var s = s0
  for i in unrollDepth..M-1:
    lz[i] = lz[i] + uint64(w)
    s = s shr 1
    case innerProductTable[s xor f]
    of 0:
      # s = s + 0
    of -1:
      s = s + SStep
    else:
      return

template search(lz: var ZeroCounter, s, f, w, i: int) =
  when i < unrollDepth:
    lz[i] = lz[i] + uint64(w)
    if i < M-1:
      let s2 = s shr 1
      case innerProductTable[s2 xor f]
      of 0:
        search(lz, s2 + 0, f, w, i+1)
      of -1:
        search(lz, s2 + SStep, f, w, i+1)
      else:
        discard
  else:
    search2(lz, s, f, w)

proc worker(base: int) {.thread.} =
  var lz: ZeroCounter
  for fi in countup(base, len(fTable)-1, numThreads):
    let (fp, w) = fTable[fi]
    let f = if (fp and (FSize div 2)) == 0: fp else: fp xor FMax
    for sp in zeroInnerProductList:
      let s = f xor sp
      search(lz, s, f, w, 0)
  acquire(lock)
  for i in 0..M-1:
    let t = lz[i].toInt128 shl (M-i).toInt128
    leadingZeros[i] = leadingZeros[i] + t
  release(lock)

proc main =
  var threads: array[numThreads, ComputeThread]
  for i in 0 .. numThreads-1:
    createThread(threads[i], worker, i)
  for i in 0 .. numThreads-1:
    joinThread(threads[i])

initLock(lock)
main()
echo(@leadingZeros)

Java ( n=22?)

Je pense que la plupart des réponses qui font mieux que d' n=16utiliser une approche similaire à celle-ci, bien qu'elles diffèrent par les symétries qu'elles exploitent et la façon dont elles divisent la tâche entre les threads.

Les vecteurs définis dans la question peuvent être remplacés par des chaînes de bits, et le produit interne avec XORing la fenêtre de chevauchement et en vérifiant qu'il existe exactement des n/2bits définis (et donc des n/2bits effacés). Il y a des n! / ((n/2)!)chaînes de nbits (à coefficient binomial central) avec n/2jeu de bits (que j'appelle des chaînes équilibrées ), donc pour toute donnée donnée, Fil y a autant de fenêtres Sdont le produit intérieur est nul. De plus, l'action de glisser le Slong de l'un et de vérifier s'il est encore possible de trouver un bit entrant qui donne un produit intérieur nul correspond à la recherche d'une arête dans un graphe dont les nœuds sont les fenêtres et dont les arêtes relient un nœud uà un nœud vdont les premiers n-1bits sont les derniersn-1des morceaux de u.

Par exemple, avec n=6et F=001001nous obtenons ce graphique:

et pour F=001011nous obtenons ce graphique:

Ensuite , nous devons compter pour chacun ide 0la nnombre de chemins de longueur , iil y a, en sommant sur les graphiques pour chaque F. Je pense que la plupart d'entre nous utilisons la recherche en profondeur d'abord.

Notez que les graphiques sont clairsemés: il est facile de prouver que chaque nœud a un degré d'au plus 1 et un degré au plus un. Cela signifie également que les seules structures possibles sont des chaînes simples et des boucles simples. Cela simplifie un peu le DFS.

J'exploite quelques symétries: les chaînes équilibrées sont fermées sous bit inverse (l' ~opération dans de nombreuses langues de la famille ALGOL), et sous rotation rotation, de sorte que nous puissions regrouper des valeurs Fqui sont liées par ces opérations et ne font que DFS une fois que.

public class CodeGolf26459v8D implements Runnable {
    private static final int NUM_THREADS = 8;

    public static void main(String[] args) {
        v8D(22);
    }

    private static void v8D(int n) {
        int[] bk = new int[1 << n];
        int off = 0;
        for (int i = 0; i < bk.length; i++) {
            bk[i] = Integer.bitCount(i) == n/2 ? off++ : -1;
        }

        int[] fwd = new int[off];
        for (int i = 0; i < bk.length; i++) {
            if (bk[i] >= 0) fwd[bk[i]] = i;
        }

        CodeGolf26459v8D[] runners = new CodeGolf26459v8D[NUM_THREADS];
        Thread[] threads = new Thread[runners.length];
        for (int i = 0; i < runners.length; i++) {
            runners[i] = new CodeGolf26459v8D(n, i, runners.length, bk, fwd);
            threads[i] = new Thread(runners[i]);
            threads[i].start();
        }

        try {
            for (int i = 0; i < threads.length; i++) threads[i].join();
        }
        catch (InterruptedException ie) {
            throw new RuntimeException("This shouldn't be reachable", ie);
        }

        long surviving = ((long)fwd.length) << (n - 1);
        for (int i = 0; i <= n; i++) {
            for (CodeGolf26459v8D runner : runners) surviving -= runner.survival[i];
            System.out.print(i == 0 ? "[" : ", ");
            java.math.BigInteger result = new java.math.BigInteger(Long.toString(surviving));
            System.out.print(result.shiftLeft(n + 1 - i));
        }
        System.out.println("]");
    }

    public final int n;
    protected final int id;
    protected final int numRunners;
    private final int[] bk;
    private final int[] fwd;

    public long[] survival;

    public CodeGolf26459v8D(int n, int id, int numRunners, int[] bk, int[] fwd) {
        this.n = n;
        this.id = id;
        this.numRunners = numRunners;

        this.bk = bk;
        this.fwd = fwd;
    }

    private int dfs2(int[] graphShape, int flip, int i) {
        if (graphShape[i] != 0) return graphShape[i];

        int succ = flip ^ (fwd[i] << 1);
        if (succ >= bk.length) succ ^= bk.length + 1;

        int j = bk[succ];
        if (j == -1) return graphShape[i] = 1;

        graphShape[i] = n + 1; // To detect cycles
        return graphShape[i] = dfs2(graphShape, flip, j) + 1;
    }

    @Override
    public void run() {
        int n = this.n;
        int[] bk = this.bk;
        int[] fwd = this.fwd;

        // NB The initial count is approx 2^(2n - 1.33 - 0.5 lg n)
        // For n=18 we overflow 32-bit
        // 64-bit is good up to n=32.
        long[] survival = new long[n + 1];
        boolean[] visited = new boolean[1 << (n - 1)];
        int th = 0;
        for (int f = 0; f < visited.length; f++) {
            if (visited[f]) continue;

            int m = 1, g = f;
            while (true) {
                visited[g] = true;
                int ng = g << 1;
                if ((ng >> (n - 1)) != 0) ng ^= (1 << n) - 1;
                if (ng == f) break;
                m++;
                g = ng;
            }

            if (th++ % numRunners != id) continue;

            int[] graphShape = new int[fwd.length];
            int flip = (f << 1) ^ f;
            for (int i = 0; i < graphShape.length; i++) {
                int life = dfs2(graphShape, flip, i);
                if (life <= n) survival[life] += m;
            }
        }

        this.survival = survival;
    }
}

Je reçois sur mon Core 2 à 2,5 GHz

# n=18
$ javac CodeGolf26459v8D.java && time java CodeGolf26459v8D
[3341140958904320, 619683355033600, 115151552380928, 21392898654208, 3982886961152, 744128512000, 141108051968, 27588886528, 5800263680, 1408761856, 438001664, 174358528, 78848000, 38050816, 18762752, 9346816, 4666496, 2333248, 1166624]

real    0m3.131s
user    0m10.133s
sys     0m0.380s

# n=20
$ javac CodeGolf26459v8D.java && time java CodeGolf26459v8D
[203141370301382656, 35792910586740736, 6316057966936064, 1114358247587840, 196906665902080, 34848574013440, 6211866460160, 1125329141760, 213330821120, 44175523840, 11014471680, 3520839680, 1431592960, 655872000, 317675520, 156820480, 78077440, 39005440, 19501440, 9750080, 4875040]

real    1m8.706s
user    4m20.980s
sys     0m0.564s

# n=22
$ javac CodeGolf26459v8D.java && time java CodeGolf26459v8D
[12410090985684467712, 2087229562810269696, 351473149499408384, 59178309967151104, 9975110458933248, 1682628717576192, 284866824372224, 48558946385920, 8416739196928, 1518499004416, 301448822784, 71620493312, 22100246528, 8676573184, 3897278464, 1860960256, 911646720, 451520512, 224785920, 112198656, 56062720, 28031360, 14015680]

real    20m10.654s
user    76m53.880s
sys     0m6.852s

Étant donné que l'ordinateur de Lembik possède 8 cœurs et exécute mon programme à thread unique précédent deux fois plus vite que le mien, je suis optimiste qu'il s'exécutera n=22dans moins de 8 minutes.

C

Il s’agit fondamentalement d’une implémentation légèrement optimisée de l’algorithme dans la question. Il peut gérer n=12dans le délai imparti.

#include <stdio.h>
#include <inttypes.h>

#define n 12
#define m (n + 1)

int main() {
    int i;
    uint64_t S, F, o[m] = {0};
    for (S = 0; S < (1LLU << (n + m - 1)); S += 2)
        for (F = 0; F < (1 << (n - 1)); F++)
            for (i = 0; i < m; i++)
                if (__builtin_popcount(((S >> i) & ((1 << n) - 1)) ^ F) == n >> 1)
                    o[i] += 4;
                else
                    break;
    for (i = 0; i < m; i++)
        printf("%" PRIu64 " ", o[i]);
    return 0;
}

Test effectué pour n=12, y compris la compilation:

$ clang -O3 -march=native -fstrict-aliasing -ftree-vectorize -Wall fast.c
$ time ./a.out 
15502147584 3497066496 792854528 179535872 41181184 9826304 2603008 883712 381952 177920 85504 42560 21280 
real    0m53.266s
user    0m53.042s
sys     0m0.068s
$

Commentaire: Je viens d’allumer mon cerveau et d’utiliser une combinatoire simple pour calculer que la première valeur sera toujours n! / ((n / 2)!)^2 * 2^(n + m - 1). Il me semble qu’il faut une solution complètement algébrique à ce problème.

Java, n=16

Pour toute valeur donnée de, Fil existe des \binom{n}{n/2}vecteurs qui ont un produit interne nul. Nous pouvons donc construire un graphe dont les sommets sont les vecteurs correspondants et dont les arêtes correspondent au décalage de S, puis il suffit de compter les chemins de longueur allant jusqu’au ngraphe.

Je n’ai pas essayé de microoptimiser cela en remplaçant les conditionnelles par des opérations au niveau des bits, mais chaque incrément double de n augmente le temps d'exécution d'environ 16 fois, donc cela ne fera pas assez de différence si je ne suis pas trop près du seuil. Sur ma machine, je ne le suis pas.

public class CodeGolf26459 {

    public static void main(String[] args) {
        v3(16);
    }

    // Order of 2^(2n-1) * n ops
    private static void v3(int n) {
        long[] counts = new long[n+1];
        int mask = (1 << n) - 1;
        for (int f = 0; f < (1 << (n-1)); f++) {
            // Find adjacencies
            long[] subcounts = new long[1 << n];
            for (int g = 0; g < (1 << n); g++) {
                subcounts[g] = Integer.bitCount(f ^ g) == n/2 ? 2 : -1;
            }

            for (int round = 0; round <= n; round++) {
                long count = 0;
                // Extend one bit.
                long[] next = new long[1 << n];
                for (int i = 0; i < (1 << n); i++) {
                    long s = subcounts[i];
                    if (s == -1) next[i] = -1;
                    else {
                        count += s;
                        int j = (i << 1) & mask;
                        if (subcounts[j] >= 0) next[j] += s;
                        if (subcounts[j + 1] >= 0) next[j + 1] += s;
                    }
                }
                counts[round] += count << (n - round);
                subcounts = next;
            }
        }

        System.out.print("[");
        for (long count : counts) System.out.print(count+", ");
        System.out.println("]");
    }
}

Je reçois sur mon Core 2 à 2,5 GHz

$ javac CodeGolf26459.java && time java -server CodeGolf26459 
[55276229099520, 10855179878400, 2137070108672, 420578918400, 83074121728, 16540581888, 3394347008, 739659776, 183838720, 57447424, 23398912, 10749184, 5223040, 2584896, 1291424, 645200, 322600, ]

real    6m2.663s
user    6m4.631s
sys     0m1.580s

RPython, N = 22 ~ 3: 23

Multi-threaded, utilisant une descente récursive sans pile. Le programme accepte deux arguments de ligne de commande: N et le nombre de threads de travail.

from time import sleep

from rpython.rlib.rthread import start_new_thread, allocate_lock
from rpython.rlib.rarithmetic import r_int64, build_int, widen
from rpython.rlib.rbigint import rbigint

r_int8 = build_int('r_char', True, 8)

class ThreadEnv:
  __slots__ = ['n', 'counts', 'num_threads',
               'v_range', 'v_num', 'running', 'lock']

  def __init__(self):
    self.n = 0
    self.counts = [rbigint.fromint(0)]
    self.num_threads = 0
    self.v_range = [0]
    self.v_num = 0
    self.running = 0
    self.lock = None

env = ThreadEnv()

bt_bits = 12
bt_mask = (1<<bt_bits)-1
# computed compile time
bit_table = [r_int8(0)]
for i in xrange(1,1<<bt_bits):
  bit_table += [((i&1)<<1) + bit_table[i>>1]]

def main(argv):
  argc = len(argv)
  if argc < 2 or argc > 3:
    print 'Usage: %s N [NUM_THREADS=2]'%argv[0]
    return 1

  if argc == 3:
    env.num_threads = int(argv[2])
  else:
    env.num_threads = 2

  env.n = int(argv[1])
  env.counts = [rbigint.fromint(0)]*env.n
  env.lock = allocate_lock()

  v_range = []
  v_max = 1<<(env.n-1)
  v_num = 0
  v = (1<<(env.n>>1))-1
  while v < v_max:
    v_num += 1
    v_range += [v]
    if v&1:
      # special case odd v
      s = (v+1)&-v
      v ^= s|(s>>1)
    else:
      s = v&-v
      r = v+s
      # s is at least 2, skip two iterations
      i = 3
      s >>= 2
      while s:
        i += 1
        s >>= 1
      v = r|((v^r)>>i)
  env.v_range = v_range
  env.v_num = v_num

  for i in xrange(env.num_threads-1):
    start_new_thread(run,())

  # use the main process as a worker
  run()

  # wait for any laggers
  while env.running:
    sleep(0.05)

  result = []
  for i in range(env.n):
    result += [env.counts[i].lshift(env.n-i+3).str()]
  result += [env.counts[env.n-1].lshift(3).str()]
  print result
  return 0

def run():
  with env.lock:
    v_start = env.running
    env.running += 1

  n = env.n
  counts = [r_int64(0)]*n
  mask = (1<<n)-1
  v_range = env.v_range
  v_num = env.v_num
  z_count = 1<<(n-2)

  for i in xrange(v_start, v_num, env.num_threads):
    v = v_range[i]
    counts[0] += z_count
    counts[1] += v_num
    r = v^(v<<1)
    for w in v_range:
      # unroll counts[2] for speed
      # ideally, we could loop over x directly,
      # rather than over all v, only to throw the majority away
      # there's a 2x-3x speed improvement to be had here...
      x = w^r
      if widen(bit_table[x>>bt_bits]) + widen(bit_table[x&bt_mask]) == n:
        counts[2] += 1
        x, y = v, x
        o, k = 2, 3
        while k < n:
          # x = F ^ S
          # y = F ^ (S<<1)
          o = k
          z = (((x^y)<<1)^y)&mask
          # z is now F ^ (S<<2), possibly xor 1
          # what S and F actually are is of no consequence

          # the compiler hint `widen` let's the translator know
          # to store the result as a native int, rather than a signed char
          bt_high = widen(bit_table[z>>bt_bits])
          if bt_high + widen(bit_table[z&bt_mask]) == n:
            counts[k] += 1
            x, y = y, z
            k += 1
          elif bt_high + widen(bit_table[(z^1)&bt_mask]) == n:
            counts[k] += 1
            x, y = y, z^1
            k += 1
          else: k = n

  with env.lock:
    for i in xrange(n):
      env.counts[i] = env.counts[i].add(rbigint.fromrarith_int(counts[i]))
    env.running -= 1

def target(*args):
  return main, None

Compiler

Créez un clone local du référentiel PyPy en utilisant mercurial, git ou ce que vous préférez. Tapez l'incantation suivante (en supposant que le script ci-dessus est nommé convolution-high.py):

$ pypy %PYPY_REPO%/rpython/bin/rpython --thread convolution-high.py

où %PYPY_REPO%représente une variable d'environnement pointant sur le référentiel que vous venez de cloner. La compilation prend environ une minute.

Exemples de timings

N = 16, 4 fils:

$ timeit convolution-high-c 16 4
[55276229099520, 10855179878400, 2137070108672, 420578918400, 83074121728, 16540581888, 3394347008, 739659776, 183838720, 57447424, 23398912, 10749184, 5223040, 2584896, 1291424, 645200, 322600]
Elapsed Time:     0:00:00.109
Process Time:     0:00:00.390

N = 18, 4 fils:

$ timeit convolution-high-c 18 4
[3341140958904320, 619683355033600, 115151552380928, 21392898654208, 3982886961152, 744128512000, 141108051968, 27588886528, 5800263680, 1408761856, 438001664, 174358528, 78848000, 38050816, 18762752, 9346816, 4666496, 2333248, 1166624]
Elapsed Time:     0:00:01.250
Process Time:     0:00:04.937

N = 20, 4 fils:

$ timeit convolution-high-c 20 4
[203141370301382656, 35792910586740736, 6316057966936064, 1114358247587840, 196906665902080, 34848574013440, 6211866460160, 1125329141760, 213330821120, 44175523840, 11014471680, 3520839680, 1431592960, 655872000, 317675520, 156820480, 78077440, 39005440, 19501440, 9750080, 4875040]
Elapsed Time:     0:00:15.531
Process Time:     0:01:01.328

N = 22, 4 fils:

$ timeit convolution-high-c 22 4
[12410090985684467712, 2087229562810269696, 351473149499408384, 59178309967151104, 9975110458933248, 1682628717576192, 284866824372224, 48558946385920, 8416739196928, 1518499004416, 301448822784, 71620493312, 22100246528, 8676573184, 3897278464, 1860960256, 911646720, 451520512, 224785920, 112198656, 56062720, 28031360, 14015680]
Elapsed Time:     0:03:23.156
Process Time:     0:13:25.437

Python 3.3, N = 20, 3.5min

Clause de non-responsabilité: mon intention n'est PAS de publier cela comme ma propre réponse, car l'algorithme que j'utilise n'est qu'un port sans vergogne de la solution RPython de primo . Mon but ici est seulement de montrer ce que vous pouvez faire en Python si vous combinez la magie de Numpy et Numba modules .

Numba a expliqué en bref:

Numba est un compilateur spécialisé juste à temps qui compile le code Python et NumPy annoté en LLVM (par le biais de décorateurs). http://numba.pydata.org/

Mise à jour 1 : Après avoir jeté les chiffres, j'ai remarqué que nous pouvions simplement ignorer complètement certains des chiffres. Alors maintenant, maxf devient (1 << n) // 2 et maxs devient maxf 2 **. Cela accélérera un peu le processus. n = 16 ne prend maintenant que ~ 48 secondes (au lieu de 4,5 minutes). J'ai aussi une autre idée que je vais essayer de voir si je peux la faire aller un peu plus vite.

Mise à jour 2: algorithme modifié (solution de primo). Bien que mon port ne supporte pas encore le multithreading, il est assez trivial d'ajouter. Il est même possible de publier CPython GIL en utilisant Numba et ctypes. Cependant, cette solution fonctionne aussi très vite sur un seul cœur!

import numpy as np
import numba as nb

bt_bits = 11
bt_mask = (1 << bt_bits) - 1
bit_table = np.zeros(1 << bt_bits, np.int32)

for i in range(0, 1 << bt_bits):
    bit_table[i] = ((i & 1) << 1) + bit_table[i >> 1]

@nb.njit("void(int32, int32, int32, int32, int64[:], int64[:])")
def run(n, m, start, re, counts, result):
    mask = (1 << n) - 1

    v_max = (1 << n) // 2
    rr = v_max // 2

    v = (1 << (n >> 1)) - 1
    while v < v_max:
        s = start

        while s < rr:
            f = v ^ s
            counts[0] += 8
            t = s << 1
            o, j = 0, 1

            while o < j and j < m:
                o = j
                w = (t ^ f) & mask
                bt_high = bit_table[w >> bt_bits]

                if bt_high + bit_table[w & bt_mask] == n:
                    counts[j] += 8
                    t <<= 1
                    j += 1
                elif bt_high + bit_table[(w ^ 1) & bt_mask] == n:
                    counts[j] += 8
                    t = (t | 1) << 1
                    j += 1
                    s += re

            s = v & -v
            r = v + s
            o = v ^ r
            o = (o >> 2) // s
            v = r | o

    for e in range(m):
        result[e] += counts[e] << (n - e)

Et enfin:

if __name__ == "__main__":
    n = 20
    m = n + 1

    result = np.zeros(m, np.int64)
    counts = np.zeros(m, np.int64)

    s1 = time.time() * 1000
    run(n, m, 0, 1, counts, result)
    s2 = time.time() * 1000

    print(result)
    print("{0}ms".format(s2 - s1))

Cela fonctionne sur ma machine en 212688ms ou ~ 3.5min.

C ++ N = 16

Je teste sur un EEEPC avec un atome .. mon temps n'a pas beaucoup de sens. :RÉ
L'atome résolvent n = 14 en 34 secondes. Et n = 16 en 20 minutes. Je veux tester n = 16 sur OP pc. Je suis optimiste.

L'idée est que chaque fois que nous trouvons une solution pour un F donné, nous avons trouvé 2 ^ i solution car nous pouvons changer la partie inférieure de S pour aboutir au même résultat.

#include <stdio.h>
#include <cinttypes>
#include <cstring>

int main()
{
   const int n = 16;
   const int m = n + 1;
   const uint64_t maxS = 1ULL << (2*n);
   const uint64_t maxF = 1ULL << n;
   const uint64_t mask = (1ULL << n)-1;
   uint64_t out[m]={0};
   uint64_t temp[m] = {0};
   for( uint64_t F = 0; F < maxF; ++F )
   {
      for( uint64_t S = 0; S < maxS; ++S )
      {
         int numSolution = 1;
         for( int i = n; i >= 0; --i )
         {
            const uint64_t window = S >> i;
            if( __builtin_popcount( mask & (window ^ F) ) == (n / 2) )
            {
               temp[i] += 1;
            } else {
               numSolution = 1 << i;
               S += numSolution - 1;
               break;
            }
         }
         for( int i = n; i >= 0; --i )
         {
            out[i] += temp[i]*numSolution;
            temp[i] = 0;
         }
      }
   }
   for( int i = n; i >= 0; --i )
   {
      uint64_t x = out[i];
      printf( "%lu ", x );
   }
   return 0;
}

compiler:

gcc 26459.cpp -std = c ++ 11 -O3 -march = natif -fstrict-aliasing -ftree-vectorize -Wall -pedantic -o 26459

JAVASCRIPT n: 12

Dans mon ordinateur, cela a pris 231,242 secondes. Dans la démo, j'utilise des webworkers pour éviter de geler le navigateur. Ceci peut être encore amélioré avec les travailleurs parallèles. Je sais que JS n'a aucune chance dans ce défi mais je l'ai fait pour le plaisir!

Cliquez pour lancer la démo en ligne

var n = 8;        
var m = n + 1;
var o = [];
var popCount = function(bits) {
  var SK5  = 0x55555555,
      SK3  = 0x33333333,
      SKF0 = 0x0f0f0f0f,
      SKFF = 0xff00ff;

  bits -= (bits >> 1) & SK5;
  bits  = (bits & SK3) + ((bits >> 2) & SK3);
  bits  = (bits & SKF0) + ((bits >> 4) & SKF0);
  bits += bits >> 8;

  return (bits + (bits >> 15)) & 63;
};
for(var S = 0; S < (1 << n + m - 1); S += 2){
  for(var F = 0; F < (1 << n - 1); F += 1){
    for (var i = 0; i < m; i++){
      var c = popCount(((S >> i) & ((1 << n) - 1)) ^ F);
      if(c == n >> 1){
        if(!o[i]) o[i] = 0;
        o[i] += 4;
      } else break;
    }
  }
}
return o;

Fortran: n = 12

Je viens de faire une version quick'n'dirty dans Fortran, aucune optimisation sauf OpenMP. Il faut presser juste au-dessous de 10 minutes pour que n = 12 sur la machine OPs, il faut 10:39 sur ma machine, ce qui est un peu plus lent.

Les entiers 64 bits ont en effet un impact négatif sur les performances. Je suppose que je devrais repenser l’ensemble de l’algorithme pour que cela soit beaucoup plus rapide. Je ne sais pas si je vais m'embêter, je pense que je vais plutôt passer un peu de temps à imaginer un bon défi qui soit plus à mon goût. Si quelqu'un d'autre veut le prendre et le suivre, allez-y :)

program golf
use iso_fortran_env
implicit none
integer, parameter ::  n=12
integer :: F(n), S(2*n)
integer(int64) :: leadingzerocounts(n+1)
integer :: k
integer(int64) :: i,j,bindec,enc

leadingzerocounts=0

!$OMP parallel do private(i,enc,j,bindec,S,F,k) reduction(+:leadingzerocounts) schedule(dynamic)
do i=0,2**(2*n)-1
  enc=i
  ! Short loop to convert i into the array S with -1s and 1s
  do j=2*n,1,-1
    bindec=2**(j-1)
    if (enc-bindec .ge. 0) then
      S(j)=1
      enc=enc-bindec
    else
      S(j)=-1
    endif
  end do
  do j=0,2**(n)-1
    ! Convert j into the array F with -1s and 1s
    enc=j
    do k=n,1,-1
      bindec=2**(k-1)
      if (enc-bindec .ge. 0) then
        F(k)=1
        enc=enc-bindec
      else
        F(k)=-1
      endif
    end do
    ! Compute dot product   
    do k=1,n+1
      if (dot_product(F,S(k:k+n-1)) /= 0) exit
      leadingzerocounts(k)=leadingzerocounts(k)+1
    end do
  end do
end do
!$OMP end parallel do

print *, leadingzerocounts

end

Lua: n = 16

Clause de non-responsabilité: mon intention n'est PAS de publier ceci comme ma propre réponse, car l'algorithme que j'utilise est volé sans vergogne de la réponse intelligente d' Anna Jokela . qui a été volé sans vergogne de la réponse intelligente de Ilmale .

En outre, ce n'est même pas valide - il y a des inexactitudes causées par des nombres à virgule flottante (il serait préférable que Lua prenne en charge les entiers 64 bits). Cependant, je suis toujours en train de le télécharger, juste pour montrer à quelle vitesse cette solution est rapide. C'est un langage de programmation dynamique, et pourtant je peux calculer n = 16 en un temps raisonnable (1 minute avec un processeur à 800 MHz).

Exécutée avec LuaJIT, l'interpréteur standard est trop lent.

local bit = require "bit"
local band = bit.band
local bor = bit.bor
local bxor = bit.bxor
local lshift = bit.lshift
local rshift = bit.rshift

-- http://stackoverflow.com/a/11283689/736054
local function pop_count(w)
    local b1 = 1431655765
    local b2 = 858993459
    local b3 = 252645135
    local b7 = 63

    w = band(rshift(w, 1), b1) + band(w, b1)
    w = band(rshift(w, 2), b2) + band(w, b2)
    w = band(w + rshift(w, 4), b3)
    return band(rshift(w, 24) + rshift(w, 16) + rshift(w, 8) + w, b7)
end

local function gen_array(n, value)
    value = value or 0
    array = {}
    for i = 1, n do
        array[i] = value
    end
    return array
end

local n = 16
local u = math.floor(n / 2)
local m = n + 1
local maxf = math.floor(lshift(1, n) / 2)
local maxs = maxf ^ 2
local mask = lshift(1, n) - 1

local out = gen_array(m, 0)
local temp = gen_array(m, 0)


for f = 0, maxf do
    local s = 0
    while s <= maxs do
        local num_solution = 1

        for i = n, 0, -1 do
            if pop_count(band(mask, bxor(rshift(s, i), f))) == u then
                temp[i + 1] = temp[i + 1] + 8
            else
                num_solution = lshift(1, i)
                s = s + num_solution - 1
                break
            end
        end

        for i = 1, m do
            out[i] = out[i] + temp[i] * num_solution
            temp[i] = 0
        end

        s = s + 1
    end
end

for i = m, 1, -1 do
    print(out[i])
end