import static choco.Choco.*;
import choco.cp.model.CPModel;
import choco.cp.solver.CPSolver;
import choco.cp.solver.search.integer.varselector.StaticVarOrder;
import choco.cp.solver.search.integer.varselector.MinDomain;
import choco.kernel.model.Model;
import choco.kernel.solver.Solver;
import choco.kernel.model.variables.integer.IntegerVariable;
import choco.kernel.solver.variables.integer.IntDomainVar;

public class BIBD {    
  
  public static void main(String[] args) {

      CPModel model = new CPModel();
      int v      = Integer.parseInt(args[0]);  // # points
      int b      = Integer.parseInt(args[1]);  // # blocks
      int r      = Integer.parseInt(args[2]);  // every point occurs r times
      int k      = Integer.parseInt(args[3]);  // size (arity) of a block
      int lambda = 1;                          // every pair of points occurs lambda times

      IntegerVariable[][] A = new IntegerVariable[v][v];
      //
      // A[i][j] = x <-> the pair of points (i,j) is in block x
      // NOTE: this forces each pair of points to occur once and once only!
      //       i.e. lambda = 1.
      //

      for (int i=0;i<v-1;i++)
	  for (int j=0;j<v;j++)
	      A[i][j] = A[j][i] = makeIntVar("A",0,b-1);
      for (int i=0;i<v;i++) A[i][i] = constant(0);

      IntegerVariable[][] block = new IntegerVariable[b][v];
      //
      // block[x][i] = 1 <-> block x contains point i
      //

      IntegerVariable[][] occurs = new IntegerVariable[v][b];
      //
      // occurs[i][x] = 1 <-> i occurs in block x
      //

      for (int x=0;x<b;x++)
	  for (int i=0;i<v;i++)
	      occurs[i][x] = block[x][i] = makeIntVar("block",0,1);
      
      for (int i=0;i<v-1;i++)
	  for (int j=i+1;j<v;j++)
	      for (int x=0;x<b;x++)
		  model.addConstraint(ifOnlyIf(eq(A[i][j],x),and(eq(block[x][i],1),eq(block[x][j],1))));
      //
      // channel between A and block
      // A[i][j] = k <-> block[k][i] && block[k][j]
      // i.e. pair (i,j) \in block k iff block[k][i] is true and block[k][j] is true
      //

      /*
      for (int i=0;i<v-1;i++)
	  for (int j=i+1;j<v;j++)
	      for (int x=0;x<b;x++)
		  model.addConstraint(ifOnlyIf(eq(A[i][j],x),and(eq(occurs[i][x],1),eq(occurs[j][x],1))));
      */

      /*
      for (int i=0;i<v-1;i++)
	  for (int j=i+1;j<v;j++){
	      model.addConstraint(domainChanneling(A[i][j],occurs[i]));
	      model.addConstraint(domainChanneling(A[i][j],occurs[j]));
	  }
      */

      for (int x=0;x<b;x++) model.addConstraint(eq(sum(block[x]),k));
      //
      // every block contains k points
      //

      for (int i=0;i<v;i++) model.addConstraint(eq(sum(occurs[i]),r));
      //
      // every point occurs in r blocks
      //

      if (args.length > 4 && Boolean.parseBoolean(args[4]))
	  for (int i=0;i<v-1;i++)
	      model.addConstraint(lex(A[i],A[i+1]));
      //
      // symmetry breaker, a redundant constraint
      //

      CPSolver solver = new CPSolver();
      solver.read(model);

      IntDomainVar[] decision = new IntDomainVar[v*(v-1)/2];
      for (int i=0,count=0;i<v-1;i++)
	  for (int j=i+1;j<v;j++,count++)
	      decision[count] = solver.getVar(A[i][j]);
      solver.setVarIntSelector(new MinDomain(solver,decision));
      // solver.setVarIntSelector(new StaticVarOrder(decision));
      //
      // decision variables are the pairs to be assigned to a block
      //

      System.out.println(solver.solve() +"   nd: "+ solver.getNodeCount() +"   ms: "+ solver.getTimeCount());
      solver.printRuntimeStatistics();
      for (int i=0;i<b;i++){
	  for (int j=0;j<v;j++)
	      if(solver.getVar(block[i][j]).getVal() == 1) System.out.print(j +" ");
	  for (int j=0;j<v;j++)
	      System.out.print(solver.getVar(block[i][j]).getVal());
	  System.out.println();
      }

      System.out.println();

      for (int i=0;i<v;i++){
	  for (int j=0;j<v;j++)
	      if (i!=j) System.out.print(solver.getVar(A[i][j]).getVal() +" ");
	      else System.out.print("- ");
	  System.out.println();
      }
  }
}