Sudoku Solver

You can find the code on my github.

• Intro
• Exact Cover
• My Implementation
• My Observations
• More Sudokus
• References

Intro

This is my solution to a piece of coursework proposed by my AI module at university, to develop a sudoku solver that gives solutions as fast as possible. In order to do this, I designed a hybrid solver that decides if the sudoku in question is difficult or easy, and utilises algorithm X or backtracking respectively.

Note Although my approach uses two separate solvers, my main focus is on Algorithm X, which will be discussed here.

Exact Cover

Exact cover is defined as a collection S of subsets of a set X such that each element in X exists in exactly one subset in S*. (Wikipedia). It might be easier to understand with an example:

Let S = {A, B, C, D, E, F} be a collection of subsets of a set X such that:

• A = {1, 4, 7}
• B = {1, 4}
• C = {4, 5, 7}
• D = {3, 5, 6}
• E = {2, 3, 6, 7}
• F = {2, 7}

The exact cover of X would be S* = {B, D, F} as:

• B ∪ D ∪ F = X
  (S* covers every element in X)
• B ∩ D ∩ F = {}
  (every element in X is covered by exactly one element in S*)

Algorithm X - how to solve exact cover problems

In order to use Algorithm X, we must first create a matrix A consisting of 0s and 1s, with a goal of selecting a subset of the rows such that the digit 1 appears in each column exactly once.

Algorithm X takes an element e in A to cover, and finds a row r that covers it. This row is added to the potential solution, and every row that also covers e is removed from A along with every column that r satisfies. It then repeats this process recursively.

The algorithm is as follows:

If A is empty, the problem is solved; terminate successfully.
Otherwise:
    choose a column, c (deterministically).
    Choose a row, r, such that A[r, c] = 1 (nondeterministically).
    Include r in the partial solution.
    For each column j such that A[r, j] = 1,
        delete column j from matrix A;
        for each row i such that A[i, j] = 1,
            delete row i from matrix A.
    Repeat this algorithm recursively on the reduced matrix A.

(Knuth, 2000) - Page 4

The above problem can be represented with the matrix:

Firstly, as matrix A is not empty, the algorithm finds the column with the lowest number of 1s. This is column 1, that has 1s in rows A and B.

The algorithm firstly selects row A (but remembers row B is a possible solution).

Row A has 1s in columns 1, 4 and 7. (This is the first for loop)

Column 1 has 1s in rows A, B. Column 4 has rows in A, B, C and column 7 has 1s in rows A, C, E and F. Therefore the only row that does not have a 1 in the same column as row A is row D (this is the second for loop.).

Therefore, row D is selected and the algorithm repeats.

As the matrix is not empty, the algorithm finds the column with the lowest number of 1s. This is column 2.

As column 2 does not contain any 1s, this branch of the algorithm terminates unsuccessfully and the algorithm moves onto the next branch, which in this case would be row B.

Continuing the algorithm, we will eventually end up with:

meaning that S* = {B, D, F} is the exact cover.

If there are no remaining unsearched branches and no solution has been found, there is no exact cover.

For a step by step version of this process, please see the Wikipedia article.

>How to solve sudoku with this method

Sudoku can be represented as an exact cover problem with a matrix A of with dimensions x and y, where:

• x represents the set of values that each cell contains, stored in the form (row, column, value)
• y represents the set of the constraints that each cell must fulfil. There are 4 constraints as outlined below:
• Every cell must contain a value.
• Every row must contain exactly one of each value.
• Every column must contain exactly one of each value.
• Every box must contain exactly one of each value.

Each value in x will satisfy a specific combination of values in y. No other value in x should satisfy these same constraints - for example, if there was a 6 in row 3, then there cannot be another 6 in row 3.

Therefore, as every value in y needs to be satisfied by exactly one value in x:

This can be represented as an exact cover problem where each column is a value in x and each row is a value in y.

My implementation

Translating sudoku to an exact cover problem

The solver firstly constructs matrix_a. The code for this is as follows:

    matrix_a = {j: set() for j in(
        [intern((f"cell {i}")) for i in product (range(9), range(9)    )] +
        [intern((f"row {i}" )) for i in product (range(9), range(1, 10))] +
        [intern((f"col {i}" )) for i in product (range(9), range(1, 10))] +
        [intern((f"box {i}" )) for i in product (range(9), range(1, 10))]
    )}

    constraints = get_constraints()

    # Populate matrix A with constraints
    # Populate matrix A with constraints
    for i, consts in constraints.items():
        for j in consts:
            matrix_a[j].add(i)

get_constraints() is a function that returns the constraints dict and caches and returns the constraints dict. Shown below.

