Sequência de jogo da vida mais longa e não repetida

16

Dado um número inteiro positivo N, determine o padrão inicial em uma grade N x N que produza a sequência não repetida mais longa sob as regras do Jogo da Vida e termina com um padrão fixo (ciclo de comprimento 1), tocado em um toro.

O objetivo não é o programa mais curto, mas o mais rápido.

Como o mundo é finito, você terminará em um ciclo, repetindo assim um estado já visitado. Se esse loop tiver o período 1, o padrão inicial será um candidato válido.

Saída: padrão inicial e número total de estados exclusivos na sequência (incluindo padrão inicial).

Agora, o toro 1x1 é especial, uma vez que uma célula pode ser considerada vizinha ou não, mas, na prática, não há problema, uma única célula viva em ambos os casos simplesmente morrerá (de superlotação ou solidão). Assim, a entrada 1 produz uma sequência de comprimento 2, sendo a sequência uma célula viva e depois morta para sempre.

A motivação para esta pergunta é que é um análogo da função de castor ocupado, mas definitivamente menos complexo, pois temos um limite na memória. Essa também será uma boa sequência para incluir no OEIS.

Para N = 3, o comprimento da sequência é 3, qualquer padrão no lado esquerdo atinge um quadrado 3x3 completamente preto e depois morre. (Todos os padrões que fazem parte de 1 ciclo são removidos).

Gráfico de estados

Per Alexandersson
fonte
11
Ah, tudo bem. O melhor código é aquele que consegue calcular o comprimento da sequência para o maior N dentro de, digamos 2 horas. A complexidade óbvia é 2 ^ (N ^ 2); portanto, se é possível melhorar isso, isso seria bom.
Por Alexandersson
11
Em tamanhos não triviais, cada padrão faz parte de uma classe de isomorfismo de padrões 8N ^ 2; portanto, se houver uma maneira rápida de canonizar, isso dará um impulso moderado.
Peter Taylor
11
Esta sequência foi adicionada ao OEIS?
mbomb007
11
Não que eu esteja ciente, ficaria feliz em vê-lo lá.
Por Alexandersson
11
Enviei esta sequência (2, 2, 3, 10, 52, 91) ao OEIS como A294241 .
Peter Kagey

Respostas:

13

C ++ até 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                                                         

Essa técnica é centrada em torno de uma função rápida do próximo estado. Cada placa é representada como uma máscara de bits de N ^ 2 bits, e truques de bits pequenos são usados ​​para calcular o próximo estado de todas as células de uma só vez. nextcompila até 200 100 instruções de montagem para N <= 8. Em seguida, fazemos apenas o teste padrão de todos os estados até que eles se repitam e escolhemos o mais longo.

Leva alguns segundos até 5x5, algumas horas para 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
fonte
11
Você pode obter uma melhoria moderada nextcontando e não classificando. #define H(x,y){x^=y;y&=(x^y);}e então algo parecido comH(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
Solução muito legal!
Por Alexandersson
@ PeterTaylor: obrigado, eu implementei (um esquema diferente para) a contagem, salvei um monte de instruções.
Keith Randall
9

Vejo as seguintes abordagens de solução possíveis:

  • Teoria pesada. Sei que há alguma literatura sobre a vida em um toro, mas ainda não li muito.
  • Força bruta para a frente: para todas as pranchas possíveis, verifique sua pontuação. Isso é basicamente o que as abordagens de Matthew e Keith fazem, embora Keith reduza o número de placas a serem verificadas por um fator de 4.
    • Otimização: representação canônica. Se pudermos verificar se um quadro está na representação canônica muito mais rápido do que é necessário para avaliar sua pontuação, obtemos uma aceleração de um fator de cerca de 8N ^ 2. (Também existem abordagens parciais com classes de equivalência menores).
    • Otimização: DP. Armazene em cache a pontuação de cada quadro, para que, em vez de percorrê-los até convergirem ou divergir, apenas andemos até encontrar um quadro que já vimos antes. Em princípio, isso aceleraria um fator da pontuação média / duração do ciclo (talvez 20 ou mais), mas, na prática, é provável que estejamos trocando bastante. Por exemplo, para N = 6, precisaríamos de capacidade para 2 ^ 36 pontuações, que em um byte por pontuação são 16 GB e precisamos de acesso aleatório, para que não possamos esperar uma boa localidade do cache.
    • Combine os dois. Para N = 6, a representação canônica completa nos permitiria reduzir o cache de DP para cerca de 60 mega-pontuações. Essa é uma abordagem promissora.
  • Força bruta para trás. Isso parece estranho no começo, mas se assumirmos que podemos encontrar facilmente naturezas-mortas e que podemos inverter facilmente a Next(board)função, vemos que ela tem alguns benefícios, embora não tantos quanto eu originalmente esperava.
    • Não nos preocupamos com placas divergentes. Não poupa muito, porque são bastante raros.
    • Não precisamos armazenar pontuações para todas as placas, portanto, deve haver menos pressão de memória do que a abordagem de DP direta.
    • Trabalhar para trás é realmente muito fácil, variando uma técnica que vi na literatura no contexto de enumerar naturezas-mortas. Ele funciona tratando cada coluna como uma letra em um alfabeto e depois observando que uma sequência de três letras determina a do meio na próxima geração. O paralelo com a enumeração de naturezas-mortas é tão próximo que eu as refatorei juntas em um método apenas um pouco estranho Prev2.
    • Pode parecer que podemos apenas canonizar as naturezas-mortas e obter algo como a aceleração de 8N ^ 2 por um custo muito baixo. No entanto, empiricamente, ainda temos uma grande redução no número de placas consideradas se canonizarmos a cada etapa.
    • Uma proporção surpreendentemente alta de placas tem uma pontuação de 2 ou 3, então ainda há pressão de memória. Eu achei necessário canonizar em tempo real, em vez de construir a geração anterior e depois canonizar. Pode ser interessante reduzir o uso de memória fazendo uma pesquisa de profundidade em vez de largura, mas sem sobrecarregar a pilha e sem fazer cálculos redundantes requer uma abordagem de rotina / continuação para enumerar as placas anteriores.

Eu não acho que a micro-otimização me permita acompanhar o código de Keith, mas por uma questão de interesse, postarei o que tenho. Isso resolve N = 5 em cerca de um minuto em uma máquina de 2 GHz usando Mono 2.4 ou .Net (sem PLINQ) e em cerca de 20 segundos usando PLINQ; N = 6 é executado por muitas horas.

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
fonte
Também estou trabalhando em outra versão para retroceder a partir de pontos fixos. Eu já enumerei os pontos fixos até N = 8 (para N = 8 existem 84396613 antes da canonização). Eu tenho um trabalho simples anterior, mas é muito lento. Parte do problema é apenas o tamanho das coisas, pois N = 6 o tabuleiro vazio tem 574384901 predecessores (antes da canonização).
Keith Randall
11
3 dias e 11 horas para confirmar que 91 é a resposta para 6x6.
31513 Peter Peter Taylor
1

Fator

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 ;

Algumas estatísticas de tempo:

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

E o teste 5 travou o REPL. Hmph. A parte mais ineficiente do programa é provavelmente a função sequência do jogo. Talvez eu consiga melhorar depois.

Matthew Rolph
fonte
Legal! Eu acho que sua saída é um muito grande, por 3x3-padrões, a mais longa sequência não-repetição tem comprimento 3, não 4 ...
Per Alexandersson