/* 

  Satisfiability Problem in Picat.

  From GLPK:s example sat.mod
  """
  SAT, Satisfiability Problem

  Written in GNU MathProg by Andrew Makhorin <mao@mai2.rcnet.ru>
  """

  There is 20160 optimal solutions (with 1 unsats).

  This Picat model was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Picat page: http://www.hakank.org/picat/

*/


import cp. % 0.016s for finding optimal solution
% import sat. % 5.78s
% import mip. % ??


main => time2(go).

go =>
  m(M),
  n(N),
  c(C),
  sat(M,N,C,optimize,_Unsat,_X,_Y),
  nl.

go2 =>
  m(M),
  n(N),
  c(C),
  time2(sat(M,N,C,optimize,Unsat,_X,_Y)),

  printf("Counting solutions where Unsat = %d\n", Unsat),
  All=findall(_, sat(M,N,C,all,Unsat,X,Y)),
  println(numsols=All.length),

  nl.


sat(M,N, C, Type, Unsat, X, Y) => 

  % """
  % main variables
  % 
  % To solve the satisfiability problem means to determine all variables
  % x[j] such that conjunction of all clauses C[1] and ... and C[m] takes
  % on the value true, i.e. all clauses are satisfied.
  %
  % Let the clause C[i] be (t or t' or ... or t''), where t, t', ..., t''
  % are either variables or their negations. The condition of satisfying
  % C[i] can be most naturally written as:
  %
  %    t + t' + ... + t'' >= 1,                                       (1)
  %
  % where t, t', t'' have to be replaced by either x[j] or (1 - x[j]).
  % The formulation (1) leads to the mip problem with no objective, i.e.
  % to a feasibility problem.
  %
  % Another, more practical way is to write the condition for C[i] as:
  %
  %    t + t' + ... + t'' + y[i] >= 1,                                (2)
  %
  % where y[i] is an auxiliary binary variable, and minimize the sum of
  % y[i]. If the sum is zero, all y[i] are also zero, and therefore all
  % clauses are satisfied. If the sum is minimal but non-zero, its value
  % shows the number of clauses which cannot be satisfied.
  % """
  X = new_list(N),
  X :: 0..1,

  % auxiliary variables
  Y = new_list(M),
  Y :: 0..1,

  % number of unsatisfied clauses 
  Unsat #= sum(Y),

  % the condition (2)
  foreach(I in 1..M)
     Len = C[I].length,
     Tmp = new_array(Len),
     Tmp :: 0..1,
     foreach(J in 1..Len) 
        CIJ = C[I,J],
        if CIJ > 0 then
          Tmp[J] #= X[CIJ]
        else  
          Tmp[J] #= 1 - X[-CIJ]
        end
     end,
     sum([T + Y[I] : T in Tmp]) #>= 1 % This is faster: 0.016s

     % This approach takes much longer: 0.21s.
     % Note: We must use X[JJ] in both branches, otherwise Picat complains about 
     % "error(type_args(out_of_bound,-1"
     % sum([cond(J #> 0, X[JJ], 1 - X[JJ]) + Y[I] : J in C[I], JJ = abs(J)]) #>= 1
  end,

  if Type = optimize then
    solve($[ff,split,min(Unsat)], X ++ Y),

    println(x=X),
    println(xsat=[I : I in 1..N, X[I] = 1]),
    println(y=Y),
    println(unsat=Unsat),
    println(unsats=[I : I in 1..M, Y[I] = 1]),
    nl
  else 
    solve($[ff,split], X ++ Y)
  end.

%
% data
% 

% """
% These data correspond to the instance hole6 (pigeon hole problem for
% 6 holes) from SATLIB, the Satisfiability Library, which is part of
% the collection at the Forschungsinstitut fuer anwendungsorientierte
% Wissensverarbeitung in Ulm Germany */
%
% The optimal solution is 1 (one clause cannot be satisfied)
% """
%
% number of clauses
m(M) => M = 133.
% number of variables
n(N) => N = 42.

% """
% clauses; each clause C[i], i = 1, ..., m, is disjunction of some
% variables or their negations; in the data section each clause is
% coded as a set of indices of corresponding variables, where negative
% indices mean negation; for example, the clause (x3 or not x7 or x11)
% is coded as the set { 3, -7, 11 }
% """
c(C) =>
C = [
[-1,-7],
[-1,-13],
[-1,-19],
[-1,-25],
[-1,-31],
[-1,-37],
[-7,-13],
[-7,-19],
[-7,-25],
[-7,-31],
[-7,-37],
[-13,-19],
[-13,-25],
[-13,-31],
[-13,-37],
[-19,-25],
[-19,-31],
[-19,-37],
[-25,-31],
[-25,-37],
[-31,-37],
[-2,-8],
[-2,-14],
[-2,-20],
[-2,-26],
[-2,-32],
[-2,-38],
[-8,-14],
[-8,-20],
[-8,-26],
[-8,-32],
[-8,-38],
[-14,-20],
[-14,-26],
[-14,-32],
[-14,-38],
[-20,-26],
[-20,-32],
[-20,-38],
[-26,-32],
[-26,-38],
[-32,-38],
[-3,-9],
[-3,-15],
[-3,-21],
[-3,-27],
[-3,-33],
[-3,-39],
[-9,-15],
[-9,-21],
[-9,-27],
[-9,-33],
[-9,-39],
[-15,-21],
[-15,-27],
[-15,-33],
[-15,-39],
[-21,-27],
[-21,-33],
[-21,-39],
[-27,-33],
[-27,-39],
[-33,-39],
[-4,-10],
[-4,-16],
[-4,-22],
[-4,-28],
[-4,-34],
[-4,-40],
[-10,-16],
[-10,-22],
[-10,-28],
[-10,-34],
[-10,-40],
[-16,-22],
[-16,-28],
[-16,-34],
[-16,-40],
[-22,-28],
[-22,-34],
[-22,-40],
[-28,-34],
[-28,-40],
[-34,-40],
[-5,-11],
[-5,-17],
[-5,-23],
[-5,-29],
[-5,-35],
[-5,-41],
[-11,-17],
[-11,-23],
[-11,-29],
[-11,-35],
[-11,-41],
[-17,-23],
[-17,-29],
[-17,-35],
[-17,-41],
[-23,-29],
[-23,-35],
[-23,-41],
[-29,-35],
[-29,-41],
[-35,-41],
[-6,-12],
[-6,-18],
[-6,-24],
[-6,-30],
[-6,-36],
[-6,-42],
[-12,-18],
[-12,-24],
[-12,-30],
[-12,-36],
[-12,-42],
[-18,-24],
[-18,-30],
[-18,-36],
[-18,-42],
[-24,-30],
[-24,-36],
[-24,-42],
[-30,-36],
[-30,-42],
[-36,-42],
[6,5,4,3,2,1],
[12,11,10,9,8,7],
[18,17,16,15,14,13],
[24,23,22,21,20,19],
[30,29,28,27,26,25],
[36,35,34,33,32,31],
[42,41,40,39,38,37]
].