Introduction
In [1] another strange Sudoku variant is presented.
The givens are just:
Obviously, more than just the standard Sudoku constraints are needed to have just one unique feasible solution. So, for this problem we have the following rules:
- Standard Sudoku rules apply. So, for each row, column, and sub-block (a.k.a. nonet [3]), we have uniqueness constraints (each cell in a region has a unique number between 1 and 9).
- Borrowing from chess: any two cells corresponding to a knight's move, must have different values.
- Similarly, any two cells that form a king's move, must also be different.
- Orthogonally adjacent cells must be non-consecutive.
Standard Sudoku setup
When we want to solve Sudokus, the easiest approach is to define the following binary decision variables [4]: \[x_{i,j,k} = \begin{cases} 1 & \text{if cell $(i,j)$ has value $k$} \\ 0 & \text{otherwise}\end{cases}\] Here \(k \in \{1,\dots,9\}\). We have 27 areas we need to check for unique values: rows, columns and nonets. We can organize this as a set: \[u_{a,i,j}\>\text{exists if and only if area $a$ contains cell $(i,j)$}\] This is data. We also have these two given cells, i.e. fixed variables for these cells. The resulting MIP model can look like [4]:
| Mixed-Integer Programming Model for standard Sudoku |
|---|
| \[\begin{align}\min\>& 0 && && \text{Dummy objective}\\ & \sum_k \color{darkred}x_{i,j,k} = 1 &&\forall i,j && \text{One $k$ in each cell}\\ & \sum_{i,j|\color{darkblue}u_{a,i,j}} \color{darkred}x_{i,j,k} = 1 && \forall a,k && \text{Unique values in each area}\\ & \color{darkred}x_{i,j,k} = 1 && \text{where } k=\color{darkblue}{\mathit{Given}}_{i,j} &&\text{Fix given values}\\ &\color{darkred}x_{i,j,k} \in \{0,1\} \end{align}\] |
During post-processing, we can calculate the completed grid using the optimal values of \(x\): \[v_{i,j} = \sum_k k \cdot x^*_{i,j,k}\] or we can include them as "accounting rows" to the model. Accounting rows are constraints that have the function to calculate some quantities for reporting purposes.
The implementation in GAMS uses a set u(a,i,j)that is populated as:
set
a 'areas' /r1*r9,c1*c9,s1*s9/ i(a) 'rows' /r1*r9/ j(a) 'columns' /c1*c9/ s(a) 'squares' /s1*s9/ u(a,i,j) 'areas with unique values' ; * populate u(a,i,j): u(a,i,j)=YES if area a contains cell (i,j) * see https://yetanothermathprogrammingconsultant.blogspot.com/2016/10/mip-modeling-from-sudoku-to-kenken.html u(i,i,j) = yes; u(j,i,j) = yes; u(s,i,j)$(ord(s)=3*ceil(ord(i)/3)+ceil(ord(j)/3)-3) = yes; |
The model equations look like:
binary variable x(i,j,k);
* v0(i,j) are the givens x.fx(i,j,k)$(v0(i,j)=ord(k)) = 1; variable z 'dummy objective'; integer variable v(i,j); v.lo(i,j) = 1; v.up(i,j) = 9; equation dummy 'objective' unique(a,k) 'all-different' values(i,j) 'one value per cell' eqv(i,j) 'calculate values' ; dummy.. z =e= 0; unique(a,k).. sum(u(a,i,j), x(i,j,k)) =e= 1; values(i,j).. sum(k, x(i,j,k)) =e= 1; eqv(i,j).. v(i,j) =e= sum(k, ord(k)*x(i,j,k)); model sudoku /all/; solve sudoku minimizing z using mip; |
We can test this approach with the above data. We will actually find a feasible solution, but of course it will not be unique. I believe there are likely millions of different solutions (or more). Counting them is not that easy. I tried to give it a shot with a solution pool approach: after setting a limit if 1,000,000 solutions, the solver stopped after reaching this limit. So, there are more than 1 million solutions for this (partial) problem.
Knight's jumps
To model how knight's moves are covering cells, we can have a look at [5].
To place knights we set up a set \(\mathit{jump}_{i,j,i',j'}\) indicating if we can jump from cell \((i,j)\) to cell \((i',j')\). If we can jump from \((i,j) \rightarrow (i',j')\) we can also jump from \((i',j') \rightarrow (i,j)\). We don't want to check both cases. To prevent this double check, we only need to look forward. So for each \((i,j)\) we need to consider just four cases:
| \(j\) | \(j+1\) | \(j+2\) | |
| \(i-2\) | \(x_{i-2,j+1}\) | ||
| \(i-1\) | \(x_{i-1,j+2}\) | ||
| \(i\) | \(x_{i,j}\) | ||
| \(i+1\) | \(x_{i+1,j+2}\) | ||
| \(i+2\) | \(x_{i+2,j+1}\) |
Note that near the border we may have fewer than four cases. In GAMS we can populate the set \(\mathit{jump}\) in a straightforward manner:
alias
(i,ii),(j,jj);
set jump(i,j,ii,jj); jump(i,j,i-2,j+1) = yes; Jump(i,j,i-1,j+2) = yes; Jump(i,j,i+1,j+2) = yes; Jump(i,j,i+2,j+1) = yes; |
There are 224 elements in the set \(\mathit{jump}\).
Now comes the more complicated part: how to make sure that the values of cell \((i,j)\) and \((i',j')\) are not the same. There are basically two approaches:
- Compare the values \(v_{i,j}\) and \(v_{i',j'}\). The constraint becomes \[|v_{i,j}-v_{i',j'}| \ge 1 \>\> \forall i,j,i',j'|\mathit{jump}(i,j,i',j')\] This can be linearized using \[v_{i,j}-v_{i',j'} \ge 1 \textbf{ or } v_{i',j'}-v_{i,j} \ge 1\] or \[\begin{align}& v_{i,j}-v_{i',j'} \ge 1 - M\cdot \delta_{i,j,i',j'} \\ & v_{i',j'}-v_{i,j} \ge 1- M (1-\delta_{i,j,i',j'}) \\ &\delta_{i,j,i',j'} \in \{0,1\}\end{align}\] Here \(M\) is a large enough constant (\(M=10\) should suffice).
- Directly work with \(x_{i,j,k}\) and \(x_{i',j',k}\). The constraint should be \[x_{i,j,k} x_{i',j',k} = 0 \>\> \forall k, \forall i,j,i',j'|\mathit{jump}(i,j,i',j') \] This non-linear constraint can be linearized as \[x_{i,j,k}+x_{i',j',k}\le 1\] This clearly is the easier route.
In GAMS this can look like:
equation
Different(i,j,ii,jj,k) 'cells forming a knights move should be different';
Different(jump(i,j,ii,jj),k).. x(i,j,k) + x(ii,jj,k) =l= 1; |
Kings's moves
Here we want to implement a similar restriction as for the knight's move. As the constraints are the same, we just have to augment the set \(\mathit{jump}\) a bit.
alias
(i,ii),(j,jj);
set jump(i,j,ii,jj) 'knights or kings move'; * knight jump(i,j,i-2,j+1) = yes; jump(i,j,i-1,j+2) = yes; jump(i,j,i+1,j+2) = yes; jump(i,j,i+2,j+1) = yes; * king jump(i,j,i+1,j+1) = yes; jump(i,j,i-1,j+1) = yes; |
After this, the set \(\mathit{jump}\) has grown from 224 to 352 elements.
Orthogonal neighbors can't be consecutive
First we need to build a set that models "orthogonal neighbors". Similar to our set \(\mathit{jump}\), we want to prevent comparing pairs twice. So given a cell \((i,j)\), we only need to consider two neighbors.
| \(j\) | \(j+1\) | |
| \(i-1\) | \(x_{i-1,j}\) | |
| \(i\) | \(x_{i,j}\) | \(x_{i,j+1}\) |
| \(i+1\) |
This is coded in GAMS as:
set nb(i,j,ii,jj) 'orthogonal neighbors';
nb(i,j,i-1,j) = yes; nb(i,j,i,j+1) = yes; |
This set has 144 elements.
The constraints are: \[\begin{align}&x_{i,j,k} + x_{i',j',k+1} \le 1\\&x_{i,j,k} + x_{i',j',k-1} \le 1\end{align}\>\> \forall i,j,i',j'|nb(i,j,i',j'), \forall k\]
The GAMS representation is:
equations
NonConsecutive1(i,j,ii,jj,k) 'orthogonal neighbors constraint' NonConsecutive2(i,j,ii,jj,k) 'orthogonal neighbors constraint' ; NonConsecutive1(nb(i,j,ii,jj),k+1).. x(i,j,k) + x(ii,jj,k+1) =l= 1; NonConsecutive2(nb(i,j,ii,jj),k-1).. x(i,j,k) + x(ii,jj,k-1) =l= 1; |
Detail: The last index of the equation definition has a \(k+1\) or \(k-1\). This tells GAMS not to generate all possile constraints for all \(k\) but rather one less. This is not strictly needed: in GAMS addressing outside the domain results in a zero. But this trick will prevent singleton constraints \(x_{i,j,k} + 0 \le 1\) to be generated.
Results
The complete model has 811 variables (808 of them binary or integer). The number of constraints is 5,878. The results look like:
---- 111 VARIABLE v.L c1 c2 c3 c4 c5 c6 c7 c8 c9 r1 4 8 3 7 2 6 1 5 9 r2 7 2 6 1 5 9 4 8 3 r3 1 5 9 4 8 3 7 2 6 r4 8 3 7 2 6 1 5 9 4 r5 2 6 1 5 9 4 8 3 7 r6 5 9 4 8 3 7 2 6 1 r7 3 7 2 6 1 5 9 4 8 r8 6 1 5 9 4 8 3 7 2 r9 9 4 8 3 7 2 6 1 5
The problem solves in the presolve phase, i.e. 0 iterations and 0 nodes.
Uniqueness of the solution
To prove the uniqueness of the solution, we add a cut that forbids the current solution \(x^*\): \[\sum_{i,j,k} x^*_{i,j,k} x_{i,j,k} \le 9^2-1\] The resulting model is infeasible. This means that the solution is indeed unique.
Conclusion
This Sudoku variant combines standard Sudoku rules with additional constraints, some of them borrowed from chess. This makes the modeling at bit trickier. In the development of the model, I have placed much of the logic in sets. This has the advantage that the constraints are fairly simple. In general, sets are easier to debug than constraints: we can display and debug sets before we have a working model.
Both the mathematical model and the GAMS representation is interesting. The mathematical formulation shows we can model this without extra (discrete) variables. The GAMS implementation details assembling sets that simplify constraints.
Sudoku models typically solve extremely fast: they are solved by the presolver. This model is no exception.
References
- The Miracle Sudoku, https://youtu.be/yKf9aUIxdb4
- Puzzled man solving 'miracle' sudoku becomes YouTube sensation,https://www.theguardian.com/lifeandstyle/2020/may/22/cracking-the-cryptic-puzzled-uk-men-become-internet-sensation-with-sudoku-channel
- Glossary of Sudoku, https://en.wikipedia.org/wiki/Glossary_of_Sudoku
- MIP Modeling: from Sudoku to KenKen via Logarithms, https://yetanothermathprogrammingconsultant.blogspot.com/2016/10/mip-modeling-from-sudoku-to-kenken.html
- Chess and solution pool, https://yetanothermathprogrammingconsultant.blogspot.com/2018/11/chess-and-solution-pool.html
No comments:
Post a Comment