The CLP Procedure

Example 3.13 Balanced Incomplete Block Design

Balanced incomplete block design (BIBD) generation is a standard combinatorial problem from design theory. The concept was originally developed in the design of statistical experiments; applications have expanded to other fields, such as coding theory, network reliability, and cryptography. A BIBD is an arrangement of $v$ distinct objects into $b$ blocks such that each block contains exactly $k$ distinct objects, each object occurs in exactly $r$ different blocks, and every two distinct objects occur together in exactly $\lambda $ blocks. A BIBD is therefore specified by its parameters $(v,b,r,k,\lambda )$. It can be proved that when a BIBD exists, its parameters must satisfy the conditions ${rv = bk}$, ${\lambda (v-1) = r(k-1)}$, and ${b \ge v}$, but these conditions are not sufficient to guarantee the existence of a BIBD (Prestwich, 2001). For instance, the parameters $(15, 21, 7, 5, 2)$ satisfy the preceding conditions, but a BIBD with these parameters does not exist. Computational methods for BIBD generation generally suffer from combinatorial explosion, in part because of the large number of symmetries: given any solution, any two objects or blocks can be exchanged to obtain another solution.

This example demonstrates how to express a BIBD problem as a CSP and how to use lexicographic ordering constraints to break symmetries. The most direct CSP model for BIBD, as described in Meseguer and Torras (2001), represents a BIBD as a ${v \times b}$ matrix $X$. Each matrix entry is a Boolean decision variable $X_{i,c}$ that satisfies $X_{i,c} = 1$ if and only if block $c$ contains object $i$. The condition that each object occurs in exactly $r$ blocks (or, equivalently, that there are $r$ 1s per row) can be expressed as $v$ linear constraints:

\[  \sum _{c=1}^ b X_{i,c} = r \quad \mr {for} \;  i=1,\dots ,v  \]

Alternatively, you can use global cardinality constraints to ensure that there are exactly $b - r$ 0s and $r$ 1s in $X_{i,1}$,…, $X_{i,b}$ for each object $i$:

\[  \mr {gcc}(X_{i,1},\dots , X_{i,b}) = ((0, 0, b-r) (1, 0, r)) \quad \mr {for} \;  i=1,\dots ,v  \]

Similarly, the condition that each block contains exactly $k$ objects (there are $k$ 1s per column) can be specified by the following constraints:

\[  \mr {gcc}(X_{1,c},\dots , X_{v,c}) = ((0, 0, v-k) (1, 0, k)) \quad \mr {for} \;  c=1,\dots ,b  \]

To enforce the final condition that every two distinct objects occur together in exactly $\lambda $ blocks (equivalently, that the scalar product of every pair of rows is equal to $\lambda $), you can introduce auxiliary variables $P_{i,j,c}$ for every ${i < j}$ that indicate whether objects $i$ and $j$ both occur in block $c$. The REIFY constraint

\[  \mr {reify} \,  P_{i,j,c}: (X_{i,c} + X_{j,c} = 2)  \]

ensures that ${P_{i,j,c} =1}$ if and only if block $c$ contains both objects $i$ and $j$. The following constraints ensure that the final condition holds:

\[  \mr {gcc}(P_{i,j,1},\dots , P_{i,j,b}) = ((0, 0, b-\lambda ) (1, 0, \lambda )) \quad \mr {for} \;  i=1,\dots ,v-1 \:  \mr {and} \:  j=i+1,\dots ,v  \]

The objects and the blocks are interchangeable, so the matrix $X$ has total row symmetry and total column symmetry. Because of the constraints on the rows, no pair of rows can be equal unless ${r = \lambda }$. To break the row symmetry, you can impose strict lexicographical ordering on the rows of $X$ as follows:

\[  (X_{i,1},\dots , X_{i,b}) <_{\mr {lex}} (X_{i-1,1},\dots , X_{i-1,b}) \quad \mr {for} \;  i=2,\dots ,v  \]

To break the column symmetry, you can impose lexicographical ordering on the columns of $X$ as follows:

\[  (X_{1,c},\dots , X_{v,c}) \le _{\mr {lex}} (X_{1,c-1},\dots , X_{v,c-1}) \quad \mr {for} \;  c=2,\dots ,b  \]

The following SAS macro incorporates all the preceding constraints. For specified parameters $(v,b,r,k,\lambda )$, the macro either finds BIBDs or proves that a BIBD does not exist.

