# Solving Sudoku as an integer programming problem

Posted 4 years ago
13982 Views
|
18 Replies
|
27 Total Likes
|

It is fairly straight forward to solve a Sudoku as an integer programming problem, by creating 9 binary variables for each cell, only one of which is one in the solution. The walk-through below and attached notebook illustrates this for the problem shown.

# Implementation

Given values, as {row, column, value}

input = {
{1,4,4},{1,5,9},{1,8,5},{2,1,6},{2,5,3},{3,1,4},
{3,2,5},{3,4,6},{3,5,2},{3,7,3},{3,9,7},{4,1,5},
{4,3,2},{4,4,7},{4,7,9},{4,8,8},{5,1,3},{5,3,6},
{5,7,2},{5,9,1},{6,2,9},{6,3,1},{6,6,2},{6,7,6},
{6,9,5},{7,1,2},{7,3,5},{7,5,1},{7,6,4},{7,8,3},
{7,9,8},{8,5,8},{8,9,9},{9,2,1},{9,5,7},{9,6,3}};


Display given values

viewmat = Table["", {9}, {9}];
Do[viewmat[[input[[i, 1]], input[[i, 2]]]] = ToString[input[[i, 3]]], {i, Length[input]}]
Grid[viewmat, Frame -> All]


Variables, as 9 x 9 x 9 matrix

varmat = Table[m[i, j, k], {i, 9}, {j, 9}, {k, 9}];


Variables as a list

vars = Flatten[varmat];


Constrain the input cells to their value

