Sudoku Dancing Links Implementation

//Federico Rampazzo, 2011
#include <iostream>
#include <vector>
#include <cstdlib>
#include <ctime>

using namespace std;

class Sudoku{
  private: // Dancing Links toroidal structures
    class header;
    class node{
    public:
      int row;
      node * up, * down, * left, * right;
      header * head;
    };
    class header : public node{
    public:
      int size;
      bool visible;
      header * left, * right;
    };
    int u, u2, u3, rowsNumber, columnsNumber, solutionsNumber, maxSolutionsNumber;
    header * root;
    vector<header> columns;
    vector<vector<node> > rows;
    vector<int> solutions;

public:
  Sudoku(int dimension) :
  u(dimension),
  u2(u*u),
  u3(u2*u2),
  rowsNumber(u3*u2),
  rows(rowsNumber, vector<node>(4, node())),
  columnsNumber(u3*4+1),
  columns(columnsNumber, header()),
  solutions(u3, 0)
  { // create the 729x324 matrix
    root = &(columns[u3*4]);
    for (int i = 0; i < columnsNumber; i++) {
      columns[i].up = &(columns[i]);
      columns[i].down = &(columns[i]);
      columns[i].left = &(columns[(i+columnsNumber-1)%columnsNumber]);
      columns[i].right = &(columns[(i+1)%columnsNumber]);
      columns[i].visible = true;
      columns[i].size = 0;
    }
    for (int i = 0; i < rowsNumber; i++) {
      int r = i/(u3), c = i%u3/u2, p = i%u2;
      for (int j = 0; j < 4; j++) {
    int columnIndex;
    // http://www.stolaf.edu/people/hansonr/sudoku/exactcover.htm
    switch (j){
      case 0: // cell constraints
        columnIndex = u3*j + r*u2+c;
        break;
      case 1: //row constraints
        columnIndex = u3*j + r*u2+p;
        break;
      case 2: //column constraints
          columnIndex = u3*j + c*u2+p;
          break;
      case 3: //block constraints
          columnIndex = u3*j + (r/u*u+c/u)*u2+p;
          break;
    }
    rows[i][j].row = i;
    rows[i][j].left = &(rows[i][(j+3)%4]);
    rows[i][j].right = &(rows[i][(j+1)%4]);
    rows[i][j].down = rows[i][j].head = &(columns[columnIndex]);
    rows[i][j].up = columns[columnIndex].up;

    columns[columnIndex].up->down = columns[columnIndex].up = &(rows[i][j]);
    columns[columnIndex].size++;
      }
    }
  };

  void print(int * s){
    for (int i=0; i<u3; i++){
      if (i!=0 && i%(u*u2)==0){
    cout << "\n";
    for (int j=0; j<u; j++){
      for (int k=0; k<u*2+2*((j==0 || j==u-1) ? 1 : 2); k++) cout << "-";
      if (j!=u-1) cout << "+";
    }
      }
      if      (i%u2==0) cout << "\n";
      else if (i%u==0) cout << "  \u007C  ";

      cout << s[i]  << " ";
    }
    cout << "\n";
  }

  int solver(int * grid, int max=1, bool fill=true){
    int i = 0;
    solutionsNumber = 0;
    maxSolutionsNumber = max;
    for (int j = 0; j<u3; j++){
      if (grid[j]) {
    int rowIndex = j*u2+grid[j]-1;
    if (rows[rowIndex][0].head->visible==false ||
      rows[rowIndex][1].head->visible==false ||
      rows[rowIndex][2].head->visible==false ||
      rows[rowIndex][3].head->visible==false) {
      return -1; // unvalid puzzle
      break;
      }
      for (int k = 0; k<4; k++)
        cover(rows[rowIndex][k].head);
      solutions[i++] = rowIndex;
      }
    }
    backtrack(i);
    if (solutionsNumber==0) return 0;
    if (fill)
      for (int j=0; j<u3; j++)
    grid[solutions[j]/u2] = solutions[j]%u2+1;

    for (int j=i-1; j>=0; j--)
      for (int k=3; k>=0; k--)
    uncover(rows[solutions[j]][k].head);
      return solutionsNumber;
  }

  int generator (int * grid, int numbers=-1, int limit=1000, int randomness = 100){
    for (int i=0; i<u3; i++)
      grid[i] = 0;
    solver(grid);

    srand((unsigned)time(0));
    for (int i=0; i<randomness; i++){
      int box = rand()%3, begin = rand()%3, end = rand()%3;
      if (rand()%2==0){
    for (int j=0; j<u2; j++){
      int beginIndex = box*u2*u+begin*u2+j,
          endIndex = box*u2*u+end*u2+j;
      int t = grid[endIndex];
      grid[endIndex] = grid[beginIndex];
      grid[beginIndex] = t;
    }
      }else{
    for (int j=0; j<u2; j++){
      int beginIndex = box*u+j*u2+begin,
          endIndex = box*u+j*u2+end;
      int t = grid[endIndex];
      grid[endIndex] = grid[beginIndex];
      grid[beginIndex] = t;
    }
      }
    }

    if (numbers==-1) numbers = u3-u3/3;
    for (int i=0; i<numbers; i++){
      int t, r, c=0;
      do{
    if (c!=0) grid[r] = t;
    r = rand()%u3;
    t = grid[r];
    grid[r] = 0;
    if (++c>limit) break;
      }while(solver(grid, 2, false)!=1);
    }
  }

private:
  void cover(header *column){ // cover a column and its nodes' rows
  column->visible = false;
  column->right->left = column->left;
  column->left->right = column->right;
  for (node *i = column->down; i!=column; i=i->down){
    for (node *j = i->right; j!=i; j=j->right) {
      j->down->up = j->up;
      j->up->down = j->down;
      j->head->size--;
    }
  }
  }

