Je promets que ce sera mon dernier défi concernant les diamants (pour un temps, en tout cas). Le bon côté des choses, ce défi n’a rien à voir avec l’art ASCII, et n’est pas non plus un code de golf, c’est donc complètement différent.

Pour rappel, chaque hexagone peut porter trois diamants différents:

Une question intéressante à poser est de savoir combien de ces pavages existent pour une taille d’hexagone donnée. Il semble que ces chiffres ont été assez étudiés et figurent dans le document OEIS A008793 .

Cependant, le problème devient plus compliqué si nous demandons combien de niveaux existent jusqu’à rotation et réflexion . Par exemple, pour une longueur de côté N = 2, il existe les 20 pavages suivants:

   ____     ____     ____     ____     ____     ____     ____     ____     ____     ____  
  /\_\_\   /\_\_\   /\_\_\   /\_\_\   /_/\_\   /_/\_\   /\_\_\   /_/\_\   /_/\_\   /_/\_\ 
 /\/\_\_\ /\/_/\_\ /\/_/_/\ /\/_/\_\ /\_\/\_\ /\_\/_/\ /\/_/_/\ /\_\/\_\ /\_\/_/\ /_/\/\_\
 \/\/_/_/ \/\_\/_/ \/\_\_\/ \/_/\/_/ \/\_\/_/ \/\_\_\/ \/_/\_\/ \/_/\/_/ \/_/\_\/ \_\/\/_/
  \/_/_/   \/_/_/   \/_/_/   \_\/_/   \/_/_/   \/_/_/   \_\/_/   \_\/_/   \_\/_/   \_\/_/ 
   ____     ____     ____     ____     ____     ____     ____     ____     ____     ____  
  /_/_/\   /\_\_\   /_/\_\   /_/_/\   /_/\_\   /_/\_\   /_/_/\   /_/_/\   /_/_/\   /_/_/\ 
 /\_\_\/\ /\/_/_/\ /_/\/_/\ /\_\_\/\ /\_\/_/\ /_/\/_/\ /_/\_\/\ /\_\_\/\ /_/\_\/\ /_/_/\/\
 \/\_\_\/ \/_/_/\/ \_\/\_\/ \/_/\_\/ \/_/_/\/ \_\/_/\/ \_\/\_\/ \/_/_/\/ \_\/_/\/ \_\_\/\/
  \/_/_/   \_\_\/   \_\/_/   \_\/_/   \_\_\/   \_\_\/   \_\/_/   \_\_\/   \_\_\/   \_\_\/ 

Mais beaucoup d'entre eux sont identiques en rotation et en réflexion. Si l'on prend en compte ces symétries, il ne reste que 6 pavages distincts:

   ____     ____     ____     ____     ____     ____  
  /\_\_\   /\_\_\   /\_\_\   /_/\_\   /_/\_\   /_/\_\ 
 /\/\_\_\ /\/_/\_\ /\/_/_/\ /\_\/_/\ /\_\/_/\ /_/\/\_\
 \/\/_/_/ \/\_\/_/ \/\_\_\/ \/\_\_\/ \/_/\_\/ \_\/\/_/
  \/_/_/   \/_/_/   \/_/_/   \/_/_/   \_\/_/   \_\/_/ 

   2        2        6        6        1        3

où les chiffres indiquent la multiplicité de chaque mosaïque. Notez que pour les grands hexagones, il y a aussi des pavages avec les multiplicités 4 et 12.

Il semble que le nombre de pavages jusqu’à la symétrie ait été étudié de manière moins approfondie. L’entrée A066931 d’ OEIS ne contient que les cinq termes:

1, 1, 6, 113, 20174

où le premier terme est pour la longueur du côté N = 0et le dernier terme pour la longueur du côté N = 4.

Je suis sûr que nous pouvons faire mieux que ça!

Votre tâche consiste à calculer le nombre de pavages pour une longueur de côté donnée.

C'est le code le plus rapide . Votre score sera la plus grande longueur Nde côté pour laquelle votre code produira le résultat correct dans les 30 minutes qui suivent sur ma machine. En cas d'égalité, je vais accepter la soumission qui produit le résultat pour que N le plus rapide.

Comme d'habitude, vous ne devez pas coder en dur les résultats que vous connaissez déjà pour gagner le match décisif. L'algorithme à résoudre N = 3doit être identique à celui qui le résout N = 5.

