Sudoku solver

[tcolgate]

This started as a port of the APL WIki sudoku solver (using DFS), but ended up deviating a fair bit. It does use the same technique outlined there, more or less. The answers for starburst and easterMonster examples are not the same as those on the APL Wiki, but seem valid. It assumes square puzzles, and I likely only works with (N^2)x(N^2) grids.

# This is zilde in APL
empty = 0 rho 0

# Pairs elements of A with elements of B, we'll use this for builings
# some arguments to calls
op A pair B = transp mix (split A) (split B)

# empty is not natively true/false in ivy
op isEmpty E = E === empty

# removes any zeros from a vector
op filterZeros M =
   (M != 0) sel M

# Common argument names
# P is a matrix containing a puzzle grid with 0s for unsolved
# R is a vector containing a row
# S is a set of possible solutions to a row or grid
# G is a matrix containing a grid (like P, but may have been transformed)
# M is any matrix, could be a regular puzzle grid, or a set of solutions

# generate an initial solution for puzzle P, with empty entries for the non-zero
# squares, and all valid numerals for zeroed, unset squares.
op possibles P =
  vs = iota max/ (rho P)
  (rho P) rho (vs empty) [1+ flatten P != 0]

# Return solutions restricted by all the values that are not already present
# in a row
op S solveNotInRow R = ,/ (not S @in R) part (not S @in R) * S

# this is a utility that zeros out all elements of R that are not present in Ns
op R solveMustBeOnceRow Ns =
  mask = Ns in@ R
  Ns * mask

# Takes solutions S, finds any solutions where one possible value only exists once
# in the set, and replaces them with that single value e.g.
# If given ((1 2 3) (1 2) (1 2) (1 2 4)), the first 3, and the last 4 must be the
# only possible values in those cells (since they *must* appear once, but only
# have one possible location, giving ((3) (1 2) (1 2) (4)).
op solveMustBeOnceRow S =
  # which items are unique in the possible answers
  0 == +/ rho flatten S : S
  Sdims = iota rho S
  Ns = (1 == +/% (flatten S) o.== Sdims) sel Sdims
  (+/ rho Ns) == 0: S
  newVals = filterZeros@ unique@ split transp mix S solveMustBeOnceRow@ Ns
  indices = filterZeros (iota rho S) * empty !==@ newVals
  S[indices] = newVals[indices]
  S

op solveRow (R S) =
  RS = S @solveNotInRow R
  RS = solveMustBeOnceRow RS
  RS

op solveGrid (G S) = mix solveRow@ G pair S

op M subgridRow (W X Y) =
  ,/ ,@ M[-1 + X + iota W][;-1 + Y + iota W]

op subgrid M =
  D = rho M
  SD = sqrt D[1]
  sgs = SD,@ 1 + SD * ,/ split (-1 +  iota SD SD)
  mix M subgridRow@ sgs

transforms = 1 2 3
# 1 - no transform
# 2 - transp rows to colums
# 3 - transp subgrids to rows

op T transform M =
  T == 1: M
  T == 2: transp M
  T == 3: subgrid M

op T solveGridWithTransform (G S) =
  TM = T transform G
  TS = T transform S
  T transform solveGrid (TM TS)

op solvePuz G =
  S = possibles G
  S = 1 solveGridWithTransform (G S)
  S = 2 solveGridWithTransform (G S)
  3 solveGridWithTransform (G S)

op P mergeFound S =
  mask = 1 == count@@ S
  P[where mask] = (first@@ S)[where mask]
  P

op isEmpty V = empty === V

op (M X Y) setP (V) =
  M[X;Y] = V
  M

op trySolve P

# Try solving P with all values of VS at position X Y
op (P X Y) tryAll VS =
  # if there are no options, we fail
  empty === VS : (P 0)
  (Pnew DONE) = trySolve (P X Y) setP 1 take VS
  # if it worked, we are done
  DONE : (Pnew DONE)
  (P X Y) tryAll 1 drop VS

# trySolve will return a vector with last found solution a true/false if solved
op trySolve P =
  not 0 in P: (P 1) # there are no unknowns, we are done
  res = solvePuz P
  Pnew = P mergeFound res
  0 == +/ +/ not@@ isEmpty@@ res : (P 0) # We have unsolved squares, but no options
  P !== Pnew : (trySolve Pnew) # We made some progress, let us go again
  # not making progress, but P has unknowns, we'll start searching
  # find result with the lowest number of options
  counts = count@@ res
  fewestCount = min/ filterZeros flatten counts
  fewestIndex = first where fewestCount == counts
  opts = first res[box fewestIndex]
  ((box P), fewestIndex) tryAll opts

op solve P =
  (ANSWER DONE) = trySolve P
  DONE : ANSWER
  print "could not solve, printing last state"
  ANSWER


s22 = mix (0 2 3 4) (3 0 0 0) (2 0 0 0) (4 0 0 0)
s33 = mix (0 0 1 6 9 0 5 0 0) (4 0 0 2 7 0 0 0 1) (0 7 0 0 0 0 0 9 0) (0 0 0 0 0 0 0 3 0) (0 0 0 4 3 0 0 0 7) (0 0 0 7 8 0 6 0 0) (0 0 6 0 0 0 8 0 0) (0 2 0 1 4 0 0 6 0) (0 1 0 3 5 0 0 4 0)
puzzle =4 4 rho 1 2 0 0 0 4 0 0 2 0 0 0 0 1 0 3
starburst = 9 9 rho 0 8 0 0 6 0 0 7 0 0 0 0 0 0 0 0 0 8 4 0 0 0 0 0 0 0 1 8 7 0 0 0 0 0 3 0 3 0 0 0 0 0 0 0 7 7 0 0 0 0 0 0 0 3 0 3 0 0 7 0 0 8 0 1 0 8 2 3 9 7 0 4
easterMonster = 9 9 rho 1 0 0 0 0 0 0 0 2 0 9 0 4 0 0 0 5 0 0 0 6 0 0 0 7 0 0 0 5 0 9 0 3 0 0 0 0 0 0 0 7 0 0 0 0 0 0 0 8 5 0 0 4 0 7 0 0 0 0 0 6 0 0 0 3 0 0 0 9 0 8 0 0 0 2 0 0 0 0 0 1