  void uncover(header *column){ // uncover a column and its nodes' rows
  column->visible = true;
  for (node *i = column->up; i!=column; i=i->up){
    for (node *j = i->left; j!=i; j=j->left) {
      j->head->size++;
      j->up->down = j;
      j->down->up = j;
    }
  }
    column->left->right = column;
  column->right->left = column;
  }

  void backtrack(int i){
    if (root->right == root) {
      solutionsNumber++;
      return;
    }
    header * column, * choosen = root->right;
  for (column = root->right; column!=root; column=column->right) {
    if (column->size < choosen->size)
      choosen = column;
    if (choosen->size==0) return;
  }
    cover(choosen);
    for (node *row = choosen->down; row!=choosen && solutionsNumber<maxSolutionsNumber; row=row->down) {
      if (solutionsNumber==0)
    solutions[i] = row->row;
      for (node *j = row->right; j!=row; j=j->right)
    cover(j->head);
      backtrack(i + 1);
      for (node *j = row->left; j!=row; j=j->left)
    uncover(j->head);
    }
    uncover(choosen);

  }

};


int main(){
  int empty[81] = {
    0, 0, 0, 0,  0, 0,  0, 0, 0,
    0, 0, 0, 0,  0, 0,  0, 0, 0,
    0, 0, 0, 0,  0, 0,  0, 0, 0,
    0, 0, 0, 0,  0, 0,  0, 0, 0,

    0, 0, 0, 0,  0, 0,  0, 0, 0,
    0, 0, 0, 0,  0, 0,  0, 0, 0,

    0, 0, 0,  0, 0, 0,  0, 0, 0,
    0, 0, 0,  0, 0, 0,  0, 0, 0,
    0, 0, 0,  0, 0, 0,  0, 0, 0,
  };
  int minusthree[81] = {
    0, 0, 0,  0, 0, 0,  0, 0, 0,
    0, 2, 3,  0, 0, 0,  7, 8, 0,
    1, 0, 0,  4, 0, 6,  0, 0, 9,

    9, 0, 0,  0, 5, 0,  0, 0, 4,
    2, 0, 0,  0, 0, 0,  0, 0, 8,
    0, 1, 0,  0, 0, 0,  0, 3, 0,

    0, 0, 8,  0, 0, 0,  3, 0, 0,
    0, 0, 0,  2, 0, 9,  0, 0, 0,
    0, 0, 0,  0, 1, 0,  0, 0, 0
  };
  int heart[81] = {
    8, 7, 2,  0, 0, 9,  0, 0, 0,
    0, 0, 0,  0, 3, 0,  4, 0, 0,
    0, 0, 0,  0, 0, 1,  2, 0, 8,

    0, 0, 9,  0, 0, 0,  0, 8, 0,
    3, 0, 7,  0, 0, 0,  6, 0, 5,
    0, 2, 0,  0, 0, 0,  1, 0, 0,

    1, 0, 8,  2, 0, 0,  0, 0, 0,
    0, 0, 5,  0, 9, 0,  0, 0, 0,
    0, 0, 0,  1, 0, 0,  9, 5, 3
  };

  clock_t begin, end;

  begin = clock();
  Sudoku sudoku(3);

  cout << "Sudoku.solver\n";
  sudoku.print(heart);
  begin = clock();
  cout << "Working.." << "\n";
  int solutions = sudoku.solver(heart, 2);
  end = clock();
  sudoku.print(heart);
  cout << solutions
       << " solutions. Solved in " << ((double)(end - begin)*1000)/CLOCKS_PER_SEC << "\n";


  cout << "\n\nSudoku.generator\n";
  int mysudoku[81];
  begin = clock();
  cout << "Working.." << "\n";
  sudoku.generator(mysudoku);
  end = clock();
  sudoku.print(mysudoku);
  cout << "Generated in " << ((double)(end - begin)*1000)/CLOCKS_PER_SEC << "\n";

  begin = clock();
  cout << "Working.." << "\n";
  sudoku.solver(mysudoku);
  end = clock();
  sudoku.print(mysudoku);
  cout << solutions
       << " solutions. Solved in " << ((double)(end - begin)*1000)/CLOCKS_PER_SEC << "\n";
}