Logic Based Puzzles (oclpe01)

/***************************************************************/
/*                                                             */
/*          S A S   S A M P L E   L I B R A R Y                */
/*                                                             */
/*    NAME: oclpe01                                            */
/*   TITLE: Logic Based Puzzles (oclpe01)                      */
/* PRODUCT: OR                                                 */
/*  SYSTEM: ALL                                                */
/*    KEYS: OR                                                 */
/*   PROCS: OPTMODEL                                           */
/*    DATA:                                                    */
/*                                                             */
/* SUPPORT:                             UPDATE:                */
/*     REF:                                                    */
/*    MISC: Example 1 from the CLP solver chapter of the       */
/*          Mathematical Programming book.                     */
/*                                                             */
/***************************************************************/

/* Given a sudoku problem */
data Indata;
   input C1-C9;
   datalines;
. . 5 . . 7 . . 1
. 7 . . 9 . . 3 .
. . . 6 . . . . .
. . 3 . . 1 . . 5
. 9 . . 8 . . 2 .
1 . . 2 . . 4 . .
. . 2 . . 6 . . 9
. . . . 4 . . 8 .
8 . . 1 . . 5 . .
;

/* Print sudoku */
%macro print_sudoku(dsn);
   goptions hsize=4in vsize=4in;
   data sdk;
      set &dsn;
      array c{9} C1-C9;
      drop C1-C9 i j;
      length function style color $8 text $1;
      retain xsys ysys hsys "3";
      i = _N_;
      do j = 1 to 9;
         function="move"; x=100*j/9; y=100*(9-i)/9; output;
         function="bar"; x=100*(j-1)/9; y=100*(10-i)/9;
            size=0.5; line=1; color="black"; output;
         function="move"; x=100*(j-1)/9+4; y=100*(10-i)/9; output;
         function="label"; text=put(c[j],1.);
         if text='.' then text='';
         position="9"; style=""; color="black"; size=5;  output;
      end;
      if _N_=9 then do;
         /* draw the out frame */
         function="move"; x=100; y=100; output;
         function="bar"; x=0; y=0; size=2; line=0; color="black"; output;
         /* Thick Vertical lines */
         function="move"; x=100/3; y=0; output;
         function="draw"; x=100/3; y=100; size=0.6; line=1; output;
         function="move"; x=200/3; y=0; output;
         function="draw"; x=200/3; y=100; output;
         /* Thick Horizontal lines */
         function="move"; x=0; y=100/3; output;
         function="draw"; x=100; y=100/3;  output;
         function="move"; x=0; y=200/3; output;
         function="draw"; x=100; y=200/3; output;
      end;
   run;
   proc ganno annotate=sdk;
   run;
%mend print_sudoku;

%print_sudoku(indata);

proc optmodel;
   /* Declare variables */
   set ROWS = 1..9;
   set COLS = ROWS; /* Use an alias for convenience and clarity */
   var X {ROWS, COLS} >= 1 <= 9 integer;

   /* Nine row constraints */
   con RowCon {i in ROWS}:
      alldiff({j in COLS} X[i,j]);

   /* Nine column constraints */
   con ColCon {j in COLS}:
      alldiff({i in ROWS} X[i,j]);

   /* Nine 3x3 block constraints */
   con BlockCon {s in 0..2, t in 0..2}:
      alldiff({i in 3*s+1..3*s+3, j in 3*t+1..3*t+3} X[i,j]);

   /* Fix variables to cell values */
   /* X[i,j] = c[i,j] if c[i,j] is not missing */
   num c {ROWS, COLS};
   read data indata into [_N_] {j in COLS} <c[_N_,j]=col('C'||j)>;
   for {i in ROWS, j in COLS: c[i,j] ne .}
      fix X[i,j] = c[i,j];

   solve;
   create data sudoku_dense from [i]=ROWS {j in COLS} <col('C'||j)=X[i,j]>;
quit;

%print_sudoku(sudoku_dense);

data Raw;
   input C1-C12;
   datalines;
3  .  .  1  5  4  .  .  1  .  9  5
.  1  .  .  3  .  .  .  .  1  3  6
.  .  4  .  .  3  .  8  .  .  2  .
5  .  .  1  .  .  9  2  5  .  .  1
.  9  .  .  5  .  .  5  .  .  .  .
5  8  1  .  .  9  .  .  3  .  6  .
.  5  .  8  .  .  2  .  .  5  5  3
.  .  .  .  5  .  .  6  .  .  1  .
2  .  .  5  1  5  .  .  5  .  .  9
.  6  .  .  4  .  1  .  .  3  .  .
1  5  1  .  .  .  .  5  .  .  5  .
5  5  .  4  .  .  3  1  6  .  .  8
;

