sysid blog

Sudoku

April 26, 2020 ☕️ 4 min read

solution

Wikipedia

a logic-based, combinatorial number-placement puzzle.

Fill a 9×9 grid with digits so that each column, each row, and each of the nine 3×3 subgrids that compose the grid contain all of the digits from 1 to 9.

The puzzle setter provides a partially completed grid, which for a well-posed puzzle has a single solution.

Challenges

  1. Understand the problem (easy)
  2. Design the index sets (medium)
  3. Translate the model into Pyomo (easy)

Tooling

  • Pyomo as LP modelling language
  • Optimizer: CBC
  • Python

Model

Indexing Schema

sudoku

Variables

yi,j;km={1,  when cell (i,j;k) contains number m 0,  else y^m_{i,j;k} = \begin{cases} 1, \ \text{ when cell $(i,j;k)$ contains number $m$ }\\ 0, \ \text{ else }\\ \end{cases}

Example: If 7 appears in (2,1;1)(2,1;1), then y2,1;17=1y^7_{2,1;1} = 1

Row Constraints

k=13j=13yi,j;km,(i=1,2,3,m=1,..,9)\sum_{k=1}^3\sum_{j=1}^3 y_{i,j;k}^m, (i=1,2,3, m=1,..,9)\\k=46j=13yi,j;km,(i=1,2,3,m=1,..,9)\sum_{k=4}^6\sum_{j=1}^3 y_{i,j;k}^m, (i=1,2,3, m=1,..,9)\\k=79j=13yi,j;km,(i=1,2,3,m=1,..,9)\sum_{k=7}^9\sum_{j=1}^3 y_{i,j;k}^m, (i=1,2,3, m=1,..,9)\\

Column Constraints

k{1,4,7}i=13yi,j;km=1,(j=1,2,3,m=1,..,9)\sum_{k \in \{1,4,7\}} \sum_{i=1}^3 y_{i,j;k}^m = 1, (j=1,2,3, m=1,..,9)\\k{2,5,8}i=13yi,j;km=1,(j=1,2,3,m=1,..,9)\sum_{k \in \{2,5,8\}} \sum_{i=1}^3 y_{i,j;k}^m = 1, (j=1,2,3, m=1,..,9)\\k{3,6,9}i=13yi,j;km=1,(j=1,2,3,m=1,..,9)\sum_{k \in \{3,6,9\}} \sum_{i=1}^3 y_{i,j;k}^m = 1, (j=1,2,3, m=1,..,9)\\

Sub-Table Constraint

Subtable contains every number mm once:

i=13j=13yi,j;km=1,(k,m=1,..,9)\sum_{i=1}^3 \sum_{j=1}^3 y_{i,j;k}^m = 1, (k,m=1,..,9)\\

Completeness Constraint

Each cell has a number in it:

m=19yi,j;km=1,(k=1,..,9 i,j=1,2,3)\sum_{m=1}^9 y_{i,j;k}^m = 1, (k=1,..,9 \ i,j=1,2,3)\\

Objective

This is not a maximization problem. Only solution feasibility is relevant, so we assume a constant objective.

Pyomo Implementation

The implementation in Pyomo is straightfoward:

# Row Constraint
model.row_c_1_3 = Constraint(model.I, model.M, rule=lambda model, i, m: sum(
    sum(model.y[i, j, k, m] for j in model.J) for k in model.K if k <= 3) == 1)

model.row_c_4_6 = Constraint(model.I, model.M, rule=lambda model, i, m: sum(
    sum(model.y[i, j, k, m] for j in model.J) for k in model.K if 4 <= k <= 6) == 1)

model.row_c_7_9 = Constraint(model.I, model.M, rule=lambda model, i, m: sum(
    sum(model.y[i, j, k, m] for j in model.J) for k in model.K if 7 <= k) == 1)

# Row Constraint
model.col_147_c = Constraint(model.J, model.M, rule=lambda model, j, m: sum(
    sum(model.y[i, j, k, m] for i in model.I) for k in model.K if k in (1, 4, 7)) == 1)

model.col_258_c = Constraint(model.J, model.M, rule=lambda model, j, m: sum(
    sum(model.y[i, j, k, m] for i in model.I) for k in model.K if k in (2, 5, 8)) == 1)

model.col_369_c = Constraint(model.J, model.M, rule=lambda model, j, m: sum(
    sum(model.y[i, j, k, m] for i in model.I) for k in model.K if k in (3, 6, 9)) == 1)

# Sub-Table constraint
model.subtable_c = Constraint(model.K, model.M, rule=lambda model, k, m: sum(
    sum(model.y[i, j, k, m] for i in model.I) for j in model.J) == 1)

# One number per cell constraint
model.number_c = Constraint(
    model.I, model.J, model.K, rule=lambda model, i, j, k: sum(model.y[i, j, k, m] for m in model.M) == 1
)

Result

Here our the input, a fairly difficult Sudoku for humans:

challenge

It takes less than a second to solve for MIP solver CBC:

---------------------------------
[ 1  2  6 ][ 4  3  7 ][ 9  5  8 ]
[ 8  9  5 ][ 6  2  1 ][ 4  7  3 ]
[ 3  7  4 ][ 9  8  5 ][ 1  2  6 ]
---------------------------------
[ 4  5  7 ][ 1  9  3 ][ 8  6  2 ]
[ 9  8  3 ][ 2  4  6 ][ 5  1  7 ]
[ 6  1  2 ][ 5  7  8 ][ 3  9  4 ]
---------------------------------
[ 2  6  9 ][ 3  1  4 ][ 7  8  5 ]
[ 5  4  8 ][ 7  6  9 ][ 2  3  1 ]
[ 7  3  1 ][ 8  5  2 ][ 6  4  9 ]
---------------------------------

Summary

The challenge here was to devise a suitable indexing schema in order to formulate the model constraints. Apart from that everything is straightforward, especially model implementation is very convenient with Pyomo.

If you are interested in the Pyomo model or the Python code contact me via mail.


always be building, by sysid.

© 2022