Graph coloring (bitwise encoding)

graphcol.cpp

#include <boost/graph/adjacency_list.hpp>
#include <boost/graph/iteration_macros.hpp>

#include <boost/range/irange.hpp>

#include <minisat/core/Dimacs.h>
#include <minisat/core/Solver.h>

#include <zlib.h>

#include <cmath>
#include <cstdio>
#include <cstdlib>

#include <fstream>
#include <iostream>
#include <sstream>

#include <string>
#include <utility>
#include <vector>

namespace dimacs {
  const std::string literal_sep(" ");
  const std::string clause_end(" 0\n");
}

namespace graphcolor {

#define NODE_TO_INT(C,V,I)                      \
  ((C) * (V) + (I))
#define SIGNAL(A,B)                             \
  (NBIT((A),(B)) ? 1 : -1)
#define NBIT(A,B)                               \
  ((A >> B) & 1)

  // File parsing based on this SO answer:
  // http://stackoverflow.com/questions/30415388/how-to-read-dimacs-vertex-coloring-graphs-in-c/30446685#30446685
  template <typename Graph>
  bool
  read_coloring_problem(std::istream& dimacs,
                        Graph& G)
  {
    size_t vertices = 0, edges = 0;
    std::string line;

    while (getline(dimacs, line)) {
      std::istringstream iss(line);
      char ch;
      if (iss >> ch) {
        size_t from, to;
        std::string format;

        switch(ch) {
        case 'c':
          // Comment.
          break;
        case 'p':
          // File specification.
          if (vertices || edges) {
            // We have already read the file specification.
            return false;
          }

          if (iss >> format >> vertices >> edges) {
            if (format != "edge") return false;
          }
          break;
        case 'e':
          if (edges-- > 0 && (iss >> from >> to) && (add_edge(from - 1, to - 1, G).second)) {
            break;
          } else {
            return false;
          }
        default:
          return false;
        }
      }
    }
    return true;
  }


  template <typename Graph>
  void
  write_different_colors_constraints(std::ostream& os,
                                     const Graph& G,
                                     const unsigned C)
  {
    using namespace boost;
    const auto index = get(vertex_index, G);

    BGL_FORALL_EDGES_T(e, G, Graph) {
      const unsigned U = index[source(e, G)];
      const unsigned V = index[target(e, G)];

      for (unsigned i : irange(1u, C + 1)) {
        os << -1 * NODE_TO_INT(C, U, i) << dimacs::literal_sep
           << -1 * NODE_TO_INT(C, V, i) << dimacs::clause_end;
      }
    }
  }


  template <typename Graph>
  void
  write_at_least_one_color_constraints(std::ostream& os,
                                       const Graph& G,
                                       const unsigned C)
  {
    using namespace boost;
    const auto index = get(vertex_index, G);

    BGL_FORALL_VERTICES_T(v, G, Graph) {
      for (unsigned i : irange(1u, C + 1)) {
        os << NODE_TO_INT(C, index[v], i) << dimacs::literal_sep;
      }
      os << dimacs::clause_end;
    }
  }


  template <typename Graph>
  void
  write_at_most_one_color_constraints_pairwise(std::ostream& os,
                                               const Graph& G,
                                               const unsigned C)
  {
    using namespace boost;
    const auto index = get(vertex_index, G);

    BGL_FORALL_VERTICES_T(v, G, Graph) {
      for (unsigned i : irange(1u, C)) {
        for (unsigned j : irange(i + 1, C + 1)) {
          os << -1 * NODE_TO_INT(C, index[v], i) << dimacs::literal_sep
             << -1 * NODE_TO_INT(C, index[v], j) << dimacs::clause_end;
        }
      }
    }
  }