%macro print_piday(dsn);
   goptions reset=title;
   goptions hsize=4in vsize=4in;
   data piday;
      set &dsn;
      array c{12} C1-C12;
      drop C1-C12 i j;
      length function style color $8 text $1;
      retain xsys ysys hsys "3";
      if _n_=1 then do;
         /* draw the color bars */
         function="move"; x=100; y=100; output; /* whole */
         function="bar"; x=0; y=0; size=1; line=0; color="cxffedaf"; style='s';
            output;
         function="move"; x=300/12; y=100; output; /* left up */
         function="bar"; x=0; y=700/12; color="cxfde28a"; output;
         function="move"; x=100; y=300/12; output; /* right down */
         function="bar"; x=600/12; y=0; output;
         function="move"; x=100; y=100; output; /* right up */
         function="bar"; x=900/12; y=700/12; color="cxffdc61"; output;
         function="move"; x=600/12; y=300/12; output; /* left down */
         function="bar"; x=0; y=0; output;
         function="move"; x=600/12; y=1000/12; output; /* purple left up */
         function="bar"; x=200/12; y=700/12; color="cxbd85ff"; output;
         function="move"; x=900/12; y=700/12; output; /* purple right down */
         function="bar"; x=700/12; y=100/12; output;
         function="move"; x=1000/12; y=1000/12; output; /* purple right up */
         function="bar"; x=600/12; y=700/12; color="cxd3aefe"; output;
         function="move"; x=500/12; y=700/12; output; /* purple left down */
         function="bar"; x=300/12; y=100/12; output;
         function="move"; x=700/12; y=300/12; color="cxffedaf"; output;
         function="bar"; x=500/12; y=100/12; output;
      end;
      i = _N_;
      do j=1 to 12;
         function="move"; x=100*j/12; y=100*(12-i)/12; output;
         function="bar"; x=100*(j-1)/12; y=100*(13-i)/12;
            size=0.5; line=1; color="black";  style='e'; output;
         function="move"; x=100*(j-1)/12+3; y=100*(13-i)/12-1; output;
         function="label"; text=put(c[j], 1.);
         if text='.' then text='';
         position="F"; style=""; color="black"; size=5;  output;
      end;
      if _n_=12 then do;
         /* draw the out frame */
         function="move"; x=100; y=100; output;
         function="bar"; x=0; y=0; size=2; line=0; position="9"; color="black";
            output;
         /* Horizontal separating lines */
         function="move"; x=1000/12; y=1000/12; output;
         function="draw"; x=200/12; y=1000/12; size=0.6; line=1; output;
         function="move"; x=100; y=700/12; output;
         function="draw"; x=0; y=700/12; output;
         function="move"; x=300/12; y=300/12; output;
         function="draw"; x=0; y=300/12; output;
         function="move"; x=100; y=300/12; output;
         function="draw"; x=900/12; y=300/12; output;
         function="move"; x=900/12; y=100/12; output;
         function="draw"; x=300/12; y=100/12; output;
         /* Vertical separating lines */
         function="move"; x=200/12; y=1000/12; output;
         function="draw"; x=200/12; y=700/12; output;
         function="move"; x=300/12; y=100; output;
         function="draw"; x=300/12; y=1000/12; output;
         function="move"; x=300/12; y=700/12; output;
         function="draw"; x=300/12; y=100/12; output;
         function="move"; x=500/12; y=700/12; output;
         function="draw"; x=500/12; y=100/12; output;
         function="move"; x=600/12; y=1000/12; output;
         function="draw"; x=600/12; y=700/12; output;
         function="move"; x=600/12; y=100/12; output;
         function="draw"; x=600/12; y=0; output;
         function="move"; x=700/12; y=700/12; output;
         function="draw"; x=700/12; y=100/12; output;
         function="move"; x=900/12; y=100; output;
         function="draw"; x=900/12; y=1000/12; output;
         function="move"; x=900/12; y=700/12; output;
         function="draw"; x=900/12; y=100/12; output;
         function="move"; x=1000/12; y=1000/12; output;
         function="draw"; x=1000/12; y=700/12; output;
      end;
   run;
   proc ganno annotate=piday;
   run;
%mend print_piday;

%print_piday(Raw);

