BranchControl.cc

Go to the documentation of this file.

/*
 *  Copyright 2007 Martin von Gagern
 *
 *
 *  This file is part of bande.
 *
 *  bande is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation; either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  bande is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */


/**
 * @file
 * Implementation of class bande::BranchControl.
 */

#include "bande.hh"

namespace bande {

  /**
   * Run search, return when all solutions have been found.
   *
   * This is the public interface to start the branching algorithm.
   */
  void BranchControl::run() {
    statistics.reset();
    branch(0);
  }

  /**
   * Recursive method to search one branch.
   *
   * This is the central method of the whole program.
   * Bit mask 1 in Settings::debug will enable debug output here.
   *
   * For every invocation, IntegerProgram::solve will be called to
   * determine whether the LP relaxation of the current problem is
   * feasible, which is a neccessary condition for any feasible
   * integral solution. If the LP is infeasible, the branch will
   * return immediately.
   *
   * If the current LP is feasible and all variables are fixed to
   * integral values, then we have a new solution, which concludes the
   * current branch as well.
   *
   * If the LP is feasible but not all variables have been fixed, one
   * variable is chosen using branchCol() and used for a case
   * destinction, leading to two subtrees which are searched by
   * recursive invocations.
   *
   * @param depth number of branching points between the root of the
   *              search tree and the current subtree.
   */
  void BranchControl::branch(int depth) {

#ifndef NDEBUG
    if (Settings::debug & 1)
      for (int i = 0; i < depth; ++i) std::cerr << '.';
#endif

    // Case 1: current LP infeasible, therefor no solutions here.
    if (!ip.solve()) { // see if current problem is feasible at all
      DBG(1, " infeasible");
      statistics.infeasible(depth);
      return;
    }

    int bcol = branchCol(); // choose one column to branch upon

    // Case 2: no branching column, which means all variables are fixed
    if (bcol == -1) {
      DBG(1, " SOLUTION");
      ip.mkIntSolution();   // round double solution to int solution
      if (!ip.checkIntSolution()) { // this should be a valid solution!
        std::cerr << "ERR!" << std::endl;
        ip.dumpSolution();
        exit(EXIT_FAILURE);
      }
      sols.solution(ip);          // remember this solution
      statistics.solution(depth); // update the statistics
      return;
    }

    // Case 3: we chose a branching column, bcol, and will branch upon this
    // Determine a suitable value at which to split, taking care of
    // rounding issues by avoiding new bounds that would coincide with
    // current bounds.
    int bval = (int)floor(ip.getColSolution(bcol));
    if (bval == ip.getColUpper(bcol)) --bval;
    if (bval < ip.getColLower(bcol)) bval = ip.getColLower(bcol);
    DBG(1, ip.getColName(bcol) << " [" << bcol << "] = " << bval
        << " (" << ip.getColSolution(bcol) << ")");
    ++depth;                            // we are about to descend one level
    UndoState state = um.save();        // create checkpoint for undoing things
    ip.setColLower(um, bcol, bval + 1); // first subtree: x_bcol >= bval+1
    branch(depth);                      // search first subtree recursively
    state.restore();                    // undo any changes from that subtree
    ip.setColUpper(um, bcol, bval);     // second subtree: x_bcol <= bval
    branch(depth);                      // search second subtree recursively
  }

  /**
   * Choose branching column.
   *
   * When there are free variables left, one of them has to be chosen
   * as the one on which the case destinction for branching will be
   * made. This method implements a heuristic to make this choice.
   *
   * The heuristic contains of two parts. First a column that
   * currently has a "most fractional" solution is chosen, where most
   * fractional means largest distance to nearest integer. If all
   * columns are integral or near enough, then the second part will
   * chose a column with minimum difference between upper and lower
   * bound, in order to increase the number of fixed columns soon.
   *
   * @return the index of the chosen column,
   *         or -1 if all columns are already fixed
   */
  int BranchControl::branchCol() {
    int bcol = -1;

    // Part 1: Try to find a column farthest from being integral.
    double limit = 1e-5; // due to rounding, treat almost int as int
    for (int col = 0; col < ip.getNumCols(); ++col) {
      double s = ip.getColSolution(col);
      double f = floor(s);
      double d = s - f;
      d *= 1 - d;
      // d now is the product of differences to the nearest integers.
      // This value will be maximal for numbers halfway between
      // integers, zero for integers, and monotone in between.
      if (limit < d) { // we have a now best choice
        limit = d;
        bcol = col;
      }
    }
    if (bcol >= 0) return bcol; // Part 1 succeeded, search no further

    // Part 2: All solutions are (almost) integral, choose one with
    // smallest difference between upper and lower bound. This will
    // help make that variable completely fixed in the near future.
    unsigned minGap = std::numeric_limits<unsigned>::max();
    for (int col = 0; col < ip.getNumCols(); ++col) {
      std::pair<int, int> bounds = ip.getColBounds(col);
      if (bounds.second == bounds.first) continue; // already fixed
      if (minGap > unsigned(bounds.second - bounds.first)) {
        minGap = bounds.second - bounds.first;
        bcol = col;
      }
    }

    // By now we either have a choice from Part 2, or all columns are
    // fixed and we return the initial value of bcol, which was -1.
    return bcol;
  }

}

---