Votre soumission ne doit pas utiliser plus de 4 Go de mémoire. Je vais laisser une marge de manœuvre à ce sujet si vous exploitez près de cette limite, mais si vous dépassez constamment cette limite ou si vous dépassez considérablement cette limite, je ne compterai pas cela Ndans votre déclaration.

Je vais tester toutes les soumissions sur mon ordinateur Windows 8, alors assurez-vous que la langue de votre choix est librement disponible sur Windows. Mathematica est la seule exception à cette règle (car j’en ai une licence). S'il vous plaît inclure des instructions sur la façon de compiler / exécuter votre code.

Bien sûr, sentez-vous libre de calculer plus de termes à votre rythme (pour la science, et pour que les autres vérifient leur nombre), mais le score de votre réponse sera déterminé par ces 30 minutes.

Algèbre, théorie des graphes, inversion de Möbius, recherche et Java

Le groupe de symétrie de l'hexagone est le groupe dièdre d'ordre 12 et est généré par une rotation de 60 degrés et un retournement de miroir sur un diamètre. Il comporte 16 sous-groupes, mais certains d'entre eux appartiennent à des groupes de conjugaison non triviaux (ceux qui n'ont que des réflexions ont 3 choix d'axe). Il y a donc 10 symétries fondamentalement différentes qu'un pavage de l'hexagone peut avoir:

Le nombre de pavages de diamant d’un sous-ensemble d’un réseau triangulaire peut être calculé comme déterminant ; mon approche initiale consistait donc à définir un déterminant pour chacune des symétries de l’hexagone, afin de calculer le nombre de pavages présentant au moins ces symétries. ; et ensuite utiliser l'inversion de Möbius dans l' algèbre d'incidence de leur poset (essentiellement une généralisation du principe d'inclusion-exclusion) pour calculer le nombre de pavages dont le groupe de symétrie correspond exactement à chacun des 10 cas. Cependant, certaines symétries ont des conditions de bord désagréables, ce qui m'a obligé à additionner de manière déterminante de nombreux déterminants. Heureusement, les valeurs obtenues pourn < 10m'a donné suffisamment de données pour pouvoir identifier les séquences pertinentes dans OEIS et reconstituer une forme fermée (pour une valeur de "fermé" qui permet des produits finis). Il y a un peu de discussion sur les séquences, et de références aux preuves, dans le résumé formel que j'ai préparé pour justifier les mises à jour de séquence OEIS.

Une fois le double comptage pris en charge, il s'avère que quatre des dix valeurs s'annulent proprement. Il ne reste donc qu'à calculer les six autres valeurs, puis à effectuer une somme pondérée.

Ce code prend moins de 30 secondes N=1000sur ma machine.

import java.math.BigInteger;

public class OptimisedCounter {
    private static int[] minp = new int[2];

    public static void main(String[] args) {
        if (args.length > 0) {
            for (String arg : args) System.out.println(count(Integer.parseInt(arg)));
        }
        else {
            for (int n = 0; n < 16; n++) {
                System.out.format("%d\t%s\n", n, count(n));
            }
        }
    }

    private static BigInteger count(int n) {
        if (n == 0) return BigInteger.ONE;

        if (minp.length < 3*n) {
            int[] wider = new int[3*n];
            System.arraycopy(minp, 0, wider, 0, minp.length);
            for (int x = minp.length; x < wider.length; x++) {
                // Find the smallest prime which divides x
                for (wider[x] = 2; x % wider[x] != 0; wider[x]++) { /* Do nothing */ }
            }
            minp = wider;
        }

        BigInteger E = countE(n), R2 = countR2(n), F = countF(n), R3 = countR3(n), R = countR(n), FR = countFR(n);
        BigInteger sum = E.add(R3);
        sum = sum.add(R2.add(R).multiply(BigInteger.valueOf(2)));
        sum = sum.add(F.add(FR).multiply(BigInteger.valueOf(3)));
        return sum.divide(BigInteger.valueOf(12));
    }

    private static BigInteger countE(int n) {
        int[] w = new int[3*n];
        for (int i = 0; i < n; i++) {
            for (int j = i + 1; j <= i + n; j++) w[j]--;
            for (int j = i + n + 1; j <= i + 2*n; j++) w[j]++;
        }
        return powerProd(w);
    }