proc optmodel;
   set ROWS = 1..12;
   /* These declarations are inexpensive and improve clarity: */
   set COLS = ROWS, REGIONS = ROWS, CELLS = ROWS cross COLS;

   /* specify a 12x12 array of region identifiers.
      The spacing is just to make the regions easier to visualize. */
   num region{CELLS} = [
       1  1  1    2  2  2  2  2  2    3  3  3
       1  1  1    2  2  2  2  2  2    3  3  3
       1  1   4  4  4  4    5  5  5  5   3  3
       1  1   4  4  4  4    5  5  5  5   3  3
       1  1   4  4  4  4    5  5  5  5   3  3
       6  6  6   7  7   8  8   9  9  10 10 10
       6  6  6   7  7   8  8   9  9  10 10 10
       6  6  6   7  7   8  8   9  9  10 10 10
       6  6  6   7  7   8  8   9  9  10 10 10
      11 11 11   7  7   8  8   9  9  12 12 12
      11 11 11   7  7   8  8   9  9  12 12 12
      11 11 11 11 11 11     12 12 12 12 12 12 ];

   /* Each area must contain two 1's, two 3's, three 5's, no 7's,
      and one for each of other values from 1 to 9. */
                    /* 1 2 3 4 5 6 7 8 9 */
   num nTimes{1..9} = [2 1 2 1 3 1 0 1 1];
   /* For convenience, create a triplet set version of nTimes.
      In this model, GCC's lower and upper bounds are the same. */
   set N_TIMES = setof{ni in 1..9} <ni,nTimes[ni],nTimes[ni]>;

   /* The number assigned to the ith row and jth column. */
   var X {CELLS} >= 1 <= 9 integer;

   /* X[i,j] = c[i,j] if c[i,j] is not missing */
   num c {CELLS};
   read data raw into [_N_] {j in COLS} <c[_N_,j]=col('C'||j)>;
   for {<i,j> in CELLS: c[i,j] ne .}
      fix X[i,j] = c[i,j];

   con RowCon {i in ROWS}:
      gcc({j in COLS} X[i,j], N_TIMES);

   con ColCon {j in COLS}:
      gcc({i in ROWS} X[i,j], N_TIMES);

   con RegionCon {ri in REGIONS}:
      gcc({<i,j> in CELLS: region[i,j] = ri} X[i,j], N_TIMES);

   solve;
   /* Replicate typical PROC CLP output from PROC OPTMODEL arrays */
   create data pdsout from
      {<i,j> in ROWS cross COLS}<col('X_'||i||'_'||j)=X[i,j]>;
quit;

/* convert solution to matrix in dense format */
%macro pds_out;
   data pdsoutsq;
   set pdsout;
   array C{12};
      %do i = 1 %to 12;
         %do j = 1 %to 12;
            C[&j.] = X_&i._&j.;
         %end;

         output;
      %end;
      drop X:;
   run;

   proc print data=pdsoutsq;
      title "Pi Day Sudoku 2008";
   run;

%mend pds_out;

%pds_out;

%print_piday(pdsoutsq);

%macro magic(n);
   proc optmodel;
      num n = &n;
      /* magic constant */
      num sum = n*(n^2+1)/2;
      set ROWS = 1..n;
      set COLS = 1..n;

      /* X[i,j] = entry (i,j) */
      var X {ROWS, COLS} >= 1 <= n^2 integer;

      /* row sums */
      con RowCon {i in ROWS}:
         sum {j in COLS} X[i,j] = sum;

      /* column sums */
      con ColCon {j in COLS}:
         sum {i in ROWS} X[i,j] = sum;

      /* diagonal: upper left to lower right */
      con DiagCon:
         sum {i in ROWS} X[i,i] = sum;

      /* diagonal: upper right to lower left */
      con AntidiagCon:
         sum {i in ROWS} X[n+1-i,i] = sum;

      /* symmetry-breaking */
      con BreakRowSymmetry:
         X[1,1] + 1 <= X[n,1];
      con BreakDiagSymmetry:
         X[1,1] + 1 <= X[n,n];
      con BreakAntidiagSymmetry:
         X[1,n] + 1 <= X[n,1];

      con alldiff(X);

      solve with CLP / varselect=minrmaxc maxtime=3;
      create data magic&n from
         {<i,j> in ROWS cross COLS}<col('X_'||i||'_'||j)=X[i,j]>;
   quit;
%mend magic;

%magic(7)

%macro convert_to_dense(n);
   /* convert solution to matrix in dense format */
   data magic7_dense;
      set magic7;
      array C{7};
      %do i = 1 %to &n;
         %do j = 1 %to &n;
            C[&j] = X_&i._&j;
         %end;
         output;
      %end;
      drop X:;
   run;
%mend convert_to_dense;
%convert_to_dense(7);

/* Print Magic Square */
%macro print_msq(dsn);
   goptions hsize=3in vsize=3in;
   data msq;
      set &dsn;
      array c{7} C1-C7;
      drop C1-C7 i j;
      length function style color $8 text $2;
      retain xsys ysys hsys "3";
      i = _N_;
      do j=1 to 7;
         function="move"; x=100*j/7; y=100*(7-i)/7; output;
         function="bar"; x=100*(j-1)/7; y=100*(8-i)/7;
            size=0.5; line=0; color="black"; output;
         function="move"; x=100*(j-1)/7+4; y=100*(8-i)/7; output;
         function="label"; text=put(c[j], 2.);
         position="9"; style=""; color="black"; size=5; output;
      end;
   run;
   proc ganno annotate=msq;
   run;
%mend print_msq;
%print_msq(magic7_dense);