gagern · September 29, 2024 00:33
diff --git a/mx4976109a.py b/mx4976109a.py
 from collections import deque
 from itertools import chain

 # Counting the number of ways how a 4×N rectangle can be filled with terominoes,
 # counting symmetric versions of the same tiling only once.
 #
 # https://math.stackexchange.com/a/4976754/35416
 # Martin von Gagern

 # We store tiles and unfinished rows as the bits of a single number:
 #
 # Bit index:    Bit value:
 # 12 13 14 15   0x1000 0x2000 0x4000 0x8000
 # 08 09 10 11   0x0100 0x0200 0x0400 0x0800
 # 04 05 06 07   0x0010 0x0020 0x0040 0x0080
 # 00 01 02 03   0x0001 0x0002 0x0004 0x0008

 # Some such bit patterns get one extra bit to denote that we are dealing with a
 # pair of tiles instead of a single tile. See "half-pair states" in the post.
 PAIRMARK = 0x10000


 def tetrominoes():
    """Construct bit patterns representing possible positions for tetrominoes.

    This function constructs all possible rotations and all possible horizontal
    positions, but vertically all tiles will have their first square on row 0
    (lowest bit value). The returned value is a quadruple of tuples. These are
    the bit patterns grouped by the leftmost (i.e. lowest) bit set in them.
    """

    # I'm too lazy to generate these programmatically, so writing them out
    # instead:
    all_tiles = """\
 ##
 ##

 ####

 #
 #
 #
 #

 ###
 #

 ###
  #

 #
 ###

  #
 ###

 #
 #
 ##

 ##
 #
 #

 #
 #
 ##

 ##
 #
 #

 ##
 ##

 ##
 ##

 #
 ##
 #

 #
 ##
 #

 ###
 #

 #
 ##
 #

 #
 ##
 #

 #
 ###
 """
    each_tile = all_tiles.split("\n\n")
    assert len(each_tile) == 19, "Wikipedia says there should be 19"
    t = [[], [], [], []]
    for shape in each_tile:
        cnt = 0
        num = 0
        # Note: In the input here we use the top row as row 0, while later code
        # puts row 0 at the bottom. So the input in string a has mirror images
        # from how the other code interprets these, but since we include all the
        # mirror images anyway, nothing changes.
        for y, row in enumerate(shape.split("\n")):
            for x, ch in enumerate(row):
                if ch != " ":
                    cnt += 1
                    num |= 1 << (4 * y + x)
        assert cnt == 4, f"Four # required:\n{shape}"
        assert num & 0xf != 0, f"# required in first row:\n{shape}"
        assert num & 0x1111 != 0, f"# required in first col:\n{shape}"
        # Find the lowest set bit so we can put it in the right group.
        i = 0
        while num & (1 << i) == 0:
            i += 1
        t[i].append(num)
        # While the rightmost column is empty, we can move the tile right.
        while num & 0x8888 == 0:
            i += 1
            num <<= 1
            t[i].append(num)
    return tuple(map(tuple, t))

 class Automaton:
    """Represents a finite state automaton."""

    def __init__(self, tiles, initial_states, print_transitions=False):
        """Constructs an automaton for given set of tiles.

        The automaton will be able to process all the given initial states.
        """
        # Breadth-first graph traversal
        seen = set(initial_states)
        q = deque(sorted(seen))
        edges = []
        while q:
            src = q.popleft()
            if src & PAIRMARK:
                # This is a "half-pair state". The only transition is to the
                # corresponding regular state.
                dst = src ^ PAIRMARK
                edges.append((src, dst))
                if dst not in seen:
                    seen.add(dst)
                    q.append(dst)
                continue
            # Find lowest unset bit.
            i = 0
            while src & (1 << i) != 0:
                i += 1
            for tile in tiles[i]:
                if src & tile != 0:
                    continue
                # dst is like src but with bits from tile set as well.
                dst = src | (tile & 0xffff)
                # But then we drop any completed rows.
                while dst & 0xf == 0xf:
                    dst >>= 4
                # If tile was a pair, then dst is marked a "half-pair state".
                dst |= tile & PAIRMARK
                edges.append((src, dst))
                if dst not in seen:
                    seen.add(dst)
                    q.append(dst)

                if print_transitions:
                    print(f"{src:04x} → {dst:04x}")
                    for y in range(3, -1, -1):
                        row = []
                        for x in range(4):
                            bit = 1 << (4*y + x)
                            if src & bit:
                                row.append("@")
                            elif tile & bit:
                                row.append("#")
                            else:
                                row.append(" ")
                        row.append("   ")
                        for x in range(4):
                            bit = 1 << (4*y + x)
                            if dst & bit:
                                row.append("@")
                            else:
                                row.append(" ")
                        print("".join(row))
                    print("")
        print(f"{len(seen):3d} states, {len(edges):3d} transitions before dead end elimination")

        # Dead end eliminiation.
        alive = set([0])
        n = 0
        while n != len(alive):
            n = len(alive)
            for src, dst in edges:
                if dst in alive:
                    alive.add(src)
        edges = [(src, dst) for src, dst in edges if src in alive and dst in alive]
        edges.sort()
        print(f"{len(alive):3d} states, {len(edges):3d} transitions  after dead end elimination")
        self.edges = edges

        # Mapping between state bit patterns (which can be large) and indices
        # (much smaller) for efficiency.
        state2idx = {}
        idx2state = []
        for s in sorted(alive):
            if s not in state2idx:
                state2idx[s] = len(state2idx)
                idx2state.append(s)
        for s in seen:
            # Mark dead end states as known but not represented. This
            # distinguishes them from states that were not even explored.
            state2idx.setdefault(s, None)
        self.state2idx = state2idx
        self.idx2state = idx2state

        # Build instructions for how to compute data for one additional step in
        # the graph.
        step = [[] for _ in idx2state]
        for src, dst in edges:
            step[state2idx[src]].append(state2idx[dst])
        self.step = tuple(map(tuple, step))
        # With 0 steps, there is 1 way to go from state 0 to the final state
        # (i.e. 0 as well). Any other state won't take you there in 0 steps.
        initial = [0] * len(idx2state)
        initial[0] = 1
        self.vectors = [tuple(initial)]

    def get(self, n, src=0):
        """Number of ways to get from state src to state 0 in n steps."""
        src_idx = self.state2idx[src]
        if src_idx is None:
            # Known to be a dead end.
            return 0
        while len(self.vectors) <= n:
            # Extend by one step. Semantically we follow edges to go from i to
            # j, then use vector to count ways to get from j to 0. The index i
            # is implicit in the iteration over self.step when constructing nxt.
            prv = self.vectors[-1]
            nxt = tuple(sum(prv[j] for j in js) for js in self.step)
            self.vectors.append(nxt)
        return self.vectors[n][src_idx]

    def dot(self, dot_name):
        """Write graph to a graphviz dot file, for visualization."""
        with open(dot_name, "w") as f:
            print("digraph {", file=f)
            print("  node [shape=box];", file=f)
            for node in self.idx2state:
                label = [f'  "{node:04x}" [label=<<TABLE CELLBORDER="0"']
                if node & PAIRMARK:
                    label.append(' COLOR="#00cccc"')
                label.append(' CELLPADDING="0" CELLSPACING="1">')
                for y in range(3, -1, -1):
                    label.append("<TR>")
                    for x in range(4):
                        label.append('<TD FIXEDSIZE="TRUE" WIDTH="7" HEIGHT="7"')
                        bit = 1 << (4*y + x)
                        if node & bit:
                            label.append(' BGCOLOR="#000099"')
                        label.append("></TD>")
                    label.append("</TR>")
                label.append("</TABLE>>];")
                print("".join(label), file=f)
            for src, dst in self.edges:
                print(f'  "{src:04x}" -> "{dst:04x}";', file=f)
            print("}", file=f)

 class Split:
    """Computes counts based on central arrangement and upper/lower parts."""

    def __init__(self, centers, automaton):
        self.centers = centers
        self.automaton = automaton

    def get(self, n, debug=False):
        """Number of ways to get a tiling with n tiles."""
        res = 0
        for nc, cs in enumerate(self.centers[:n+1]):
            if (n & 1) != (nc & 1):
                continue
            for c in cs:
                x = self.automaton.get((n - nc) >> 1, src=c)
                if debug:
                    print(f"{x:3d} ways to go from 0x{c:04x} to 0 in {(n - nc) >> 1} steps")
                res += x
        return res


 def reflect_one(tile):
    """Reflect a single bit pattern in the vertical axis."""
    return (((tile & 0x1111) << 3) |
            ((tile & 0x2222) << 1) |
            ((tile & 0x4444) >> 1) |
            ((tile & 0x8888) >> 3))

 def reflected(singles):
    """Reflect the given single tiles in the vertical axis.

    Each returned bit pattern will either be a single tile which is its own
    mirror image, or a pair of mirror image tiles with the PAIRMARK bit set in
    addition to the bits of the tiles. The returned format is like that of
    tetrominoes, namely a quadruple of tuples of bit patterns, where the
    quadruple groups patterns based on the lowest set bit. The last two of these
    are of course empty.
    """
    res = [[], [], [], []]
    for i, lst in enumerate(singles[:2]):
        for t in lst:
            r = reflect_one(t)
            if r == t:
                res[i].append(t)
            if r & t == 0:
                pair = t | r | PAIRMARK
                res[i].append(pair)
    return tuple(map(tuple, res))

 def reflective_centers():
    """Centeral arrangements for reflection in the horizontal axis.

    This returns a quintuple of tuples. The index in the outer tuple indicates
    how many tiles are used in the central arrangements, while the value inside
    the inner tuple is a bit pattern for the resultin shape that protrudes into
    the upper half of the tiling. Shapes may occur more than once.
    """
    # The upper halves of vertical bar and square:
    halves = [0x11, 0x22, 0x44, 0x88, 0x3, 0x6, 0xc]
    lst = []
    for i, a in enumerate(halves):
        for b in halves[i+1:]:
            if a & b == 0:
                s = a | b
                while s & 0xf == 0xf:
                    s >>= 4
                lst.append(s)
    lst.sort()
    # 0 tiles: straight line.
    # 1 tile: horizontal bar.
    # 2 tiles: as computed above.
    # 3 tiles: impossible.
    # 4 tiles: four vertical bars, straight line.
    res = ((0,), (0,), tuple(lst), (), (0,))
    print(f"reflective_centers: {tuple(map(len,res))!r}")
    return res

 def rotate_one(tile):
    """Rotate a single 4x4 bit pattern by 180°."""
    # This is effectively a bit reversal of a 16 bit integer.
    tile = ((tile & 0xff00) >> 8) | ((tile & 0x00ff) << 8)
    tile = ((tile & 0xf0f0) >> 4) | ((tile & 0x0f0f) << 4)
    tile = ((tile & 0xcccc) >> 2) | ((tile & 0x3333) << 2)
    tile = ((tile & 0xaaaa) >> 1) | ((tile & 0x5555) << 1)
    return tile

 def rotated_centers(singles, print_transitions=False):
    """Centeral arrangements for 180° rotation.

    This returns a quintuple of tuples. The index in the outer tuple indicates
    how many tiles are used in the central arrangements, while the value inside
    the inner tuple is a bit pattern for the resultin shape that protrudes into
    the upper half of the tiling. Shapes may occur more than once.
    """
    singles = list(chain.from_iterable(singles))

    # We model even pairs as 6 rows high and odd pairs as 7 rows high.
    even_pairs = []
    for a in singles:
        c = rotate_one(a)
        b = c << 8
        # Drop c down to row 0, so it becomes the canonical bit pattern for that
        # tile.
        while c & 0xf == 0:
            c >>= 4
        if a > c:
            # We will generate this pair when the loop iterates over c, so by
            # continuing now we avoid adding each pair twice.
            continue
        # Move the tile up (and its image down) until it protrudes into the
        # upper half.
        while a & 0xf000 == 0:
            a <<= 4
            b >>= 4
        # Keep moving the tile further up as long as it still protrudes into the
        # lower half.
        while a & 0xfff != 0:
            if a & b == 0:
                # The tile and its mirror image don't overlap, so it's a pair.
                even_pairs.append(a | b)
            a <<= 4
            b >>= 4

    odd_pairs = []
    for a in singles:
        c = rotate_one(a)
        b = c << 12
        while c & 0xf == 0:
            c >>= 4
        if a > c:
            continue
        # Move the tile up until it protrudes into the central row.
        while a & 0xf000 == 0:
            a <<= 4
            b >>= 4
        # This time we keep moving up until the tile still intersects the
        # central row, so one more f in that bit pattern here.
        while a & 0xffff != 0:
            if a & b == 0:
                odd_pairs.append(a | b)
            a <<= 4
            b >>= 4

    # Identify tiles that are their own rotated copy.
    solo_centers = set()
    for a in singles:
        b = rotate_one(a)
        if a << 12 == b:
            solo_centers.add(b)
        if a << 4 == b:
            solo_centers.add(a << 8)
    solo_centers = sorted(solo_centers)

    shapes = [[0], solo_centers, [], [], []]
    # Outer loop iterates over pairs, not over i (the number of tiles). This
    # avoids counting a configuration multiple times if the same pairs get added
    # to it in different order. Here the outer loop fixes the order.
    for pairs, parity in ((even_pairs, 0), (odd_pairs, 1)):
        for pair in pairs:
            for i in range(parity, 3, 2): # (0, 2) for even, just (1,) for odd.
                for src in shapes[i]:
                    if src & pair != 0:
                        # pair overlaps with what's already there.
                        continue
                    dst = src | pair
                    shapes[i + 2].append(dst)
                    if print_transitions:
                        s = dst >> 12
                        while s & 0xf == 0xf:
                            s >>= 4
                        print(f"{src:07x} → {dst:07x} ({i + 2} tiles, state {s:04x})")
                        for y in range(5 + parity, -1, -1):
                            row = []
                            for x in range(4):
                                bit = 1 << (4*y + x)
                                if src & bit:
                                    row.append("@")
                                elif pair & bit:
                                    row.append("#")
                                else:
                                    row.append(" ")
                            print("".join(row))

    # Take shapes of upper half of each central arrangement and store them in
    # the result.
    res = [[] for _ in shapes]
    for i, ss in enumerate(shapes):
        for s in ss:
            if i & 1 and s & 0xf000 != 0xf000:
                # Omit odd center with incomplete central row.
                continue
            s >>= 12
            while s & 0xf == 0xf:
                # Drop complete rows from state representation.
                s >>= 4
            res[i].append(s)
    print(f"rotated_centers: {tuple(map(len,res))!r}")
    return tuple(map(tuple, res))

 def fully_symmetric_centers(reflective_centers):
    """Centeral arrangements for full rectangle symmetry.

    This is like reflective_centers, but counting only the shapes that are
    symmetric around the vertical axis as well.
    """
    # Since squares and vertical bars are the only tiles that may get used
    # together, and they leave different heights in the resulting shape
    # protrusion, there is no danger that a symmetric shape could result from an
    # unsymmetric placement of tiles.
    return tuple(tuple(j for j in i if j == reflect_one(j)) for i in reflective_centers)

 def main():
    t = tetrominoes()
    rt = reflected(t)
    rc = reflective_centers()
    rot = rotated_centers(t)
    sym = fully_symmetric_centers(rc)

    b = Automaton(t, chain(chain.from_iterable(rc), chain.from_iterable(rot)))
    c = Automaton(rt, [0])
    d = Split(rc, b)
    e = Split(rot, b)
    f = Split(sym, c)

    print("| {:>3s} | {:>32s} | {:>32s} | {:>20s} | {:>20s} | {:>20s} | {:>10s} |".format("$N$", "$a_N$", "$b_N$", "$c_N$", "$d_N$", "$e_N$", "$f_N$"))
    sep = "-"*50
    print(f"|-{sep:.3s}:|-{sep:.32s}:|-{sep:.32s}:|-{sep:.20s}:|-{sep:.20s}:|-{sep:.20s}:|-{sep:.10s}:|")
    for i in range(51):
        bi = b.get(i)
        ci = c.get(i)
        di = d.get(i)
        ei = e.get(i)
        fi = f.get(i)
        ai = bi + ci + di + ei
        assert ai & 3 == 0, f"{ai} not divisible by 4"
        ai >>= 2
        print(f"| {i:3d} | {ai:32d} | {bi:32d} | {ci:20d} | {di:20d} | {ei:20d} | {fi:10d} |")

 if __name__ == "__main__":
    main()
	from collections import deque
	from itertools import chain

	# Counting the number of ways how a 4×N rectangle can be filled with terominoes,
	# counting symmetric versions of the same tiling only once.
	#
	# https://math.stackexchange.com/a/4976754/35416
	# Martin von Gagern

	# We store tiles and unfinished rows as the bits of a single number:
	#
	# Bit index: Bit value:
	# 12 13 14 15 0x1000 0x2000 0x4000 0x8000
	# 08 09 10 11 0x0100 0x0200 0x0400 0x0800
	# 04 05 06 07 0x0010 0x0020 0x0040 0x0080
	# 00 01 02 03 0x0001 0x0002 0x0004 0x0008

	# Some such bit patterns get one extra bit to denote that we are dealing with a
	# pair of tiles instead of a single tile. See "half-pair states" in the post.
	PAIRMARK = 0x10000


	def tetrominoes():
	"""Construct bit patterns representing possible positions for tetrominoes.

	This function constructs all possible rotations and all possible horizontal
	positions, but vertically all tiles will have their first square on row 0
	(lowest bit value). The returned value is a quadruple of tuples. These are
	the bit patterns grouped by the leftmost (i.e. lowest) bit set in them.
	"""

	# I'm too lazy to generate these programmatically, so writing them out
	# instead:
	all_tiles = """\
	##
	##

	####

	#
	#
	#
	#

	###
	#

	###
	#

	#
	###

	#
	###

	#
	#
	##

	##
	#
	#

	#
	#
	##

	##
	#
	#

	##
	##

	##
	##

	#
	##
	#

	#
	##
	#

	###
	#

	#
	##
	#

	#
	##
	#

	#
	###
	"""
	each_tile = all_tiles.split("\n\n")
	assert len(each_tile) == 19, "Wikipedia says there should be 19"
	t = [[], [], [], []]
	for shape in each_tile:
	cnt = 0
	num = 0
	# Note: In the input here we use the top row as row 0, while later code
	# puts row 0 at the bottom. So the input in string a has mirror images
	# from how the other code interprets these, but since we include all the
	# mirror images anyway, nothing changes.
	for y, row in enumerate(shape.split("\n")):
	for x, ch in enumerate(row):
	if ch != " ":
	cnt += 1
	num \|= 1 << (4 * y + x)
	assert cnt == 4, f"Four # required:\n{shape}"
	assert num & 0xf != 0, f"# required in first row:\n{shape}"
	assert num & 0x1111 != 0, f"# required in first col:\n{shape}"
	# Find the lowest set bit so we can put it in the right group.
	i = 0
	while num & (1 << i) == 0:
	i += 1
	t[i].append(num)
	# While the rightmost column is empty, we can move the tile right.
	while num & 0x8888 == 0:
	i += 1
	num <<= 1
	t[i].append(num)
	return tuple(map(tuple, t))

	class Automaton:
	"""Represents a finite state automaton."""

	def __init__(self, tiles, initial_states, print_transitions=False):
	"""Constructs an automaton for given set of tiles.

	The automaton will be able to process all the given initial states.
	"""
	# Breadth-first graph traversal
	seen = set(initial_states)
	q = deque(sorted(seen))
	edges = []
	while q:
	src = q.popleft()
	if src & PAIRMARK:
	# This is a "half-pair state". The only transition is to the
	# corresponding regular state.
	dst = src ^ PAIRMARK
	edges.append((src, dst))
	if dst not in seen:
	seen.add(dst)
	q.append(dst)
	continue
	# Find lowest unset bit.
	i = 0
	while src & (1 << i) != 0:
	i += 1
	for tile in tiles[i]:
	if src & tile != 0:
	continue
	# dst is like src but with bits from tile set as well.
	dst = src \| (tile & 0xffff)
	# But then we drop any completed rows.
	while dst & 0xf == 0xf:
	dst >>= 4
	# If tile was a pair, then dst is marked a "half-pair state".
	dst \|= tile & PAIRMARK
	edges.append((src, dst))
	if dst not in seen:
	seen.add(dst)
	q.append(dst)

	if print_transitions:
	print(f"{src:04x} → {dst:04x}")
	for y in range(3, -1, -1):
	row = []
	for x in range(4):
	bit = 1 << (4*y + x)
	if src & bit:
	row.append("@")
	elif tile & bit:
	row.append("#")
	else:
	row.append(" ")
	row.append(" ")
	for x in range(4):
	bit = 1 << (4*y + x)
	if dst & bit:
	row.append("@")
	else:
	row.append(" ")
	print("".join(row))
	print("")
	print(f"{len(seen):3d} states, {len(edges):3d} transitions before dead end elimination")

	# Dead end eliminiation.
	alive = set([0])
	n = 0
	while n != len(alive):
	n = len(alive)
	for src, dst in edges:
	if dst in alive:
	alive.add(src)
	edges = [(src, dst) for src, dst in edges if src in alive and dst in alive]
	edges.sort()
	print(f"{len(alive):3d} states, {len(edges):3d} transitions after dead end elimination")
	self.edges = edges

	# Mapping between state bit patterns (which can be large) and indices
	# (much smaller) for efficiency.
	state2idx = {}
	idx2state = []
	for s in sorted(alive):
	if s not in state2idx:
	state2idx[s] = len(state2idx)
	idx2state.append(s)
	for s in seen:
	# Mark dead end states as known but not represented. This
	# distinguishes them from states that were not even explored.
	state2idx.setdefault(s, None)
	self.state2idx = state2idx
	self.idx2state = idx2state

	# Build instructions for how to compute data for one additional step in
	# the graph.
	step = [[] for _ in idx2state]
	for src, dst in edges:
	step[state2idx[src]].append(state2idx[dst])
	self.step = tuple(map(tuple, step))
	# With 0 steps, there is 1 way to go from state 0 to the final state
	# (i.e. 0 as well). Any other state won't take you there in 0 steps.
	initial = [0] * len(idx2state)
	initial[0] = 1
	self.vectors = [tuple(initial)]

	def get(self, n, src=0):
	"""Number of ways to get from state src to state 0 in n steps."""
	src_idx = self.state2idx[src]
	if src_idx is None:
	# Known to be a dead end.
	return 0
	while len(self.vectors) <= n:
	# Extend by one step. Semantically we follow edges to go from i to
	# j, then use vector to count ways to get from j to 0. The index i
	# is implicit in the iteration over self.step when constructing nxt.
	prv = self.vectors[-1]
	nxt = tuple(sum(prv[j] for j in js) for js in self.step)
	self.vectors.append(nxt)
	return self.vectors[n][src_idx]

	def dot(self, dot_name):
	"""Write graph to a graphviz dot file, for visualization."""
	with open(dot_name, "w") as f:
	print("digraph {", file=f)
	print(" node [shape=box];", file=f)
	for node in self.idx2state:
	label = [f' "{node:04x}" [label=<<TABLE CELLBORDER="0"']
	if node & PAIRMARK:
	label.append(' COLOR="#00cccc"')
	label.append(' CELLPADDING="0" CELLSPACING="1">')
	for y in range(3, -1, -1):
	label.append("<TR>")
	for x in range(4):
	label.append('<TD FIXEDSIZE="TRUE" WIDTH="7" HEIGHT="7"')
	bit = 1 << (4*y + x)
	if node & bit:
	label.append(' BGCOLOR="#000099"')
	label.append("></TD>")
	label.append("</TR>")
	label.append("</TABLE>>];")
	print("".join(label), file=f)
	for src, dst in self.edges:
	print(f' "{src:04x}" -> "{dst:04x}";', file=f)
	print("}", file=f)

	class Split:
	"""Computes counts based on central arrangement and upper/lower parts."""

	def __init__(self, centers, automaton):
	self.centers = centers
	self.automaton = automaton

	def get(self, n, debug=False):
	"""Number of ways to get a tiling with n tiles."""
	res = 0
	for nc, cs in enumerate(self.centers[:n+1]):
	if (n & 1) != (nc & 1):
	continue
	for c in cs:
	x = self.automaton.get((n - nc) >> 1, src=c)
	if debug:
	print(f"{x:3d} ways to go from 0x{c:04x} to 0 in {(n - nc) >> 1} steps")
	res += x
	return res


	def reflect_one(tile):
	"""Reflect a single bit pattern in the vertical axis."""
	return (((tile & 0x1111) << 3) \|
	((tile & 0x2222) << 1) \|
	((tile & 0x4444) >> 1) \|
	((tile & 0x8888) >> 3))

	def reflected(singles):
	"""Reflect the given single tiles in the vertical axis.

	Each returned bit pattern will either be a single tile which is its own
	mirror image, or a pair of mirror image tiles with the PAIRMARK bit set in
	addition to the bits of the tiles. The returned format is like that of
	tetrominoes, namely a quadruple of tuples of bit patterns, where the
	quadruple groups patterns based on the lowest set bit. The last two of these
	are of course empty.
	"""
	res = [[], [], [], []]
	for i, lst in enumerate(singles[:2]):
	for t in lst:
	r = reflect_one(t)
	if r == t:
	res[i].append(t)
	if r & t == 0:
	pair = t \| r \| PAIRMARK
	res[i].append(pair)
	return tuple(map(tuple, res))

	def reflective_centers():
	"""Centeral arrangements for reflection in the horizontal axis.

	This returns a quintuple of tuples. The index in the outer tuple indicates
	how many tiles are used in the central arrangements, while the value inside
	the inner tuple is a bit pattern for the resultin shape that protrudes into
	the upper half of the tiling. Shapes may occur more than once.
	"""
	# The upper halves of vertical bar and square:
	halves = [0x11, 0x22, 0x44, 0x88, 0x3, 0x6, 0xc]
	lst = []
	for i, a in enumerate(halves):
	for b in halves[i+1:]:
	if a & b == 0:
	s = a \| b
	while s & 0xf == 0xf:
	s >>= 4
	lst.append(s)
	lst.sort()
	# 0 tiles: straight line.
	# 1 tile: horizontal bar.
	# 2 tiles: as computed above.
	# 3 tiles: impossible.
	# 4 tiles: four vertical bars, straight line.
	res = ((0,), (0,), tuple(lst), (), (0,))
	print(f"reflective_centers: {tuple(map(len,res))!r}")
	return res

	def rotate_one(tile):
	"""Rotate a single 4x4 bit pattern by 180°."""
	# This is effectively a bit reversal of a 16 bit integer.
	tile = ((tile & 0xff00) >> 8) \| ((tile & 0x00ff) << 8)
	tile = ((tile & 0xf0f0) >> 4) \| ((tile & 0x0f0f) << 4)
	tile = ((tile & 0xcccc) >> 2) \| ((tile & 0x3333) << 2)
	tile = ((tile & 0xaaaa) >> 1) \| ((tile & 0x5555) << 1)
	return tile

	def rotated_centers(singles, print_transitions=False):
	"""Centeral arrangements for 180° rotation.

	This returns a quintuple of tuples. The index in the outer tuple indicates
	how many tiles are used in the central arrangements, while the value inside
	the inner tuple is a bit pattern for the resultin shape that protrudes into
	the upper half of the tiling. Shapes may occur more than once.
	"""
	singles = list(chain.from_iterable(singles))

	# We model even pairs as 6 rows high and odd pairs as 7 rows high.
	even_pairs = []
	for a in singles:
	c = rotate_one(a)
	b = c << 8
	# Drop c down to row 0, so it becomes the canonical bit pattern for that
	# tile.
	while c & 0xf == 0:
	c >>= 4
	if a > c:
	# We will generate this pair when the loop iterates over c, so by
	# continuing now we avoid adding each pair twice.
	continue
	# Move the tile up (and its image down) until it protrudes into the
	# upper half.
	while a & 0xf000 == 0:
	a <<= 4
	b >>= 4
	# Keep moving the tile further up as long as it still protrudes into the
	# lower half.
	while a & 0xfff != 0:
	if a & b == 0:
	# The tile and its mirror image don't overlap, so it's a pair.
	even_pairs.append(a \| b)
	a <<= 4
	b >>= 4

	odd_pairs = []
	for a in singles:
	c = rotate_one(a)
	b = c << 12
	while c & 0xf == 0:
	c >>= 4
	if a > c:
	continue
	# Move the tile up until it protrudes into the central row.
	while a & 0xf000 == 0:
	a <<= 4
	b >>= 4
	# This time we keep moving up until the tile still intersects the
	# central row, so one more f in that bit pattern here.
	while a & 0xffff != 0:
	if a & b == 0:
	odd_pairs.append(a \| b)
	a <<= 4
	b >>= 4

	# Identify tiles that are their own rotated copy.
	solo_centers = set()
	for a in singles:
	b = rotate_one(a)
	if a << 12 == b:
	solo_centers.add(b)
	if a << 4 == b:
	solo_centers.add(a << 8)
	solo_centers = sorted(solo_centers)

	shapes = [[0], solo_centers, [], [], []]
	# Outer loop iterates over pairs, not over i (the number of tiles). This
	# avoids counting a configuration multiple times if the same pairs get added
	# to it in different order. Here the outer loop fixes the order.
	for pairs, parity in ((even_pairs, 0), (odd_pairs, 1)):
	for pair in pairs:
	for i in range(parity, 3, 2): # (0, 2) for even, just (1,) for odd.
	for src in shapes[i]:
	if src & pair != 0:
	# pair overlaps with what's already there.
	continue
	dst = src \| pair
	shapes[i + 2].append(dst)
	if print_transitions:
	s = dst >> 12
	while s & 0xf == 0xf:
	s >>= 4
	print(f"{src:07x} → {dst:07x} ({i + 2} tiles, state {s:04x})")
	for y in range(5 + parity, -1, -1):
	row = []
	for x in range(4):
	bit = 1 << (4*y + x)
	if src & bit:
	row.append("@")
	elif pair & bit:
	row.append("#")
	else:
	row.append(" ")
	print("".join(row))

	# Take shapes of upper half of each central arrangement and store them in
	# the result.
	res = [[] for _ in shapes]
	for i, ss in enumerate(shapes):
	for s in ss:
	if i & 1 and s & 0xf000 != 0xf000:
	# Omit odd center with incomplete central row.
	continue
	s >>= 12
	while s & 0xf == 0xf:
	# Drop complete rows from state representation.
	s >>= 4
	res[i].append(s)
	print(f"rotated_centers: {tuple(map(len,res))!r}")
	return tuple(map(tuple, res))

	def fully_symmetric_centers(reflective_centers):
	"""Centeral arrangements for full rectangle symmetry.

	This is like reflective_centers, but counting only the shapes that are
	symmetric around the vertical axis as well.
	"""
	# Since squares and vertical bars are the only tiles that may get used
	# together, and they leave different heights in the resulting shape
	# protrusion, there is no danger that a symmetric shape could result from an
	# unsymmetric placement of tiles.
	return tuple(tuple(j for j in i if j == reflect_one(j)) for i in reflective_centers)

	def main():
	t = tetrominoes()
	rt = reflected(t)
	rc = reflective_centers()
	rot = rotated_centers(t)
	sym = fully_symmetric_centers(rc)

	b = Automaton(t, chain(chain.from_iterable(rc), chain.from_iterable(rot)))
	c = Automaton(rt, [0])
	d = Split(rc, b)
	e = Split(rot, b)
	f = Split(sym, c)

	print("\| {:>3s} \| {:>32s} \| {:>32s} \| {:>20s} \| {:>20s} \| {:>20s} \| {:>10s} \|".format("$N$", "$a_N$", "$b_N$", "$c_N$", "$d_N$", "$e_N$", "$f_N$"))
	sep = "-"*50
	print(f"\|-{sep:.3s}:\|-{sep:.32s}:\|-{sep:.32s}:\|-{sep:.20s}:\|-{sep:.20s}:\|-{sep:.20s}:\|-{sep:.10s}:\|")
	for i in range(51):
	bi = b.get(i)
	ci = c.get(i)
	di = d.get(i)
	ei = e.get(i)
	fi = f.get(i)
	ai = bi + ci + di + ei
	assert ai & 3 == 0, f"{ai} not divisible by 4"
	ai >>= 2
	print(f"\| {i:3d} \| {ai:32d} \| {bi:32d} \| {ci:20d} \| {di:20d} \| {ei:20d} \| {fi:10d} \|")

	if __name__ == "__main__":
	main()
No results found