Monday, November 10, 2008

Sudoku: Constraint Programming vs MIP

The Sudoku puzzle can be coded easily using constraint programming thanks to built-in constructs that handle the all-different constraint. Here is an example using the Microsoft Solver Foundation CSP solver:


(click to enlarge)

The only difficulty was to model: for all i,j such that d[i,j]>0 require that x[i,j]=d[i,j]. I used the following formulation: Foreach[{i,I},{j,J},x[i,j]==d[i,j] | d[i,j]==0] where | is an "or" operator.

The complete model looks like:





// Alternative 1: (does not work in the beta)
// Alternative 2:
Foreach[{i,I},{j,J},x[i,j]==d[i,j] | d[i,j]==0],



The above sudoku can also be coded in GAMS:

Sudoku solver
September 2005

Code written by:
Paul van der Eijk
GAMS Development Corporation
1217 Potomac Street NW
Washington DC, 20007
Tel: (202) 342-0180 Fax: (202) 342-0181

Based on Latin Square example in GAMS model library.

Adapted to solve a Sudoku problem by:
Sherman Robinson
Economics Department
University of Sussex
Brighton BN1 9SN


r rows / r1 * r9/
c columns / c1 * c9/
v values / v1 * v9/
sc sub-col /sc1 * sc3/
sr sub-row /sr1 * sr3/
scs(sc, c) /sc1.(c1*c3), sc2.(c4*c6), sc3.(c7*c9)/
srs(sr, r) /sr1.(r1*r3), sr2.(r4*r6), sr3.(r7*r9)/

table problem(r,c)
c1 c2 c3 c4 c5 c6 c7 c8 c9
r1 5 3 7
r2 6 1 9 5
r3 9 8 6
r4 8 6 3
r5 4 8 3 1
r6 7 2 6
r7 6 2 8
r8 4 1 9 5
r9 8 7 9
display problem;

Parameter value(v) "Values";
value(v) = ord(v) ;

binary variable x(r,c,v);
variable w;

eq1(r,c) "exactly one value for each cell"
eq2(c, v) "columns entries have to be unique"
eq3(r, v) "row entries have to be unique"
eq4(sr, sc, v) "sub-squares have to be unique"
nobj "definition of objective - anything"

eq1(r, c).. sum(v, x(r, c, v)) =E= 1;
eq2(c, v).. sum(r, x(r, c, v)) =E= 1;
eq3(r, v).. sum(c, x(r, c, v)) =E= 1;
eq4(sr, sc, v).. sum((c,r)$(scs(sc, c) and srs(sr, r)), x(r, c, v)) =E= 1;

nobj.. w =E= sum((r, c, v), x(r, c, v));

*SR Fix values in the problem
x.l(r,c,v)$(problem(r,c) = value(v)) = 1 ;
x.fx(r,c,v)$(x.l(r,c,v) = 1) = 1 ;

option limrow=1000, limcol=0;
option mip=cplex;
model sudoku /eq1, eq2, eq3, eq4, nobj/;
solve sudoku minimizing w using mip;

parameter solution(r,c) "Display of results";
solution(r, c)$(x.l(r, c, v)) = Ord(v);
option solution:0:1:1;
display solution;

I wrote an Excel based GUI for this model as an example, see