Advent of Code, in Erlang: Day 22

Published Wednesday, December 22, 2021 by Bryan

Twenty-two days into Advent of Code, it finally happened: I predicted what Part 2 would require! Not that I let the prediction guide me away from an overly simple Part 1 implementation, that I knew I would throw away. But I saw it coming.

Day 22's problem is basically calculating volume of a simplified version of Constructive Solid Geometry. The simplifications are that the only operators are "union" (called "on" in this puzzle) and "difference" (called "off") - no "intersection" - and space is made of integer cubes. Part 1 limits the problem further, to a small 101x101x101-cube space.

pointmap_apply([{Dir, {XMin,XMax}, {YMin,YMax}, {ZMin,ZMax}}|Rest],
                   Cubes) ->
    NewCubes = lists:foldl(
                 case Dir of
                     on ->
                         fun turn_on/2;
                     off ->
                         fun turn_off/2
                 end,
                 Cubes,
                 [{X,Y,Z}
                  || X <- lists:seq(XMin, XMax),
                     Y <- lists:seq(YMin, YMax),
                     Z <- lists:seq(ZMin, ZMax)]),
    pointmap_apply(Rest, NewCubes);
pointmap_apply([], Cubes) ->
    Cubes.

pointmap_on(Point, Cubes) ->
    Cubes#{Point => on}.

pointmap_off(Point, Cubes) ->
    maps:remove(Point, Cubes).

pointmap_count(Cubes) ->
    maps:size(Cubes).

