Counting types of neighbour - in Python.

(I tried using the Spyder IDE to write a more literal Python program. Output samples needed to be pasted in to the source, but I didn't particularly want many of Jupyters strengths, so tried Spyder together with http://hilite.me/)

Best viewed on larger than 7" screens.

# -*- coding: utf-8 -*-
"""
Created on Sat Jan 30 10:24:54 2021

@author: Paddy3118
"""

#from conway_cubes_2_plus_dims import *
#from conway_cubes_2_plus_dims import neighbours
from itertools import product


# %% #Full neighbourhood
"""
In my [previous post][1] I showed how the number of neighbours to a point in
N-Dimensional space increases exponentially:

[1]: https://paddy3118.blogspot.com/2021/01/conways-game-of-life-in-7-dimensions.html

"""

# %% ##Function

def neighbours(point):
    return {n
            for n in product(*((n-1, n, n+1)
                               for n in point))
            if n != point}

for dim in range(1,10):
    n = neighbours(tuple([0] * dim))
    print(f"A {dim}-D grid point has {len(n):_} neighbours.")
    assert len(n) == 3**dim - 1

# %% ##Output
"""
A 1-D grid point has 2 neighbours.
A 2-D grid point has 8 neighbours.
A 3-D grid point has 26 neighbours.
A 4-D grid point has 80 neighbours.
A 5-D grid point has 242 neighbours.
A 6-D grid point has 728 neighbours.
A 7-D grid point has 2_186 neighbours.
A 8-D grid point has 6_560 neighbours.
A 9-D grid point has 19_682 neighbours.
"""

# %% ##Visuals
"""
Here we define neighbours as all cells, n whose indices differ by at most 1,
i.e. by -1, 0, or +1 from the point X itself; apart from the point X.

In 1-D:

    .nXn.

In 2-D:

    .....
    .nnn.
    .nXn.
    .nnn.
    .....

# 'Axial' neighbours

 There is another definition of neighbour that "cuts out the diagonals" in 2-D
 to form:

     .....
     ..n..
     .nXn.
     ..n..
     .....

In 3-D this would add two extra cells, one directly above and one below X.

In 1-D this would be the two cells one either side of X.

I call this axial because if X is translated to the origin then axial
neighbours are the 2N cells in which only one coordinate can differ by either
+/- 1.

# Let's calclate them.
"""

# %% #Axial neighbourhood
from typing import Tuple, Set, Dict

Point = Tuple[int, int, int]
PointSet = Set[Point]


def axial_neighbours(point: Point) -> PointSet:
    dim, pt = len(point), list(point)
    return {tuple(pt[:d] + [pt[d] + delta] + pt[d+1:])
            for d in range(dim) for delta in (-1, +1)}

print("\n" + '='*60 + '\n')
for dim in range(1,10):
    n = axial_neighbours(tuple([0] * dim))
    text = ':\n  ' + str(sorted(n)) if dim <= 2 else ''
    print(f"A {dim}-D grid point has {len(n):_} axial neighbours{text}")
    assert len(n) == 2*dim

# %% ##Output
"""
============================================================

A 1-D grid point has 2 axial neighbours:
  [(-1,), (1,)]
A 2-D grid point has 4 axial neighbours:
  [(-1, 0), (0, -1), (0, 1), (1, 0)]
A 3-D grid point has 6 axial neighbours
A 4-D grid point has 8 axial neighbours
A 5-D grid point has 10 axial neighbours
A 6-D grid point has 12 axial neighbours
A 7-D grid point has 14 axial neighbours
A 8-D grid point has 16 axial neighbours
A 9-D grid point has 18 axial neighbours
"""

# %% #Another way of classifying...
"""
The 3-D axial neighbourhood around X can be described as:
    Think of X as a cube with six sides oriented with the axes of the 3
    dimensions. The axial neighbourhood is six similar cubes with one face
    aligned with each of X's faces. A kind of 3-D cross of cubes.

In 3-D, the "full" neighbourhood around point X describes a 3x3x3 units cube
centered on X.

In 3-D: Thinking of that 3x3x3 cube around X:

    What if we excluded the six corners of that cube?
    * Those excluded corner points have all their three coordinates different
      from those of X,
      i.e. if excluded e = Point(e_x, e_y, e_z) and X = Point(x_x, x_y, x_z):
          e_x != x_x and e_y != x_y and e_z != x_z
      Also:
          |e_x - x_x| == 1 and |e_y - x_y| == 1 and |e_z - x_z| == 1

    * We _Keep_ all points in the cube that have 1 and 2 different coords to
      those of X

    Another definition of the _axial_ neighbourhood case is
    * We _Keep_ all points in the cube that have only 1 different coord to
      those of X

This can be generalised to thinking in N dimensions of neighbourhood types
keeping from 1 to _up to_ N differences in coordinates (DinC).
"""
def axial_neighbours(point: Point) -> PointSet:
    dim, pt = len(point), list(point)
    return {tuple(pt[:d] + [pt[d] + delta] + pt[d+1:])
            for d in range(dim) for delta in (-1, +1)}

#%% #Neighbourhood by Differences in Coordinates

from collections import defaultdict
from pprint import pformat
from math import comb, factorial as fact


def d_in_c_neighbourhood(dim: int) -> Dict[int, PointSet]:
    """
    Split neighbourhood cube around origin point in dim-dimensions mapping
    count of coords not-equal to origin -> those neighbours
    """
    origin = tuple([0] * dim)
    cube = {n
            for n in product(*((n-1, n, n+1)
                               for n in origin))
            if n != origin}
    d_in_c = defaultdict(set)
    for n in cube:
        d_in_c[dim - n.count(0)].add(n)
    return dict(d_in_c)

def _check_counts(d: int, c: int, n: int) -> None:
    """
    Some checks on counts

    Parameters
    ----------
    d : int
        Dimensionality.
    c : int
        Count of neighbour coords UNequal to origin.
    n : int
        Number of neighbour points with exactly c off-origin coords.

    Returns
    -------
    None.

    """
    # # From OEIS
    # if   c == 1: assert n == d* 2
    # elif c == 2: assert n == d*(d-1)* 2
    # elif c == 3: assert n == d*(d-1)*(d-2)* 4 / 3
    # elif c == 4: assert n == d*(d-1)*(d-2)*(d-3)* 2 / 3
    # elif c == 5: assert n == d*(d-1)*(d-2)*(d-3)*(d-4)* 4 / 15
    # # Getting the hang of it
    # elif c == 6: assert n == d*(d-1)*(d-2)*(d-3)*(d-4)*(d-5)* 4 / 45

    # Noticed some powers of 2 leading to
    # if c == 6: assert n == d*(d-1)*(d-2)*(d-3)*(d-4)*(d-5)* 2**6 / fact(6)
    # if c == 6: assert n == fact(d) / fact(c) / fact(d - c) * 2**c

    # Finally:
    assert n == fact(d) / fact(c) / fact(d - c) * 2**c


print("\n" + '='*60 + '\n')
for dim in range(1,10):
    d = d_in_c_neighbourhood(dim)
    tot = sum(len(n_set) for n_set in d.values())
    print(f"A {dim}-D point has a total of {tot:_} full neighbours of which:")
    for diff_count in sorted(d):
        n_count = len(d[diff_count])
        print(f"  {n_count:4_} have exactly {diff_count:2} different coords.")
        _check_counts(dim, diff_count, n_count)

# %% ##Output
"""
============================================================

A 1-D point has a total of 2 full neighbours of which:
     2 have exactly  1 different coords.
A 2-D point has a total of 8 full neighbours of which:
     4 have exactly  1 different coords.
     4 have exactly  2 different coords.
A 3-D point has a total of 26 full neighbours of which:
     6 have exactly  1 different coords.
    12 have exactly  2 different coords.
     8 have exactly  3 different coords.
A 4-D point has a total of 80 full neighbours of which:
     8 have exactly  1 different coords.
    24 have exactly  2 different coords.
    32 have exactly  3 different coords.
    16 have exactly  4 different coords.
A 5-D point has a total of 242 full neighbours of which:
    10 have exactly  1 different coords.
    40 have exactly  2 different coords.
    80 have exactly  3 different coords.
    80 have exactly  4 different coords.
    32 have exactly  5 different coords.
A 6-D point has a total of 728 full neighbours of which:
    12 have exactly  1 different coords.
    60 have exactly  2 different coords.
   160 have exactly  3 different coords.
   240 have exactly  4 different coords.
   192 have exactly  5 different coords.
    64 have exactly  6 different coords.
A 7-D point has a total of 2_186 full neighbours of which:
    14 have exactly  1 different coords.
    84 have exactly  2 different coords.
   280 have exactly  3 different coords.
   560 have exactly  4 different coords.
   672 have exactly  5 different coords.
   448 have exactly  6 different coords.
   128 have exactly  7 different coords.
A 8-D point has a total of 6_560 full neighbours of which:
    16 have exactly  1 different coords.
   112 have exactly  2 different coords.
   448 have exactly  3 different coords.
  1_120 have exactly  4 different coords.
  1_792 have exactly  5 different coords.
  1_792 have exactly  6 different coords.
  1_024 have exactly  7 different coords.
   256 have exactly  8 different coords.
A 9-D point has a total of 19_682 full neighbours of which:
    18 have exactly  1 different coords.
   144 have exactly  2 different coords.
   672 have exactly  3 different coords.
  2_016 have exactly  4 different coords.
  4_032 have exactly  5 different coords.
  5_376 have exactly  6 different coords.
  4_608 have exactly  7 different coords.
  2_304 have exactly  8 different coords.
   512 have exactly  9 different coords.
"""

# %% #Manhattan
"""
What, up until now, I have called the Differences in Coordinates is better
known as the Manhattan distance, or [Taxicab geometry][2].


[2]: https://en.wikipedia.org/wiki/Taxicab_geometry
"""

# %% #Formula for Taxicab counts of cell neighbours
"""
Function _check_counts was originally set up as a crude check of taxicab
distance of 1 as the formula was obvious.

The formulea for taxcab distances 2, 3, and 4 were got from a search on [OEIS][3]

The need for factorials is obvious. The count for a taxicab distance of N in
N-D space is 2**N which allowed me to work out a final formula:

If:
    d is the number of spatial dimiensions.
    t is the taxicab distance of a neighbouring point to the origin.
Then:
    n, the count of neighbours with exactly this taxicab distance is

        n = f(d, t)
          = d! / t! / (d-t)! * 2**t


We can use the assertion from the exercising of function neighbours() at the
beginning to state that:

        sum(f(d, t)) for t from 0 .. d
          = g(d)
          = 3**d - 1


[3]: https://oeis.org/
"""

# %% ##Taxicab visualisation
"""
Lets create a visualisation of the difference in coords for neighbours in <= 4-D.
(The function is general, but I'm lost after 3-D)!

The origin will show as zero, 0; and neighbours surround it as digits which are
the taxicab distance from 0.
"""

def to_str(taxi: Dict[int, PointSet], indent: int=4) -> str:
    """


    Parameters
    ----------
    taxi : Dict[int, PointSet]
        Map taxicab distance from origin
        -> set off neighbours with that difference.
    indent : int
        indent output with spaces

    Returns
    -------
    str
        Visusalisation of region.

    """

    if not taxi:
        return ''

    ap = set()              # all points
    neighbour2taxi = {}
    for taxi_distance, nbrs in taxi.items():
        ap |= nbrs
        for n in nbrs:
            neighbour2taxi[n] = taxi_distance

    # Dimensionality
    dims = len(n)

    # Add in origin showing as 0
    origin = tuple([0] * dims)
    ap.add(origin)
    neighbour2taxi[origin] = 0

    # Plots neighbourhood of origin (plus extra space)
    space = 1
    minmax = [range(-(1 + space), (1 + space) + 1)
              for dim in range(dims)]

    txt = []
    indent_txt = ' ' * indent
    for plane_coords in product(*minmax[2:]):
        ptxt = ['\n' + indent_txt
                + ', '.join(f"dim{dim}={val}"
                            for dim, val in enumerate(plane_coords, 2))]

        ptxt += [''.join((str(neighbour2taxi[tuple([x, y] + list(plane_coords))])
                          if tuple([x, y] + list(plane_coords)) in ap
                          else '.')
                         for x in minmax[0])
                 for y in minmax[1]]

        # Don't plot planes with no neighbours (due to extra space).
        if ''.join(ptxt).count('.') < (3 + space*2)**2:
            txt += ptxt

    return '\n'.join(indent_txt + t for t in txt)


print("\n" + '='*60 + '\n')
for dim in range(2,5):
    d = d_in_c_neighbourhood(dim)
    tot = sum(len(n_set) for n_set in d.values())
    print(f"\nA {dim}-D point has a total of {tot:_} full neighbours of which:")
    for diff_count in sorted(d):
        n_count = len(d[diff_count])
        print(f"  {n_count:4_} have taxicab distance {diff_count:2} from the origin.")
        _check_counts(dim, diff_count, n_count)
    print(to_str(d))

# %% ##Output
"""
============================================================


A 2-D point has a total of 8 full neighbours of which:
     4 have taxicab distance  1 from the origin.
     4 have taxicab distance  2 from the origin.


    .....
    .212.
    .101.
    .212.
    .....

A 3-D point has a total of 26 full neighbours of which:
     6 have taxicab distance  1 from the origin.
    12 have taxicab distance  2 from the origin.
     8 have taxicab distance  3 from the origin.

    dim2=-1
    .....
    .323.
    .212.
    .323.
    .....

    dim2=0
    .....
    .212.
    .101.
    .212.
    .....

    dim2=1
    .....
    .323.
    .212.
    .323.
    .....

A 4-D point has a total of 80 full neighbours of which:
     8 have taxicab distance  1 from the origin.
    24 have taxicab distance  2 from the origin.
    32 have taxicab distance  3 from the origin.
    16 have taxicab distance  4 from the origin.

    dim2=-1, dim3=-1
    .....
    .434.
    .323.
    .434.
    .....

    dim2=-1, dim3=0
    .....
    .323.
    .212.
    .323.
    .....

    dim2=-1, dim3=1
    .....
    .434.
    .323.
    .434.
    .....

    dim2=0, dim3=-1
    .....
    .323.
    .212.
    .323.
    .....

    dim2=0, dim3=0
    .....
    .212.
    .101.
    .212.
    .....

    dim2=0, dim3=1
    .....
    .323.
    .212.
    .323.
    .....

    dim2=1, dim3=-1
    .....
    .434.
    .323.
    .434.
    .....

    dim2=1, dim3=0
    .....
    .323.
    .212.
    .323.
    .....

    dim2=1, dim3=1
    .....
    .434.
    .323.
    .434.
    .....

"""

Conways Game of Life in 7 Dimensions!

 

I stumbled across the site "Advent of Code" and its puzzle Conways cube by accident. I really liked the idea of the problem, even though it took me some time to understand the given examples, but then I really got into it.

I wrote a 3D answer; updated it to a 4D answer, then set about changing the Python source to handle any dimension >= 2.

The problem 

This is the gist of what is requested

• Conways Game of life, (GOL), is usually played on a 2-D rectangular grid of cells that are in one of two states: inactive or active at the end of any cycle. 
• We map a coordinate system on the grid where grid intersections are atall whole  integer coordinates. (Negative too).
  
• The neighbours of a cell on the grid are defined as any cell with all coordinates that differs from the target cell by either -1, 0, or +1 independently - apart from the target cell in the center of course.
  For the original 2-D GOL, each cell has 8 neighbours.
  
• The grid, (called a  Pocket Domain in the puzzle), is first seeded/initialised with an arrangement of activated  cells in the 2-D plane coresponding to cycle(0).
• The configuration of active cells in the next cycle cycle(n+1) is calculated from the configuration in the current cycle cycle(n) by applying the following two simple rules to all cells:
1. If a cell is active in cycle(n) it will stay active if it has 2 or 3 active neighbours; otherwise it becomes inactive in cycle(n+1).
2. If a cell is inactive in cycle(n) it will become active if it has exactly 3 active neighbours; otherwise it becomes inactive in cycle(n+1).
•  Note that although the pocket domain/grid is infinite, the rules need only be applied in the neighbourhood of any active cedlls as all cells outside that are inactive and can only remain inactive in the next cycle.
• The puzzle defines a starting configuration/format and asks how many cells are active in some future cycle.
  

The textual format

For 2-D it is a rectangular plane of full stops '.' and hashes '#' where hashes represent active cells on the plane, and full-stops inactive cells. Lets call this the x-y coordinate plane

• The rectangle shown is only large enough to show all active cells.

 

.#.
..#
###

A 3-D pocket domain is visualised as one or more x-y slices through the cubic x,y,z extent of the domain,The z coordinate for each slice in the cycle is different and stated above each slice like so:

z=-1
#..
..#
.#.

z=0
#.#
.##
.#.

z=1
#..
..#
.#.

The range of x and y coords shown, of all x-y slices for a cycle stay the same; but may be different between cycles as the extent of all active cells changes between cycles.

Their 4-D textual format introduces w as the forth dimension and now x-y planes through the encompassing hyper-rectangle of all active cells need both the z and w coords specified for each plane.

z=-1, w=-1
#..
..#
.#.

z=0, w=-1
#..
..#
.#.

z=1, w=-1
#..
..#
.#.

z=-1, w=0
#..
..#
.#.

N-Dimensional Solution

My N-Dimensional extension to the Text Format

My format is a slight change in that it will read the above, but names the dimensions greater than two, (assuming the x-y plane  as dim0 and dim1), dim2, dim3, ... for example:

dim2=0, dim3=0, dim4=0
.#.
..#
###

The internal grid data-structure

I started with no idea if there was to be a limit to the extent of the 3-D pocket dimension, and new that 2-D GOL looked pretty sparse w.r.t. active cells.

 

So I opted for a sparse grid.

A bit of thought and I confirmed I only needed to keep active cells, from which I could calculate all neighbours that could be affected then calculate what the active cells of the next cycle become.

If I keep the coordinate of cells as a tuple then I will be doing a lot of lookup on whether a coord is active or not. I made the pocket domain/grid/pd a set of active points where points were tuples. Checking if a neighbour is active amounts to checking if it is in the set.

Reading the textual format

My initial 3-D reader only reads the 3-d format. It uses a regexp to extract the Z and XY sections and turns them into a set of the active points defined.

from itertools import product
import re

#%% to/From textual format

# Help: https://regex101.com/codegen?language=python
pocket_plane_re = r"""
	# NOTE: add an extra \n\n to end of string being matched

	^z=(?P<Z>[-+0-9]+?)    # Single z-plane coord
	\s*?\n
	(?P<XY>.*?\n)\n        # Following xy-plane description
	"""

def _plane_from_str(text, z):
    """single pocket dim z-plane to set of active cell coords"""
    pd = {(x, y, z)         # Pocket dimension
          for y, row in enumerate(text.strip().split('\n'))
          for x, state in enumerate(row)
          if state == '#'}
    return pd

def from_str(text):
    """Pocket dim text to set of active cell coords, the internal format"""
    pd = set()              # pocket dim active cells accumulator
    matches = re.finditer(pocket_plane_re, text + '\n\n',
                          re.MULTILINE | re.DOTALL | re.VERBOSE)
    for match in matches:
        z = int(match.groupdict()['Z'])
        xy =    match.groupdict()['XY']
        pd |= _plane_from_str(xy, z)
    return pd

function from_str separates and extracts the z coord and the xy plane from the text then passes them to function _plane_from_str to return a set of active cell coordinates to add to the pocket domain.

The N-Dimensional reader 

from itertools import product
import re

#%% to/From textual format

# Help: https://regex101.com/codegen?language=python

# Parses textaul version of state.
pocket_plane_re = r"""
	# Note: add an extra \n\n to end of string being matched,
    #       and a \n to the beginning.

	# Coordinates of plane, (optional for 2D).
	^(?P<PLANE>(?: [a-z][a-z0-9]*=(?P<NUM1>-?\d+))
	    (?:,\s+ (?: [a-z][a-z0-9]*=(?P<NUM2>-?\d+)))*)?

	\s*?\n
	(?P<XY>[.#\n]+?\n)\n        # xy-plane description
	"""

def _plane_from_str(text, plane_coords):
    """single pocket dim z-plane to set of active cell coords"""
    pcoords = list(plane_coords)
    pd = {tuple([x, y] + pcoords)         # Pocket dimension
          for y, row in enumerate(text.strip().split('\n'))
          for x, state in enumerate(row)
          if state == '#'}
    return pd

def from_str(text: str) -> set:
    """Pocket dim text to set of active cell coords, the internal format"""
    pd = set()              # pocket dim active cells accumulator
    matches = re.finditer(pocket_plane_re, '\n' + text + '\n\n',
                          re.MULTILINE | re.DOTALL | re.VERBOSE)
    for match in matches:
        # Extra coords for the xy plane
        plane_coords_txt = match.groupdict()['PLANE'] or ''
        plane_coords = [int(matchp.group()[1:])
                        for matchp in
                        re.finditer(r"=-?\d+", plane_coords_txt)]
        # xy plane
        xy = match.groupdict()['XY']
        # Accumulate active cells in plane
        pd |= _plane_from_str(xy, plane_coords)
    return pd

 

An extended regexp handles their being no extra coord needed for the x-y plane, when working in 2-D, as well as the unbounded list of extra dimension coords in more than 2-D.

 

Writing the textual format.

Creates a textual representation of the internal format of a cycle.

 

In the 3-D case I first finds the minimum/maximum extents of activity on all three axis in minmax. Stepping through the range of the z axis it then finds and prints all the activity in the x-y plane for that value of z

I had had difficulty in understanding that the extend of x and y changes between cycles so argument print_extent was used in debugging.

def to_str(pd, print_extent=False):
    "From set of active cell coords to output textual format"
    # Extent of pocket dimenson
    minmax = [range(min(pd, key=lambda pt:pt[dim])[dim],
                    max(pd, key=lambda pt:pt[dim])[dim] + 1)
              for dim in range(3)]

    txt = [f"\n// FROM: {tuple(r.start for r in minmax)}"
           f" TO: {tuple(r.stop -1 for r in minmax)}"] if print_extent else []

    append = txt.append
    for z in minmax[2]:
        append(f"\nz={z}")
        for y in minmax[1]:
            append(''.join(('#' if (x, y, z) in pd else '.')
                           for x in minmax[0]))

    return '\n'.join(txt)

The N-Dimensional writer

def to_str(pd, print_extent=False) -> str:
    "From set of active cell coords to output textual format"
    if not pd:
        return ''
    # Dimensionality of pocket dimension
    point = pd.pop()
    dims = len(point)
    pd |= {point}           # Put that point back

    # Extent of pocket dimenson
    minmax = [range(min(pd, key=lambda pt:pt[dim])[dim],
                    max(pd, key=lambda pt:pt[dim])[dim] + 1)
              for dim in range(dims)]

    txt = [f"\n// FROM: {tuple(r.start for r in minmax)}"
           f" TO: {tuple(r.stop -1 for r in minmax)}"] if print_extent else []

    for plane_coords in product(*minmax[2:]):

        ptxt = ['\n' + ', '.join(f"dim{dim}={val}"
                                 for dim, val in enumerate(plane_coords, 2))]
        ptxt += [''.join(('#'
                          if tuple([x, y] + list(plane_coords)) in pd
                          else '.')
                         for x in minmax[0])
                 for y in minmax[1]]
        if '#' in ''.join(ptxt):
            txt += ptxt

    return '\n'.join(txt)

To find out how many dimensions there are, a point is poped from, (and later added back to), the set of active points, pd, and stored in dims. minmax and assignments are similar,

The creation of the text for each x-y plane is separated out. many of the x-y planes can be of all inactive cells - this modification will not display x-y planes if they have no set cells.

Generating all neighbours to a cell

An important part of the program that is repeated many times, as well as many times for the same cell coordinates.

3-D neighbour function exploration

Knowing that it is calculated many times for the same cell I decided to first have a cached function. I skipped using lru_cache for a simple, fast homebrew of using a dict _neighbour_cache.

_neighbour_cache = {}
def neighbours(point):
    "return the points 26 neighbour coords"
    if point not in _neighbour_cache:
        x, y, z = point
        _neighbour_cache[point] = {xyz
                                   for xyz in product((x-1, x, x+1),
                                                      (y-1, y, y+1),
                                                      (z-1, z, z+1))
                                   if xyz != point}
    return _neighbour_cache[point]

That cache was large and in higher dimensions would get huge or inefficient, so I decided to time some alternate functions for the 3-D case.

Alternate implementations and timings

def neighboursX(point):
    x, y, z = point
    return {
         (x-1, y-1, z-1), (x  , y-1, z-1), (x+1, y-1, z-1),
         (x-1, y  , z-1), (x  , y  , z-1), (x+1, y  , z-1),
         (x-1, y+1, z-1), (x  , y+1, z-1), (x+1, y+1, z-1),

         (x-1, y-1, z  ), (x  , y-1, z  ), (x+1, y-1, z  ),
         (x-1, y  , z  ),                  (x+1, y  , z  ),
         (x-1, y+1, z  ), (x  , y+1, z  ), (x+1, y+1, z  ),

         (x-1, y-1, z+1), (x  , y-1, z+1), (x+1, y-1, z+1),
         (x-1, y  , z+1), (x  , y  , z+1), (x+1, y  , z+1),
         (x-1, y+1, z+1), (x  , y+1, z+1), (x+1, y+1, z+1),
        }

def neighboursY(point):
    x, y, z = point
    return {xyz
            for xyz in product((x-1, x, x+1),
                               (y-1, y, y+1),
                               (z-1, z, z+1))
            if xyz != point}

"""

In [80]: p = (1, 1, 1)

In [81]: neighbours(p) == neighboursX(p) == neighboursY(p)
Out[81]: True

In [82]: %timeit neighbours(p)
719 ns ± 20.2 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)

In [83]: %timeit neighboursX(p)
9.75 µs ± 82 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

In [84]: %timeit neighboursY(p)
11.2 µs ± 265 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

In [85]:

"""

Caching using a simple global dict in neighbours is an order of magnitude
faster than calculations, but calculations take no extra memory

neighboursY seems to hit the sweet spot of maintainability speed and memory use and was adopted for use in the 4-D version.

4-D neighbour function

def neighbours(point):
    x, y, z, w = point
    return {xyzw
            for xyzw in product((x-1, x, x+1),
                                (y-1, y, y+1),
                                (z-1, z, z+1),
                                (w-1, w, w+1))
            if xyzw != point}

It is just a 4-D extension of the 3-D neighboursY function.

The N-D neighbour function

def neighbours(point):
    return {n
            for n in product(*((n-1, n, n+1)
                               for n in point))
            if n != point}

 

Working that out was sweet :-)

 