    private static BigInteger countR2(int n) {
        int[] w = new int[3*n];
        for (int i = 0; i < n; i++) {
            w[3*i+2]++;
            for (int j = 3*i + 1; j <= 2*i + n + 1; j++) w[j]--;
            for (int j = 2*i + n + 1; j <= i + n + n; j++) w[j]++;
        }
        return powerProd(w);
    }

    private static BigInteger countF(int n) {
        int[] w = new int[3*n];
        for (int i = 0; i < n; i++) {
            for (int j = 2*i + 1; j <= 2*i + n; j++) w[j]--;
            for (int j = i + n + 1; j <= i + 2*n; j++) w[j]++;
        }
        return powerProd(w);
    }

    private static BigInteger countR3(int n) {
        if ((n & 1) == 1) return BigInteger.ZERO;
        return countE(n / 2).pow(2);
    }

    private static BigInteger countR(int n) {
        if ((n & 1) == 1) return BigInteger.ZERO;
        int m = n / 2;
        int[] w = new int[3*m-1];
        for (int i = 0; i < m; i++) {
            for (int j = 1; j <= 3*i+1; j++) w[j] += 2;
            for (int j = 1; j <= i + m; j++) w[j] -= 2;
        }
        return powerProd(w);
    }

    private static BigInteger countFR(int n) {
        if ((n & 1) == 1) return BigInteger.ZERO;
        int m = n / 2;
        int[] w = new int[3*n-2];
        for (int j = 1; j <= m; j++) w[j]--;
        for (int j = 2*m; j <= 3*m-1; j++) w[j]++;
        for (int i = 0; i <= 2*m-3; i++) {
            for (int j = i + 2*m + 1; j <= i + 4*m; j++) w[j]++;
            for (int j = 2*i + 3; j <= 2*i + 2*m + 2; j++) w[j]--;
        }
        return powerProd(w);
    }

    private static BigInteger powerProd(int[] w) {
        BigInteger result = BigInteger.ONE;
        for (int x = w.length - 1; x > 1; x--) {
            if (w[x] == 0) continue;

            int p = minp[x];
            if (p == x) result = result.multiply(BigInteger.valueOf(p).pow(w[p]));
            else {
                // Redistribute it. This should ensure we avoid negatives.
                w[p] += w[x];
                w[x / p] += w[x];
            }
        }

        return result;
    }
}

C

introduction

Comme l'a commenté David Carraher, le moyen le plus simple d'analyser le pavage hexagonal semblait être de tirer parti de son isomorphisme avec le diagramme de Young tridimensionnel, essentiellement un carré x, y rempli de barres de hauteur entière dont la hauteur z doit rester la même ou augmenter. lorsque l’axe z est approché.

J'ai commencé avec un algorithme permettant de trouver les totaux, plus facile à adapter pour le comptage de symétrie que l'algorithme publié, basé sur un biais par rapport à l'un des trois axes cartésiens.

Algorithme

Je commence par remplir les cellules des plans x, y et z avec des 1, tandis que le reste de la zone contient des zéros. Une fois cela fait, je construis le motif couche par couche, chaque couche contenant les cellules qui ont une distance de manhattan 3D commune à l'origine. Une cellule ne peut contenir un 1 que si les trois cellules en dessous contiennent également un 1. Si l'une d'entre elles contient un 0, la cellule doit être un 0.

L'avantage de construire le motif de cette manière est que chaque couche est symétrique autour de la ligne x = y = z. Cela signifie que chaque couche peut être vérifiée indépendamment pour sa symétrie.

Vérification de la symétrie

Les symétries du solide sont les suivantes: Rotation 3 fois autour de la ligne x = y = z -> Rotation 3 fois autour du centre de l'hexagone; et 3 réflexions x autour des 3 plans contenant la ligne x = y = z et chacun des axes x, y, z -> réflexion autour des lignes passant par les angles hexagonaux.

Cela ne fait qu'ajouter une symétrie 6 fois. Pour obtenir la symétrie complète de l'hexagone, un autre type de symétrie doit être envisagé. Chaque solide (construit à partir de 1) a un solide complémentaire (construit à partir de 0). Lorsque N est impair, le solide complémentaire doit être différent du solide d'origine (car il ne leur est pas possible d'avoir le même nombre de cubes). Pourtant, lorsque le solide complémentaire est retourné, on constate que sa représentation 2D sous forme de mosaïque de diamant est identique (à l’exception d’une opération à symétrie double) au solide initial. Lorsque N est pair, il est possible que le solide soit auto-inverse.

