Example: Magic Squares

The model will eliminate symmetries by imposing order constraints. Propagation will be drastically improved by exploiting a redundant constraint.

Problem Specification

The Magic Square Puzzle consists in finding for given N an N x N - matrix such that: A magic square for N = 3 looks as follows:
2 7 6
9 5 1
4 3 8
The Magic Square Puzzle is extremely hard for large N. Even for N = 4, our script will have to explore more than 10000 nodes to find a solution.

Viewpoint and Constraints

We model the problem by having a variable Fij for every field (i,j) of the matrix. Moreover, we have one additional variable S and require that the sum of every row, column, and main diagonal equals S.

Without loss of generality, we can impose the following order constraints eliminating symmetries:

F11 < FNN FN1 < F1N F11 < FN1
Since the sum of the sums of the rows must equal the sum of all fields, we have the redundant constraint
$ {\frac{{N^2}}{{2}}}$ * (N2 + 1) = N * S
To see this, note that sum of all fields equals
1 + 2 + 3 + ... + N2 = $ {\frac{{N^2}}{{2}}}$ * (N2 + 1)
and that the sum of each of the N rows must be S.

Branching Strategy

We branch on the variables F11,..., FNN with a first-fail strategy splitting the domain of the selected variable and trying the lower half first.


val maxValue = (valOf (Int.maxInt))div 2

fun magicsquares n space =
        val nn = n * n
        val matrix = Vector.tabulate(n,fn x => 
        val matrixv =Vector.concat(Vector.toList(matrix))
        val sum = FD.range(space,(1,maxValue))
        val sumn = FD.range(space,(1,maxValue))
        val diagonal1 = Vector.tabulate(n,fn x =>(x,x))
        val diagonal2 = Vector.tabulate(n,fn x =>(n-1-x,x))
        val diagonal1vars =,y) => 
        val diagonal2vars =,y) => 
        val field11 = Vector.sub(Vector.sub(matrix,0),0) 
        val field1N = Vector.sub(Vector.sub(matrix,0),n-1)
        val fieldN1 = Vector.sub(Vector.sub(matrix,n-1),0)
        val fieldNN = Vector.sub(Vector.sub(matrix,n-1),n-1)                 
        (* diagonals *)
        post(space,SUMV(diagonal1vars) `= FD(sum),FD.DOM);
        post(space,SUMV(diagonal2vars) `= FD(sum),FD.DOM);
        (* columns *)
            val colmn = Vector.tabulate(n,fn x =>
            post(space,SUMV(colmn) `= FD(sum),FD.DOM)
        (* rows *)
        Vector.appi(fn(x,y) =>
                     post(space,SUMV(y) `= FD(sum),FD.DOM))matrix;
        (* Eliminate symmetries *)
        (*redundant constraints *)
        post(space,FD(sum) `* `n `= FD(sumn),FD.DOM);
        FD.relI(space,sumn,FD.EQ,(nn *(nn + 1))div 2);

With the Explorer you can find out that for N = 3 there is exactly one magic square satisfying the ordering constraints of our model. Without the ordering constraints there are 8 different solutions. Omitting the propagator for the redundant constraint will increase the search tree by one order of magnitude.

Andreas Rossberg 2006-08-28