  template <typename Graph>
  void
  write_at_most_one_color_constraints_bitwise(std::ostream& os,
                                              const Graph& G,
                                              const unsigned C)
  {
    using namespace boost;
    const auto index = get(vertex_index, G);

    const unsigned V = num_vertices(G);
    const unsigned K = std::ceil(std::log2(C)); // Number of new variables

    // for (tie(vi, vi_end) = vertices(G); vi != vi_end; ++vi) {
    BGL_FORALL_VERTICES_T(v, G, Graph) {
      const unsigned I      = index[v];
      const unsigned OFFSET = (V * C) + (I * K) + 1;

      // Usual bitwise constraints.
      for (unsigned j : irange(1u, C + 1)) {
        for (unsigned b : irange(0u, K)) {
          os << -1 * NODE_TO_INT(C, I, j)       << dimacs::literal_sep
             << SIGNAL(j - 1, b) * (OFFSET + b) << dimacs::clause_end;
        }
      }

      // Domain restriction.
      for (unsigned j : irange(C, (unsigned) std::pow(2, K))) {
        for (unsigned b : irange(0u, K)) {
          os << -1 * SIGNAL(j, b) * (OFFSET + b)
             << dimacs::literal_sep;
        }
        os << dimacs::clause_end;
      }
    }
  }

  template <typename Graph>
  void
  write_symmetry_breaking_constraints(std::ostream& os,
                                      const Graph& G,
                                      const unsigned C)
  {
    const auto index = get(boost::vertex_index, G);
    typename boost::graph_traits<Graph>::vertex_iterator vi;

    // We would gain the biggest benefit out of this if the vertex whose color is
    // fixed was the one with the greatest degree.
    os << NODE_TO_INT(C, index[*vi], 1) << dimacs::clause_end;
  }


  template <typename Graph>
  void
  write_encoding(std::ostream& os,
                 const Graph& G,
                 const unsigned C)
  {
    std::ostringstream aux_oss;

    const unsigned V = num_vertices(G);
    const unsigned E = num_edges(G);
    const unsigned K = std::ceil(std::log2(C));

    // Bitwise encoding with symmetry breaking
    const unsigned NUM_CLAUSES   = (E * C) + V + (V * (C * (K - 1) + std::pow(2, K))) + 1;
    const unsigned NUM_VARIABLES = V * (C + K);
    write_at_most_one_color_constraints_bitwise(aux_oss, G, C);
    write_symmetry_breaking_constraints(aux_oss, G, C);

    // Pairwise encoding
    // const unsigned NUM_CLAUSES   = (E * C) + V + V * (C * (C - 1) / 2);
    // const unsigned NUM_VARIABLES = V * C;
    // write_at_most_one_color_constraints_pairwise(aux_oss, G, C);

    write_different_colors_constraints(aux_oss, G, C);
    write_at_least_one_color_constraints(aux_oss, G, C);

    os << "p cnf " << NUM_VARIABLES << " " << NUM_CLAUSES << "\n"
       << aux_oss.str() << std::flush;
  }


  template <typename Graph>
  bool
  solve_coloring_problem(const Graph& g,
                         const unsigned C)
  {
    const char *dimacs_filename = tmpnam(nullptr);
    std::ofstream dimacs(dimacs_filename);

    write_encoding(dimacs, g, C);
    dimacs.close();

    gzFile in = gzopen(dimacs_filename, "rb");

    Minisat::Solver solver;
    parse_DIMACS(in, solver);

    gzclose(in);
    std::remove(dimacs_filename);

    return solver.solve();
  }

} // namespace graphcolor


template <typename Graph>
unsigned search(const Graph& G,
                const unsigned unsat,
                const unsigned sat)
{
  if (unsat == sat || unsat + 1 == sat)
    return sat;

  const unsigned middle = sat - (sat - unsat) / 2;
  if (graphcolor::solve_coloring_problem(G, middle)) {
    return search(G, unsat, middle);
  }
  return search(G, middle, sat);
};

int main(int argc, char* argv[])
{
  if (argc < 2 || (argc == 2 && !strcmp(argv[1], "--help"))) {
    std::cout << "USAGE: graphcolor <input-file>" << std::endl;
    return 0;
  }

  const char* filename = argv[1];
  std::ifstream infile(filename);

  using namespace boost;
  adjacency_list<vecS, vecS, undirectedS> G;

  if (graphcolor::read_coloring_problem(infile, G)) {
    std::cout << search(G, 0, num_vertices(G)) << std::endl;
  } else {
    std::cout << "** Error: bad DIMACS graph configuration." << std::endl;
  }

  return 0;
}