Cela se voit dans les exemples pour N = 2 dans la question. Vu de gauche, le premier hexagone ressemble à un cube solide avec 8 petits cubes, tandis que le dernier hexagone ressemble à une coque vide avec 0 petit cube. Si vu de la droite, l'inverse est vrai. Les 3ème, 4ème et 5ème hexagones et les 16ème, 17ème et 18ème hexagones semblent contenir 2 ou 6 cubes et se complètent donc en 3 dimensions. Ils sont reliés les uns aux autres en 2 dimensions par une opération de symétrie 2 (rotation 2 fois, ou réflexion autour d'un axe passant par les arêtes de l'hexagone). D'autre part, les 9ème, 10ème, 11ème et 12ème hexagones montrent des motifs 3D qui sont leurs propres compléments et ont donc une symétrie plus élevée (ce sont donc les seuls modèles avec une multiplicité impaire).

Notez qu'avoir (N ^ 3) / 2 cubes est une condition nécessaire pour être auto-complémentaire, mais en général, il ne suffit pas si N> 2. Le résultat de tout cela est que pour les N impairs, les pavages apparaissent toujours par paires (N ^ 3) / 2 cubes doivent être soigneusement inspectés.

Code actuel (génère le total correct pour N = 1,2,3,5. Erreur comme indiqué pour N = 4.)

int n;                     //side length

char t[11][11][11];        //grid sized for N up to 10

int q[29][192], r[29];     //tables of coordinates for up to 10*3-2=28 layers 

int c[9];                  //counts arrangements found by symmetry class. c[8] contains total.


//recursive layer counting function. m= manhattan distance, e= number of cells in previous layers, s=symmetry class.
void f(int m,int e,int s){

  int u[64], v[64], w[64]; //shortlists for x,y,z coordinates of cells in this layer
  int j=0;                 
  int x,y,z;

  for (int i=r[m]*3; i; i-=3){
    // get a set of coordinates for a cell in the current layer.
    x=q[m][i-3]; y= q[m][i-2]; z= q[m][i-1];
    // if the three cells in the previous layer are filled, add it to the shortlist u[],v[],w[]. j indicates the length of the shortlist.
    if (t[x][y][z-1] && t[x][y-1][z] && t[x-1][y][z]) u[j]=x, v[j]=y, w[j++]=z ;
  }


  // there are 1<<j possible arrangements for this layer.   
  for (int i = 1 << j; i--;) {

    int d = 0;

    // for each value of i, set the 1's bits of t[] to the 1's bits of i. Count the number of 1's into d as we go.
    for (int k = j; k--;) d+=(t[u[k]][v[k]][w[k]]=(i>>k)&1);

    // we have no interest in i=0 as it is the empty layer and therefore the same as the previous recursion step. 
    // Still we loop through it to ensure t[] is properly cleared.      

    if(i>0){
      int s1=s;    //local copy of symmetry class. 1's bit for 3 fold rotation, 2's bit for reflection in y axis.
      int sc=0;    //symmetry of self-complement.

      //if previous layers were symmetrical, test if the symmetry has been reduced by the current layer 
      if (s1) for (int k = j; k--;) s1 &= (t[u[k]][v[k]][w[k]]==t[w[k]][u[k]][v[k]]) | (t[u[k]][v[k]][w[k]]==t[w[k]][v[k]][u[k]])<<1;

      //if exactly half the cells are filled, test for self complement
      if ((e+d)*2==n*n*n){
        sc=1;
        for(int A=1; A<=(n>>1); A++)for(int B=1; B<=n; B++)for(int C=1; C<=n; C++) sc&=t[A][B][C]^t[n+1-A][n+1-B][n+1-C];
      }

      //increment counters for total and for symmetry class.
      c[8]++; c[s1+(sc<<2)]++;

      //uncomment for graphic display of each block stacking with metadata. not recommended for n>3.
      //printf("m=%d  j=%d  i=%d c1=%d-2*%d=%d c3=%d cy=%d(cs=%d) c3v=%d ctot=%d\n",m,j,i,c[0],c[2],c[0]-2*c[2],c[1],c[2],c[2]*3,c[3],c[8]);
      //printf("m=%d  j=%d  i=%d C1=%d-2*%d=%d C3=%d CY=%d(CS=%d) C3V=%d ctot=%d\n",m,j,i,c[4],c[6],c[4]-2*c[6],c[5],c[6],c[6]*3,c[7],c[8]);
      //for (int A = 0; A<4; A++, puts(""))for (int B = 0; B<4; B++, printf(" "))for (int C = 0; C<4; C++) printf("%c",34+t[A][B][C]);

      //recurse to next level.
      if(m<n*3-2)f(m + 1,e+d,s1);

    }
  } 
}

