044: Steiner triple systems

/*

  Steiner triplets in Comet.

  http://www.probp.com/examples/clpset/steiner.pl
  """
  The ternary Steiner problem of order n is to find n(n-1)/6 sets of elements 
  in {1,2,...,n} such that each set contains three elements and any two 
  sets have at most one element in common.

  For example, the following shows a solution for size n=7:

      {1,2,3}, {1,4,5}, {1,6,7}, {2,4,6}, {2,5,7}, {3,4,7}, {3,5,6}

  Problem taken from:
  C. Gervet: Interval Propagation to Reason about Sets: Definition and 
             Implementation of a PracticalLanguage,  
             Constraints, An International Journal, vol.1, pp.191-246, 1997.
  """

  Also see:
  - http://mathworld.wolfram.com/SteinerTripleSystem.html
  - http://en.wikipedia.org/wiki/Steiner_system


  Note: This model uses arrays of booleans as an representation of sets.


  This Comet model was created by Hakan Kjellerstrand (hakank@gmail.com)
  Also, see my Comet page: http://www.hakank.org/comet 

 */

% Licenced under CC-BY-4.0 : http://creativecommons.org/licenses/by/4.0/

import cotfd;

int t0 = System.getCPUTime();

int n = 7;
int nb = (n * (n-1)) / 6;

cout << "n: " << n << " nb: " << nb << endl;

Solver m();

var{bool} sets[1..nb, 1..n](m);

Integer num_solutions(0);

exploreall {
// explore {

  forall(i in 1..nb) 
    m.post(sum(j in 1..n) sets[i,j] == 3);

  forall(i in 1..nb, j in i+1..nb) {
    m.post( sum(k in 1..n) ( sets[i,k] && sets[j,k]) <= 1);    
  }

  // symmetry breaking
  forall(i in 1..nb, j in i+1..nb)
    m.post(lexleq(all(k in 1..n) sets[i,k], all(k in 1..n) sets[j,k]) );

} using {

  label(m);

  num_solutions := num_solutions + 1;


  forall(i in 1..nb) {
    forall(j in 1..n) {
      if (sets[i,j]) 
        cout << j << " ";
    }
    cout << endl;
  }
  cout << endl;

}

cout << "\nnum_solutions: " << num_solutions << endl;
    
int t1 = System.getCPUTime();
cout << "time:      " << (t1-t0) << endl;
cout << "#choices = " << m.getNChoice() << endl;
cout << "#fail    = " << m.getNFail() << endl;
cout << "#propag  = " << m.getNPropag() << endl;