Checking if a cell becomes active in the next cycle

 Take a point and the set of all active points in the current cycle then checks the points neighbours to see if this point will be active or not in the next cycle

def is_now_active(point, pd):
    "Is point active w.r.t. pd"
    ncount = 0              # active neighbours count
    for npoint in neighbours(point):
        ncount += npoint in pd
        if ncount >3:
            break
    if point in pd:         # Currently active
        return 2 <= ncount <= 3
    else:                   # Currently inactive
        return ncount == 3

 

The above didn't need changing from being written for the 3-D case and being reused in the 4-D and N-D versions.

 

in higher dimensions points can have many neighbours and many possible active neighbours but we only need a count of up to 4 for the following checks to work.

 

 Generating the next cycle from the current.

This again did not need changing from the 3-D case

def next_cycle(pd: set) -> set:
    "generate next cycle from current Pocket Dimension pd"
    possibles = set()           # cells that could become active
    for active in pd:
        possibles |= neighbours(active)
    possibles |= pd

    return {point for point in possibles if is_now_active(point, pd)}

it checks all current points and their neighbours to see if they will be active in the next cycle.

Examples of N-Dimensional results

if __name__ == '__main__':
    cycle0 = """
.#.
..#
###
"""

    pd2 =  from_str(cycle0)
    start_pd_for_dim = {d: {tuple(list(pt) + [0]*(d-2)) for pt in pd2}
                        for d in range(2, 8)}
    for dim in range(2, 8):
        print('='*32)
        n, pd = 0, start_pd_for_dim[dim]
        lines, active = 0, 1
        while True:
            text = to_str(pd)
            active = text.count('#')
            active_planes = ('\n' + text).count('\n\n')
            lines = text.count('\n')
            print(f"\n{dim}-D: After {n} cycles, {active:_} cells on "
                  f"{active_planes:_} active plane(s):\n"
                  f"{text if lines < 60 else ''}")
            if not (lines <2**12 and 0 < active < 2_000 and n < 10):
                break
            n, pd = n + 1, next_cycle(pd)