main()
{
  scanf("%d",&n);

  int x,y,z;

  // Fill x,y and z planes of t[] with 1's
  for (int a=0; a<9; a++) for (int b=0; b<9; b++) t[a][b][0]= t[0][a][b]= t[b][0][a]= 1;

  // Build table of coordinates for each manhattan layer
  for (int m=1; m < n*3-1; m++){
    printf("m=%d : ",m);
    int j=0;
    for (x = 1; x <= n; x++) for (y = 1; y <= n; y++) {
      z=m+2-x-y;
      if (z>0 && z <= n) q[m][j++] = x, q[m][j++] = y, q[m][j++]=z, printf(" %d%d%d ",x,y,z);
      r[m]=j/3;
    }
    printf(" : r=%d\n",r[m]);
  }

  // Set count to 1 representing the empty box (symmetry c3v)
  c[8]=1; c[3]=1; 

  // Start searching at f=1, with 0 cells occupied and symmetry 3=c3v
  f(1,0,3); 

  // c[2 and 6] only contain reflections in y axis, therefore must be multiplied by 3.
  // Similarly the reflections in x and z axis must be subtracted from c[0] and c[4].
  c[0]-=c[2]*2; c[2]*=3; 
  c[4]-=c[6]*2; c[6]*=3;



  int cr[9];cr[8]=0;
  printf("non self-complement                   self-complement\n");
  printf("c1  %9d/12=%9d           C1  %9d/6=%9d\n",   c[0], cr[0]=c[0]/12,     c[4], cr[4]=c[4]/6);
  if(cr[0]*12!=c[0])puts("c1 division error");if(cr[4]*6!=c[4])puts("C1 division error");

  printf("c3  %9d/4 =%9d           C3  %9d/2=%9d\n",   c[1], cr[1]=c[1]/4,      c[5], cr[5]=c[5]/2);
  if(cr[1]*4!=c[1])puts("c3 division error");if(cr[5]*2!=c[5])puts("C3 division error");

  printf("cs  %9d/6 =%9d           CS  %9d/3=%9d\n",   c[2], cr[2]=c[2]/6,      c[6], cr[6]=c[6]/3);
  if(cr[2]*6!=c[2])puts("cs division error");if(cr[6]*3!=c[6])puts("CS division error");

  printf("c3v %9d/2 =%9d           C3V %9d/1=%9d\n",   c[3], cr[3]=c[3]/2,      c[7], cr[7]=c[7]);
  if(cr[3]*2!=c[3])puts("c3v division error");  

  for(int i=8;i--;)cr[8]+=cr[i]; 
  printf("total =%d unique =%d",c[8],cr[8]);    
}

Sortie

Le programme génère une table de sortie de 8 entrées, conformément aux 8 symétries du solide. Le solide peut avoir l’une des 4 symétries suivantes (notation Schoenflies)

c1: no symmetry
c3: 3-fold axis of rotation (produces 3-fold axis of rotation in hexagon tiling)
cs: plane of reflection (produces line of reflection in hexagon tiling)
c3v both of the above (produces 3-fold axis of rotation and three lines of reflection through the hexagon corners)

De plus, lorsque le solide a exactement la moitié des cellules avec 1 et l'autre avec 0, il est possible d'inverser tous les 1 et les 0, puis d'inverser les coordonnées par le centre de l'espace du cube. C'est ce que j'appelle l'auto-complément, mais un terme plus mathématique serait "antisymétrique par rapport à un centre d'inversion".

Cette opération de symétrie donne un axe de rotation de 2 fois dans le pavage hexagonal.

Les motifs qui ont cette symétrie sont listés dans une colonne séparée. Ils ne se produisent que lorsque N est pair.

