prolog sudoku solver running out of global stack

Question

I tried writing a binary sudoku solver in swi-prolog. (a binary sudoku is explained here)

The problem is that I am now running out of global stack. I am giving it 2 gb which should be more than enough. Am I using a flawed algorithm? Is there something I could be doing better to avoid running in to the lack of global stack error with such small puzzles?

A bit more information : I am already running out of stack on 4X4 puzzles which have with the first constraint applied only 6^4 possibilities. You can query this problem with :

problems(2,Field),binary_sudoku(Field).

code here :

:-use_module(library(clpfd)).

valid_row(Row) :-
    Row ins 0..1,
    length(Row,L),
    sum(Row,#=,L/2).

matrixNth1(Matr,X,Y,El) :-
    nth1(Y,Matr,CurRow),
    nth1(X,CurRow,El).

all_diff([]).
all_diff([X|Y]) :-
    maplist(dif(X),Y),
    all_diff(Y).


valid(_,1,1).
valid(Rows,1,Y) :-
    length(Rows,Y).
valid(Rows,X,1) :-
    length(Rows,X).
valid(Rows,X,X) :-
    length(Rows,X).

valid(Rows,X,Y) :-
    matrixNth1(Rows,X,Y,0).
valid(Rows,X,Y):-
    AboveY is Y-1,
    matrixNth1(Rows,X,AboveY,0).
valid(Rows,X,Y):-
    BelowY is Y+1,
    matrixNth1(Rows,X,BelowY,0).
valid(Rows,X,Y):-
    LeftX is X-1,
    matrixNth1(Rows,LeftX,Y,0).
valid(Rows,X,Y):-
    RightX is X+1,
    matrixNth1(Rows,RightX,Y,0).

binary_sudoku(Rows) :-
    length(Rows,Height),
    transpose(Rows,Cols),
    length(Cols,Height),
    maplist(valid_row,Rows),
    foreach(between(1,Height,X),foreach(between(1,Height,Y),valid(Rows,X,Y))),
    all_diff(Rows),all_diff(Cols).


problems(1,[[_,_],[_,_]]).

problems(2,[[_,_,_,_],[_,_,_,_],[_,_,_,_],[_,_,_,_]]).

Answer 1

There are a couple of problems with your code. If you are really a Prolog/clpfd beginner, you might consider starting with something a bit easier.

The idea of CLP programs is to first set up the constraints (which are deterministic), and then do the search for solutions (which is nondeterministic). In your code, valid_row/1 and all_diff/1 can be considered constraints, but valid/3 is nondeterministic, and makes lots of choices.

So the first change you should make is to move your call to valid/3 to the very end.

Then you should change valid/3 so that it systematically enumerates all possible variable assignments. With your current code, and a 3x3 matrix, the first few solutions that

foreach(between(1,Height,X),foreach(between(1,Height,Y),valid(Rows,X,Y)))

generates are

[[A, 0, C], [0, 0, 0], [G, 0, I]]
[[A, 0, C], [0, 0, 0], [G, 0, 0]]
[[A, 0, C], [0, 0, 0], [G, 0, I]]
[[A, 0, C], [0, 0, 0], [G, 0, I]]
[[A, 0, 0], [0, 0, F], [G, 0, I]]
...

That's bad, because they still contain variables, and they are not even mutually exclusive. Both of which means that you are visiting the same assignments over and over again. Rewrite it such that every time it succeeds, all variables are set to 0/1, and all solutions are different. Then you will be sure to traverse the search space only once, and you might have a chance to find a solution in reasonable time.

Answer 2

Here is a compact solution in ECLiPSe (a Prolog with constraint and modelling extensions, http://eclipseclp.org). It uses sum-constraints for the number of 0/1s per row/column, sequence/4 constraints for the no-three-1s condition, and lex_ne/2 to enforce difference between rows. The solution search is done by the labeling/1 call at the end. Also, Matrix-notation is used, which is more convenient than lists in such settings.

:- lib(gfd).

solve(Name, Mat) :-
    problem(Name, Mat),
    dim(Mat, [N,N]),
    Mat #:: 0..1,
    N #= 2*K,
    ( for(I,1,N), param(Mat,K,N) do
        sum(Mat[I,1..N]) #= K,
        sum(Mat[1..N,I]) #= K,
        sequence(1, 2, 3, Mat[I,1..N]),
        sequence(1, 2, 3, Mat[1..N,I]),
        ( for(J,I+1,N), param(Mat,I,N) do
            lex_ne(Mat[I,1..N], Mat[J,1..N]),
            lex_ne(Mat[1..N,I], Mat[1..N,J])
        )
    ),
    labeling(Mat).

problem(2, [](
    [](_,1,0,_,_,_,_,0,_,_,0,_),
    [](_,1,1,_,_,1,_,_,_,_,_,_),
    [](_,_,_,_,_,_,_,_,1,_,_,0),
    [](_,_,0,0,_,_,_,_,_,_,_,0),
    [](_,_,_,_,_,_,1,1,_,0,_,_),
    [](_,1,_,0,_,1,1,_,_,_,1,_),
    [](_,_,_,_,_,_,_,_,1,_,_,_),
    [](1,_,_,1,_,_,_,_,_,_,0,_),
    [](_,1,_,_,_,_,_,_,0,_,_,_),
    [](_,_,_,_,_,_,_,0,_,_,_,_),
    [](1,_,_,_,_,_,_,_,_,_,_,1),
    [](_,1,_,1,_,_,_,_,_,0,0,_))).

This brings up the (unique) solution quickly:

?- solve(2, M).
M = []([](1, 1, 0, 0, 1, 0, 1, 0, 1, 1, 0, 0),
       [](0, 1, 1, 0, 0, 1, 0, 1, 0, 0, 1, 1),
       [](1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0),
       [](1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0),
       [](0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1),
       [](0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0),
       [](1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1),
       [](1, 0, 1, 1, 0, 0, 1, 0, 0, 1, 0, 1),
       [](0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 1, 0),
       [](0, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0),
       [](1, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1),
       [](0, 1, 0, 1, 0, 1, 0, 1, 1, 0, 0, 1))
Yes (0.03s cpu, solution 1, maybe more)