Ignore:
Timestamp:
Apr 1, 2013 1:09:22 PM (6 years ago)
Author:
forrest
Message:

make it easier to use slow exotic cuts
more on cutoff as constraint and multiple root solver fixes
general fixing of bugs found on MIQP etc

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Cbc/src/CbcFollowOn.cpp

    r1854 r1880  
    2323#include "CbcFollowOn.hpp"
    2424#include "CbcBranchActual.hpp"
     25#include "CbcBranchCut.hpp"
    2526#include "CoinSort.hpp"
    2627#include "CoinError.hpp"
     
    525526//##############################################################################
    526527
     528//##############################################################################
     529
     530// Default Constructor
     531CbcIdiotBranch::CbcIdiotBranch ()
     532        : CbcObject()
     533{
     534  id_ = 1000000000 +  CutBranchingObj;
     535}
     536
     537// Useful constructor
     538CbcIdiotBranch::CbcIdiotBranch (CbcModel * model)
     539        : CbcObject(model)
     540{
     541    assert (model);
     542    id_ = 1000000000 +  CutBranchingObj;
     543}
     544
     545// Copy constructor
     546CbcIdiotBranch::CbcIdiotBranch ( const CbcIdiotBranch & rhs)
     547  : CbcObject(rhs),
     548    randomNumberGenerator_(rhs.randomNumberGenerator_),
     549    savedRandomNumberGenerator_(rhs.savedRandomNumberGenerator_)
     550{
     551}
     552
     553// Clone
     554CbcObject *
     555CbcIdiotBranch::clone() const
     556{
     557    return new CbcIdiotBranch(*this);
     558}
     559
     560// Assignment operator
     561CbcIdiotBranch &
     562CbcIdiotBranch::operator=( const CbcIdiotBranch & rhs)
     563{
     564    if (this != &rhs) {
     565        CbcObject::operator=(rhs);
     566        randomNumberGenerator_ = rhs.randomNumberGenerator_;
     567        savedRandomNumberGenerator_ = rhs.savedRandomNumberGenerator_;
     568    }
     569    return *this;
     570}
     571
     572// Destructor
     573CbcIdiotBranch::~CbcIdiotBranch ()
     574{
     575}
     576double
     577CbcIdiotBranch::infeasibility(const OsiBranchingInformation * info,
     578                           int &preferredWay) const
     579{
     580  randomNumberGenerator_ = savedRandomNumberGenerator_;
     581  double rhs = buildCut(info,0,preferredWay).ub();
     582  double fraction = rhs-floor(rhs);
     583  if (fraction>0.5)
     584    fraction=1.0-fraction;
     585  return fraction;
     586}
     587
     588// This looks at solution and sets bounds to contain solution
     589void
     590CbcIdiotBranch::feasibleRegion()
     591{
     592}
     593
     594CbcBranchingObject *
     595CbcIdiotBranch::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * info, int way)
     596{
     597  // down and up
     598  randomNumberGenerator_ = savedRandomNumberGenerator_;
     599  int preferredWay;
     600  OsiRowCut downCut = buildCut(info,0,preferredWay);
     601  double rhs = downCut.ub();
     602  assert(rhs == downCut.lb());
     603  OsiRowCut upCut =downCut;
     604  downCut.setUb(floor(rhs));
     605  downCut.setLb(-COIN_DBL_MAX);
     606  upCut.setLb(ceil(rhs));
     607  upCut.setUb(COIN_DBL_MAX);
     608  CbcBranchingObject * branch
     609    = new CbcCutBranchingObject(model_, downCut,upCut,true);
     610  return branch;
     611}
     612// Initialize for branching
     613void
     614CbcIdiotBranch::initializeForBranching(CbcModel * )
     615{
     616  savedRandomNumberGenerator_ = randomNumberGenerator_;
     617}
     618// Build "cut"
     619OsiRowCut
     620CbcIdiotBranch::buildCut(const OsiBranchingInformation * info,int /*type*/,int & preferredWay) const
     621{
     622  int numberIntegers = model_->numberIntegers();
     623  const int * integerVariable = model_->integerVariable();
     624  int * which = new int [numberIntegers];
     625  double * away = new double [numberIntegers];
     626  const double * solution = info->solution_;
     627  const double * lower = info->lower_;
     628  const double * upper = info->upper_;
     629  double integerTolerance = model_->getIntegerTolerance();
     630  //int nMax=CoinMin(4,numberIntegers/2);
     631  int n=0;
     632  for (int i=0;i<numberIntegers;i++) {
     633    int iColumn=integerVariable[i];
     634    double value = solution[iColumn];
     635    value = CoinMax(value, lower[iColumn]);
     636    value = CoinMin(value, upper[iColumn]);
     637    //#define JUST_SMALL
     638#ifndef JUST_SMALL
     639    double nearest = floor(value + 0.5);
     640    if (fabs(value-nearest)>integerTolerance) {
     641      double random = 0.0; //0.25*randomNumberGenerator_.randomDouble();
     642      away[n]=-fabs(value-nearest)*(1.0+random);
     643      which[n++]=iColumn;
     644    }
     645#else
     646    double distanceDown = value-floor(value);
     647    if (distanceDown>10.0*integerTolerance&&distanceDown<0.1) {
     648      double random = 0.0; //0.25*randomNumberGenerator_.randomDouble();
     649      away[n]=distanceDown*(1.0+random);
     650      which[n++]=iColumn;
     651    }
     652#endif
     653  }
     654  CoinSort_2(away,away+n,which);
     655  OsiRowCut possibleCut;
     656  possibleCut.setUb(0.0);
     657  if (n>1) {
     658    int nUse=0;
     659    double useRhs=0.0;
     660    double best=0.0;
     661    double rhs = 0.0;
     662    double scaleFactor=1.0;
     663    /* should be able to be clever to find best cut
     664       maybe don't do if away >0.45
     665       maybe don't be random
     666       maybe do complete enumeration on biggest k (where sum of k >= 1.0)
     667       maybe allow +-1
     668     */
     669    for (int i=0;i<n;i++) {
     670      int iColumn=which[i];
     671      double value = solution[iColumn];
     672      value = CoinMax(value, lower[iColumn]);
     673      value = CoinMin(value, upper[iColumn]);
     674#define PLUS_MINUS
     675#ifndef PLUS_MINUS
     676      away[i]=1.0;
     677      rhs += value;
     678#else
     679      if (value-floor(value)<=0.5) {
     680        away[i]=1.0;
     681        rhs += value;
     682      } else {
     683        away[i]=-1.0;
     684        rhs -= value;
     685      }
     686#endif
     687      double nearest = floor(rhs + 0.5);
     688      double distance=fabs(rhs-nearest);
     689      // scale back
     690      distance *= scaleFactor;
     691      scaleFactor *= 0.95;
     692      if (distance>best) {
     693        nUse=i+1;
     694        best=distance;
     695        useRhs=rhs;
     696      }
     697    }
     698    // create cut
     699    if (nUse>1) {
     700      possibleCut.setRow(nUse,which,away);
     701      possibleCut.setLb(useRhs);
     702      possibleCut.setUb(useRhs);
     703    }
     704  }
     705  delete [] which;
     706  delete [] away;
     707  return possibleCut;
     708}
     709#if 0
     710//##############################################################################
     711
     712// Default Constructor
     713CbcBoundingBranchingObject::CbcBoundingBranchingObject()
     714        : CbcBranchingObject()
     715{
     716  branchCut_.setUb(0.0);
     717}
     718
     719// Useful constructor
     720CbcBoundingBranchingObject::CbcBoundingBranchingObject (CbcModel * model,
     721                                                        int way , const OsiRowCut * cut)
     722        : CbcBranchingObject(model, 0, way, 0.5)
     723{
     724  branchCut_ = *cut;
     725}
     726
     727// Copy constructor
     728CbcBoundingBranchingObject::CbcBoundingBranchingObject ( const CbcBoundingBranchingObject & rhs) : CbcBranchingObject(rhs)
     729{
     730  branchCut_ = rhs.branchCut_;
     731}
     732
     733// Assignment operator
     734CbcBoundingBranchingObject &
     735CbcBoundingBranchingObject::operator=( const CbcBoundingBranchingObject & rhs)
     736{
     737    if (this != &rhs) {
     738        CbcBranchingObject::operator=(rhs);
     739        branchCut_ = rhs.branchCut_;
     740    }
     741    return *this;
     742}
     743CbcBranchingObject *
     744CbcBoundingBranchingObject::clone() const
     745{
     746    return (new CbcBoundingBranchingObject(*this));
     747}
     748
     749
     750// Destructor
     751CbcBoundingBranchingObject::~CbcBoundingBranchingObject ()
     752{
     753}
     754double
     755CbcBoundingBranchingObject::branch()
     756{
     757    decrementNumberBranchesLeft();
     758    OsiSolverInterface * solver = model_->solver();
     759    const double * columnLower = solver->getColLower();
     760    int i;
     761    // *** for way - up means bound all those in up section
     762    if (way_ < 0) {
     763#ifdef FULL_PRINT
     764        printf("Down Bound ");
     765#endif
     766        //printf("Down Bound %d\n",numberDown_);
     767        for (i = 0; i < numberDown_; i++) {
     768            int iColumn = downList_[i];
     769            model_->solver()->setColUpper(iColumn, columnLower[iColumn]);
     770#ifdef FULL_PRINT
     771            printf("Setting bound on %d to lower bound\n", iColumn);
     772#endif
     773        }
     774        way_ = 1;         // Swap direction
     775    } else {
     776#ifdef FULL_PRINT
     777        printf("Up Bound ");
     778#endif
     779        //printf("Up Bound %d\n",numberUp_);
     780        for (i = 0; i < numberUp_; i++) {
     781            int iColumn = upList_[i];
     782            model_->solver()->setColUpper(iColumn, columnLower[iColumn]);
     783#ifdef FULL_PRINT
     784            printf("Setting bound on %d to lower bound\n", iColumn);
     785#endif
     786        }
     787        way_ = -1;        // Swap direction
     788    }
     789#ifdef FULL_PRINT
     790    printf("\n");
     791#endif
     792    return 0.0;
     793}
     794void
     795CbcBoundingBranchingObject::print()
     796{
     797  OsiRowCut cut = branchCut_;
     798  if (way_ < 0) {
     799    printf("Down Fix ");
     800    cut.setUb(floor(branchCut_.ub()));
     801    cut.setLb(-COIN_DBL_MAX);
     802  } else {
     803    printf("Up Fix ");
     804    cut.setLb(ceil(branchCut_.lb()));
     805    cut.setUb(COIN_DBL_MAX);
     806  }
     807  printf("\n");
     808  cut.print();
     809}
     810
     811/** Compare the original object of \c this with the original object of \c
     812    brObj. Assumes that there is an ordering of the original objects.
     813    This method should be invoked only if \c this and brObj are of the same
     814    type.
     815    Return negative/0/positive depending on whether \c this is
     816    smaller/same/larger than the argument.
     817*/
     818int
     819CbcBoundingBranchingObject::compareOriginalObject
     820(const CbcBranchingObject* /*brObj*/) const
     821{
     822    throw("must implement");
     823}
     824
     825/** Compare the \c this with \c brObj. \c this and \c brObj must be os the
     826    same type and must have the same original object, but they may have
     827    different feasible regions.
     828    Return the appropriate CbcRangeCompare value (first argument being the
     829    sub/superset if that's the case). In case of overlap (and if \c
     830    replaceIfOverlap is true) replace the current branching object with one
     831    whose feasible region is the overlap.
     832   */
     833CbcRangeCompare
     834CbcBoundingBranchingObject::compareBranchingObject
     835(const CbcBranchingObject* /*brObj*/, const bool /*replaceIfOverlap*/)
     836{
     837#ifdef JJF_ZERO //ndef NDEBUG
     838    const CbcBoundingBranchingObject* br =
     839        dynamic_cast<const CbcBoundingBranchingObject*>(brObj);
     840    assert(br);
     841#endif
     842    // If two BoundingBranchingObject's have the same base object then it's pretty
     843    // much guaranteed
     844    throw("must implement");
     845}
     846
     847//##############################################################################
     848
     849#endif
Note: See TracChangeset for help on using the changeset viewer.