Just using the one 2-D initial cycle to create the same initial x-y plane in 2 to seven pocket domains in dict start_pd_for_dim.

Then for increasing pocket domain dimensions run the Conway cycles. cycle printing is stopped for large outputs. and whole dimensions are truncated when they get "large"

Output

A sample from the end of the output for the __main__ above

================================

6-D: After 0 cycles, 5 cells on 1 active plane(s):

dim2=0, dim3=0, dim4=0, dim5=0
.#.
..#
###

6-D: After 1 cycles, 245 cells on 81 active plane(s):


6-D: After 2 cycles, 464 cells on 48 active plane(s):


6-D: After 3 cycles, 15_744 cells on 1_504 active plane(s):

================================

7-D: After 0 cycles, 5 cells on 1 active plane(s):

dim2=0, dim3=0, dim4=0, dim5=0, dim6=0
.#.
..#
###

7-D: After 1 cycles, 731 cells on 243 active plane(s):


7-D: After 2 cycles, 1_152 cells on 112 active plane(s):


7-D: After 3 cycles, 106_400 cells on 10_320 active plane(s):

Extra: Just how many neighbours does a point have?

>>>  for d in range(1,15):
         n = len(neighbours(tuple([0] * d)))
         print(f"A {d}-D grid point has {n:_} neighbours.")
         assert n == 3**d - 1
A 1-D grid point has 2 neighbours.
A 2-D grid point has 8 neighbours.
A 3-D grid point has 26 neighbours.
A 4-D grid point has 80 neighbours.
A 5-D grid point has 242 neighbours.
A 6-D grid point has 728 neighbours.
A 7-D grid point has 2_186 neighbours.
A 8-D grid point has 6_560 neighbours.
A 9-D grid point has 19_682 neighbours.
A 10-D grid point has 59_048 neighbours.
A 11-D grid point has 177_146 neighbours.
A 12-D grid point has 531_440 neighbours.
A 13-D grid point has 1_594_322 neighbours.
A 14-D grid point has 4_782_968 neighbours.

In [3]: 

 

There's a lot of neighbours in higher dimensions!

 

 END.