La plus longue séquence Game-of-Life non répétitive

16

Étant donné un entier N positif, déterminez le motif de départ sur une grille N x N qui donne la séquence non répétitive la plus longue selon les règles du jeu de la vie, et se termine par un motif fixe (cycle de longueur 1), joué sur un tore.

L'objectif n'est pas le programme le plus court, mais le plus rapide.

Puisque le monde est fini, vous finirez par finir dans une boucle, répétant ainsi un état déjà visité. Si cette boucle a la période 1, alors le modèle de départ est un candidat valide.

Sortie: modèle de départ et nombre total d'états uniques dans la séquence (y compris le modèle de départ).

Maintenant, le tore 1x1 est spécial, car une cellule peut être considérée comme voisine d'elle-même ou non, mais en pratique, il n'y a pas de problème, une seule cellule vivante mourra dans les deux cas (de surpopulation ou de solitude). Ainsi, l'entrée 1 donne une séquence de longueur 2, la séquence étant une cellule vivante, puis morte à jamais.

La motivation de cette question est qu'il s'agit d'un analogue de la fonction de castor occupé, mais certainement moins complexe, car nous avons une limite sur la mémoire. Ce sera une belle séquence à inclure également sur OEIS.

Pour N = 3, la longueur de la séquence est de 3, tout motif sur le côté gauche atteint un carré 3x3 complètement noir, puis meurt. (Tous les motifs faisant partie d'un cycle ont été supprimés).

Graphique des états

Per Alexandersson
la source
1
Ah, d'accord. Le meilleur code est celui qui parvient à calculer la longueur de séquence pour le plus grand N dans, disons 2 heures. La complexité évidente est 2 ^ (N ^ 2), donc s'il est possible d'améliorer cela, ce serait bien.
Per Alexandersson
1
À des tailles non triviales, chaque motif fait partie d'une classe d'isomorphisme de 8N ^ 2 motifs, donc s'il existe un moyen rapide de canoniser, cela donne un coup de pouce modéré.
Peter Taylor
1
Cette séquence a-t-elle été ajoutée à OEIS?
mbomb007
1
Pas que je sache, serait heureux de le voir là-bas.
Per Alexandersson
1
J'ai soumis cette séquence (2, 2, 3, 10, 52, 91) à l'OEIS sous le numéro A294241 .
Peter Kagey

Réponses:

13

C ++ jusqu'à N = 6

3x3 answer 3: 111 000 000                                                                                        
4x4 answer 10: 1110 0010 1100 0000                                                                               
5x5 answer 52: 11010 10000 11011 10100 00000                                                                     
6x6 answer 91: 100011 010100 110011 110100 101000 100000                                                         

Cette technique est centrée sur une fonction d'état suivant rapide. Chaque carte est représentée sous la forme d'un masque de bits de N ^ 2 bits, et des astuces bit-twiddly sont utilisées pour calculer l'état suivant pour toutes les cellules à la fois. nextcompile jusqu'à environ 200 100 instructions d'assemblage pour N <= 8. Ensuite, nous faisons simplement l'essai standard tous les états jusqu'à ce qu'ils se répètent et choisissons le plus long.

Prend quelques secondes jusqu'à 5x5, quelques heures pour 6x6.

#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <algorithm>
using namespace std;

#define N 6

typedef uint64_t board;

// board is N bits by N bits, with indexes like this (N=4):                                                        
//  0  1  2  3                                                                                                     
//  4  5  6  7                                                                                                     
//  8  9 10 11                                                                                                     
// 12 13 14 15                                                                                                     

#if N==3
#define LEFT_COL (1 + (1<<3) + (1<<6))
#define RIGHT_COL ((1<<2) + (1<<5) + (1<<8))
#define ALL 0x1ffULL
#elif N==4
#define LEFT_COL 0x1111ULL
#define RIGHT_COL 0x8888ULL
#define ALL 0xffffULL
#elif N==5
#define LEFT_COL (1ULL + (1<<5) + (1<<10) + (1<<15) + (1<<20))
#define RIGHT_COL ((1ULL<<4) + (1<<9) + (1<<14) + (1<<19) + (1<<24))
#define ALL 0x1ffffffULL
#elif N==6
#define LEFT_COL (1 + (1<<6) + (1<<12) + (1<<18) + (1<<24) + (1ULL<<30))
#define RIGHT_COL ((1<<5) + (1<<11) + (1<<17) + (1<<23) + (1<<29) + (1ULL<<35))
#define ALL 0xfffffffffULL
#else
#error "bad N"
#endif

static inline board north(board b) {
  return (b >> N) + (b << N*N-N);
}
static inline board south(board b) {
  return (b << N) + (b >> N*N-N);
}
static inline board west(board b) {
  return ((b & ~LEFT_COL) >> 1) + ((b & LEFT_COL) << N-1);
}
static inline board east(board b) {
  return ((b & ~RIGHT_COL) << 1) + ((b & RIGHT_COL) >> N-1);
}

board next(board b) {
  board n1 = north(b);
  board n2 = south(b);
  board n3 = west(b);
  board n4 = east(b);
  board n5 = north(n3);
  board n6 = north(n4);
  board n7 = south(n3);
  board n8 = south(n4);

  // add all the bits bitparallel-y to a 2-bit accumulator with overflow
  board a0 = 0;
  board a1 = 0;
  board overflow = 0;
#define ADD(x) overflow |= a0 & a1 & x; a1 ^= a0 & x; a0 ^= x;

  a0 = n1; // no carry yet
  a1 ^= a0 & n2; a0 ^= n2; // no overflow yet
  a1 ^= a0 & n3; a0 ^= n3; // no overflow yet
  ADD(n4);
  ADD(n5);
  ADD(n6);
  ADD(n7);
  ADD(n8);
  return (a1 & (b | a0)) & ~overflow & ALL;
}
void print(board b) {
  for (int i = 0; i < N; i++) {
    for (int j = 0; j < N; j++) {
      printf("%d", (int)(b >> i*N+j & 1));
    }
    printf(" ");
  }
  if (b >> N*N) printf("*");
  printf("\n");
}

int main(int argc, char *argv[]) {
  // Somewhere in the starting pattern there are a 1 and 0 together.  Using translational                          
  // and rotational symmetry, we can put these in the first two bits.  So we need only consider                    
  // 1 mod 4 boards.                                                                                               

  board steps[10000];
  int maxsteps = -1;
  for (board b = 1; b < 1ULL << N*N; b += 4) {
    int nsteps = 0;
    board x = b;
    while (true) {
      int repeat = steps + nsteps - find(steps, steps + nsteps, x);
      if (repeat > 0) {
        if (repeat == 1 && nsteps > maxsteps) {
          printf("%d: ", nsteps);
          print(b);
          maxsteps = nsteps;
        }
        break;
      }
      steps[nsteps++] = x;
      x = next(x);
    }
  }
}
Keith Randall
la source
1
Vous pourriez obtenir une amélioration modérée nexten comptant plutôt qu'en triant. #define H(x,y){x^=y;y&=(x^y);}puis quelque chose commeH(n1,n2);H(n3,n4);H(n5,n6);H(n7,n8);H(n1,n3);H(n5,n7);H(n2,n4);H(n6,n8);H(n1,n5);H(n3,n7);H(n2,n6);H(n2,n3);H(n2,n5); return n2 & (b | n1) & ~(n3|n4|n5|n6|n7|n8) & ALL;
Peter Taylor
Solution vraiment cool!
Per Alexandersson
@PeterTaylor: merci, j'ai implémenté (un schéma différent pour) le comptage, enregistré un tas d'instructions.
Keith Randall
9

Je vois les approches de solution possibles suivantes:

  • Théorie lourde. Je sais qu'il existe de la littérature sur la vie sur un tore, mais je n'en ai pas beaucoup lu.
  • Force brute en avant: pour chaque planche possible, vérifiez son score. C'est essentiellement ce que font les approches de Matthew et Keith, bien que Keith réduit le nombre de planches à vérifier par un facteur de 4.
    • Optimisation: représentation canonique. Si nous pouvons vérifier si une planche est en représentation canonique beaucoup plus rapidement qu'il n'en faut pour évaluer son score, nous obtenons une accélération d'un facteur d'environ 8N ^ 2. (Il existe également des approches partielles avec des classes d'équivalence plus petites).
    • Optimisation: DP. Mettez en cache le score de chaque tableau, de sorte qu'au lieu de les parcourir jusqu'à ce qu'ils convergent ou divergent, nous marchons jusqu'à ce que nous trouvions un tableau que nous avons vu auparavant. En principe, cela donnerait une accélération d'un facteur de la durée moyenne du score / cycle (peut-être 20 ou plus), mais dans la pratique, nous allons probablement échanger fortement. Par exemple, pour N = 6, nous aurions besoin d'une capacité de 2 ^ 36 scores, qui à un octet par score est de 16 Go, et nous avons besoin d'un accès aléatoire, donc nous ne pouvons pas nous attendre à une bonne localité de cache.
    • Combinez les deux. Pour N = 6, la représentation canonique complète nous permettrait de réduire le cache DP à environ 60 méga-scores. Il s'agit d'une approche prometteuse.
  • Force brute en arrière. Cela semble étrange au début, mais si nous supposons que nous pouvons facilement trouver des natures mortes et que nous pouvons facilement inverser la Next(board)fonction, nous constatons qu'elle présente certains avantages, bien que pas autant que je l'espérais à l'origine.
    • Nous ne nous soucions pas du tout des planches divergentes. Pas beaucoup d'économies, car elles sont assez rares.
    • Nous n'avons pas besoin de stocker les scores pour toutes les cartes, il devrait donc y avoir moins de pression sur la mémoire que l'approche DP vers l'avant.
    • Le travail en arrière est en fait assez facile en faisant varier une technique que j'ai vue dans la littérature dans le contexte de l'énumération des natures mortes. Il fonctionne en traitant chaque colonne comme une lettre dans un alphabet, puis en observant qu'une séquence de trois lettres détermine celle du milieu dans la prochaine génération. Le parallèle avec l' énumération des natures mortes est si proche que je l' ai refactorisé - les ensemble dans une méthode peu maladroite, Prev2.
    • Il peut sembler que nous pouvons simplement canoniser les natures mortes et obtenir quelque chose comme l'accélération 8N ^ 2 pour très peu de frais. Cependant, empiriquement, nous obtenons toujours une grande réduction du nombre de cartes prises en compte si nous canonisons à chaque étape.
    • Une proportion étonnamment élevée de planches ont un score de 2 ou 3, il y a donc toujours une pression mémoire. J'ai trouvé nécessaire de canoniser à la volée plutôt que de construire la génération précédente, puis de canoniser. Il peut être intéressant de réduire l'utilisation de la mémoire en effectuant une recherche en profondeur d'abord plutôt qu'en largeur, mais le faire sans déborder la pile et sans effectuer de calculs redondants nécessite une approche de co-routine / continuation pour énumérer les cartes précédentes.

Je ne pense pas que la micro-optimisation me permettra de rattraper le code de Keith, mais pour l'intérêt, je posterai ce que j'ai. Cela résout N = 5 en environ une minute sur une machine à 2 GHz en utilisant Mono 2.4 ou .Net (sans PLINQ) et en environ 20 secondes en utilisant PLINQ; N = 6 s'exécute pendant plusieurs heures.

using System;
using System.Collections.Generic;
using System.Linq;

namespace Sandbox {
    class Codegolf9393 {
        internal static void Main() {
            new Codegolf9393(4).Solve();
        }

        private readonly int _Size;
        private readonly uint _AlphabetSize;
        private readonly uint[] _Transitions;
        private readonly uint[][] _PrevData1;
        private readonly uint[][] _PrevData2;
        private readonly uint[,,] _CanonicalData;

        private Codegolf9393(int size) {
            if (size > 8) throw new NotImplementedException("We need to fit the bits in a ulong");

            _Size = size;
            _AlphabetSize = 1u << _Size;

            _Transitions = new uint[_AlphabetSize * _AlphabetSize * _AlphabetSize];
            _PrevData1 = new uint[_AlphabetSize * _AlphabetSize][];
            _PrevData2 = new uint[_AlphabetSize * _AlphabetSize * _AlphabetSize][];
            _CanonicalData = new uint[_Size, 2, _AlphabetSize];

            InitTransitions();
        }

        private void InitTransitions() {
            HashSet<uint>[] tmpPrev1 = new HashSet<uint>[_AlphabetSize * _AlphabetSize];
            HashSet<uint>[] tmpPrev2 = new HashSet<uint>[_AlphabetSize * _AlphabetSize * _AlphabetSize];
            for (int i = 0; i < tmpPrev1.Length; i++) tmpPrev1[i] = new HashSet<uint>();
            for (int i = 0; i < tmpPrev2.Length; i++) tmpPrev2[i] = new HashSet<uint>();

            for (uint i = 0; i < _AlphabetSize; i++) {
                for (uint j = 0; j < _AlphabetSize; j++) {
                    uint prefix = Pack(i, j);
                    for (uint k = 0; k < _AlphabetSize; k++) {
                        // Build table for forwards checking
                        uint jprime = 0;
                        for (int l = 0; l < _Size; l++) {
                            uint count = GetBit(i, l-1) + GetBit(i, l) + GetBit(i, l+1) + GetBit(j, l-1) + GetBit(j, l+1) + GetBit(k, l-1) + GetBit(k, l) + GetBit(k, l+1);
                            uint alive = GetBit(j, l);
                            jprime = SetBit(jprime, l, (count == 3 || (alive + count == 3)) ? 1u : 0u);
                        }
                        _Transitions[Pack(prefix, k)] = jprime;

                        // Build tables for backwards possibilities
                        tmpPrev1[Pack(jprime, j)].Add(k);
                        tmpPrev2[Pack(jprime, i, j)].Add(k);
                    }
                }
            }

            for (int i = 0; i < tmpPrev1.Length; i++) _PrevData1[i] = tmpPrev1[i].ToArray();
            for (int i = 0; i < tmpPrev2.Length; i++) _PrevData2[i] = tmpPrev2[i].ToArray();

            for (uint col = 0; col < _AlphabetSize; col++) {
                _CanonicalData[0, 0, col] = col;
                _CanonicalData[0, 1, col] = VFlip(col);
                for (int rot = 1; rot < _Size; rot++) {
                    _CanonicalData[rot, 0, col] = VRotate(_CanonicalData[rot - 1, 0, col]);
                    _CanonicalData[rot, 1, col] = VRotate(_CanonicalData[rot - 1, 1, col]);
                }
            }
        }

        private ICollection<ulong> Prev2(bool stillLife, ulong next, ulong prev, int idx, ICollection<ulong> accum) {
            if (stillLife) next = prev;

            if (idx == 0) {
                for (uint a = 0; a < _AlphabetSize; a++) Prev2(stillLife, next, SetColumn(0, idx, a), idx + 1, accum);
            }
            else if (idx < _Size) {
                uint i = GetColumn(prev, idx - 2), j = GetColumn(prev, idx - 1);
                uint jprime = GetColumn(next, idx - 1);
                uint[] succ = idx == 1 ? _PrevData1[Pack(jprime, j)] : _PrevData2[Pack(jprime, i, j)];
                foreach (uint b in succ) Prev2(stillLife, next, SetColumn(prev, idx, b), idx + 1, accum);
            }
            else {
                // Final checks: does the loop round work?
                uint a0 = GetColumn(prev, 0), a1 = GetColumn(prev, 1);
                uint am = GetColumn(prev, _Size - 2), an = GetColumn(prev, _Size - 1);
                if (_Transitions[Pack(am, an, a0)] == GetColumn(next, _Size - 1) &&
                    _Transitions[Pack(an, a0, a1)] == GetColumn(next, 0)) {
                    accum.Add(Canonicalise(prev));
                }
            }

            return accum;
        }

        internal void Solve() {
            DateTime start = DateTime.UtcNow;
            ICollection<ulong> gen = Prev2(true, 0, 0, 0, new HashSet<ulong>());
            for (int depth = 1; gen.Count > 0; depth++) {
                Console.WriteLine("Length {0}: {1}", depth, gen.Count);
                ICollection<ulong> nextGen;

                #if NET_40
                nextGen = new HashSet<ulong>(gen.AsParallel().SelectMany(board => Prev2(false, board, 0, 0, new HashSet<ulong>())));
                #else
                nextGen = new HashSet<ulong>();
                foreach (ulong board in gen) Prev2(false, board, 0, 0, nextGen);
                #endif

                // We don't want the still lifes to persist or we'll loop for ever
                if (depth == 1) {
                    foreach (ulong stilllife in gen) nextGen.Remove(stilllife);
                }

                gen = nextGen;
            }
            Console.WriteLine("Time taken: {0}", DateTime.UtcNow - start);
        }

        private ulong Canonicalise(ulong board)
        {
            // Find the minimum board under rotation and reflection using something akin to radix sort.
            Isomorphism canonical = new Isomorphism(0, 1, 0, 1);
            for (int xoff = 0; xoff < _Size; xoff++) {
                for (int yoff = 0; yoff < _Size; yoff++) {
                    for (int xdir = -1; xdir <= 1; xdir += 2) {
                        for (int ydir = 0; ydir <= 1; ydir++) {
                            Isomorphism candidate = new Isomorphism(xoff, xdir, yoff, ydir);

                            for (int col = 0; col < _Size; col++) {
                                uint a = canonical.Column(this, board, col);
                                uint b = candidate.Column(this, board, col);

                                if (b < a) canonical = candidate;
                                if (a != b) break;
                            }
                        }
                    }
                }
            }

            ulong canonicalValue = 0;
            for (int i = 0; i < _Size; i++) canonicalValue = SetColumn(canonicalValue, i, canonical.Column(this, board, i));
            return canonicalValue;
        }

        struct Isomorphism {
            int xoff, xdir, yoff, ydir;

            internal Isomorphism(int xoff, int xdir, int yoff, int ydir) {
                this.xoff = xoff;
                this.xdir = xdir;
                this.yoff = yoff;
                this.ydir = ydir;
            }

            internal uint Column(Codegolf9393 _this, ulong board, int col) {
                uint basic = _this.GetColumn(board, xoff + col * xdir);
                return _this._CanonicalData[yoff, ydir, basic];
            }
        }

        private uint VRotate(uint col) {
            return ((col << 1) | (col >> (_Size - 1))) & (_AlphabetSize - 1);
        }

        private uint VFlip(uint col) {
            uint replacement = 0;
            for (int row = 0; row < _Size; row++)
                replacement = SetBit(replacement, row, GetBit(col, _Size - row - 1));
            return replacement;
        }

        private uint GetBit(uint n, int bit) {
            bit %= _Size;
            if (bit < 0) bit += _Size;

            return (n >> bit) & 1;
        }

        private uint SetBit(uint n, int bit, uint value) {
            bit %= _Size;
            if (bit < 0) bit += _Size;

            uint mask = 1u << bit;
            return (n & ~mask) | (value == 0 ? 0 : mask);
        }

        private uint Pack(uint a, uint b) { return (a << _Size) | b; }
        private uint Pack(uint a, uint b, uint c) {
            return (((a << _Size) | b) << _Size) | c;
        }

        private uint GetColumn(ulong n, int col) {
            col %= _Size;
            if (col < 0) col += _Size;
            return (_AlphabetSize - 1) & (uint)(n >> (col * _Size));
        }

        private ulong SetColumn(ulong n, int col, uint value) {
            col %= _Size;
            if (col < 0) col += _Size;

            ulong mask = (_AlphabetSize - 1) << (col * _Size);
            return (n & ~mask) | (((ulong)value) << (col * _Size));
        }
    }
}
Peter Taylor
la source
Je travaille également sur une autre version pour reculer à partir de points fixes. J'ai déjà énuméré les points fixes jusqu'à N = 8 (pour N = 8 il y en a 84396613 avant la canonisation). J'ai un simple aperçu qui fonctionne, mais il est trop lent. Une partie du problème vient de la taille des choses, pour N = 6, le plateau vide a 574384901 prédécesseurs (avant la canonisation).
Keith Randall
1
3 jours et 11 heures pour confirmer que 91 est la réponse pour 6x6.
Peter Taylor
1

Facteur

USING: arrays grouping kernel locals math math.functions math.parser math.order math.ranges math.vectors sequences sequences.extras ;
IN: longest-gof-pattern

:: neighbors ( x y game -- neighbors )
game length :> len 
x y game -rot 2array {
    { -1 -1 }
    { -1 0 }
    { -1 1 }
    { 0 -1 }
    { 0 1 }
    { 1 -1 }
    { 1 0 }
    { 1 1 }
} [
    v+ [
        dup 0 <
        [ dup abs len mod - abs len mod ] [ abs len mod ]
        if
    ] map
] with map [ swap [ first2 ] dip nth nth ] with map ;

: next ( game -- next )
dup [
    [
        neighbors sum
        [ [ 1 = ] [ 2 3 between? ] bi* and ]
        [ [ 0 = ] [ 3 = ] bi* and ] 2bi or 1 0 ?
    ] curry curry map-index
] curry map-index ;

: suffixes ( seq -- suffixes )
{ }
[ [ [ suffix ] curry map ] [ 1array 1array ] bi append ]
reduce ;

! find largest repeating pattern
: LRP ( seq -- pattern )
dup length iota
[ 1 + [ reverse ] dip group [ reverse ] map reverse ] with
map dup [ dup last [ = ] curry map ] map
[ suffixes [ t [ and ] reduce ] map [ ] count ] map
dup supremum [ = ] curry find drop swap nth last ;

: game-sequence ( game -- seq )
1array [
    dup [
        dup length 2 >
        [ 2 tail-slice* [ first ] [ last ] bi = not ]
        [ drop t ] if
    ] [ LRP length 1 > not ] bi and
] [ dup last next suffix ] while ;

: pad-to-with ( str len padstr -- rstr )
[ swap dup length swapd - ] dip [ ] curry replicate ""
[ append ] reduce prepend ;

:: all-NxN-games ( n -- games )
2 n sq ^ iota [
    >bin n sq "0" pad-to-with n group
    [ [ 48 = 0 1 ? ] { } map-as ] map
] map ;

: longest-gof-pattern ( n -- game )
all-NxN-games [ game-sequence ] map [ length ] supremum-by but-last ;

Quelques statistiques de temps:

IN: longest-gof-pattern [ 3 longest-gof-pattern ] time dup length . . 
Running time: 0.08850873500000001 seconds

3
{
   { { 1 1 1 } { 0 0 0 } { 0 0 0 } }
   { { 1 1 1 } { 1 1 1 } { 1 1 1 } }
   { { 0 0 0 } { 0 0 0 } { 0 0 0 } }
}

IN: longest-gof-pattern [ 4 longest-gof-pattern ] time dup length . . 
Running time: 49.667698828 seconds

10
{
  { { 0 1 1 0 } { 0 1 0 0 } { 0 1 0 0 } { 1 1 0 1 } }
  { { 0 1 1 0 } { 0 1 0 0 } { 0 1 0 0 } { 0 0 0 1 } }
  { { 0 1 1 0 } { 0 1 0 0 } { 0 0 1 0 } { 1 1 0 0 } }
  { { 0 1 1 0 } { 0 1 0 0 } { 0 0 1 0 } { 0 0 0 1 } }
  { { 0 1 1 0 } { 0 1 0 0 } { 0 0 1 0 } { 1 1 0 1 } }
  { { 0 1 1 0 } { 0 1 0 0 } { 0 0 1 1 } { 0 0 0 1 } }
  { { 0 1 0 1 } { 0 1 0 1 } { 0 0 1 1 } { 1 1 0 1 } }
  { { 1 1 0 1 } { 1 1 0 1 } { 0 0 0 0 } { 1 1 0 0 } }
  { { 1 1 0 1 } { 1 1 0 1 } { 0 0 1 1 } { 1 1 1 1 } }
  { { 0 0 0 0 } { 0 0 0 0 } { 0 0 0 0 } { 0 0 0 0 } }
}

Et le test 5 a fait planter le REPL. Hmph. La partie la plus inefficace du programme est probablement la fonction jeu-séquence. Je pourrais peut-être l'améliorer plus tard.

Matthew Rolph
la source
Cool! Je pense que votre sortie est 1 trop grande, pour les motifs 3x3, la séquence non répétitive la plus longue a une longueur 3, pas 4 ...
Par Alexandersson