I don't think it's optimal in any way. I seem to spend a lot of time folding and unfolding the matricies. I've just realised that f@@ for ranks is a thing, that helps in a couple of places. Also updated to used the recent where and indexed assignment improvements

[tcolgate]

I had a go at implementing the Position Vector approach shown here. I've editte this to include the optimisation to pass the contention map from the top level. It seems to work, but only really for simpler puzzles. It has the "downside" of calculating all possible solutions, which is pretty catastrophic for the starburst and easterMonster puzzles above.
I may be missing something int he original position that makes the calculation exit early somehow (but if so, I really don't understand how that can be the case)

# sbox calculate the subgrid boxes for an NxN puzzle
op sbox N = transp N sel transp N sel (N N) rho iota N ** 2

# rcb creates a grid with (ROW COLUMN SUBGRID) tupples for each position
# in a puzzle
op rcb DIM = flatten@@ ,@@ transp mix (iota DIM) (sbox first DIM ** 1/2)

# cmap calculates the contention map for a puzzle, each entry has a 1 in each square that contend with the item and the row/column position in the map
op cmap M = box@@ 1 in@@@@ M o.== M

# avl calculates all the possible values for position X Y in puzzle P. CMAP and POSSIBLES
# are the precalculated contention map and all possible numeric values in a grid
op (P CMAP POSSIBLES)  avl (X Y) =
  (not POSSIBLES in (P * (first CMAP[X ; Y]))) sel POSSIBLES

# set sets position X Y in matrix M to value V
op (X Y V) set M =
  M + (rho M) take (-1 * X Y) take V

# pvec gives a vector of all possible grid P for a given ampty position X Y
op (X Y) pvec (CMAP POSSIBLES P) =
  ((X Y) ,@ (((P CMAP POSSIBLES) avl X Y))) @set P

# Not sure if this is quite equivalent to the bind/curry approach in the video
# ,/ 4 3 pvec@ ,/ 4 2 pvec@ ,/ 3 4 pvec@ ,/ 3 2 pvec@ ,/ 2 4 pvec@ ,/ 2 3 pvec@ ,/ 4 4 pvec@ ,/ 3 3 pvec@ ,/ 2 2 pvec@ 1 1  pvec s22

# pvex applies pvec for every empty position in the EV, the vector of empty positions
op EV pvex (CMAP POSSIBLES P) =
  (count EV) == 1 : (flatten -1 take EV) pvec (CMAP POSSIBLES P)
  res = ((-1 drop EV) pvex (CMAP POSSIBLES P))
  args = (CMAP POSSIBLES) ,@ box@ res
  ,/ (flatten -1 take EV) pvec@ args

# emt calculates a vector of empty (zeroed) position in puzzle P
op emt P =
  (,P == 0)  sel , iota rho P

# solves puzzle P
op solve P =
  POSSIBLES = iota first rho P
  CMAP = cmap rcb rho P
  (emt P) pvex (CMAP POSSIBLES P)

[tcolgate]

A couple of thoughts after having gone through these:

• )debug trace is nice, but often does quite get you "line level" debug. debug/break/step would be great
• native zilde / () support would be nice. A native empty that can count as false, or even have 0 1 0 0 sel v smartly fill with empties if v is a vector of vectors, if that makes sense.
• can the V[where x] = V[where x] assignment work x and v are vectors of vectors? I can't quite see why it wouldn't be possible.
• first class funcs, currying and function/operator literals would be fantastic!

[tcolgate]

Things didn't seem complete without a printer. This is a bit crude but:

op padOuter (R TR) = 
  TR[2] , (first R) , TR[4]

op SGW padInner (R TR) = 
  pRow = (SGW sel iota SGW) part R
  flatten (pRow[iota SGW - 1] @, TR[3]), box pRow[SGW]

op SGW padLine (R TR) =
  themedInner = SGW padInner (R TR)
  args = (box themedInner) TR
  padOuter args

op show_sudoku P =
  theme = 4 4 rho ' ┃│┃━┏┯┓─┠┼┨━┗┷┛'
  cols = first rho P
  sgw = sqrt cols
  grid = text@@ P
  empties = first@@ grid == "0"
  grid[where empties] = ((rho P) rho theme[1;1])[where empties]
  lines = split grid
  args = (box@ lines) @, box theme[1]
  mid = sgw padLine@ args
  midPad = box sgw padLine ((cols rho theme[3;1]) theme[3])
  midBoxes = (sgw sel iota sgw) part mid
  mid = (,/ (midBoxes[iota sgw - 1] @, midPad)), midBoxes[sgw]
  top = box sgw padLine ((cols rho theme[2;1]) theme[2])
  bot = box sgw padLine ((cols rho theme[4;1]) theme[4])
  mix top, mid , bot

giving

show_sudoku easterMonster
┏━━━┯━━━┯━━━┓
┃1  │   │  2┃
┃ 9 │4  │ 5 ┃
┃  6│   │7  ┃
┠───┼───┼───┨
┃ 5 │9 3│   ┃
┃   │ 7 │   ┃
┃   │85 │ 4 ┃
┠───┼───┼───┨
┃7  │   │6  ┃
┃ 3 │  9│ 8 ┃
┃  2│   │  1┃
┗━━━┷━━━┷━━━┛