cons1 = (varmat[[Sequence @@ #]] == 1 &) /@ input


The sum of the binary variables for each cell is 1

cons2 = Flatten @ Table[ (Sum[varmat[[i, j, k]], {k, 9}] == 1), {i, 9}, {j, 9}];


All different constraint for the rows

cons3 = Flatten @ Table[ (Sum[varmat[[i, j, k]], {i, 9}] == 1), {j, 9}, {k, 9}];


All different constraint for the columns

cons4 = Flatten @ Table[ (Sum[varmat[[i, j, k]], {j, 9}] == 1), {i, 9}, {k, 9}];


All different constraint for the submatrices

sm[di_, dj_] := Flatten [Table[{i, j}, {i, 1 + 3*(di - 1), 3*di}, {j, 1 + 3*(dj - 1), 3*dj}],1]
cons5 = Flatten @ Table[(Total[m[Sequence @@ #, k] & /@ sm[i, j]] == 1), {i, 3}, {j, 3}, {k, 9}];


Confine the variables to the range 0 to 1

cons6 = Thread[0 <= vars <= 1];


Combine the constraints

Length[allcons = Join[cons1, cons2, cons3, cons4, cons5, cons6]]


1089

Solve the problem, specifying that the variables are integers.

AbsoluteTiming[sol = FindMinimum[{0, allcons, Element[vars, Integers]}, vars];]


{0.0946335, Null}

Find the values for each cell

resmat = Table[Sum[k*m[i, j, k], {k, 9}], {i, 9}, {j, 9}] /. sol[[2]]


Display the input and result

{Grid[viewmat, Frame -> All], Grid[resmat, Frame -> All]}


Check the result

And @@ Table[Unequal[Sequence @@ resmat[[i]]], {i, 9}]


True

And @@ Table[Unequal[Sequence @@ Transpose[resmat][[i]]], {i, 9}]


True

And @@ Flatten @ Table[Unequal[resmat[[Sequence @@ #]] & /@ sm[i, j]], {i, 3}, {j, 3}]


True

Attachments:
18 Replies
Sort By:
Posted 4 years ago
 That's a neat way of doing it! I actually wrote a solved on my last flight (a week ago or so) from USA back to Europe. I was having doubts with the in-flight entertainment system version of Sudoku; it generated non-unique Sudokus...so the puzzles had many many solutions... big shame. I'm still considering sending Delta airlines an email about it...
Posted 4 years ago
 Some years ago I wrote a "brute force" Sudoku solver which goes row by row, finding all possible solutions. It's in the Library Archive. The Integer programming method runs much faster. I also tried using Reduce with Unequal constraints as well but ran out of patience waiting for it to return.
Posted 4 years ago
Posted 4 years ago
 It would be nice if Reduce applied this technique automatically.
Posted 4 years ago
 Reduce uses exact methods. Which might be somewhat slower than FindMinimum, but still should work for purposes of solving a Sudoku.
Posted 4 years ago
 using Unequal constraints or Integer Programming?
Posted 4 years ago
 Using integer programming. Unequal constraints will very possibly get blown up into a logical disjunction, which is a combinatorial explosion in terms of the computational complexity. ILP has exponential complexity also, on a bad day. But sudoku puzzles rarely give it quite that bad a day.
Posted 4 years ago
 FindMinimum takes about 0.1 sec. Reduce does not return a solution after 10 minutes.
Posted 4 years ago
 Difficult for me to comment in the absence of actual code. If it was something like this Reduce[allcons, vars, Integers]; I will agree: exact ILP can be slow. That's why I often write branch-and-prune loops using NMinimize.
Posted 4 years ago
 yes, I tried Reduce[ allcons, vars, Integers ].
Posted 4 years ago
 - you have earned "Featured Contributor" badge, congratulations! This is a great post and it has been selected for the curated Staff Picks group. Your profile is now distinguished by a "Featured Contributor" badge and displayed on the "Featured Contributor" board.
Posted 4 years ago
 Very cool. I see this post is trending now on Reddit-Programming: https://redd.it/5h1xgq There is a comment there mentioning MS-OML. Just for the sake of comparison i give their code that solves the problem. I wonder if we can learn anything form it. Model [ Parameters[Sets,I,J], Parameters[Integers,d[I,J]], Decisions[Integers[1,9],x[I,J]], Constraints [ FilteredForeach[{i,I},{j,J},d[i,j]>0,x[i,j]==d[i,j]], Foreach [{i,I},Unequal[Foreach[{j,J},x[i,j]]]], Foreach [{j,J},Unequal[Foreach[{i,I},x[i,j]]]], Foreach [{ib,3}, Foreach [{jb,3}, Unequal[Foreach[{i,ib*3,ib*3+3},{j,jb*3,jb*3+3},x[i,j] ] ] ] ] ] ] 
Posted 4 years ago
 I don't understand why the Foreach statements are nested in the MS-OML code.
Posted 4 years ago
 Don't shoot the messenger. I just reposted it here verbatim for easy access by Community folks. I do not understand MS-OML language at all. Perhaps some knowledgable folks would be able to comment.
Posted 4 years ago
 I am not familiar with MS-OML but it seems clear enough what this particular nesting does. We'll look at one such."Foreach [{i,I},Unequal[Foreach[{j,J},x[i,j]]]]"Taking this apart, for each i in the "row" subscript set, we set up an unequal involving ALL of {x[i,1],x[i,2],...x[i,9]}.
Posted 3 years ago
 Much clever thought went into your solution, thanks for sharing! I did not understand at first, what and why was happening in your code and i had to refresh my memory on numerical optimization. Understanding now the principle of your algorithm (maybe it cannot be called an algorithm because it is just setting up the maths to be solved by the solver in 1 step), i actually cannot think of any other alternative sensible way of solving the sudoku puzzle on a computer. I believe that there are many different ways (actual algorithms!) which solve the sudoku puzzle on a computer, e.g. working with permutations and such. But i always wanted to fully understand just one way of getting there, and if it is a mathematical way, the better. And yours is a purely mathematical way, solving a huge system of equations!!To not lose my train of thought, i've edited your *.nb-file (minor code rearrangements and improvements, added explanations) and am sharing it here as well, please see attached.Hope this helps, God bless! Attachments:
 Single line Sudoku solver: https://rosettacode.org/wiki/Sudoku#Mathematica solve[sudoku_] := NestWhile[ Join @@ Table[ Table[ReplacePart[s, #1 -> n], {n, #2}] & @@ First@SortBy[{#, Complement[Range@9, s[[First@#]], s[[;; , Last@#]], Catenate@ Extract[Partition[s, {3, 3}], Quotient[#, 3, -2]]]} & /@ Position[s, 0, {2}], Length@Last@# &], {s, #}] &, {sudoku}, ! FreeQ[#, 0] &]