%macro bibd(v, b, r, k, lambda, out=bibdout);
   /* Arrange v objects into b blocks such that:
         (i) each object occurs in exactly r blocks,
         (ii) each block contains exactly k objects,
         (iii) every pair of objects occur together in exactly lambda blocks.

      Equivalently, create a binary matrix with v rows and b columns,
      with r 1s per row, k 1s per column,
      and scalar product lambda between any pair of distinct rows.
   */

   /* Check necessary conditions */
   %if (%eval(&r * &v) ne %eval(&b * &k)) or
      (%eval(&lambda * (&v - 1)) ne %eval(&r * (&k - 1))) or
      (&v > &b) %then %do;
      %put BIBD necessary conditions are not met.;
      %goto EXIT;
   %end;

   proc clp out=&out(keep=x:) domain=[0,1] varselect=FIFO;
      /* Decision variables: */
      /* Decision variable X_i_c = 1 iff object i occurs in block c. */
      var (
           %do i=1 %to &v;
              x&i._1-x&i._&b.
           %end;
          ) = [0,1];

      /* Mandatory constraints: */
      /* (i) Each object occurs in exactly r blocks. */
      %let q = %eval(&b.-&r.);  /* each row has &q 0s and &r 1s */
      %do i=1 %to &v;
         gcc( x&i._1-x&i._&b. ) = ((0,0,&q.) (1,0,&r.));
      %end;

      /* (ii) Each block contains exactly k objects. */
      %let h = %eval(&v.-&k.);  /* each column has &h 0s and &k 1s */
      %do c=1 %to &b;
         gcc(
             %do i=1 %to &v;
                x&i._&c.
             %end;
            ) = ((0,0,&h.) (1,0,&k.));
      %end;

      /* (iii) Every pair of objects occurs in exactly lambda blocks. */
      %let t = %eval(&b.-&lambda.);
      %do i=1 %to %eval(&v.-1);
         %do j=%eval(&i.+1) %to &v;
            /* auxiliary variable p_i_j_c =1 iff both i and j occur in c */
            var ( p&i._&j._1-p&i._&j._&b. ) = [0,1];
            %do c=1 %to &b;
               reify p&i._&j._&c.: (x&i._&c. + x&j._&c. = 2);
            %end;

            gcc(p&i._&j._1-p&i._&j._&b.) = ((0,0,&t.) (1,0,&lambda.));
         %end;
      %end;

      /* Symmetry breaking constraints: */
      /* Break row symmetry via lexicographic ordering constraints. */
      %do i = 2 %to &v.;
         %let i1 = %eval(&i.-1);
         lexico( (x&i._1-x&i._&b.) LEX_LT (x&i1._1-x&i1._&b.) );
      %end;

      /* Break column symmetry via lexicographic ordering constraints. */
      %do c = 2 %to &b.;
         %let c1 = %eval(&c.-1);
         lexico( ( %do i = 1 %to &v.;
                      x&i._&c.
                   %end; )
                 LEX_LE
                 ( %do i = 1 %to &v.;
                      x&i._&c1.
                   %end; ) );
      %end;
   run;
   %put &_orclp_;
%EXIT:
%mend bibd;

The following statement invokes the macro to find a BIBD design for the parameters $(15,15,7,7,3)$:


%bibd(15,15,7,7,3);

The output is displayed in Output 3.13.1.

Output 3.13.1: Balanced Incomplete Block Design for (15,15,7,7,3)

Balanced Incomplete Block Design Problem
(15, 15, 7, 7, 3)

Obs Block1 Block2 Block3 Block4 Block5 Block6 Block7 Block8 Block9 Block10 Block11 Block12 Block13 Block14 Block15
1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0
2 1 1 1 0 0 0 0 1 1 1 1 0 0 0 0
3 1 1 0 1 0 0 0 1 0 0 0 1 1 1 0
4 1 0 1 0 1 0 0 0 1 0 0 1 1 0 1
5 1 0 0 1 0 1 0 0 0 1 1 1 0 0 1
6 1 0 0 0 1 0 1 0 0 1 1 0 1 1 0
7 1 0 0 0 0 1 1 1 1 0 0 0 0 1 1
8 0 1 1 0 0 0 1 0 0 1 0 1 0 1 1
9 0 1 0 1 0 0 1 0 1 0 1 0 1 0 1
10 0 1 0 0 1 1 0 1 0 1 0 0 1 0 1
11 0 1 0 0 1 1 0 0 1 0 1 1 0 1 0
12 0 0 1 1 1 0 0 1 0 0 1 0 0 1 1
13 0 0 1 1 0 1 0 0 1 1 0 0 1 1 0
14 0 0 1 0 0 1 1 1 0 0 1 1 1 0 0
15 0 0 0 1 1 0 1 1 1 1 0 1 0 0 0