# Cache the constraints dict; this is not updated and can be re-used.
@memoize
def get_constraints() -> dict:
    '''
    Gets the general constraint dict for all sudokus.

    @returns:
        constraints (dict) : Constraints dictionary
    '''
    constraints = {}
    for row, col, cell in product(range(9), range(9), range(1, 10)):
        box = (row // 3) * 3 + (col // 3)
        # Boxes are labelled like this:
        #   0 1 2
        #   3 4 5
        #   6 7 8
        constraints[(row, col, cell)] = [
            # Each cell must have a value
            intern((f"cell ({row}, {col})")),
            # Each row must have each value
            intern((f"row ({row}, {cell})")),
            # Each column must have each value
            intern((f"col ({col}, {cell})")),
            # Each box must have each value
            intern((f"box ({box}, {cell})"))
        ]
    return constraints

This creates a matrix as follows:

Once the initial matrix has been generated, we can update the constraints to reflect the initial state - to remove the rows and columns that we know are incorrect. If our original sudoku contains the value (0,5,3) for example, then we can remove all columns that imply a different value at that location, such as (0,5,2). If the value of any given cell does not exist, the sudoku is not solvable.

This is completed with the code below:

    # Update constraints to reflect initial state
    for (row, col), cell in ndenumerate(sudoku):
        if cell != 0:
            try:
                select(matrix_a, constraints, (row, col, cell))
            except KeyError:
                # Sudoku is not solvable
                return full((9, 9), fill_value=-1)

We can then proceed onto Algorithm X .

Algorithm X

def find_solution(matrix_a, constraints, solution):
    '''
    Recursively attempts to find a solution to a given exact cover problem

    @args:
        matrix_a: The search space matrix
        constraints: the constraints dict
        solution: The state to find (or not find) soltion for

    @returns:
        list: The solution
    '''
    if not matrix_a: yield solution
    else:
        col = choose_col(matrix_a)
        for row in list(matrix_a[col]):
            solution.append(row)
            cols = select(matrix_a, constraints, row)

            # Keep trying to find solution recursively
            for i in find_solution(matrix_a, constraints, solution): yield i

            # This branch does not have a solution, deselect it
            deselect(matrix_a, constraints, row, cols)

            # Remove value from this possible solution
            solution.pop()

This is the recursive backtracking algorithm that attempts to find the solution to the given sudoku.

If matrix_a is empty, we have fulfilled all the constraints and thus solved the sudoku.

Otherwise, find the column with the fewest values in matrix_a using the algorithm below:

def choose_col(matrix_a, constraints):
    """
    Returns col with fewest possible values.

    @args
        matrix_a: the search space matrix
        constraints: The constraints dict

    @returns
        col : the column with the fewest possible values
    """

    best_col_val = float("inf")
    best_col = None

    for col in matrix_a:
        cur_col_val = len(matrix_a[col])

        if best_col_val > cur_col_val:

            best_col = col
            best_col_val = cur_col_val

            # Do not waste time if we have already found a column with only
            # one value
            if cur_col_val == 1:
                return best_col

    return best_col

This returns the selected column. Append it to the partial solution.

The associated rows and columns are then removed from the matrix:

def select(matrix_a, constraints, row) -> list:
    '''
    removes associated rows, cols from matrix
    @args:
        matrix_a: the search space matrix
        constraints: the constraints dict
        row: The row to be selected

    @returns:
        cols (list): List of removed columns
    '''
    # Keep removed columns so they can be added back into sudoku
    cols = []

    # For each constraint this row satisfies:
    for i in constraints[row]:

        # For all other constraints that also satisfy i:
        for j in matrix_a[i]:

            # For all other constraints that j satisfies:
            for k in constraints[j]:

                # Remove all other constraints except i
                if k != i:
                    matrix_a[k].remove(j)

        cols.append(matrix_a.pop(i))

    return cols

And the program recursively tries again.

If the program has exhausted all the possible solutions on a given branch, it will then deselect it and return the removed values back into matrix A:

def deselect(matrix_a, constraints, row, cols) -> None:
    '''
    Restores a branch with a no solutions back into matrix_a

    @args:
        matrix_a: The search space matrix
        constraints: the constraints dict
        cols: Columns to restore into matrix_a
    '''
    for i in reversed(constraints[row]):

        # Get top column from list of removed columns
        matrix_a[i] = cols.pop()

        # For each other value that satisfies constraint:
        for j in matrix_a[i]:

            # For other constraints that value satisfies:
            for k in constraints[j]:

                # Add value back into matrix_a
                matrix_a[k].add(j)

Translating the solved exact cover problem back to a sudoku

When find_solution has found and returned a solution, it will be in the form of a list of tuples containing the remaining values to put back into the sudoku. The solver then fills those values into the original sudoku, avoiding making a copy in order to save processing time.

    # find solution and update initial state with it
    for solution in find_solution(matrix_a, constraints, []):
        for (row, col, val) in solution:
            # Fill original values directly into input sudoku to save time
            sudoku[row][col] = val
        return sudoku

If find_solution has exhausted all possible branches, then there is no possible solution and hence the solver will return the error grid:

return full((9, 9), fill_value=-1)

My Observations

While running some tests on the solver, I noticed some rather weird behaviour, which I will do my best to document here.

When given a blank sudoku, the solver always returns the same value.

When the sudoku solver is given an array of zeros, it will consistently return the same solution, as it is the first solution it reaches. Why it is this specific solution, I am unsure. I have not found any explanation online, although I think it would be interesting to see if other implementations of Algorithm X reach the same solution when given an empty initial state.

+-------+-------+-------+
| 4 7 1 | 3 8 6 | 5 9 2 |
| 9 3 2 | 5 4 7 | 6 1 8 |
| 8 5 6 | 2 1 9 | 7 4 3 |
+-------+-------+-------+
| 2 9 3 | 1 6 8 | 4 5 7 |
| 6 8 7 | 9 5 4 | 3 2 1 |
| 1 4 5 | 7 3 2 | 8 6 9 |
+-------+-------+-------+
| 7 6 9 | 8 2 5 | 1 3 4 |
| 3 2 4 | 6 7 1 | 9 8 5 |
| 5 1 8 | 4 9 3 | 2 7 6 |
+-------+-------+-------+

I believe that it has something to do with how the compiler treats data stored in lists and dictionaries, as when the code is run using pypy3, it produces a different result when given the same test:

+-------+-------+-------+
| 1 2 3 | 4 5 6 | 7 8 9 |
| 7 8 9 | 1 2 3 | 4 5 6 |
| 4 5 6 | 7 8 9 | 1 2 3 |
+-------+-------+-------+
| 3 1 2 | 8 4 5 | 9 6 7 |
| 6 9 7 | 3 1 2 | 8 4 5 |
| 8 4 5 | 6 9 7 | 3 1 2 |
+-------+-------+-------+
| 2 3 1 | 5 7 4 | 6 9 8 |
| 9 6 8 | 2 3 1 | 5 7 4 |
| 5 7 4 | 9 6 8 | 2 3 1 |
+-------+-------+-------+

More Sudokus

These are a few resources I used to test the solver. I used the 1 million sudoku dataset in order to get an average solve time.

1 million sudoku games: https://www.kaggle.com/datasets/bryanpark/sudoku

9 million sudoku games: https://www.kaggle.com/datasets/rohanrao/sudoku

You can also find a weekly "Unsolvable" sudoku posted here.

References

Knuth, D. 2000. Dancing Links. https://www.ocf.berkeley.edu/~jchu/publicportal/sudoku/0011047.pdf accessed Dec 2022

Wikipedia, 2022. Dancing Links. https://wikipedia.org/wiki/Dancing_Links?lang=en accessed Dec 2022

Wikipedia, 2022. Exact cover. https://wikipedia.org/wiki/Exact_cover?lang=en accessed Dec 2022

Wikipedia, 2022. Knuth's Algorithm X. https://wikiless.org/wiki/Knuth%27s_Algorithm_X?lang=en accessed Dec 2022

Gillespie, P. (n.d). Text to Ascii Art Generator (TAAG): https://patorjk.com/software/taag/#p=display&f=Big accessed Dec 2022

Inkala, A. 2012. Arto Inkala's Monster sudoku: https://www.sudokuwiki.org/Arto_Inkala_Sudoku accessed Dec 2022

Inkala, A. 2006. Escargot: https://www.sudokuwiki.org/Escargot accessed Dec 2022

pyutils, 2022. line_profiler https://github.com/pyutils/line_profiler accessed Dec 2022

Park, K. 2014. 1 million Sudoku games: https://www.kaggle.com/datasets/bryanpark/sudoku accessed Dec 2022

Rohanrao, 2019. 9 million Sudoku Puzzles and Solutions: https://www.kaggle.com/datasets/rohanrao/sudoku accessed Dec 2022

Sudoku Wiki. (n.d.). Weekly Unsolvable puzzle: https://www.sudokuwiki.org/Weekly_Sudoku.asp?puz=28 accessed Dec 2022

Blender, 2013. How to import all modules in dir: https://stackoverflow.com/questions/16852811/python-how-to-import-from-all-modules-in-dir accessed Dec 2022