Add 2023/21 py

Author: erikw

.vscode/settings.json

@@ -19,18 +19,24 @@
     },
 
     "cSpell.words": [
+        "adjm",
         "adventofcode",
         "BSLASH",
         "cloc",
         "csym",
         "dcount",
         "multitime",
+        "ncol",
         "ndir",
         "ndist",
+        "nidx",
         "npos",
         "nprev",
+        "nrow",
         "pdir",
         "ppos",
-        "pprev"
+        "pprev",
+        "starte",
+        "Vandermonde"
     ],
 }

2023/21/README.md

@@ -0,0 +1,12 @@
+# Advent of Code - 2023 Day 21
+Here are my solutions to this puzzle.
+
+* Problem instructions: [adventofcode.com/2023/day/21](https://adventofcode.com/2023/day/21)
+* Input: [adventofcode.com/2023/day/21/input](https://adventofcode.com/2023/day/21/input)
+
+Fetch input by exporting `$AOC_SESSION` in your shell and then:
+```bash
+curl -OJLsb session=$AOC_SESSION adventofcode.com/2023/day/21/input
+```
+
+or run `fetch_input.sh` in this directory.

2023/21/fetch_input.sh

@@ -0,0 +1,5 @@
+#!/usr/bin/env sh
+# Make sure $AOC_SESSION is exported before running this script.
+
+curl --remote-name --remote-header-name --silent --fail -A 'https://erikw.me/contact' --cookie "session=$AOC_SESSION" "https://adventofcode.com/2023/day/21/input"
+test "$?" -eq 0 && echo "Fetched input" && exit 0 || echo "Failed to fetch input" && exit 1

2023/21/input1.0

@@ -0,0 +1,11 @@
+...........
+.....###.#.
+.###.##..#.
+..#.#...#..
+....#.#....
+.##..S####.
+.##..#...#.
+.......##..
+.##.#.####.
+.##..##.##.
+...........

2023/21/instructions.url

@@ -0,0 +1,2 @@
+[InternetShortcut]
+URL = https://adventofcode.com/2023/day/21

2023/21/output1.0

@@ -0,0 +1 @@
+16

2023/21/output2.0

@@ -0,0 +1,73 @@
+#!/usr/bin/env python3
+import fileinput
+from pprint import pprint
+from collections import deque
+
+# STEPS = 1
+# STEPS = 6
+STEPS = 64
+
+SYM_START = "S"
+SYM_PLOT = "."
+SYM_ROCK = "#"
+
+NEIGHBOR_DELTAS = [
+     (-1, 0),
+     (1, 0),
+     (0, -1),
+     (0, 1),
+]
+
+def read_map():
+    map = []
+    pos_start = None
+    for i, line in enumerate(fileinput.input()):
+        row = list(line.rstrip("\n"))
+        try:
+            col_start = row.index(SYM_START)
+        except ValueError:
+            pass
+        else:
+             row[col_start] = SYM_PLOT
+             pos_start = (i, col_start)
+        map.append(row)
+    return map, pos_start
+
+def print_map(map):
+     for row in map:
+          print("".join(row))
+
+def positions_steps_away(map, pos_start, steps_target):
+    que = deque([(*pos_start, 0)])
+    visited = set()
+    pos_targets = set()
+
+    while que:
+        elem = que.popleft()
+        if elem in visited:
+            continue
+        visited.add(elem)
+        row, col, steps = elem
+
+        if steps == steps_target:
+            pos_targets.add((row, col))
+            continue
+
+        for dr, dc in NEIGHBOR_DELTAS:
+            nrow, ncol = row + dr, col + dc
+            if not (0 <= nrow < len(map) and 0 <= ncol < len(map[0]) and map[nrow][ncol] == SYM_PLOT):
+                continue
+            que.append((nrow, ncol, steps + 1))
+
+    return len(pos_targets)
+
+
+
+
+def main():
+    map, pos_start = read_map()
+    garden_positions = positions_steps_away(map, pos_start, STEPS)
+    print(garden_positions)
+
+if __name__ == '__main__':
+	main()

2023/21/part1.py

@@ -0,0 +1,76 @@
+#!/usr/bin/env python3
+import fileinput
+from pprint import pprint
+import scipy.sparse as sparse
+
+# STEPS = 1
+# STEPS = 6
+STEPS = 64
+
+SYM_START = "S"
+SYM_PLOT = "."
+SYM_ROCK = "#"
+
+NEIGHBOR_DELTAS = [
+     (-1, 0),
+     (1, 0),
+     (0, -1),
+     (0, 1),
+]
+
+def read_map():
+    map = []
+    pos_start = None
+    for i, line in enumerate(fileinput.input()):
+        row = list(line.rstrip("\n"))
+        try:
+            col_start = row.index(SYM_START)
+        except ValueError:
+            pass
+        else:
+             row[col_start] = SYM_PLOT
+             pos_start = (i, col_start)
+        map.append(row)
+    return map, pos_start
+
+def print_map(map):
+     for row in map:
+          print("".join(row))
+
+def adjacency_matrix(map):
+    num_poss = len(map) * len(map[0])
+    # adj = sparse.csr_matrix((num_poss, num_poss))
+    adj = sparse.dok_matrix((num_poss, num_poss))
+    for row in range(len(map)):
+        for col in range(len(map[0])):
+            if map[row][col] != SYM_PLOT:
+                continue
+            for dr, dc in NEIGHBOR_DELTAS:
+                nrow, ncol = row + dr, col + dc
+                if 0 <= nrow < len(map) and 0 <= ncol < len(map[0]) and map[nrow][ncol] == SYM_PLOT:
+                   adj_idx =  row * len(map[0]) + col
+                   adj_nidx = nrow * len(map[0]) + ncol
+                   adj[adj_idx, adj_nidx] = 1
+                   adj[adj_nidx, adj_idx] = 1
+    return adj
+
+
+
+
+def main():
+    map, pos_start = read_map()
+
+    # With an adjacency matrix A and then A^k:
+    # cell (x,y) in A^k means the number of ways we can end up in y starting form x with k steps.
+    # Thus we are looking for sum(A^k[start_plot]), meaning the sum of all positions we can end up at after k steps starting from the start plot.
+    # Ref: https://math.stackexchange.com/a/207835/147723
+
+    map_adjm = adjacency_matrix(map)
+    map_adjm_pow = sparse.linalg.matrix_power(map_adjm, STEPS)
+    adj_start_idx = pos_start[0] * len(map[0]) + pos_start[1]
+
+    garden_positions = map_adjm_pow.count_nonzero(axis=1)[adj_start_idx]
+    print(garden_positions)
+
+if __name__ == '__main__':
+	main()

2023/21/part1_matrix.py

@@ -0,0 +1,109 @@
+#!/usr/bin/env python3
+import fileinput
+from collections import deque
+from math import ceil
+import numpy as np
+
+STEPS = 26501365
+
+SYM_START = "S"
+SYM_PLOT = "."
+SYM_ROCK = "#"
+
+NEIGHBOR_DELTAS = [
+     (-1, 0),
+     (1, 0),
+     (0, -1),
+     (0, 1),
+]
+
+
+def read_map():
+    map = []
+    pos_start = None
+    for i, line in enumerate(fileinput.input()):
+        row = list(line.rstrip("\n"))
+        try:
+            col_start = row.index(SYM_START)
+        except ValueError:
+            pass
+        else:
+             row[col_start] = SYM_PLOT
+             pos_start = (i, col_start)
+        map.append(row)
+    return map, pos_start
+
+
+def positions_steps_away(map, pos_start, steps_target):
+    que = deque([(*pos_start, 0)])
+    visited = set()
+    pos_targets = set()
+
+    while que:
+        elem = que.popleft()
+        if elem in visited:
+            continue
+        visited.add(elem)
+        row, col, steps = elem
+
+        if steps == steps_target:
+            pos_targets.add((row, col))
+            continue
+
+        for dr, dc in NEIGHBOR_DELTAS:
+            nrow, ncol = row + dr, col + dc
+            if not (0 <= nrow < len(map) and 0 <= ncol < len(map[0]) and map[nrow][ncol] == SYM_PLOT):
+                continue
+            que.append((nrow, ncol, steps + 1))
+
+    return len(pos_targets)
+
+
+# Expand map to cover the max number steps in one straight line from the start position (assuming start position in the middle of the grid).
+def expand_map(map, pos_start, max_steps):
+    times = ceil(max_steps / (len(map)/2)) 
+    mape = []
+    for _ in range(times):
+        for row in map:
+            mape.append(times * row)
+    return mape, (pos_start[0] + len(map) * (times//2), pos_start[1] + len(map[0]) * (times//2))
+
+def main():
+    map, pos_start = read_map()
+
+    # h/t https://www.reddit.com/r/adventofcode/comments/18nevo3/comment/keao6r4/?utm_source=share&utm_medium=web3x&utm_name=web3xcss&utm_term=1&utm_content=share_button
+    # As STEPS = 202300 * len(map) + 65
+    # or: 26501365 = 202300 * 131 + 65
+    n = STEPS // len(map) # 202300
+    r = STEPS % len(map) # 65
+    # it suggests a recurrence on period 131 * n + 65
+
+    # As of the visual shape of the input data, it looks like we could fit a quadratic function to it y = ax^2 + bx + c
+    # To fit a function, we can solve with linear matrix equations with  A * x = b <=> x = A^-1 * b where A is coefficient matrix.
+
+    # First get some data points for x=[0, 1, 2]
+    ans = []
+    for x in range(3):
+        steps = x * len(map) + r
+        mape, pos_starte = expand_map(map, pos_start, steps)
+        ans.append(positions_steps_away(mape, pos_starte, steps))
+    b = np.array(ans)
+
+    # A is the Vandermonde matrix with each row the geometric progression (regression?) [x^2 x^1 x^0] for x:s (0, 1, 2) (same x's used to calculate y1, y2, y3 above)
+    # The vandermonde matrix usually goes in the other direction with each row [x^0 x^1 x^2]. Having it reverse just makes the solution vector x's elements appear in the same order we want to use it.
+    # See 'polynomial interpolation' at https://en.wikipedia.org/wiki/Vandermonde_matrix#Applications
+    # [ 0^2 0^1 0^0 ]   [ 0 0 1 ]
+    # [ 1^2 1^1 1^0 ] = [ 1 1 1 ]
+    # [ 2^2 2^1 2^0 ]   [ 4 2 1 ]
+    A = np.matrix([[0, 0, 1], [1, 1, 1], [4, 2, 1]])
+
+    # Solve for solution vector x!
+    x = np.linalg.solve(A, b).astype(np.int64)
+
+    # Now we have found the coefficients for the quadratic function, namely x.
+    # Thus we can calculate the number of garden_positions for the give period n we want (according to earlier reasoning)
+    garden_positions = x[0] * n**2 + x[1] * n + x[2]
+    print(garden_positions)
+
+if __name__ == '__main__':
+	main()

2023/21/part2.py

@@ -1,7 +1,9 @@
 # Runtime Dependencies
 graphviz
+numpy
 python-mermaid
 termcolor
+scipy
 z3-solver
 
 # Development Tools