bbsolver.cpp

/********************************************************************************************
 Name:       DiscretizationSolver
 Discrete optimization for graph discretization
 Author:     J.Omer
 Sources:    C++
 License:    GNU General Public License v.2
 History:
 *********************************************************************************************/

#include "bbsolver.hpp"

#include <algorithm>


ILOSTLBEGIN
//

ILOLAZYCONSTRAINTCALLBACK2(CyclesLazyConstraintsBB,
                           Instance &, inst,
                           DiscretizationSolver&, solver) {
#ifdef VERBOSE
  std::cout << "revorder: check lazy dicycle constraints" << std::endl;
#endif
  // Search for cycles in the current solution digraph
  //
  // generate the adjacency lists of the vertices corresponding to the current
  // cplex solution
  float tolerance = 1e-06;
  std::map<Vertex*, std::vector<Vertex*> > adjlist;
  std::map<Vertex*, int > nbrefs;
  for (Vertex* v : inst.vertices_) {
    adjlist[v] = std::vector<Vertex*>();
    nbrefs[v] = 0;
  }
  for (Vertex* u : inst.vertices_) {
    if (solver.modeltype_ == WITNESS) {
      if (getValue(solver.isvertexinclique(u) >= 1 - tolerance)) continue;
    }
    for (Vertex* v : u->neighbors_) {
      if (getValue(solver.isbefore(u, v)) >= 1 - tolerance) {
        adjlist[u].push_back(v);
        nbrefs[v]++;
      }
    }
  }
  //
  // search for the root of the digraph
  Vertex* root = nullptr;
  int minnbrefs = inst.U();
  if (solver.modeltype_ == WITNESS) {
    for (Vertex* v : inst.vertices_) {
      if (getValue(solver.isvertexinclique(v) >= 1 - tolerance)) {
        root = v;
        break;
      }
    }
  } else {
    for (Vertex* v : inst.vertices_) {
      if (nbrefs[v] < minnbrefs) {
        minnbrefs = nbrefs[v];
        root = v;
        break;
      }
    }
  }
  //
  // call the enumeration of cycles, returns false if there are none
  std::vector<std::vector<Vertex*> > cycles;
  bool iscyclic = solver.enumeratecycles(adjlist, root, cycles);

  // Look for dicycle inequalities : we restrict the search to the dicycles
  // that contain at least one edge in the distance graph
  // besides that, the search is performed using brute force
  int nbaddedcuts = 0;
  IloEnv env = getEnv();
  if (iscyclic) {
    for (std::vector<Vertex*> path : cycles) {
      IloExpr sumedges(env);
      int cyclelength = path.size();
      Vertex* u = path.back();
      for (int i = 0; i < cyclelength; i++) {
        Vertex* v = path[i];
        sumedges += solver.isbefore(u, v);
        u = v;
      }
      if ((solver.modeltype_ == WITNESS) && (cyclelength <= inst.L() + 1)) {
        add(sumedges - solver.isvertexinclique(path.front())
                <= cyclelength - 1).end();
        nbaddedcuts++;
      } else {
        add(sumedges <= cyclelength - 1).end();
        nbaddedcuts++;
      }
    }
#ifdef VERBOSE
    std::cout << "revorder: number of violated dicycle constraints : " << nbaddedcuts << std::endl;
#endif
  } else {
#ifdef VERBOSE
    std::cout << "revorder: the graph is acyclic " << std::endl;
#endif
  }
  return;
}
/*************************************************************************************************************************************************************/
//
// sort the nodes by ascending number of ordered vertices
bool bbnodesnborderedcompare(BBNode* n1, BBNode* n2) {
  return (n1->nbordered() < n2->nbordered());
}
bool bbnodesbestcompare(BBNode*n1, BBNode* n2) {
  return (static_cast<float>(n1->cost_) / n1->nbordered()
      > static_cast<float>(n2->cost_) / n2->nbordered());
}

// Public methods
//
// constructor and destructor
//
// constructor used only for the root
BBNode::BBNode(Instance& inst):
    depth_(0),
    inst_(&inst),
    nbordered_(0),
    nbfullref_(0),
    istreated_(true) {}
//
// constructor used only for the initial cliques
BBNode::BBNode(Instance& inst,
               BBNode* root,
               Clique& c,
               bool branchonsmallerindex):
    depth_(1),
    father_(root),
    inst_(&inst),
    nbordered_(0),
    nbfullref_(0),
    istreated_(false) {
  this->orderedvertices_.clear();
  for (Vertex* v : inst.vertices_) {
    this->rank_[v] = -1;
    this->nbrefs_[v] = 0;
    this->isfullref_[v] = false;
    this->isordered_[v] = false;
    this->ispotentialchild_[v] = false;
    this->mustbefullref_[v] = false;
    this->maxrefs_[v] = v->degree_;
  }
  int currentrank = 1;
  for (Vertex* u : c.vertices_) {
    this->rank_[u] = currentrank++;
    this->orderedvertices_.push_back(u);
    this->switchtoordered(u);
    this->nbordered_++;
    //
    // update the number of references of the neighbors of the vertex
    for (Vertex* v : u->neighbors_) {
      if (!this->isordered_[v]) {
        this->nbrefs_[v]++;
      }
    }
  }

  // record the list of partially and fully referenced (unordered) vertices
  for (Vertex* v : inst.vertices_) {
    if (!this->isordered_[v]) {
      // record the list of unordered fully-referenced vertices
      if (this->nbrefs_[v] >= inst.U()) {
        this->isfullref_[v] = true;
        this->fullrefs_.push_back(v);
      } else if (this->nbrefs_[v] >= inst.L()) {
        // add partially referenced vertices to the list of potential choices
        // for branching
        if (branchonsmallerindex) {
          bool largerindex = true;
          Vertex* u = c.vertices_.back();
          if (!u->isneighbor(v)) {
            if (u->id_ >= v->id_) largerindex = false;
          }
          //
          // to break symmetries, make sure that the chosen potential child
          // is always that with biggest index
          if (largerindex) {
            this->ispotentialchild_[v] = true;
            this->potentialchildren_.push_back(v);
            //
            // partially referenced vertices with less than U neighbors can
            // be ordered right away since they will never be fully referenced
            if (this->maxrefs_[v] <= inst.U() - 1) {
              this->readytoorder_.push_back(v);
#ifdef VERBOSE
              std::cout << "revorder: bb: order vertex " << v->id_ << " without branching, maximum nb of references = " << this->maxrefs_[v] << std::endl;
#endif
            }
          } else {
            this->mustbefullref_[v] = true;
            // if a vertex must be fully-referenced, but it has exactly U
            // neighbors, it will for sure come after all its neighbors in
            // the order
            if (this->maxrefs_[v] <= inst.U()) {