Mon compte semble être légèrement en retrait pour N = 4. En discutant avec Peter Taylor, il apparaît que je ne détecte pas de pavages présentant uniquement une symétrie d'une ligne sur les arêtes de l'hexagone. Ceci est probablement dû au fait que je n'ai pas testé l'auto-complément (antisymétrie) pour des opérations autres que (inversion) x (identité.). ) peut découvrir les symétries manquantes. Je m'attendrais alors à ce que la première ligne des données pour N = 4 ressemble à ceci (16 de moins en c1 et 32 ​​de plus en C1):

c1   224064/12=18672          C1  534/6=89

Cela alignerait les totaux sur la réponse de Peter et sur https://oeis.org/A066931/a066931.txt.

la sortie actuelle est la suivante.

N=1
non self-complement     self-complement
c1      0/12= 0           C1  0/6= 0
c3      0/4 = 0           C3  0/2= 0
cs      0/6 = 0           CS  0/3= 0
c3v     2/2 = 1           C3V 0/1= 0
total =2 unique =1

non self-complement     self-complement
N=2
c1      0/12= 0           C1  0/6= 0
c3      0/4 = 0           C3  0/2= 0
cs     12/6 = 2           CS  3/3= 1
c3v     4/2 = 2           C3V 1/1= 1
total =20 unique =6

N=3
non self-complement     self-complement
c1    672/12=56           C1  0/6= 0
c3      4/4 = 1           C3  0/2= 0
cs    288/6 =48           CS  0/3= 0
c3v    16/2 = 8           C3V 0/1= 0
total =980 unique =113

N=4 (errors as discussed)
non self-complement     self-complement
c1   224256/12=18688          C1  342/6=57
c3       64/4 =16             C3  2/2= 1
cs     8064/6 =1344           CS  54/3=18
c3v      64/2 =32             C3V 2/1= 2
total =232848 unique =20158

N=5
non self-complement     self-complement
c1  266774112/12=22231176        C1  0/6= 0
c3       1100/4 =275             C3  0/2= 0
cs     451968/6 =75328           CS  0/3= 0
c3v       352/2 =176             C3V 0/1= 0
total =267227532 unique =22306955

Liste de tâches (mise à jour)

Ranger le code actuel.

Fait plus ou moins

Implémentez la vérification de symétrie pour la couche en cours et transmettez un paramètre pour la symétrie de la couche précédente (cela ne sert à rien de vérifier si la dernière couche était asymétrique.)

Fait, les résultats pour N impair correspondent aux données publiées.

Ajoute une option pour supprimer le comptage des figures asymétriques (devrait courir beaucoup plus vite)

Cela peut être fait en ajoutant une autre condition à l'appel récursif: if(s1 && m<n*3-2)f(m + 1,e+d,s1)Cela réduit le temps d'exécution pour N = 5 de 5 minutes à environ une seconde. En conséquence, la première ligne de la sortie devient le total des déchets (tout comme les totaux globaux), mais si le total est déjà connu d’OEIS, le nombre de pavages asymétriques peut être reconstitué, au moins pour les N. impairs.

Mais même pour N, le nombre de solides asymétriques (selon les symétries c3v) qui s'auto-complétaient serait perdu. Dans ce cas, un programme séparé dédié aux solides avec exactement (N ** 3) / 2 cellules avec un 1 peut être utile. Avec cette option disponible (et en comptant correctement), il peut être possible d’essayer N = 6, mais l’exécution prendra beaucoup de temps.

Implémentez le comptage des cellules pour réduire la recherche à (N ^ 3) / 2 cubes.

Non fait, les économies devraient être marginales

Implémentez la vérification de symétrie (solide complémentaire) pour les motifs contenant exactement (N ^ 3) / 2 cubes.

Fait, mais semble avoir des omissions, voir N = 4.

Trouvez un moyen de choisir le chiffre lexical le plus bas parmi un chiffre asymétrique.

Les économies ne devraient pas être aussi bonnes. Supprimer les figures asymétriques élimine la plupart de cela. La seule réflexion vérifiée est le plan passant par l’axe des y (x et z sont calculés ultérieurement en les multipliant par 3.). Peut-être que cela courrait presque deux fois plus vite si on n'en comptait qu'un.

Pour faciliter ceci, améliorez éventuellement la façon dont les coordonnées dans chaque couche sont listées (elles forment des groupes dégénérés de 6 ou 3, avec éventuellement un groupe de 1 au centre exact de la couche.)

Intéressant, mais il y a probablement d'autres questions sur le site à explorer.