Solving Advent of Code 2020 day 17 by not solving it
20210122
... in which we take a look at the day 17 problem from Advent of Code 2020, Conway Cubes. We will end up solving it eventually, but in a very roundabout way; you will see why shortly.
Along the way, we will:
 test and refactor,
 avoid hard problems, and
 write some idiomatic Python, from scratch.
There are both shorter and faster ways of solving this, but we'll take the scenic route, in a way that helps you translate those skills to your regular, nonpuzzle coding.
You can read this without having solved the problem, but to get the most from it, I recommend you do it on your own first, and then follow along as we do it again here.
Contents
The problem #
The problem might look familiar (if you didn't read it yet, go do that first):
 we have a grid of cells that can be either active or inactive
 cells change state at each step in time depending on their neighbors
These things are called cellular automata, the most famous of which is probably Conway's Game of Life. Our problem seems to have the same rules as Life:
 any active cell with 2 or 3 active neighbors remains active
 any inactive cell with 3 active neighbors becomes active
 all other cells become/remain inactive
... except the grid is 3dimensional.
Minor spoiler alert:
In part two of the problem, the grid is 4dimensional.
You might notice that this is relatively hard to visualize: the example output has 5 layers after 3 cycles in 3 dimensions, and 25 layers after only 2 cycles in 4 dimensions! Understanding what's happening in 4D seems difficult even assuming the number of planes remains 25 for future iterations.^{1}
The normal flow of solving the problem would be something like this:
 implement the 3D version; then,
 after seeing the second part, either
 do 4D in a similar fashion, or
 generalize 3D to an arbitrary number of dimensions.
We'll do the opposite: implement 2D, which is trivial to check by hand, and once we have some tests, generalize that.
Parsing the input #
First, let's parse the input.
We'll transform the input string into a list of lists of integers, where 0 and 1 represent inactive and active cubes, respectively. For the test input, it should look like this:
[[0, 1, 0], [0, 0, 1], [1, 1, 1]]
Using numbers separates concerns, but also allows for nicer code:
we can write if active
instead of if active == '#'
, and
sum(cubes)
instead of sum(active == '#' for active in cubes)
.


Later on, we'll need to turn that back into a string; let's take care of it now:


We're not bothering to minimize the output like the example does (hide empty rows/columns), since we won't be looking at more than one layer per cycle anyway.
Having both functions, we can write a basic roundtrip test:


The world #
There's more than one way to represent the world, but most straightforward should be the same we did the input, as a nested list of integers.
The grid is supposed to be infinite; there are no infinite lists in Python, so we have to cheat: we'll make a big enough world, and check if we've reached the edges.
For now, 8 cubes ought to be enough for anyone:


[[0] * x_size] * y_size
won't work, because the outer list
will contain the (same) inner list y_size times:
>>> l = [[0] * 3] * 3
>>> l
[[0, 0, 0], [0, 0, 0], [0, 0, 0]]
>>> l[0][0] = 1
>>> l
[[1, 0, 0], [1, 0, 0], [1, 0, 0]]
We'll place the input in the middle, to give it room to grow in all directions:


At this point, the script should output this:
........
........
...#....
....#...
..###...
........
........
........
Simulation #
Let's simulate one cycle.
There are many optimizations to this, but again, we'll do the straightforward thing: go through each cell, count its neighbors, and change its state according to the rules.
The cube state changes simultaneously for all cubes, so we cannot update the world inplace: given two neighboring cubes updated in order, the first cube's old state will be gone when we get to the second one. To avoid this, we will use two worlds, always reading from the old world and changing the new one.
We'll put the logic in a function that takes the two worlds as arguments:


Going through each cell:


The cell's neighbors are all the cells in a 3x3 square centered on it, except the cell itself:


old[ny][nx]
should now get us the neighbor's state;
however, there are a couple of (literal) edge cases:
 The current cell has one of the coordinates 0 (is on the top or left edge). Due to how list indexing works, the neighbors to the top/left (coordinate 1) will be cells from the opposite edge of the world. We don't want this: the world is supposed to be infinite, not wrap around.
 The current cell has one of the coordinates
size  1
(is on the bottom or right edge). We'll get an IndexError when trying to get the bottom/right neighbor.
In both cases, we want to know when the current cube is active, since it's possible one of its out of bounds neighbors should become active:


For now, the exceptions will remain unhandled, and we'll make the world bigger by hand when required. We could do something else, like catch the exception, scale the world by some constant factor, and continue simulating.
Having that, we set the cell's new state according to the rules:


Let's call that 6 times. Instead instead of creating a new world on each iteration, we can reuse the old one from the previous iteration, since we'd be throwing it away anyway.


Running the script should output this:
........
........
...#....
....#...
..###...
........
........
........
after cycle #1
........
........
........
..#.#...
...##...
...#....
........
........
after cycle #2
........
........
........
....#...
..#.#...
...##...
........
........
after cycle #3
........
........
........
...#....
....##..
...##...
........
........
after cycle #4
........
........
........
....#...
.....#..
...###..
........
........
after cycle #5
........
........
........
........
...#.#..
....##..
....#...
........
after cycle #6
........
........
........
........
.....#..
...#.#..
....##..
........
A quick visual inspection confirms our simulation is correct.
You may notice that after 4 cycles, we have the initial pattern, but moved one cell to the bottom and right; this pattern is called a glider.
To wrap up, we'll calculate the puzzle answer (the number of cubes left in the active state after the sixth cycle), and add a test to make sure we won't break anything later on:


Cleaning up #
Now that our 2D implementation is correct, we can make it more general.
We won't do that straight away, though; to make it possible to call in multiple ways and generalize incrementally, we'll split the script into functions.
Throughout these changes, our test suite should keep passing.
While working on the script, I used entr to rerun it automatically on every save:
echo conway_cubes.py  entr rc python3 conway_cubes.py
Your editor may have have a "save and run" function; if it doesn't, entr is a nifty editorindependent way of achieving the same thing.
First, let's extract the neighborcounting logic:




Then, the construction of worlds:






Then, the input copying:




Finally, let's wrap the whole solution in a function.
Again, to separate logic from presentation, we'll split this into a generator that yields new states, forever; it yields the initial state first for convenience:


... and a function that drives it for 6 cycles, printing the resulting worlds.
Since we want to use the function for more dimensions, but don't want to implement formatting for anything except 2D, we'll make printing optional through the magic of dependency injection.
Also, to keep the tests working, we'll return the final state and the number of active cubes.


Going multidimensional #
OK, we're ready now:
let's follow simulate_forever
and make it work with N dimensions.
First, we'll fix make_world
.
Instead of one size argument per dimension,
it will take a single argument, a tuple of sizes.
Using size
or sizes
as the name might be confusing
(since they are different kinds of sizes);
how about we call it shape
?
To handle the arbitrary number of dimensions, make_world
will be recursive:


We now need to bubble this change up to run
.
For convenience, it will keep its size
, but get a new dimensions
argument,
so we can have an ndimensional "square" world.






Next, input copying. For 2D, the destination was the world itself; for higher dimensions, the destination is the (2D) plane in the middle of each higher dimension:


Continuing with simulate
.
We need a generic version of the nested loop that goes over all the cells.
We'll write a version of enumerate that works with nested lists,
and instead of i, value pairs, yields (..., k, j, i), value pairs:


You may notice we stole the function name from numpy.
In the first draft of this article it was called coord_enumerate
and it generated ..., k, j, i, value tuples,
but I changed it to make things easier for people already familiar with numpy.
Unlike with numpy arrays, we have to infer the number of dimensions, since our world doesn't have an explicit shape.
Also unlike the numpy version, ours only works with lists; to make it work with other sequences, we could use something like:
isinstance(temp, collections.abc.Sequence)
We'll leave get_active_neighbors
last, and finish simulate
first:


run
also goes through all the cubes to count the active ones,
so we'll change it to use ndenumerate
as well:


Finally, let's do get_active_neighbors
.
First, we need to flatten the loop that generates neighbor coordinates.
We could use recursion, like we did with make_world
and ndenumerate
,
but the standard library already has a tool for this:
itertools.product generates the cartesian product of some iterables;
even better, it can generate the product of an iterable with itself.




Even if itertools.product did not have a repeat
argument,
we could still emulate it with argument unpacking:
product(* [(1, 0, 1)] * dimensions)
After using make_directions
,
and making the neighbor counting parts generic as well,
get_active_neighbors
becomes:


The last remaining thing is to fix its call site in simulate
:


We're done!
3D #
After 6 cycles with the test input, there should be 112 active cubes:


It turns out 8 cubes is not enough:
after cycle #0: ...
after cycle #1: ...
after cycle #2: ...
after cycle #3: ...
Traceback (most recent call last):
File "conway_cubes.py", line 161, in <module>
_, active_cubes = run(8, 3, parse_input(TEST_INPUT))
File "conway_cubes.py", line 131, in run
for cycle, world in zip(range(6 + 1), worlds):
File "conway_cubes.py", line 123, in simulate_forever
simulate(old, new)
File "conway_cubes.py", line 84, in simulate
active_neighbors = get_active_neighbors(old, active, coords)
File "conway_cubes.py", line 47, in get_active_neighbors
raise RuntimeError(f"active on edge: {coords}")
RuntimeError: active on edge: (2, 3, 0)
Good thing we're checking for that. Bumping the world size to 16 fixes it:
after cycle #0: ...
after cycle #1: ...
after cycle #2: ...
after cycle #3: ...
after cycle #4: ...
after cycle #5: ...
after cycle #6: ...
the result is 112
Time to plug the puzzle input; I assume they are different for everyone, so feel free to use yours instead.
At first I got another active on edge
error;
since running with the test input already takes about a second,
I increased the world size by smaller increments instead of doubling it;
20 is enough for my input:


the result is 386
, which unlocks the second part of the problem.
4D #
Moving on to even higher dimensions:
The 4D cube count for the test input should be 848:


... and after 10 seconds the first cycle hasn't ended yet.
Let's add some timings in run
before we let this run any further.




We use time.perf_counter because it is the "clock with the highest available resolution to measure a short duration"; since the timescales we're dealing with don't probably need it, we could have used time.monotonic (on Unix, they're actually the same).
In general, you should avoid time.time for measuring durations, because it can skip backward or forward (e.g. daylight saving time), or move slower or faster than real time (e.g. leap smearing).
After about 1 minute, I get the (correct) result for the test input.
After 2 more minutes, the result is 2276
; this is also correct.
Conclusions #
What have we learned?
Sometimes, it's easier to start with a simpler problem, and iterate towards the real one: we did most of the coding for 2D, which was easier to understand and debug; once we got that general enough, 3D and 4D just worked.
Tests make refactoring quicker. Even for "exploratory" coding like this, when you don't know exactly what the final code will look like, or even what you're trying to achieve, some parts "settle" earlier than others; once that happens, a simple highlevel test allows you to focus on the next thing without having to constantly check that stuff that was working is still working.
The 4D simulation is really slow ... exponentially slow. We'll see at least one way of improving that in a future post.
Bonus: CLI and tests #
When run, the script goes through all the tests and then prints the results for the two parts of the problem (3 and 4 dimensions).
Let's make it a bit easier to run a single simulation with a specific world size, dimension, and number of cycles; this will be useful when profiling.
First, we'll move the formatting test (line 22) and the other tests
to a new file, test_conway_cubes.py
, and into test functions;
we leave TEST_INPUT
in the module since we'll use it from the CLI as well.


We want to be able to run fewer cycles, so we add a new argument to run
:




And then remove the rest of the top level code in favor of a rudimentary CLI:


We now can run the tests using pytest.
The 4D one takes a while; we use k
to skip it:
$ pytest test_conway_cubes.py q k 'not 4'
... [100%]
3 passed, 1 deselected in 1.14s
The script is invoked like this:
$ python3 conway_cubes.py real 20 4 6
after cycle #0 (0.01s): ...
after cycle #1 (24.15s): ...
after cycle #2 (23.45s): ...
after cycle #3 (23.79s): ...
after cycle #4 (24.31s): ...
after cycle #5 (24.15s): ...
after cycle #6 (24.07s): ...
the result is 2276 (143.94s)
Bonus: more tests #
Our tests cover the most common cases, but there's something we neglected entirely: edge cases.
At the moment, that's not necessarily an issue: the code dealing with it is straightforward, and the puzzle itself is testing it for us – if we go over the edge and the script doesn't fail, we won't get the right result.
However, the code may become less straightforward after optimizing it, and we don't want to rely on the website for validation every time we make a change.
Let's fix that.
The first test checks we get errors when there are active cells on the edges:




Here, the parametrize decorator defines
4 different input
and 3 different dimensions
arguments,
so that the test_edge_errors
function will run using them in turn,
in all possible combinations.
The second test checks that the new world can have active cells on the edges, as long as the old one didn't:


Because of parametrization, we get 14 actual tests:
12 for test_edge_errors
, and 2 for test_edge_ok
:
$ pytest test_conway_cubes.py q k edge
.............. [100%]
14 passed, 4 deselected in 0.06s
The final version of the script and test file.
If you think that's a very conservative assumption, trust your instincts. ^{[return]}