[Math] Gröbner basis for Sudoku

ag.algebraic-geometryalgorithmsgroebner-basesrecreational-mathematicssudoku

I'm trying to write a program that solves Sudokus using Gröbner bases. I introduced $81$ variables, $x_1$ to $x_{81}$, this is a linearisation of the Sudoku board.

The space of valid Sudokus is defined by:

  • for $i=1,\ldots,81$: $$F_i = (x_i – 1)(x_i – 2)\cdots(x_i – 9)$$
    This represents the fact that all squares have integer values between 1 and 9.

  • for all $x_i$ and $x_j$ which are not equal but in the same row, column or block: $$G_{ij} = (F_i – F_j)/(x_i – x_j)$$ This represents that the variables $x_i$ and $x_j$ can not be equal.

All these $F_i$ and $G_{ij}$ together define the space of valid Sudokus. This consists of $891$ polynomials.

Now to solve a Sudoku we can add the clues to the space, so by example if the clue of a Sudoku is the first square is a $5$, then we add $(x_1 – 5)$ to the space. If we now take the Gröbner basis of this space we can directly see the solution for it.

I understand what I am doing this far. But I have trouble finding a computable manner for finding the gröbner bases. I have succesfully done everything for $4 \times 4$ Sudokus (or so-called Shidokus). But neither Maple nor Singular are giving me a result for the Gröbner basis of the $9 \times 9$ Sudoku space.
You can see the commands I gave to Maple here. First I define the $891$ polynomials, then I ask for a basis of it.

I read papers saying it's feasible although non-performant to do what I strive for but I don't see how to find the solution, as they don't include many implementation details. Can anyone point me to a direction, making this problem easier for Maple or other software?

Best Answer

Here is a Singular Code that works quite well:

ring A = 0,(t,x(1..9)),lp;
/* Characteristic 0 works suprisingly well for this problem. */
/* We choose a lexicographic ordering since we will compute an
   elimination ideal. */

poly p = (t-x(1))*(t-x(2))*(t-x(3))*(t-x(4))*(t-x(5))*(t-x(6))*(t-x(7))*(t-x(8))*(t-x(9))-(t-1)*(t-2)*(t-3)*(t-4)*(t-5)*(t-6)*(t-7)*(t-8)*(t-9);
/* p(x)=0 in Q[t] implies that x is a permutation of the numbers 1 to 9. */

matrix c = coeffs(p,t);
ideal J = (c[1..9,1]);
/* J expresses that x is a permutation of the numbers 1 to 9. However,
   surprisingly, it is better to use only constraints saying that x(8)
   and x(9) are distinct integers between 1 and 9. This is done by
   computing an elimination ideal. */
ideal JG = groebner(J);
ideal J2 = (JG[1],JG[2]);
/* J2 is the ideal expressing that x(1) and x(2) are distinct integers
   between 1 and 9. */

ring R=0,(x(1..81)),dp;
ideal I;
map psi;
proc f(k,l,m,n,o,p,q,r,s)
{intvec v = k,l,m,n,o,p,q,r,s;
 int i,j;
 for (i=1; i<=8; i++) {for (j=i+1; j<=9; j++)
 {psi = A,0,1,2,3,4,5,6,7,x(v[i]),x(v[j]); I = I + psi(J2);}}}

/* Code the rules into the ideal. */
f(1,2,3,4,5,6,7,8,9);
f(10,11,12,13,14,15,16,17,18);
f(19,20,21,22,23,24,25,26,27);
f(28,29,30,31,32,33,34,35,36);
f(37,38,39,40,41,42,43,44,45);
f(46,47,48,49,50,51,52,53,54);
f(55,56,57,58,59,60,61,62,63);
f(64,65,66,67,68,69,70,71,72);
f(73,74,75,76,77,78,79,80,81);
f(1,10,19,28,37,46,55,64,73);
f(2,11,20,29,38,47,56,65,74);
f(3,12,21,30,39,48,57,66,75);
f(4,13,22,31,40,49,58,67,76);
f(5,14,23,32,41,50,59,68,77);
f(6,15,24,33,42,51,60,69,78);
f(7,16,25,34,43,52,61,70,79);
f(8,17,26,35,44,53,62,71,80);
f(9,18,27,36,45,54,63,72,81);
f(1,2,3,10,11,12,19,20,21);
f(4,5,6,13,14,15,22,23,24);
f(7,8,9,16,17,18,25,26,27);
f(28,29,30,37,38,39,46,47,48);
f(31,32,33,40,41,42,49,50,51);
f(34,35,36,43,44,45,52,53,54);
f(55,56,57,64,65,66,73,74,75);
f(58,59,60,67,68,69,76,77,78);
f(61,62,63,70,71,72,79,80,81);

/* Code a uniquely solvable Sudoku problem into the ideal.  */
I=I+(x(3)-4);
I=I+(x(6)-3);
I=I+(x(7)-6);
I=I+(x(9)-9);
I=I+(x(12)-8);
I=I+(x(13)-9);
I=I+(x(16)-2);
I=I+(x(18)-1);
I=I+(x(22)-8);
I=I+(x(23)-1);
I=I+(x(27)-7);
I=I+(x(28)-6);
I=I+(x(33)-7);
I=I+(x(37)-8);
I=I+(x(40)-3);
I=I+(x(42)-9);
I=I+(x(45)-5);
I=I+(x(49)-6);
I=I+(x(54)-3);
I=I+(x(55)-4);
I=I+(x(59)-7);
I=I+(x(60)-6);
I=I+(x(64)-2);
I=I+(x(66)-7);
I=I+(x(69)-5);
I=I+(x(70)-1);
I=I+(x(73)-9);
I=I+(x(75)-1);
I=I+(x(76)-2);
I=I+(x(79)-4);

option(redSB);
groebner(I,30);
/* You get the solution quickly. */
Related Question