I'm trying to search for still lifes in Conway's Game of Life.
An object in Conway's Game of Life can be represented by an array of booleans, say, Array[b, {w, h}]. A still life is a stable object, so it needs to satisfy the following conditions:
- For each living cell, the number of its living neighbors must be 2 or 3,
- For each dead cell, the number of its living neighbors must NOT be 3.
Written in Mathematica code, these conditions are:
NeighborCount[k_, {i_, j_}] :=
BooleanCountingFunction[{k},
Delete[Catenate[Array[b, {3, 3}, {i, j} - 1]], 5]];
StillLifeCondition[i_, j_] :=
(b[i, j] && NeighborCount[{2, 3}, {i, j}]) ||
(! b[i, j] && ! NeighborCount[{3}, {i, j}]);
Moreover, all cells outside the boundary are dead. So the condition for the whole object can be written as:
Array[StillLifeCondition[##] /.
b[i_, j_] /; i < 1 || i > w || j < 1 || j > h :> False &,
{w, h} + 2, 0, And]
Now this is a Boolean satisfiability problem, so we can apply SatisfiabilityInstances:
SearchStillLife[w_, h_] :=
ArrayReshape[
SatisfiabilityInstances[
Array[StillLifeCondition[##] /.
b[i_, j_] /; i < 1 || i > w || j < 1 || j > h :> False &,
{w, h} + 2, 0, And],
Catenate[Array[b, {w, h}]]][[1]], {w, h}];
But it is slow. Searching for a still life of size 10x10 takes 6 seconds:
ArrayPlot[Boole@SearchStillLife[10, 10], Mesh -> All] // AbsoluteTiming
We can specify the return form of BooleanCountingFunction. I tested all the possible forms. Among them, "BFF" is the fastest. So we can rewrite the NeighborCount function:
NeighborCount[k_, {i_, j_}] :=
BooleanCountingFunction[{k},
Delete[Catenate[Array[b, {3, 3}, {i, j} - 1]], 5], "BFF"];
Now searching for a still life of size 10x10 takes ~0.5 seconds.
But it is still very slow for large sizes. Searching for a still life of size 30x30 takes 1246 seconds:
ArrayPlot[Boole@SearchStillLife[30, 30], Mesh -> All] // AbsoluteTiming
How can I speed up the code?
Another problem is that there is no randomness in the code, so the result is completely determined by the input. The simplest way to solve it is to Xor it with a random array:
SearchRandomStillLife[w_, h_] :=
Block[{r = RandomChoice[{True, False}, {w, h}]},
MapThread[Xor, {r,
ArrayReshape[
SatisfiabilityInstances[
Array[StillLifeCondition[##] /.
b[i_, j_] /; i < 1 || i > w || j < 1 || j > h :> False /.
b[i_, j_] :> Xor[b[i, j], r[[i, j]]] &,
{w, h} + 2, 0, And],
Catenate[Array[b, {w, h}]]][[1]], {w, h}]}, 2]]
Now the timing depends on this random array. Usually it is slower. The following takes 1.5 seconds.
SeedRandom[233];
ArrayPlot[Boole@SearchRandomStillLife[10, 10], Mesh -> All] // AbsoluteTiming
Is there a better way to generate random result?
Another algorithm:
Here is the code I wrote in the first version of this question. It is extremely slow, but works for oscillators, and generates random results:
(* An association that, given the first two rows of an 3x3 grid, and \
the next generation of the central cell, return a list of all \
possible third rows. *)
a = Merge[
k[#[[;; 2]], CellularAutomaton["GameOfLife", #][[2, 2]]] -> #[[
3]] & /@ Tuples[Tuples[{0, 1}, 3], 3], Identity];
(* Given two rows of cells and the next generation of the second row, \
return a list of all possible third rows. *)
NextRow[x_, y_, z_] :=
FindPath[Flatten@
MapIndexed[Outer[v[#2 - 1, Most@#] -> v[#2, Rest@#1] &, ##, 1] &,
Lookup[a,
MapThread[
k, {First[
Partition[{x, y}, {2, 3}, {2, 1}, {{1, -2}, {-1, 2}}, 0]],
z}], {}]], v[0, {0, 0}], v[Length@y, {0, 0}], Length@y + 2,
All][[;; , 2 ;;, 2, 1]];
(* Search for an oscillator with width w, height h, and period p. *)
SearchOscillators[w_, h_, p_] := Module[{o = Table[0, w + 2], i, f},
i = NextRow[o, o, o];
f[l : {___, x_, y_}] :=
If[Length@l == h + 2,
If[SameQ[##, o] & @@ y && SubsetQ[i, x], l, Missing["NotFound"]],
Catch[Do[
If[MissingQ[#], , Throw[#]] &[f[Append[l, z]]], {z,
Tuples[RandomSample /@
MapThread[NextRow, {x, y, RotateLeft[y]}]]}];
Missing["NotFound"]]];
ArrayPlot[#, Mesh -> All] & /@
Transpose[
NestWhile[f[{Table[o, p], RandomChoice[i, p]}] &,
Missing["NotFound"], MissingQ]]]
