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"


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 ((float) n1->cost_ / n1->nbordered() > (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);
            }
            //
            // add partially referenced vertices to the list of potential choices for branching
            else if (this->nbrefs_[v] >= inst.L()) {
                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()) {
                            for (Vertex* u: v->neighbors_) {
                                if (this->mustbefullref_[u]) {
                                    if (this->maxrefs_[u] >= inst.U()+1) {
                                        this->maxrefs_[u]--;
                                    }
                                }
                                else this->maxrefs_[u]--;
                            }
                        }
                    }
                }
                else {
                    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;