This implementation uses a map to track which cubes are on. The key in the map is the coordinates of the cube. If the key has a value in the map, the key is on. (I've also used the atom 'on' as the value.) If the key does not have a value in the map, the cube is off. So, the number of keys in the map is equivalent to the number of cubes that are turned on (the answer to the puzzle). It just barely handles the Part 1 space, requiring a million entries if the entire space is lit. The reason I wrote it anyway is that I recognized the opportunity to use list-comprehension permutation. That middle part, where {X,Y,Z} comes from an X <- generator, and a Y <- generator, and a Z <- generator, automatically creates every combination of the values of the three. It has never been what I have wanted any time I've tried to use it, but this time, where I want the 27 individual coordinates in x=10..12,y=10..12,z=10..12, it sure is neat.

Lim2 = puzzle22:limit_instructions(Inst2, {-50,50}, {-50,50}, {-50,50})
590784 = puzzle22:pointmap_count(puzzle22:pointmap_apply(Lim2)).

"For now, ignore cubes outside this region," said Part 1. My code did, but my eyes didn't. The instructions for the cubes outside of the reduces space weren't just far outside, they were also large. The entire initial space held just over a million cubes. The first example instruction outside that space was for 18.7 trillion cubes. And, as predicted, Part 2's task was to turn on (and off) those regions as well.

So the one-element-per-coordinate solutions isn't going to work. We need to summarize a region of cubes as being "on". That's just the region from the instruction itself, until we turn part of it off. Then we need to slice out a region. In theory, this should come naturally to me, since slicing cubes into smaller cubes is a big part of what I do when I'm not coding.

region_apply([{on, X, Y, Z}|Rest], Cubes) ->
    region_apply(Rest, [{X, Y, Z}|Cubes]);
region_apply([{off, X, Y, Z}|Rest], Cubes) ->
    region_apply(Rest, subtract({X, Y, Z}, Cubes));
region_apply([], Cubes) ->
    unique(Cubes).

region_count(Cubes) ->
    lists:foldl(fun({X, Y, Z}, Acc) ->
                        axis_size(X) * axis_size(Y) * axis_size(Z) + Acc
                end,
                0,
                Cubes).

axis_size({Min, Max}) ->
    (Max - Min) + 1.

unique([C|Rest]) ->
    [C|unique(subtract(C, Rest))];
unique([]) ->
    [].

Above is the high-level code. You can think of my storage as basically being nothing but "on" instructions, for only the regions that would still be turned on after any "off" instructions are processed. The number of "on" cubes is then the volume of all remaining regions.

I've chosen to handle on instructions by simply adding the region to the accumulated list, and then reducing to a unique set of non-overlapping regions at the end, once all instructions have been processed. This might be a poor choice, if the instructions happen to create a lot of small on-regions. That last subtraction will compare a lot of regions to each other, possibly for no reason. But, it could also be that some regions are entirely turned off by some later instruction. In that case, we didn't need to deduplicate that region anyway.

Space subtraction was a hard thing for me to get right. The woodworker in me wanted to make as few cuts as possible, to keep regions as large as possible. This would be a good idea for efficiency as well - the point of this rewrite was to reduce the number of things we needed to keep track of. Alas, the logic for planning large cuts kept getting too gnarly for me to feel confident in, so I came up with the following.

subtract2({{X1Min, X1Max}, {Y1Min, Y1Max}, {Z1Min, Z1Max}},
          {{X2Min, X2Max}, {Y2Min, Y2Max}, {Z2Min, Z2Max}}=C)
  when X1Max < X2Min; X1Min > X2Max;
       Y1Max < Y2Min; Y1Min > Y2Max;
       Z1Max < Z2Min; Z1Min > Z2Max ->
    %% not overlapping
    [C];
subtract2({{X1Min, X1Max}, {Y1Min, Y1Max}, {Z1Min, Z1Max}},
          {{X2Min, X2Max}, {Y2Min, Y2Max}, {Z2Min, Z2Max}})
  when X1Min =< X2Min, X1Max >= X2Max,
       Y1Min =< Y2Min, Y1Max >= Y2Max,
       Z1Min =< Z2Min, Z1Max >= Z2Max ->
    %% 100% overlap
    [];
subtract2({X1, Y1, Z1}=C, {X2, Y2, Z2}) ->
    Xs = split_axis(X1, X2),
    Ys = split_axis(Y1, Y2),
    Zs = split_axis(Z1, Z2),
    Splits = [ {X,Y,Z} || X <- Xs, Y <- Ys, Z <- Zs ],
    subtract(C, Splits).

split_axis({A1Min, A1Max}, {A2Min, A2Max}) ->
    As = min_compare(A1Min, A2Min) ++ max_compare(A1Max, A2Max),
    group_split(As).

group_split([AMin, AMax | Rest]) ->
    [{AMin, AMax}|group_split(Rest)];
group_split([]) ->
    [].

min_compare(P1, P2) when P1 =< P2 ->
    [P2];
min_compare(P1, P2) ->
    [P2, P1-1, P1].

max_compare(P1, P2) when P1 >= P2 ->
    [P2];
max_compare(P1, P2) ->
    [P1, P1+1, P2].

If the existing region is completely outside of the subtracting region, leave the existing region unmodified. If the existing region is completely within the subtracting region, remove the existing region completely. If the intersection is more complicated than either of those, split the existing region, and then reapply these rules.

How to split the region was the hard part to figure out. After figuring out how to choose the points on each axis, I realized what I had was this: any intersection of two rectangular solids ("cuboids" in the puzzle description) can be modeled by breaking one of the solids into at most 27 pieces: each corner (8), each edge between the corners (12), each face between the edges (6), and the center (1). The split_axis function, and the functions it calls, chooses the planes in which to dice a cube this way. It may choose fewer pieces, but 27 is the simplified model. The recursive call to subtract at the end of subtract2 finds the correct splits to dispose of. This function never recusively calls itself more than once, because the splits are guaranteed to match either the first or second clause.

Sidenote: I got to use the permutation list comprehension again!

590784 = puzzle22:region_count(puzzle22:region_apply(Lim2, [])).
460 = length(puzzle22:region_apply(Lim2, [])).

39769202357779 = puzzle22:region_count(puzzle22:region_apply(Inst2, [])).
462 = length(puzzle22:region_apply(Inst2, [])).

2758514936282235 = puzzle22:region_count(puzzle22:region_apply(Inst3, [])).
1311 = length(puzzle22:region_apply(Inst3, [])).

This implementation can track the 590,784 on cubes from the limited region of Part 1 with only 460 elements - less than 1/1000 of the map implementation. It handles the two large regions that we ignored before, with just two more elements. Part 2 included a larger example, with more regions outside the "reboot" zone, and this implementation needs only 1311 elements to track 100x more cubes.

12173 = length(puzzle22:region_apply(Inst, [])).
timer:tc(fun() -> puzzle22:region_count(puzzle22:region_apply(Inst, [])) end).
% {2786539,no spoilers}

My full-puzzle input produces 12173 distinct regions covering about the same number of cubes turned on. Just under 3sec to do it. I'll bet there's a great implementation of this in any of the common CAD packages. I'll look in OpenSCAD's source code if I need to know it some day. For now, my own source is on github.