diff --git a/.envrc b/.envrc
new file mode 100644
index 0000000..3550a30
--- /dev/null
+++ b/.envrc
@@ -0,0 +1 @@
+use flake
diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000..c6f8b46
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,3 @@
+.direnv
+*.smt
+__pycache__
diff --git a/README b/README
new file mode 100644
index 0000000..61ab2d2
--- /dev/null
+++ b/README
@@ -0,0 +1 @@
+https://www.hempuli.com/blogblog/archives/2213
diff --git a/flake.lock b/flake.lock
new file mode 100644
index 0000000..124fe6c
--- /dev/null
+++ b/flake.lock
@@ -0,0 +1,127 @@
+{
+ "nodes": {
+ "flake-utils": {
+ "locked": {
+ "lastModified": 1644229661,
+ "narHash": "sha256-1YdnJAsNy69bpcjuoKdOYQX0YxZBiCYZo4Twxerqv7k=",
+ "owner": "numtide",
+ "repo": "flake-utils",
+ "rev": "3cecb5b042f7f209c56ffd8371b2711a290ec797",
+ "type": "github"
+ },
+ "original": {
+ "owner": "numtide",
+ "repo": "flake-utils",
+ "type": "github"
+ }
+ },
+ "mach-nix": {
+ "inputs": {
+ "flake-utils": [
+ "flake-utils"
+ ],
+ "nixpkgs": [
+ "nixpkgs"
+ ],
+ "pypi-deps-db": [
+ "pypi-deps-db"
+ ]
+ },
+ "locked": {
+ "lastModified": 1654084003,
+ "narHash": "sha256-j/XrVVistvM+Ua+0tNFvO5z83isL+LBgmBi9XppxuKA=",
+ "owner": "DavHau",
+ "repo": "mach-nix",
+ "rev": "7e14360bde07dcae32e5e24f366c83272f52923f",
+ "type": "github"
+ },
+ "original": {
+ "owner": "DavHau",
+ "ref": "3.5.0",
+ "repo": "mach-nix",
+ "type": "github"
+ }
+ },
+ "nixpkgs": {
+ "locked": {
+ "lastModified": 1645630103,
+ "narHash": "sha256-9fAKdW8ijvx1LIobMDnoyHwbE6T/cqQ/b4OTZ+idOHE=",
+ "owner": "nixos",
+ "repo": "nixpkgs",
+ "rev": "3e7847061787ddcb8c33abc035f6effb801f05c9",
+ "type": "github"
+ },
+ "original": {
+ "owner": "nixos",
+ "repo": "nixpkgs",
+ "type": "github"
+ }
+ },
+ "nixpkgsPy36": {
+ "locked": {
+ "lastModified": 1601475821,
+ "narHash": "sha256-7AI8j/xq5slauMGwC3Dp2K9TKDyDtBXBebeyWsE9euE=",
+ "owner": "NixOS",
+ "repo": "nixpkgs",
+ "rev": "b4db68ff563895eea6aab4ff24fa04ef403dfe14",
+ "type": "github"
+ },
+ "original": {
+ "id": "nixpkgs",
+ "rev": "b4db68ff563895eea6aab4ff24fa04ef403dfe14",
+ "type": "indirect"
+ }
+ },
+ "pypi-deps-db": {
+ "inputs": {
+ "mach-nix": [
+ "mach-nix"
+ ],
+ "nixpkgs": [
+ "nixpkgs"
+ ],
+ "nixpkgsPy36": "nixpkgsPy36",
+ "pypiIndex": "pypiIndex"
+ },
+ "locked": {
+ "lastModified": 1648888082,
+ "narHash": "sha256-CIh1Qal1zHxqZHdbICk7bWzDj+k1//AZU+OUvvDssuw=",
+ "owner": "DavHau",
+ "repo": "pypi-deps-db",
+ "rev": "edfdd998cd2ab0ff4901e47eca89ac0363b9c440",
+ "type": "github"
+ },
+ "original": {
+ "owner": "DavHau",
+ "repo": "pypi-deps-db",
+ "type": "github"
+ }
+ },
+ "pypiIndex": {
+ "flake": false,
+ "locked": {
+ "lastModified": 1648875522,
+ "narHash": "sha256-smpS9pU1Z7dxYIcsTXj3BuNKTeoPW6RNtRS+m3TbDlw=",
+ "owner": "davhau",
+ "repo": "nix-pypi-fetcher",
+ "rev": "58ec5b020f8996d713297a937d31d37c0d5782d7",
+ "type": "github"
+ },
+ "original": {
+ "owner": "davhau",
+ "repo": "nix-pypi-fetcher",
+ "type": "github"
+ }
+ },
+ "root": {
+ "inputs": {
+ "flake-utils": "flake-utils",
+ "mach-nix": "mach-nix",
+ "nixpkgs": "nixpkgs",
+ "pypi-deps-db": "pypi-deps-db"
+ }
+ }
+ },
+ "root": "root",
+ "version": 7
+}
diff --git a/flake.nix b/flake.nix
new file mode 100644
index 0000000..78f8367
--- /dev/null
+++ b/flake.nix
@@ -0,0 +1,53 @@
+{
+ description = "A very basic flake";
+
+ inputs = {
+ nixpkgs.url = "github:nixos/nixpkgs";
+ flake-utils.url = "github:numtide/flake-utils";
+
+ # This section will allow us to create a python environment
+ # with specific predefined python packages from PyPi
+ pypi-deps-db = {
+ url = "github:DavHau/pypi-deps-db";
+ inputs.nixpkgs.follows = "nixpkgs";
+ inputs.mach-nix.follows = "mach-nix";
+ };
+ mach-nix = {
+ url = "github:DavHau/mach-nix/3.5.0";
+ inputs.nixpkgs.follows = "nixpkgs";
+ inputs.flake-utils.follows = "flake-utils";
+ inputs.pypi-deps-db.follows = "pypi-deps-db";
+ };
+ };
+
+ outputs = { self, nixpkgs, flake-utils, mach-nix, ... }@attr:
+ flake-utils.lib.eachDefaultSystem (system:
+ let
+ pkgs = import nixpkgs { inherit system; };
+
+ # create a custom python environment
+ myPython = mach-nix.lib.${system}.mkPython {
+ # specify the base version of python you with to use
+ python = "python310";
+
+ requirements = ''
+ setuptools
+ numpy
+ z3-solver
+ '';
+ _.z3-solver = {
+ nativeBuildInputs.add = with pkgs; [ cmake ];
+ dontUseCmakeConfigure = true;
+ };
+ };
+ in {
+ devShell = pkgs.mkShell {
+ nativeBuildInputs = [
+ # Now you can use your custom python environemt!
+ # This should also work for `buildInputs` and so on!
+ myPython
+ ];
+ };
+ }
+ );
+}
diff --git a/solve.py b/solve.py
new file mode 100644
index 0000000..02b7cba
--- /dev/null
+++ b/solve.py
@@ -0,0 +1,204 @@
+import logging
+from z3 import *
+from tools import *
+from time import time
+from itertools import islice, accumulate
+import numpy as np
+
+logging.basicConfig(level=logging.DEBUG)
+
+SOLUTION_COUNT = None # None for all solutions
+WIDTH, HEIGHT = 5, 5
+SIZE = WIDTH * HEIGHT
+print(f"WIDTH = {WIDTH}, HEIGHT = {HEIGHT}")
+coordinate, cell_number = coordinate_l(WIDTH), cell_number_l(WIDTH)
+
+CONSTANTS = {
+ # (0, 0): False,
+ # (WIDTH - 1, HEIGHT - 1): True,
+ # (0, HEIGHT - 1): True,
+ # (WIDTH - 1, 0): True,
+ # (4, 1): 6,
+ # (2, 2): 5,
+ # (1, 3): 3,
+ # (3, 4): 3,
+ # (0, 4): 3,
+ # (2, 4): 1,
+ # (4, 4): 4,
+}
+
+s = Solver()
+ts = time()
+
+logging.debug("defining value spaces")
+t = time()
+grid = np.empty((WIDTH, HEIGHT), dtype=ExprRef)
+for index in np.ndindex(*grid.shape):
+ grid[index] = Bool(f"cell_{'_'.join(str(v) for v in index)}")
+logging.debug(f"{time() - t:f} seconds")
+
+logging.debug("constraining: no 2x2 shaded")
+t = time()
+# 2x2 shaded cells are not allowed
+for x in range(0, WIDTH - 1):
+ for y in range(0, HEIGHT - 1):
+ s.add(Not(And([grid[x][y], grid[x][y + 1], grid[x + 1][y], grid[x + 1][y + 1]])))
+logging.debug(f"{time() - t:f} seconds")
+
+logging.debug("constructing: view matrix")
+t = time()
+right, down = np.array([1, 0]), np.array([0, 1])
+cardinals = [right, down, -right, -down]
+edge = np.asarray(grid.shape)
+view = np.empty(grid.shape, dtype=ExprRef)
+for index in np.ndindex(*view.shape):
+ start = np.asarray(index)
+ visible = []
+ for direction in cardinals:
+ ray = []
+ a = 1
+ p = start + a * direction
+ while min(p) >= 0 and min(edge - p) > 0:
+ ray.append(p)
+ a = a + 1
+ p = start + a * direction
+ if len(ray) == 0:
+ continue
+ cells = [grid[tuple(p)] for p in ray]
+ visible_in_direction = list(accumulate(cells, And))
+ visible.extend(visible_in_direction)
+ view[index] = Sum([If(cell, 1, 0) for cell in visible])
+logging.debug(f"{time() - t:f} seconds")
+
+logging.debug("constraining: CONSTANTS")
+# hard set these coordinates
+for key, value in CONSTANTS.items():
+ t = time()
+ if key < (0, 0) or key >= (WIDTH, HEIGHT):
+ raise ValueError(f"Constant key '{key}' is outside of board range")
+ if type(value) is bool:
+ logging.debug(f"constraining: grid at {key} is {value}")
+ s.add(grid[key] == value)
+ elif type(value) is int:
+ if 0 > value or value > WIDTH + HEIGHT - 2:
+ raise ValueError(f"Invalid constant number for {value} in CONSTANTS at key {key}")
+ logging.debug(f"constraining: view at {key} is {value}, grid at {key} is {False}")
+ s.add(view[key] == value)
+ s.add(Not(grid[key]))
+ else:
+ raise ValueError(f"Invalid constant value type '{type(value)}'")
+ logging.debug(f"{time() - t:f} seconds")
+
+ logging.debug("checking sat")
+ t = time()
+ sat_result = s.check()
+ logging.debug(f"{time() - t:f} seconds")
+
+ if sat_result == unsat:
+ print("The latest constraint caused an unsat D:")
+ exit(1)
+
+logging.debug("constructing: adjacency matrix")
+t0 = time()
+adjacency = np.empty((SIZE - 2, SIZE, SIZE), dtype=ExprRef)
+# the adjacency matrix adjacency[0][i][j] equals 1 when cell#i and cell#j in grid are shaded and connected, otherwise 0
+for index in np.ndindex(*adjacency[0].shape):
+ i, j = index
+ difference = abs(i - j)
+ cardinal_neighbors = difference == WIDTH or (
+ difference == 1 and not abs(i % WIDTH - j % WIDTH) != 1)
+ if cardinal_neighbors:
+ adjacency[0][i][j] = And(grid[coordinate(i)], grid[coordinate(j)])
+ else:
+ adjacency[0][i][j] = False
+logging.debug(f"{time() - t:f} seconds")
+
+logging.debug("constructing: adjacency_k")
+t0 = time()
+
+# powers of The Adjacency Matrix
+for k in range(1, adjacency.shape[0]):
+ mat_mul = z3_bool_mat_mul(adjacency[0], adjacency[k - 1])
+ for index in np.ndindex(*adjacency[k].shape):
+ adjacency[k][index] = mat_mul(*index)
+logging.debug(f"{time() - t:f} seconds")
+
+logging.debug("constructing: sum_adjacency_k")
+t = time()
+
+# matrix for the sum of all adjacency^k
+mat_sum = z3_bool_mat_sum(adjacency)
+adjacency_k_sum = np.empty(adjacency[0].shape, dtype=ExprRef)
+for index in np.ndindex(*adjacency_k_sum.shape):
+ adjacency_k_sum[index] = mat_sum(*index)
+logging.debug(f"{time() - t:f} seconds")
+
+logging.debug("constraining: sum of adjacency^k is nonzero for shaded cell pairs")
+t = time()
+# constrain the sum of adjacency^k for k in [1..SIZE-1], is positive for all shaded cells
+for index in np.ndindex(*adjacency_k_sum.shape):
+ i, j = index
+ nonzero = adjacency_k_sum[index]
+ shaded = And(grid[coordinate(i)], grid[coordinate(j)])
+ s.add(shaded == nonzero)
+t1 = time()
+logging.debug(f"{t1 - t:f} seconds")
+logging.debug(f"constructing constraints total time: {t1 - ts:f} seconds")
+
+logging.debug("finished constraining puzzle rules")
+logging.debug("exporting solver assertions to file:")
+t = time()
+file_name = f"{WIDTH}_{HEIGHT}.smt"
+with open(file_name, "w") as f:
+ f.write(s.to_smt2())
+logging.debug(f"{time() - t:f} seconds")
+
+logging.debug("checking sat")
+t = time()
+# are we SAT?
+sat_result = s.check()
+t1 = time()
+logging.debug(f"{t1 - t:f} seconds")
+logging.debug(f"Total time: {t1 - ts:f} seconds")
+
+if sat_result == unsat:
+ print("We are not SAT D:")
+ exit(1)
+
+for key, value in CONSTANTS.items():
+ if type(value) is bool:
+ value = ('' if value else 'un') + 'shaded'
+ print(f"{key} will be {value}")
+
+free_terms = [grid[index] for index in np.ndindex(*grid.shape) if index not in CONSTANTS]
+solutions = all_smt(s, free_terms)
+# for i, m in enumerate(solutions, start=1):
+# print(i)
+# exit(0)
+sl = islice(solutions, SOLUTION_COUNT)
+
+t = time()
+for i, m in enumerate(sl, start=1):
+ logging.debug(f"{time() - t:f} seconds")
+ m: ModelRef
+ eval_bool_func = np.vectorize(lambda expr: is_true(m.eval(expr, model_completion=True)))
+ eval_int_func = np.vectorize(lambda expr: m.eval(expr, model_completion=True).as_long())
+ shading = eval_bool_func(grid)
+ numbers = eval_int_func(view)
+ cell_display_func = np.vectorize(cell_display_l(shading, numbers))
+
+ constant_numbers = np.empty(grid.shape, dtype=str)
+ for index in np.ndindex(*constant_numbers.shape):
+ constant_numbers[index] = ' '
+ for key, value in CONSTANTS.items():
+ if type(value) is int:
+ constant_numbers[key] = str(value)
+
+ shaded_or_constant_number_display = np.vectorize(shaded_or_display_l(shading, constant_numbers))
+ cells = np.fromfunction(shaded_or_constant_number_display, grid.shape, dtype=int)
+ print(f"Solution #{i}:")
+ mat_display(cells)
+ print()
+ t = time()
+if not next(solutions, None):
+ print("No more solutions")
diff --git a/tools.py b/tools.py
new file mode 100644
index 0000000..11e962b
--- /dev/null
+++ b/tools.py
@@ -0,0 +1,84 @@
+from z3 import And, Or, sat
+
+
+# implement python operators for some z3 objects
+# z3.BoolRef.__radd__ = lambda self, other: self + other
+# z3.BoolRef.__add__ = lambda self, other: And(self, other)
+# z3.BoolRef.__rmul__ = lambda self, other: self * other
+# z3.BoolRef.__mul__ = lambda self, other: Or(self, other)
+
+
+def coordinate_l(width):
+ return lambda cell_number: (cell_number % width, cell_number // width)
+
+
+def cell_number_l(width):
+ return lambda x, y: x + y * width
+
+
+def z3_bool_mat_mul(mat1, mat2):
+ return lambda i, j: simplify_or([simplify_and(mat1[i][k], mat2[k][j]) for k in range(len(mat1))])
+
+
+def simplify_and(*terms):
+ if any(t is False for t in terms):
+ return False
+ return And(terms)
+
+
+def simplify_or(terms):
+ a = [t for t in terms if t is not False]
+ if len(a) == 0:
+ return False
+ return Or(a)
+
+
+def z3_bool_mat_sum(mats):
+ return lambda *index: simplify_or([mat[index] for mat in mats])
+
+
+def mat_display(mat):
+ height = len(mat[0])
+ width = len(mat)
+ print('x', ' '.join([str(x) for x in range(width)]))
+ for y in range(height):
+ print(y, end=" ")
+ for x in range(width):
+ print(mat[x][y], end=" ")
+ print()
+
+
+def bool_display(v):
+ return '#' if v else ' '
+
+
+def cell_display_l(shading, numbers):
+ return lambda *index: "#" if shading[index] else str(numbers[index])
+ # return lambda *index: "#" if shading[index] else ' '
+
+
+def shaded_or_display_l(shading, other):
+ return lambda *index: "#" if shading[index] else other[index]
+ # return lambda *index: "#" if shading[index] else ' '
+
+
+def all_smt(s, initial_terms):
+ def block_term(s, m, t):
+ s.add(t != m.eval(t, model_completion=True))
+
+ def fix_term(s, m, t):
+ s.add(t == m.eval(t, model_completion=True))
+
+ def all_smt_rec(terms):
+ if sat == s.check():
+ m = s.model()
+ yield m
+ for i in range(len(terms)):
+ s.push()
+ block_term(s, m, terms[i])
+ for j in range(i):
+ fix_term(s, m, terms[j])
+ yield from all_smt_rec(terms[i:])
+ s.pop()
+
+ yield from all_smt_rec(list(initial_terms))
|
