Changeset 1880


Ignore:
Timestamp:
Apr 1, 2013 1:09:22 PM (7 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

Location:
trunk/Cbc/src
Files:
23 edited

Legend:

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

    r1854 r1880  
    138138                return equalityTest(x, y); // so ties will be broken in consistent manner
    139139        }
     140    }
     141    if (!weight_) {
     142      double testX =  x->objectiveValue() + 1.0e-9 * x->numberUnsatisfied();
     143      double testY = y->objectiveValue() + 1.0e-9 * y->numberUnsatisfied();
     144      if (testX != testY)
     145        return testX > testY;
     146      else
     147        return equalityTest(x, y); // so ties will be broken in consistent manner
    140148    }
    141149    //weight_=0.0;
  • trunk/Cbc/src/CbcCutGenerator.cpp

    r1839 r1880  
    5050        numberActiveCutsAtRoot_(0),
    5151        numberShortCutsAtRoot_(0),
    52         switches_(1)
     52        switches_(1),
     53        maximumTries_(-1)
    5354{
    5455}
     
    7374        numberActiveCutsAtRoot_(0),
    7475        numberShortCutsAtRoot_(0),
    75         switches_(1)
     76        switches_(1),
     77        maximumTries_(-1)
    7678{
    7779    if (howOften < -1900) {
     
    111113    generatorName_ = CoinStrdup(rhs.generatorName_);
    112114    switches_ = rhs.switches_;
     115    maximumTries_ = rhs.maximumTries_;
    113116    timeInCutGenerator_ = rhs.timeInCutGenerator_;
    114117    savedCuts_ = rhs.savedCuts_;
     
    141144        generatorName_ = CoinStrdup(rhs.generatorName_);
    142145        switches_ = rhs.switches_;
     146        maximumTries_ = rhs.maximumTries_;
    143147        timeInCutGenerator_ = rhs.timeInCutGenerator_;
    144148        savedCuts_ = rhs.savedCuts_;
     
    202206    if (howOften == -100)
    203207        return false;
     208    int pass = model_->getCurrentPassNumber() - 1;
     209    if (maximumTries_>0) {
     210      // howOften means what it says
     211      if ((pass%howOften)!=0||depth)
     212        return false;
     213      else
     214        howOften=1;
     215    }
    204216    if (howOften > 0)
    205217        howOften = howOften % 1000000;
     
    210222    bool returnCode = false;
    211223    //OsiSolverInterface * solver = model_->solver();
    212     int pass = model_->getCurrentPassNumber() - 1;
    213224    // Reset cuts on first pass
    214225    if (!pass)
     
    287298        CglProbing* generator =
    288299            dynamic_cast<CglProbing*>(generator_);
     300        //if (!depth&&!pass)
     301        //printf("Cut generator %s when %d\n",generatorName_,whenCutGenerator_);
    289302        if (!generator) {
    290303            // Pass across model information in case it could be useful
     
    11191132        }
    11201133#endif
    1121         {
    1122             int numberRowCutsAfter = cs.sizeRowCuts() ;
    1123             if (numberRowCutsBefore < numberRowCutsAfter) {
    1124                 for (int k = numberRowCutsBefore; k < numberRowCutsAfter; k++) {
    1125                     OsiRowCut thisCut = cs.rowCut(k) ;
    1126                     int n = thisCut.row().getNumElements();
    1127                     numberElements_ += n;
    1128                 }
     1134        int numberRowCutsAfter = cs.sizeRowCuts() ;
     1135        int numberColumnCutsAfter = cs.sizeColCuts() ;
     1136        if (numberRowCutsBefore < numberRowCutsAfter) {
     1137            for (int k = numberRowCutsBefore; k < numberRowCutsAfter; k++) {
     1138                OsiRowCut thisCut = cs.rowCut(k) ;
     1139                int n = thisCut.row().getNumElements();
     1140                numberElements_ += n;
     1141            }
    11291142#ifdef JJF_ZERO
    1130                 printf("generator %s generated %d row cuts\n",
    1131                        generatorName_, numberRowCutsAfter - numberRowCutsBefore);
    1132 #endif
    1133                 numberCuts_ += numberRowCutsAfter - numberRowCutsBefore;
    1134             }
    1135             int numberColumnCutsAfter = cs.sizeColCuts() ;
    1136             if (numberColumnCutsBefore < numberColumnCutsAfter) {
     1143            printf("generator %s generated %d row cuts\n",
     1144                   generatorName_, numberRowCutsAfter - numberRowCutsBefore);
     1145#endif
     1146            numberCuts_ += numberRowCutsAfter - numberRowCutsBefore;
     1147        }
     1148        if (numberColumnCutsBefore < numberColumnCutsAfter) {
    11371149#ifdef JJF_ZERO
    1138                 printf("generator %s generated %d column cuts\n",
    1139                        generatorName_, numberColumnCutsAfter - numberColumnCutsBefore);
    1140 #endif
    1141                 numberColumnCuts_ += numberColumnCutsAfter - numberColumnCutsBefore;
    1142             }
     1150            printf("generator %s generated %d column cuts\n",
     1151                   generatorName_, numberColumnCutsAfter - numberColumnCutsBefore);
     1152#endif
     1153            numberColumnCuts_ += numberColumnCutsAfter - numberColumnCutsBefore;
    11431154        }
    11441155        if (timing())
    11451156            timeInCutGenerator_ += CoinCpuTime() - time1;
    1146 #ifdef JJF_ZERO
    11471157        // switch off if first time and no good
    1148         if (node == NULL && !pass) {
    1149             if (cs.sizeCuts() - cutsBefore < CoinAbs(switchOffIfLessThan_)) {
    1150                 whenCutGenerator_ = -99;
    1151                 whenCutGeneratorInSub_ = -200;
     1158        if (node == NULL && !pass ) {
     1159            if (numberRowCutsAfter - numberRowCutsBefore
     1160                < switchOffIfLessThan_ /*&& numberCuts_ < switchOffIfLessThan_*/) {
     1161              // switch off
     1162              maximumTries_ = 0;
     1163              whenCutGenerator_=-100;
     1164              //whenCutGenerator_ = -100;
     1165              //whenCutGeneratorInSub_ = -200;
    11521166            }
    11531167        }
    1154 #endif
     1168        if (maximumTries_>0) {
     1169          maximumTries_--;
     1170          if (!maximumTries_)
     1171            whenCutGenerator_=-100;
     1172        }
    11551173    }
    11561174    return returnCode;
  • trunk/Cbc/src/CbcCutGenerator.hpp

    r1839 r1880  
    166166        return depthCutGeneratorInSub_;
    167167    }
     168    /// Set maximum number of times to enter
     169    inline void setMaximumTries(int value)
     170    { maximumTries_ = value;}
     171    /// Get maximum number of times to enter
     172    inline int maximumTries() const
     173    { return maximumTries_;}
    168174
    169175    /// Get switches (for debug)
     
    451457    /// Switches - see gets and sets
    452458    int switches_;
     459    /// Maximum number of times to enter
     460    int maximumTries_;
    453461};
    454462
  • 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
  • trunk/Cbc/src/CbcFollowOn.hpp

    r1854 r1880  
    149149};
    150150
     151/** Define an idiotic idea class.
     152    The idea of this is that we take some integer variables away from integer and
     153    sum them with some randomness to get signed sum close to 0.5.  We then can
     154    branch to exclude that gap.
     155
     156    This branching rule should be in addition to normal rules and have a high priority.
     157*/
     158
     159class CbcIdiotBranch : public CbcObject {
     160
     161public:
     162
     163    // Default Constructor
     164    CbcIdiotBranch ();
     165
     166    /** Useful constructor
     167    */
     168    CbcIdiotBranch (CbcModel * model);
     169
     170    // Copy constructor
     171    CbcIdiotBranch ( const CbcIdiotBranch &);
     172
     173    /// Clone
     174    virtual CbcObject * clone() const;
     175
     176    // Assignment operator
     177    CbcIdiotBranch & operator=( const CbcIdiotBranch& rhs);
     178
     179    // Destructor
     180    ~CbcIdiotBranch ();
     181
     182    /// Infeasibility - large is 0.5
     183    virtual double infeasibility(const OsiBranchingInformation * info,
     184                                 int &preferredWay) const;
     185
     186    using CbcObject::feasibleRegion ;
     187    /// This looks at solution and sets bounds to contain solution
     188    virtual void feasibleRegion();
     189
     190    /// Creates a branching object
     191    virtual CbcBranchingObject * createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * info, int way) ;
     192    /// Initialize for branching
     193    virtual void initializeForBranching(CbcModel * );
     194protected:
     195    /// Build "cut"
     196    OsiRowCut buildCut(const OsiBranchingInformation * info,int type,int & preferredWay) const;
     197    /// data
     198    /// Thread specific random number generator
     199    mutable CoinThreadRandom randomNumberGenerator_;
     200    /// Saved version of thread specific random number generator
     201    mutable CoinThreadRandom savedRandomNumberGenerator_;
     202};
     203
    151204#endif
    152205
  • trunk/Cbc/src/CbcGeneralDepth.cpp

    r1854 r1880  
    2323#include "CbcGeneralDepth.hpp"
    2424#include "CbcBranchActual.hpp"
     25#include "CbcHeuristicDive.hpp"
    2526#include "CoinSort.hpp"
    2627#include "CoinError.hpp"
     
    160161        = dynamic_cast<OsiClpSolverInterface *> (solver);
    161162        if (clpSolver) {
     163          if ((model_->moreSpecialOptions()&33554432)==0) {
    162164            ClpNodeStuff * info = nodeInfo_;
    163165            info->integerTolerance_ = model_->getIntegerTolerance();
     
    246248            simplex->setLogLevel(saveLevel);
    247249            numberNodes_ = info->nNodes_;
    248             int numberDo = numberNodes_;
    249             if (whichSolution_ >= 0)
    250                 numberDo--;
    251             if (numberDo > 0) {
    252                 return 0.5;
    253             } else {
    254                 // no solution
    255                 return COIN_DBL_MAX; // say infeasible
    256             }
     250          } else {
     251            // Try diving
     252            // See if any diving heuristics set to do dive+save
     253            CbcHeuristicDive * dive=NULL;
     254            for (int i = 0; i < model_->numberHeuristics(); i++) {
     255              CbcHeuristicDive * possible = dynamic_cast<CbcHeuristicDive *>(model_->heuristic(i));
     256              if (possible&&possible->maxSimplexIterations()==COIN_INT_MAX) {
     257                // if more than one then rotate later?
     258                //if (possible->canHeuristicRun()) {
     259                dive=possible;
     260                break;
     261              }
     262            }
     263            assert (dive); // otherwise moreSpecial should have been turned off
     264            CbcSubProblem ** nodes=NULL;
     265            int branchState=dive->fathom(model_,numberNodes_,nodes);
     266            if (branchState) {
     267              printf("new solution\n");
     268              whichSolution_=numberNodes_-1;
     269            } else {
     270              whichSolution_=-1;
     271            }
     272#if 0
     273            if (0) {
     274              for (int iNode=0;iNode<numberNodes;iNode++) {
     275                //tree_->push(nodes[iNode]) ;
     276              }
     277              assert (node->nodeInfo());
     278              if (node->nodeInfo()->numberBranchesLeft()) {
     279                tree_->push(node) ;
     280              } else {
     281                node->setActive(false);
     282              }
     283            }
     284#endif
     285            //delete [] nodes;
     286            model_->setTemporaryPointer(reinterpret_cast<void *>(nodes));
     287            // end try diving
     288          }
     289          int numberDo = numberNodes_;
     290          if (numberDo > 0 || whichSolution_ >= 0) {
     291            return 0.5;
     292          } else {
     293            // no solution
     294            return COIN_DBL_MAX; // say infeasible
     295          }
    257296        } else {
    258297            return -1.0;
     
    284323#endif
    285324CbcBranchingObject *
    286 CbcGeneralDepth::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int /*way*/)
     325CbcGeneralDepth::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * info, int /*way*/)
    287326{
    288327    int numberDo = numberNodes_;
    289     if (whichSolution_ >= 0)
     328    if (whichSolution_ >= 0 && (model_->moreSpecialOptions()&33554432)==0)
    290329        numberDo--;
    291330    assert (numberDo > 0);
     
    308347    ClpSimplex * simplex = clpSolver->getModelPtr();
    309348    int numberColumns = simplex->numberColumns();
    310     double * lowerBefore = CoinCopyOfArray(simplex->getColLower(),
    311                                            numberColumns);
    312     double * upperBefore = CoinCopyOfArray(simplex->getColUpper(),
    313                                            numberColumns);
    314     ClpNodeStuff * info = nodeInfo_;
    315     double * weight = new double[numberNodes_];
    316     int * whichNode = new int [numberNodes_];
    317     // Sort
    318     for (iNode = 0; iNode < numberNodes_; iNode++) {
     349    if ((model_->moreSpecialOptions()&33554432)==0) {
     350      double * lowerBefore = CoinCopyOfArray(simplex->getColLower(),
     351                                             numberColumns);
     352      double * upperBefore = CoinCopyOfArray(simplex->getColUpper(),
     353                                             numberColumns);
     354      ClpNodeStuff * info = nodeInfo_;
     355      double * weight = new double[numberNodes_];
     356      int * whichNode = new int [numberNodes_];
     357      // Sort
     358      for (iNode = 0; iNode < numberNodes_; iNode++) {
    319359        if (iNode != whichSolution_) {
    320             double objectiveValue = info->nodeInfo_[iNode]->objectiveValue();
    321             double sumInfeasibilities = info->nodeInfo_[iNode]->sumInfeasibilities();
    322             int numberInfeasibilities = info->nodeInfo_[iNode]->numberInfeasibilities();
    323             double thisWeight = 0.0;
     360          double objectiveValue = info->nodeInfo_[iNode]->objectiveValue();
     361          double sumInfeasibilities = info->nodeInfo_[iNode]->sumInfeasibilities();
     362          int numberInfeasibilities = info->nodeInfo_[iNode]->numberInfeasibilities();
     363          double thisWeight = 0.0;
    324364#if 1
    325             // just closest
    326             thisWeight = 1.0e9 * numberInfeasibilities;
    327             thisWeight += sumInfeasibilities;
    328             thisWeight += 1.0e-7 * objectiveValue;
    329             // Try estimate
    330             thisWeight = info->nodeInfo_[iNode]->estimatedSolution();
     365          // just closest
     366          thisWeight = 1.0e9 * numberInfeasibilities;
     367          thisWeight += sumInfeasibilities;
     368          thisWeight += 1.0e-7 * objectiveValue;
     369          // Try estimate
     370          thisWeight = info->nodeInfo_[iNode]->estimatedSolution();
    331371#else
    332             thisWeight = 1.0e-3 * numberInfeasibilities;
    333             thisWeight += 1.0e-5 * sumInfeasibilities;
    334             thisWeight += objectiveValue;
     372          thisWeight = 1.0e-3 * numberInfeasibilities;
     373          thisWeight += 1.0e-5 * sumInfeasibilities;
     374          thisWeight += objectiveValue;
    335375#endif
    336             whichNode[iProb] = iNode;
    337             weight[iProb++] = thisWeight;
    338         }
    339     }
    340     assert (iProb == numberDo);
    341     CoinSort_2(weight, weight + numberDo, whichNode);
    342     for (iProb = 0; iProb < numberDo; iProb++) {
     376          whichNode[iProb] = iNode;
     377          weight[iProb++] = thisWeight;
     378        }
     379      }
     380      assert (iProb == numberDo);
     381      CoinSort_2(weight, weight + numberDo, whichNode);
     382      for (iProb = 0; iProb < numberDo; iProb++) {
    343383        iNode = whichNode[iProb];
    344384        ClpNode * node = info->nodeInfo_[iNode];
     
    353393#ifdef CHECK_PATH
    354394        if (simplex->numberColumns() == numberColumns_Z) {
    355             bool onOptimal = true;
    356             const double * columnLower = simplex->columnLower();
    357             const double * columnUpper = simplex->columnUpper();
    358             for (int i = 0; i < numberColumns_Z; i++) {
    359                 if (iNode == gotGoodNode_Z)
    360                     printf("good %d %d %g %g\n", iNode, i, columnLower[i], columnUpper[i]);
    361                 if (columnUpper[i] < debuggerSolution_Z[i] || columnLower[i] > debuggerSolution_Z[i] && simplex->isInteger(i)) {
    362                     onOptimal = false;
    363                     break;
    364                 }
    365             }
    366             if (onOptimal) {
    367                 printf("adding to node %x as %d - objs\n", this, iProb);
    368                 for (int j = 0; j <= iProb; j++)
    369                     printf("%d %g\n", j, sub[j].objectiveValue_);
    370             }
     395          bool onOptimal = true;
     396          const double * columnLower = simplex->columnLower();
     397          const double * columnUpper = simplex->columnUpper();
     398          for (int i = 0; i < numberColumns_Z; i++) {
     399            if (iNode == gotGoodNode_Z)
     400              printf("good %d %d %g %g\n", iNode, i, columnLower[i], columnUpper[i]);
     401            if (columnUpper[i] < debuggerSolution_Z[i] || columnLower[i] > debuggerSolution_Z[i] && simplex->isInteger(i)) {
     402              onOptimal = false;
     403              break;
     404            }
     405          }
     406          if (onOptimal) {
     407            printf("adding to node %x as %d - objs\n", this, iProb);
     408            for (int j = 0; j <= iProb; j++)
     409              printf("%d %g\n", j, sub[j].objectiveValue_);
     410          }
    371411        }
    372412#endif
    373     }
    374     delete [] weight;
    375     delete [] whichNode;
    376     const double * lower = solver->getColLower();
    377     const double * upper = solver->getColUpper();
    378     // restore bounds
    379     for ( int j = 0; j < numberColumns; j++) {
     413      }
     414      delete [] weight;
     415      delete [] whichNode;
     416      const double * lower = solver->getColLower();
     417      const double * upper = solver->getColUpper();
     418      // restore bounds
     419      for ( int j = 0; j < numberColumns; j++) {
    380420        if (lowerBefore[j] != lower[j])
    381             solver->setColLower(j, lowerBefore[j]);
     421          solver->setColLower(j, lowerBefore[j]);
    382422        if (upperBefore[j] != upper[j])
    383             solver->setColUpper(j, upperBefore[j]);
    384     }
    385     delete [] upperBefore;
    386     delete [] lowerBefore;
     423          solver->setColUpper(j, upperBefore[j]);
     424      }
     425      delete [] upperBefore;
     426      delete [] lowerBefore;
     427    } else {
     428      // from diving
     429      CbcSubProblem ** nodes = reinterpret_cast<CbcSubProblem **>
     430        (model_->temporaryPointer());
     431      assert (nodes);
     432      int adjustDepth=info->depth_;
     433      assert (numberDo);
     434      numberNodes_=0;
     435      for (iProb = 0; iProb < numberDo; iProb++) {
     436        if ((nodes[iProb]->problemStatus_&2)==0) {
     437          // create subproblem (and swap way and/or make inactive)
     438          sub[numberNodes_].takeOver(*nodes[iProb],true);
     439          // but adjust depth
     440          sub[numberNodes_].depth_+=adjustDepth;
     441          numberNodes_++;
     442        }
     443        delete nodes[iProb];
     444      }
     445      branch->numberSubProblems_ = numberNodes_;
     446      branch->numberSubLeft_ = numberNodes_;
     447      branch->setNumberBranches(numberNodes_);
     448      if (!numberNodes_) {
     449        // infeasible
     450        delete branch;
     451        branch=NULL;
     452      }
     453      delete [] nodes;
     454    }
    387455    return branch;
    388456}
  • trunk/Cbc/src/CbcGeneralDepth.hpp

    r1854 r1880  
    7070    inline void setMaximumDepth(int value) {
    7171        maximumDepth_ = value;
     72    }
     73    /// Return number of nodes
     74    inline int numberNodes() const {
     75        return numberNodes_;
    7276    }
    7377    /// Get which solution
  • trunk/Cbc/src/CbcHeuristic.cpp

    r1864 r1880  
    470470    }
    471471    randomNumberGenerator_.setSeed(value);
     472}
     473// Get seed
     474int
     475CbcHeuristic::getSeed() const
     476{
     477  return randomNumberGenerator_.getSeed();
    472478}
    473479
     
    648654    // valueNow,numberRowsNow,numberColumnsNow,
    649655    // valueStart,numberRowsStart,numberColumnsStart);
    650     if (10*numberRowsNow < 8*numberRowsStart)
     656    if (10*numberRowsNow < 8*numberRowsStart || 10*numberColumnsNow < 7*numberColumnsStart)
    651657        return valueNow / valueStart;
    652658    else if (10*numberRowsNow < 9*numberRowsStart)
     
    826832                    }
    827833                }
    828             } else if (ratio > fractionSmall && after > 300) {
     834            } else if (ratio > fractionSmall && after > 300 && numberNodes >=0) {
    829835                returnCode = -1;
    830836            }
  • trunk/Cbc/src/CbcHeuristic.hpp

    r1822 r1880  
    237237    /// Set random number generator seed
    238238    void setSeed(int value);
     239    /// Get random number generator seed
     240    int getSeed() const;
    239241    /// Sets decay factor (for howOften) on failure
    240242    inline void setDecayFactor(double value) {
  • trunk/Cbc/src/CbcHeuristicDive.cpp

    r1854 r1880  
    1111#include "CbcHeuristicDive.hpp"
    1212#include "CbcStrategy.hpp"
     13#include "CbcModel.hpp"
     14#include "CbcSubProblem.hpp"
    1315#include "OsiAuxInfo.hpp"
    1416#include  "CoinTime.hpp"
     
    190192}
    191193
    192 // See if diving will give better solution
    193 // Sets value of solution
    194 // Returns 1 if solution, 0 if not
    195 int
    196 CbcHeuristicDive::solution(double & solutionValue,
    197                            double * betterSolution)
    198 {
    199     int nodeCount = model_->getNodeCount();
    200     if (feasibilityPumpOptions_>0 && (nodeCount % feasibilityPumpOptions_) != 0)
    201         return 0;
    202 #ifdef DIVE_DEBUG
    203     std::cout << "solutionValue = " << solutionValue << std::endl;
    204 #endif
    205     ++numCouldRun_;
    206 
    207     // test if the heuristic can run
    208     if (!canHeuristicRun())
    209         return 0;
    210 
    211 #ifdef JJF_ZERO
    212     // See if to do
    213     if (!when() || (when() % 10 == 1 && model_->phase() != 1) ||
    214             (when() % 10 == 2 && (model_->phase() != 2 && model_->phase() != 3)))
    215         return 0; // switched off
    216 #endif
    217 
     194// inner part of dive
     195int
     196CbcHeuristicDive::solution(double & solutionValue, int & numberNodes,
     197                           int & numberCuts, OsiRowCut ** cuts,
     198                           CbcSubProblem ** & nodes,
     199                           double * newSolution)
     200{
    218201#ifdef DIVE_DEBUG
    219202    int nRoundInfeasible = 0;
    220203    int nRoundFeasible = 0;
     204#endif
    221205    int reasonToStop = 0;
    222 #endif
    223206    double time1 = CoinCpuTime();
    224207    int numberSimplexIterations = 0;
    225208    int maxSimplexIterations = (model_->getNodeCount()) ? maxSimplexIterations_
    226209                               : maxSimplexIterationsAtRoot_;
    227 
     210    // but can't be exactly coin_int_max
     211    maxSimplexIterations = CoinMin(maxSimplexIterations,COIN_INT_MAX>>3);
    228212    OsiSolverInterface * solver = cloneBut(6); // was model_->solver()->clone();
    229213# ifdef COIN_HAS_CLP
     
    231215    = dynamic_cast<OsiClpSolverInterface *> (solver);
    232216    if (clpSolver) {
     217      ClpSimplex * clpSimplex = clpSolver->getModelPtr();
     218      int oneSolveIts = clpSimplex->maximumIterations();
     219      oneSolveIts = CoinMin(1000+2*(clpSimplex->numberRows()+clpSimplex->numberColumns()),oneSolveIts);
     220      clpSimplex->setMaximumIterations(oneSolveIts);
     221      if (!nodes) {
    233222        // say give up easily
    234         ClpSimplex * clpSimplex = clpSolver->getModelPtr();
    235223        clpSimplex->setMoreSpecialOptions(clpSimplex->moreSpecialOptions() | 64);
     224      } else {
     225        // get ray
     226        int specialOptions = clpSimplex->specialOptions();
     227        specialOptions &= ~0x3100000;
     228        specialOptions |= 32;
     229        clpSimplex->setSpecialOptions(specialOptions);
     230        clpSolver->setSpecialOptions(clpSolver->specialOptions() | 1048576);
     231        if ((model_->moreSpecialOptions()&16777216)!=0) {
     232          // cutoff is constraint
     233          clpSolver->setDblParam(OsiDualObjectiveLimit, COIN_DBL_MAX);
     234        }
     235      }
    236236    }
    237237# endif
     
    268268    // Get solution array for heuristic solution
    269269    int numberColumns = solver->getNumCols();
    270     double * newSolution = new double [numberColumns];
    271270    memcpy(newSolution, solution, numberColumns*sizeof(double));
    272271
    273272    // vectors to store the latest variables fixed at their bounds
    274273    int* columnFixed = new int [numberIntegers];
    275     double* originalBound = new double [numberIntegers];
     274    double* originalBound = new double [numberIntegers+2*numberColumns];
     275    double * lowerBefore = originalBound+numberIntegers;
     276    double * upperBefore = lowerBefore+numberColumns;
     277    memcpy(lowerBefore,lower,numberColumns*sizeof(double));
     278    memcpy(upperBefore,upper,numberColumns*sizeof(double));
     279    double * lastDjs=newSolution+numberColumns;
    276280    bool * fixedAtLowerBound = new bool [numberIntegers];
    277281    PseudoReducedCost * candidate = new PseudoReducedCost [numberIntegers];
     
    279283
    280284    int maxNumberAtBoundToFix = static_cast<int> (floor(percentageToFix_ * numberIntegers));
     285    assert (!maxNumberAtBoundToFix||!nodes);
    281286
    282287    // count how many fractional variables
     
    308313        bool canRound = selectVariableToBranch(solver, newSolution,
    309314                                               bestColumn, bestRound);
    310 
    311315        // if the solution is not trivially roundable, we don't try to round;
    312316        // if the solution is trivially roundable, we try to round. However,
     
    340344                nRoundFeasible++;
    341345#endif
    342                 // Round all the fractional variables
    343                 for (int i = 0; i < numberIntegers; i++) {
     346                if (!nodes||bestColumn<0) {
     347                  // Round all the fractional variables
     348                  for (int i = 0; i < numberIntegers; i++) {
    344349                    int iColumn = integerVariable[i];
    345350                    double value = newSolution[iColumn];
    346351                    if (fabs(floor(value + 0.5) - value) > integerTolerance) {
    347                         assert(downLocks_[i] == 0 || upLocks_[i] == 0);
    348                         if (downLocks_[i] == 0 && upLocks_[i] == 0) {
    349                             if (direction * objective[iColumn] >= 0.0)
    350                                 newSolution[iColumn] = floor(value);
    351                             else
    352                                 newSolution[iColumn] = ceil(value);
    353                         } else if (downLocks_[i] == 0)
    354                             newSolution[iColumn] = floor(value);
    355                         else
    356                             newSolution[iColumn] = ceil(value);
     352                      assert(downLocks_[i] == 0 || upLocks_[i] == 0);
     353                      if (downLocks_[i] == 0 && upLocks_[i] == 0) {
     354                        if (direction * objective[iColumn] >= 0.0)
     355                          newSolution[iColumn] = floor(value);
     356                        else
     357                          newSolution[iColumn] = ceil(value);
     358                      } else if (downLocks_[i] == 0)
     359                        newSolution[iColumn] = floor(value);
     360                      else
     361                        newSolution[iColumn] = ceil(value);
    357362                    }
    358                 }
    359                 break;
    360             }
     363                  }
     364                  break;
     365                } else {
     366                  // can't round if going to use in branching
     367                  int i;
     368                  for (i = 0; i < numberIntegers; i++) {
     369                    int iColumn = integerVariable[i];
     370                    double value = newSolution[bestColumn];
     371                    if (fabs(floor(value + 0.5) - value) > integerTolerance) {
     372                      if (iColumn==bestColumn) {
     373                        assert(downLocks_[i] == 0 || upLocks_[i] == 0);
     374                        double obj = objective[bestColumn];
     375                        if (downLocks_[i] == 0 && upLocks_[i] == 0) {
     376                          if (direction * obj >= 0.0)
     377                            bestRound=-1;
     378                          else
     379                            bestRound=1;
     380                        } else if (downLocks_[i] == 0)
     381                          bestRound=-1;
     382                        else
     383                          bestRound=1;
     384                        break;
     385                      }
     386                    }
     387                  }
     388                }
     389            }
    361390#ifdef DIVE_DEBUG
    362391            else
     
    549578
    550579        double originalBoundBestColumn;
     580        double bestColumnValue;
     581        int whichWay;
    551582        if (bestColumn >= 0) {
     583            bestColumnValue = newSolution[bestColumn];
    552584            if (bestRound < 0) {
    553585                originalBoundBestColumn = upper[bestColumn];
    554                 solver->setColUpper(bestColumn, floor(newSolution[bestColumn]));
     586                solver->setColUpper(bestColumn, floor(bestColumnValue));
     587                whichWay=0;
    555588            } else {
    556589                originalBoundBestColumn = lower[bestColumn];
    557                 solver->setColLower(bestColumn, ceil(newSolution[bestColumn]));
     590                solver->setColLower(bestColumn, ceil(bestColumnValue));
     591                whichWay=1;
    558592            }
    559593        } else {
     
    562596        int originalBestRound = bestRound;
    563597        int saveModelOptions = model_->specialOptions();
     598       
    564599        while (1) {
    565600
     
    567602            solver->resolve();
    568603            model_->setSpecialOptions(saveModelOptions);
    569             if (!solver->isAbandoned()) {
     604            if (!solver->isAbandoned()&&!solver->isIterationLimitReached()) {
    570605                numberSimplexIterations += solver->getIterationCount();
    571606            } else {
    572607                numberSimplexIterations = maxSimplexIterations + 1;
     608                reasonToStop += 100;
    573609                break;
    574610            }
    575611
    576612            if (!solver->isProvenOptimal()) {
     613                if (nodes) {
     614                  if (solver->isProvenPrimalInfeasible()) {
     615                    if (maxSimplexIterationsAtRoot_!=COIN_INT_MAX) {
     616                      // stop now
     617                      printf("stopping on first infeasibility\n");
     618                      break;
     619                    } else if (cuts) {
     620                      // can do conflict cut
     621                      printf("could do intermediate conflict cut\n");
     622                      bool localCut;
     623                      OsiRowCut * cut = model_->conflictCut(solver,localCut);
     624                      if (cut) {
     625                        if (!localCut) {
     626                          model_->makePartialCut(cut,solver);
     627                          cuts[numberCuts++]=cut;
     628                        } else {
     629                          delete cut;
     630                        }
     631                      }
     632                    }
     633                  } else {
     634                    reasonToStop += 10;
     635                    break;
     636                  }
     637                }
    577638                if (numberAtBoundFixed > 0) {
    578639                    // Remove the bound fix for variables that were at bounds
     
    587648                } else if (bestRound == originalBestRound) {
    588649                    bestRound *= (-1);
     650                    whichWay |=2;
    589651                    if (bestRound < 0) {
    590652                        solver->setColLower(bestColumn, originalBoundBestColumn);
    591                         solver->setColUpper(bestColumn, floor(newSolution[bestColumn]));
     653                        solver->setColUpper(bestColumn, floor(bestColumnValue));
    592654                    } else {
    593                         solver->setColLower(bestColumn, ceil(newSolution[bestColumn]));
     655                        solver->setColLower(bestColumn, ceil(bestColumnValue));
    594656                        solver->setColUpper(bestColumn, originalBoundBestColumn);
    595657                    }
     
    602664        if (!solver->isProvenOptimal() ||
    603665                direction*solver->getObjValue() >= solutionValue) {
    604 #ifdef DIVE_DEBUG
    605             reasonToStop = 1;
    606 #endif
    607             break;
    608         }
    609 
    610         if (iteration > maxIterations_) {
    611 #ifdef DIVE_DEBUG
    612             reasonToStop = 2;
    613 #endif
    614             break;
    615         }
    616 
    617         if (CoinCpuTime() - time1 > maxTime_) {
    618 #ifdef DIVE_DEBUG
    619             reasonToStop = 3;
    620 #endif
    621             break;
    622         }
    623 
    624         if (numberSimplexIterations > maxSimplexIterations) {
    625 #ifdef DIVE_DEBUG
    626             reasonToStop = 4;
    627 #endif
     666            reasonToStop += 1;
     667        } else if (iteration > maxIterations_) {
     668            reasonToStop += 2;
     669        } else if (CoinCpuTime() - time1 > maxTime_) {
     670            reasonToStop += 3;
     671        } else if (numberSimplexIterations > maxSimplexIterations) {
     672            reasonToStop += 4;
    628673            // also switch off
    629674#ifdef CLP_INVESTIGATE
     
    632677#endif
    633678            when_ = 0;
    634             break;
    635         }
    636 
    637         if (solver->getIterationCount() > 1000 && iteration > 3) {
    638 #ifdef DIVE_DEBUG
    639             reasonToStop = 5;
    640 #endif
     679        } else if (solver->getIterationCount() > 1000 && iteration > 3 && !nodes) {
     680            reasonToStop += 5;
    641681            // also switch off
    642682#ifdef CLP_INVESTIGATE
     
    645685#endif
    646686            when_ = 0;
    647             break;
    648687        }
    649688
    650689        memcpy(newSolution, solution, numberColumns*sizeof(double));
    651690        numberFractionalVariables = 0;
     691        double sumFractionalVariables=0.0;
    652692        for (int i = 0; i < numberIntegers; i++) {
    653693            int iColumn = integerVariable[i];
    654694            double value = newSolution[iColumn];
    655             if (fabs(floor(value + 0.5) - value) > integerTolerance) {
     695            double away = fabs(floor(value + 0.5) - value);
     696            if (away > integerTolerance) {
    656697                numberFractionalVariables++;
    657             }
    658         }
    659 
    660     }
    661 
    662 
    663     double * rowActivity = new double[numberRows];
    664     memset(rowActivity, 0, numberRows*sizeof(double));
    665 
     698                sumFractionalVariables += away;
     699            }
     700        }
     701        if (nodes) {
     702          // save information
     703          //branchValues[numberNodes]=bestColumnValue;
     704          //statuses[numberNodes]=whichWay+(bestColumn<<2);
     705          //bases[numberNodes]=solver->getWarmStart();
     706          ClpSimplex * simplex = clpSolver->getModelPtr();
     707          CbcSubProblem * sub =
     708            new CbcSubProblem(clpSolver,lowerBefore,upperBefore,
     709                          simplex->statusArray(),numberNodes);
     710          nodes[numberNodes]=sub;
     711          // other stuff
     712          sub->branchValue_=bestColumnValue;
     713          sub->problemStatus_=whichWay;
     714          sub->branchVariable_=bestColumn;
     715          sub->objectiveValue_ = simplex->objectiveValue();
     716          sub->sumInfeasibilities_ = sumFractionalVariables;
     717          sub->numberInfeasibilities_ = numberFractionalVariables;
     718          printf("DiveNode %d column %d way %d bvalue %g obj %g\n",
     719                 numberNodes,sub->branchVariable_,sub->problemStatus_,
     720                 sub->branchValue_,sub->objectiveValue_);
     721          numberNodes++;
     722          if (solver->isProvenOptimal()) {
     723            memcpy(lastDjs,solver->getReducedCost(),numberColumns*sizeof(double));
     724            memcpy(lowerBefore,lower,numberColumns*sizeof(double));
     725            memcpy(upperBefore,upper,numberColumns*sizeof(double));
     726          }
     727        }
     728        if (!numberFractionalVariables||reasonToStop)
     729          break;
     730    }
     731    if (nodes) {
     732      printf("Exiting dive for reason %d\n",reasonToStop);
     733      if (reasonToStop>1) {
     734        printf("problems in diving\n");
     735        int whichWay=nodes[numberNodes-1]->problemStatus_;
     736        CbcSubProblem * sub;
     737        if ((whichWay&2)==0) {
     738          // leave both ways
     739          sub = new CbcSubProblem(*nodes[numberNodes-1]);
     740          nodes[numberNodes++]=sub;
     741        } else {
     742          sub = nodes[numberNodes-1];
     743        }
     744        if ((whichWay&1)==0)
     745          sub->problemStatus_=whichWay|1;
     746        else
     747          sub->problemStatus_=whichWay&~1;
     748      }
     749      if (!numberNodes) {
     750        // was good at start! - create fake
     751        clpSolver->resolve();
     752        ClpSimplex * simplex = clpSolver->getModelPtr();
     753        CbcSubProblem * sub =
     754          new CbcSubProblem(clpSolver,lowerBefore,upperBefore,
     755                            simplex->statusArray(),numberNodes);
     756        nodes[numberNodes]=sub;
     757        // other stuff
     758        sub->branchValue_=0.0;
     759        sub->problemStatus_=0;
     760        sub->branchVariable_=-1;
     761        sub->objectiveValue_ = simplex->objectiveValue();
     762        sub->sumInfeasibilities_ = 0.0;
     763        sub->numberInfeasibilities_ = 0;
     764        printf("DiveNode %d column %d way %d bvalue %g obj %g\n",
     765               numberNodes,sub->branchVariable_,sub->problemStatus_,
     766               sub->branchValue_,sub->objectiveValue_);
     767        numberNodes++;
     768        assert (solver->isProvenOptimal());
     769      }
     770      nodes[numberNodes-1]->problemStatus_ |= 256*reasonToStop;
     771      // use djs as well
     772      if (solver->isProvenPrimalInfeasible()&&cuts) {
     773        // can do conflict cut and re-order
     774        printf("could do final conflict cut\n");
     775        bool localCut;
     776        OsiRowCut * cut = model_->conflictCut(solver,localCut);
     777        if (cut) {
     778          printf("cut - need to use conflict and previous djs\n");
     779          if (!localCut) {
     780            model_->makePartialCut(cut,solver);
     781            cuts[numberCuts++]=cut;
     782          } else {
     783            delete cut;
     784          }
     785        } else {
     786          printf("bad conflict - just use previous djs\n");
     787        }
     788      }
     789    }
     790   
    666791    // re-compute new solution value
    667792    double objOffset = 0.0;
     
    669794    newSolutionValue = -objOffset;
    670795    for (int i = 0 ; i < numberColumns ; i++ )
    671         newSolutionValue += objective[i] * newSolution[i];
     796      newSolutionValue += objective[i] * newSolution[i];
    672797    newSolutionValue *= direction;
    673798    //printf("new solution value %g %g\n",newSolutionValue,solutionValue);
    674     if (newSolutionValue < solutionValue) {
    675         // paranoid check
    676         memset(rowActivity, 0, numberRows*sizeof(double));
    677         for (int i = 0; i < numberColumns; i++) {
    678             int j;
    679             double value = newSolution[i];
    680             if (value) {
    681                 for (j = columnStart[i];
    682                         j < columnStart[i] + columnLength[i]; j++) {
    683                     int iRow = row[j];
    684                     rowActivity[iRow] += value * element[j];
    685                 }
    686             }
    687         }
    688         // check was approximately feasible
    689         bool feasible = true;
    690         for (int i = 0; i < numberRows; i++) {
    691             if (rowActivity[i] < rowLower[i]) {
    692                 if (rowActivity[i] < rowLower[i] - 1000.0*primalTolerance)
    693                     feasible = false;
    694             } else if (rowActivity[i] > rowUpper[i]) {
    695                 if (rowActivity[i] > rowUpper[i] + 1000.0*primalTolerance)
    696                     feasible = false;
    697             }
    698         }
    699         for (int i = 0; i < numberIntegers; i++) {
    700             int iColumn = integerVariable[i];
    701             double value = newSolution[iColumn];
    702             if (fabs(floor(value + 0.5) - value) > integerTolerance) {
    703                 feasible = false;
    704                 break;
    705             }
    706         }
    707         if (feasible) {
    708             // new solution
    709             memcpy(betterSolution, newSolution, numberColumns*sizeof(double));
    710             solutionValue = newSolutionValue;
    711             //printf("** Solution of %g found by CbcHeuristicDive\n",newSolutionValue);
    712             returnCode = 1;
    713         } else {
    714             // Can easily happen
    715             //printf("Debug CbcHeuristicDive giving bad solution\n");
    716         }
     799    if (newSolutionValue < solutionValue && !reasonToStop) {
     800      double * rowActivity = new double[numberRows];
     801      memset(rowActivity, 0, numberRows*sizeof(double));
     802      // paranoid check
     803      memset(rowActivity, 0, numberRows*sizeof(double));
     804      for (int i = 0; i < numberColumns; i++) {
     805        int j;
     806        double value = newSolution[i];
     807        if (value) {
     808          for (j = columnStart[i];
     809               j < columnStart[i] + columnLength[i]; j++) {
     810            int iRow = row[j];
     811            rowActivity[iRow] += value * element[j];
     812          }
     813        }
     814      }
     815      // check was approximately feasible
     816      bool feasible = true;
     817      for (int i = 0; i < numberRows; i++) {
     818        if (rowActivity[i] < rowLower[i]) {
     819          if (rowActivity[i] < rowLower[i] - 1000.0*primalTolerance)
     820            feasible = false;
     821        } else if (rowActivity[i] > rowUpper[i]) {
     822          if (rowActivity[i] > rowUpper[i] + 1000.0*primalTolerance)
     823            feasible = false;
     824        }
     825      }
     826      for (int i = 0; i < numberIntegers; i++) {
     827        int iColumn = integerVariable[i];
     828        double value = newSolution[iColumn];
     829        if (fabs(floor(value + 0.5) - value) > integerTolerance) {
     830          feasible = false;
     831          break;
     832        }
     833      }
     834      if (feasible) {
     835        // new solution
     836        solutionValue = newSolutionValue;
     837        //printf("** Solution of %g found by CbcHeuristicDive\n",newSolutionValue);
     838        //if (cuts)
     839        //clpSolver->getModelPtr()->writeMps("good8.mps", 2);
     840        returnCode = 1;
     841      } else {
     842        // Can easily happen
     843        printf("Debug CbcHeuristicDive giving bad solution\n");
     844      }
     845      delete [] rowActivity;
    717846    }
    718847
     
    726855#endif
    727856
    728     delete [] newSolution;
    729857    delete [] columnFixed;
    730858    delete [] originalBound;
    731859    delete [] fixedAtLowerBound;
    732860    delete [] candidate;
    733     delete [] rowActivity;
    734861    delete [] random;
    735862    delete [] downArray_;
     
    739866    delete solver;
    740867    return returnCode;
     868}
     869// See if diving will give better solution
     870// Sets value of solution
     871// Returns 1 if solution, 0 if not
     872int
     873CbcHeuristicDive::solution(double & solutionValue,
     874                           double * betterSolution)
     875{
     876    int nodeCount = model_->getNodeCount();
     877    if (feasibilityPumpOptions_>0 && (nodeCount % feasibilityPumpOptions_) != 0)
     878        return 0;
     879#ifdef DIVE_DEBUG
     880    std::cout << "solutionValue = " << solutionValue << std::endl;
     881#endif
     882    ++numCouldRun_;
     883
     884    // test if the heuristic can run
     885    if (!canHeuristicRun())
     886        return 0;
     887
     888#ifdef JJF_ZERO
     889    // See if to do
     890    if (!when() || (when() % 10 == 1 && model_->phase() != 1) ||
     891            (when() % 10 == 2 && (model_->phase() != 2 && model_->phase() != 3)))
     892        return 0; // switched off
     893#endif
     894    // Get solution array for heuristic solution
     895    int numberColumns = model_->solver()->getNumCols();
     896    double * newSolution = new double [numberColumns];
     897    int numberCuts=0;
     898    int numberNodes=-1;
     899    CbcSubProblem ** nodes=NULL;
     900    int returnCode=solution(solutionValue,numberNodes,numberCuts,
     901                            NULL,nodes,
     902                            newSolution);
     903  if (returnCode==1)
     904    memcpy(betterSolution, newSolution, numberColumns*sizeof(double));
     905
     906    delete [] newSolution;
     907    return returnCode;
     908}
     909/* returns 0 if no solution, 1 if valid solution
     910   with better objective value than one passed in
     911   also returns list of nodes
     912   This does Fractional Diving
     913*/
     914int
     915CbcHeuristicDive::fathom(CbcModel * model, int & numberNodes,
     916                         CbcSubProblem ** & nodes)
     917{
     918  double solutionValue = model->getCutoff();
     919  numberNodes=0;
     920  // Get solution array for heuristic solution
     921  int numberColumns = model_->solver()->getNumCols();
     922  double * newSolution = new double [4*numberColumns];
     923  double * lastDjs = newSolution+numberColumns;
     924  double * originalLower = lastDjs+numberColumns;
     925  double * originalUpper = originalLower+numberColumns;
     926  memcpy(originalLower,model_->solver()->getColLower(),
     927         numberColumns*sizeof(double));
     928  memcpy(originalUpper,model_->solver()->getColUpper(),
     929         numberColumns*sizeof(double));
     930  int numberCuts=0;
     931  OsiRowCut ** cuts = NULL; //new OsiRowCut * [maxIterations_];
     932  nodes=new CbcSubProblem * [maxIterations_+2];
     933  int returnCode=solution(solutionValue,numberNodes,numberCuts,
     934                          cuts,nodes,
     935                          newSolution);
     936
     937  if (returnCode==1) {
     938    // copy to best solution ? or put in solver
     939    printf("Solution from heuristic fathom\n");
     940  }
     941  int numberFeasibleNodes=numberNodes;
     942  if (returnCode!=1)
     943    numberFeasibleNodes--;
     944  if (numberFeasibleNodes>0) {
     945    CoinWarmStartBasis * basis = nodes[numberFeasibleNodes-1]->status_;
     946    //double * sort = new double [numberFeasibleNodes];
     947    //int * whichNode = new int [numberFeasibleNodes];
     948    //int numberNodesNew=0;
     949    // use djs on previous unless feasible
     950    for (int iNode=0;iNode<numberFeasibleNodes;iNode++) {
     951      CbcSubProblem * sub = nodes[iNode];
     952      double branchValue = sub->branchValue_;
     953      int iStatus=sub->problemStatus_;
     954      int iColumn = sub->branchVariable_;
     955      bool secondBranch = (iStatus&2)!=0;
     956      bool branchUp;
     957      if (!secondBranch)
     958        branchUp = (iStatus&1)!=0;
     959      else
     960        branchUp = (iStatus&1)==0;
     961      double djValue=lastDjs[iColumn];
     962      sub->djValue_=fabs(djValue);
     963      if (!branchUp&&floor(branchValue)==originalLower[iColumn]
     964          &&basis->getStructStatus(iColumn) == CoinWarmStartBasis::atLowerBound) {
     965        if (djValue>0.0) {
     966          // naturally goes to LB
     967          printf("ignoring branch down on %d (node %d) from value of %g - branch was %s - dj %g\n",
     968                 iColumn,iNode,branchValue,secondBranch ? "second" : "first",
     969                 djValue);
     970          sub->problemStatus_ |= 4;
     971          //} else {
     972          // put on list
     973          //sort[numberNodesNew]=djValue;
     974          //whichNode[numberNodesNew++]=iNode;
     975        }
     976      } else if (branchUp&&ceil(branchValue)==originalUpper[iColumn]
     977          &&basis->getStructStatus(iColumn) == CoinWarmStartBasis::atUpperBound) {
     978        if (djValue<0.0) {
     979          // naturally goes to UB
     980          printf("ignoring branch up on %d (node %d) from value of %g - branch was %s - dj %g\n",
     981                 iColumn,iNode,branchValue,secondBranch ? "second" : "first",
     982                 djValue);
     983          sub->problemStatus_ |= 4;
     984          //} else {
     985          // put on list
     986          //sort[numberNodesNew]=-djValue;
     987          //whichNode[numberNodesNew++]=iNode;
     988        }
     989      }
     990    }
     991    // use conflict to order nodes
     992    for (int iCut=0;iCut<numberCuts;iCut++) {
     993    }
     994    //CoinSort_2(sort,sort+numberNodesNew,whichNode);
     995    // create nodes
     996    // last node will have one way already done
     997  }
     998  for (int iCut=0;iCut<numberCuts;iCut++) {
     999    delete cuts[iCut];
     1000  }
     1001  delete [] cuts;
     1002  delete [] newSolution;
     1003  return returnCode;
    7411004}
    7421005
  • trunk/Cbc/src/CbcHeuristicDive.hpp

    r1854 r1880  
    88
    99#include "CbcHeuristic.hpp"
     10class CbcSubProblem;
     11class OsiRowCut;
    1012struct PseudoReducedCost {
    1113    int var;
     
    5961    virtual int solution(double & objectiveValue,
    6062                         double * newSolution);
     63    /// inner part of dive
     64  int solution(double & objectiveValue, int & numberNodes,
     65                 int & numberCuts, OsiRowCut ** cuts,
     66                 CbcSubProblem ** & nodes,
     67                 double * newSolution);
     68    /** returns 0 if no solution, 1 if valid solution
     69        with better objective value than one passed in
     70        also returns list of nodes
     71        This does Fractional Diving
     72    */
     73    int fathom(CbcModel * model, int & numberNodes,CbcSubProblem ** & nodes);
    6174
    6275    /// Validate model i.e. sets when_ to 0 if necessary (may be NULL)
     
    7992    void setMaxSimplexIterations(int value) {
    8093        maxSimplexIterations_ = value;
     94    }
     95    /// Get maximum number of simplex iterations
     96    inline int maxSimplexIterations() const {
     97        return maxSimplexIterations_;
    8198    }
    8299
  • trunk/Cbc/src/CbcHeuristicFPump.cpp

    r1871 r1880  
    496496            maximumPasses = 100; // feasibility problem?
    497497        }
     498        randomNumberGenerator_.randomize();
    498499        if (model_->getRandomSeed()!=-1)
    499500          clpSolver->getModelPtr()->setRandomSeed(randomNumberGenerator_.getSeed());
     501        clpSolver->getModelPtr()->randomNumberGenerator()->randomize();
    500502      }
    501503    }
     
    10231025                    offset += costValue * newSolution[iColumn];
    10241026                }
     1027                if (numberPasses==1 && !totalNumberPasses && (model_->specialOptions()&8388608)!=0) {
     1028                  // doing multiple solvers - make a real difference - flip 5%
     1029                  for (i = 0; i < numberIntegers; i++) {
     1030                    int iColumn = integerVariable[i];
     1031                    double value = floor(newSolution[iColumn]+0.5);
     1032                    if (fabs(value-solution[iColumn])>primalTolerance) {
     1033                      value = randomNumberGenerator_.randomDouble();
     1034                      if(value<0.05) {
     1035                        //printf("Flipping %d - random %g\n",iColumn,value);
     1036                        solver->setObjCoeff(iColumn,-solver->getObjCoefficients()[iColumn]);
     1037                      }
     1038                    }
     1039                  }
     1040                }
    10251041                solver->setDblParam(OsiObjOffset, -offset);
    10261042                if (!general && false) {
  • trunk/Cbc/src/CbcModel.cpp

    r1876 r1880  
    2424#include <cmath>
    2525#include <cfloat>
    26 
    2726#ifdef COIN_HAS_CLP
    2827// include Presolve from Clp
     
    447446        if (!numberContinuousObj && numberIntegerObj <= 5 && numberIntegerWeight <= 100 &&
    448447                numberIntegerObj*3 < numberObjects_ && !parentModel_ && solver_->getNumRows() > 100)
    449             iType = 3 + 4;
     448          iType = 1 + 4 + ((moreSpecialOptions_&536870912)==0) ? 2 : 0;
    450449        else if (!numberContinuousObj && numberIntegerObj <= 100 &&
    451450                 numberIntegerObj*5 < numberObjects_ && numberIntegerWeight <= 100 &&
    452451                 !parentModel_ &&
    453452                 solver_->getNumRows() > 100 && cost != -COIN_DBL_MAX)
    454             iType = 2 + 4;
     453          iType = 4 + ((moreSpecialOptions_&536870912)==0) ? 2 : 0;
    455454        else if (!numberContinuousObj && numberIntegerObj <= 100 &&
    456455                 numberIntegerObj*5 < numberObjects_ &&
     
    511510                << general << CoinMessageEol ;
    512511                // switch off clp type branching
    513                 fastNodeDepth_ = -1;
     512                // no ? fastNodeDepth_ = -1;
    514513                int highPriority = (branchOnSatisfied) ? -999 : 100;
    515514                for (int i = 0; i < numberObjects_; i++) {
     
    10721071                                          messages())
    10731072                << value << CoinMessageEol ;
    1074                 setDblParam(CbcModel::CbcCutoffIncrement, value*0.999) ;
     1073                setDblParam(CbcModel::CbcCutoffIncrement, CoinMax(value*0.999,value-1.0e-4)) ;
    10751074            }
    10761075        }
     
    15961595
    15971596{
     1597  if (!parentModel_)
    15981598    /*
    15991599      Capture a time stamp before we start (unless set).
     
    21592159            // switch off fast nodes for now
    21602160            fastNodeDepth_ = -1;
     2161            moreSpecialOptions_ &= ~33554432; // no diving
    21612162        }
    21622163        if (numberThreads_ > 0) {
     
    22272228
    22282229    // add cutoff as constraint if wanted
    2229     if ((moreSpecialOptions_&16777216)!=0) {
     2230    if (cutoffRowNumber_==-2) {
    22302231      if (!parentModel_) {
    22312232        int numberColumns=solver_->getNumCols();
     
    22412242        if (n) {
    22422243          double cutoff=getCutoff();
     2244          // relax a little bit
     2245          cutoff += 1.0e-4;
    22432246          double offset;
    22442247          solver_->getDblParam(OsiObjOffset, offset);
     2248          cutoffRowNumber_ = solver_->getNumRows();
    22452249          solver_->addRow(n,indices,obj,-COIN_DBL_MAX,CoinMin(cutoff,1.0e25)+offset);
    22462250        } else {
    22472251          // no objective!
    2248           moreSpecialOptions_ &= ~16777216;
     2252          cutoffRowNumber_ = -1;
    22492253        }
    22502254        delete [] indices;
     
    22522256      } else {
    22532257        // switch off
    2254         moreSpecialOptions_ &= ~16777216;
     2258        cutoffRowNumber_ = -1;
    22552259      }
    22562260    }
     
    24892493      }
    24902494      CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis *> (solver_->getEmptyWarmStart());
     2495      int numberIterations=0;
    24912496      for (int i=0;i<numberModels;i++) {
    24922497        rootModels[i]=new CbcModel(*this);
     
    24982503        // use seed
    24992504        rootModels[i]->setSpecialOptions(specialOptions_ |(4194304|8388608));
     2505        rootModels[i]->setMoreSpecialOptions(moreSpecialOptions_ &
     2506                                             (~134217728));
    25002507        rootModels[i]->solver_->setWarmStart(basis);
    25012508#ifdef COIN_HAS_CLP
     
    25032510          = dynamic_cast<OsiClpSolverInterface *> (rootModels[i]->solver_);
    25042511        if (clpSolver) {
    2505           clpSolver->getModelPtr()->setRandomSeed(newSeed+20000000*i);
    2506           clpSolver->getModelPtr()->allSlackBasis();
     2512          ClpSimplex * simplex = clpSolver->getModelPtr();
     2513          if (defaultHandler_)
     2514            simplex->setDefaultMessageHandler();
     2515          simplex->setRandomSeed(newSeed+20000000*i);
     2516          simplex->allSlackBasis();
     2517          int logLevel=simplex->logLevel();
     2518          if (logLevel==1)
     2519            simplex->setLogLevel(0);
     2520          if (i!=0) {
     2521            double random = simplex->randomNumberGenerator()->randomDouble();
     2522            int bias = static_cast<int>(random*(numberIterations/4));
     2523            simplex->setMaximumIterations(numberIterations/2+bias);
     2524            simplex->primal();
     2525            simplex->setMaximumIterations(COIN_INT_MAX);
     2526            simplex->dual();
     2527          } else {
     2528            simplex->primal();
     2529            numberIterations=simplex->numberIterations();
     2530          }
     2531          simplex->setLogLevel(logLevel);
     2532          clpSolver->setWarmStart(NULL);
    25072533        }
    25082534#endif
     2535        for (int j=0;j<numberHeuristics_;j++)
     2536          rootModels[i]->heuristic_[j]->setSeed(rootModels[i]->heuristic_[j]->getSeed()+100000000*i);
     2537        for (int j=0;j<numberCutGenerators_;j++)
     2538          rootModels[i]->generator_[j]->generator()->refreshSolver(rootModels[i]->solver_);
    25092539      }
    25102540      delete basis;
     
    26632693    }
    26642694    // replace solverType
     2695    double * tightBounds = NULL;
    26652696    if (solverCharacteristics_->tryCuts())  {
    26662697
     
    26692700            if (!eventHappened_ && feasible) {
    26702701                if (rootModels) {
     2702                  // for fixings
     2703                  int numberColumns=solver_->getNumCols();
     2704                  tightBounds = new double [2*numberColumns];
     2705                  {
     2706                    const double * lower = solver_->getColLower();
     2707                    const double * upper = solver_->getColUpper();
     2708                    for (int i=0;i<numberColumns;i++) {
     2709                      tightBounds[2*i+0]=lower[i];
     2710                      tightBounds[2*i+1]=upper[i];
     2711                    }
     2712                  }
    26712713                  int numberModels = multipleRootTries_%100;
    26722714                  const OsiSolverInterface ** solvers = new
     
    26762718                  for (int i=0;i<numberModels;i++) {
    26772719                    solvers[i]=rootModels[i]->solver();
     2720                    const double * lower = solvers[i]->getColLower();
     2721                    const double * upper = solvers[i]->getColUpper();
     2722                    for (int j=0;j<numberColumns;j++) {
     2723                      tightBounds[2*j+0]=CoinMax(lower[j],tightBounds[2*j+0]);
     2724                      tightBounds[2*j+1]=CoinMin(upper[j],tightBounds[2*j+1]);
     2725                    }
    26782726                    int numberRows2=solvers[i]->getNumRows();
    26792727                    assert (numberRows2>=numberRows);
     
    26872735                    generator_[j]->scaleBackStatistics(numberModels);
    26882736                  }
    2689                   CbcRowCuts rowCut(maxCuts);
     2737                  //CbcRowCuts rowCut(maxCuts);
     2738                  const OsiRowCutDebugger *debugger = NULL;
     2739                  if ((specialOptions_&1) != 0)
     2740                    debugger = solver_->getRowCutDebugger() ;
    26902741                  for (int iModel=0;iModel<numberModels;iModel++) {
    26912742                    int numberRows2=solvers[iModel]->getNumRows();
     
    27032754                      CoinBigIndex start = rowStart[iRow];
    27042755                      rc.setRow(rowLength[iRow],column+start,elements+start,false);
    2705                       rowCut.addCutIfNotDuplicate(rc);
     2756                      if (debugger)
     2757                        CoinAssert (!debugger->invalidCut(rc));
     2758                      globalCuts_.addCutIfNotDuplicate(rc);
    27062759                    }
    2707                     //int cutsAdded=rowCut.numberCuts()-numberCuts;
     2760                    //int cutsAdded=globalCuts_.numberCuts()-numberCuts;
    27082761                    //numberCuts += cutsAdded;
    27092762                    //printf("Model %d gave %d cuts (out of %d possible)\n",
    27102763                    //     iModel,cutsAdded,numberRows2-numberRows);
    27112764                  }
    2712                   rowCut.addCuts(cuts);
     2765                  // normally replace global cuts
     2766                  //if (!globalCuts_.())
     2767                  //globalCuts_=rowCutrowCut.addCuts(globalCuts_);
     2768                  //rowCut.addCuts(globalCuts_);
     2769                  int nTightened=0;
     2770                  bool feasible=true;
     2771                  {
     2772                    double tolerance=1.0e-5;
     2773                    const double * lower = solver_->getColLower();
     2774                    const double * upper = solver_->getColUpper();
     2775                    for (int i=0;i<numberColumns;i++) {
     2776                      if (tightBounds[2*i+0]>tightBounds[2*i+1]) {
     2777                        feasible=false;
     2778                        printf("Bad bounds on %d %g,%g was %g,%g\n",
     2779                               i,tightBounds[2*i+0],tightBounds[2*i+1],
     2780                               lower[i],upper[i]);
     2781                      }
     2782                      //int k=0;
     2783                      double oldLower=lower[i];
     2784                      double oldUpper=upper[i];
     2785                      if (tightBounds[2*i+0]>oldLower+tolerance) {
     2786                        nTightened++;
     2787                        //k++;
     2788                        solver_->setColLower(i,tightBounds[2*i+0]);
     2789                      }
     2790                      if (tightBounds[2*i+1]<oldUpper-tolerance) {
     2791                        nTightened++;
     2792                        //k++;
     2793                        solver_->setColUpper(i,tightBounds[2*i+1]);
     2794                      }
     2795                      //if (k)
     2796                      //printf("new bounds on %d %g,%g was %g,%g\n",
     2797                      //       i,tightBounds[2*i+0],tightBounds[2*i+1],
     2798                      //       oldLower,oldUpper);
     2799                    }
     2800                    if (!feasible)
     2801                      abort(); // deal with later
     2802                  }
     2803                  delete [] tightBounds;
     2804                  tightBounds=NULL;
    27132805                  char printBuffer[200];
    27142806                  sprintf(printBuffer,"%d solvers added %d different cuts out of pool of %d",
    2715                           numberModels,cuts.sizeRowCuts(),maxCuts);
     2807                          numberModels,globalCuts_.sizeRowCuts(),maxCuts);
    27162808                  messageHandler()->message(CBC_GENERAL, messages())
    27172809                    << printBuffer << CoinMessageEol ;
     2810                  if (nTightened) {
     2811                    sprintf(printBuffer,"%d bounds were tightened",
     2812                          nTightened);
     2813                    messageHandler()->message(CBC_GENERAL, messages())
     2814                      << printBuffer << CoinMessageEol ;
     2815                  }
    27182816                  delete [] solvers;
     2817                }
     2818                if (!parentModel_&&(moreSpecialOptions_&67108864) != 0) {
     2819                  // load cuts from file
     2820                  FILE * fp = fopen("global.cuts","rb");
     2821                  if (fp) {
     2822                    size_t nRead;
     2823                    int numberColumns=solver_->getNumCols();
     2824                    int numCols;
     2825                    nRead = fread(&numCols, sizeof(int), 1, fp);
     2826                    if (nRead != 1)
     2827                      throw("Error in fread");
     2828                    if (numberColumns!=numCols) {
     2829                      printf("Mismatch on columns %d %d\n",numberColumns,numCols);
     2830                      fclose(fp);
     2831                    } else {
     2832                      // If rootModel just do some
     2833                      double threshold=-1.0;
     2834                      if (!multipleRootTries_)
     2835                        threshold=0.5;
     2836                      int initialCuts=0;
     2837                      int initialGlobal = globalCuts_.sizeRowCuts();
     2838                      double * elements = new double [numberColumns+2];
     2839                      int * indices = new int [numberColumns];
     2840                      int numberEntries=1;
     2841                      while (numberEntries>0) {
     2842                        nRead = fread(&numberEntries, sizeof(int), 1, fp);
     2843                        if (nRead != 1)
     2844                          throw("Error in fread");
     2845                        double randomNumber=randomNumberGenerator_.randomDouble();
     2846                        if (numberEntries>0) {
     2847                          initialCuts++;
     2848                          nRead = fread(elements, sizeof(double), numberEntries+2, fp);
     2849                          if (nRead != static_cast<size_t>(numberEntries+2))
     2850                            throw("Error in fread");
     2851                          nRead = fread(indices, sizeof(int), numberEntries, fp);
     2852                          if (nRead != static_cast<size_t>(numberEntries))
     2853                            throw("Error in fread");
     2854                          if (randomNumber>threshold) {
     2855                            OsiRowCut rc;
     2856                            rc.setLb(elements[numberEntries]);
     2857                            rc.setUb(elements[numberEntries+1]);
     2858                            rc.setRow(numberEntries,indices,elements,
     2859                                      false);
     2860                            rc.setGloballyValidAsInteger(2);
     2861                            globalCuts_.addCutIfNotDuplicate(rc) ;
     2862                          }
     2863                        }
     2864                      }
     2865                      fclose(fp);
     2866                      // fixes
     2867                      int nTightened=0;
     2868                      fp = fopen("global.fix","rb");
     2869                      if (fp) {
     2870                        nRead = fread(indices, sizeof(int), 2, fp);
     2871                        if (nRead != 2)
     2872                          throw("Error in fread");
     2873                        if (numberColumns!=indices[0]) {
     2874                          printf("Mismatch on columns %d %d\n",numberColumns,
     2875                                 indices[0]);
     2876                        } else {
     2877                          indices[0]=1;
     2878                          while (indices[0]>=0) {
     2879                            nRead = fread(indices, sizeof(int), 2, fp);
     2880                            if (nRead != 2)
     2881                              throw("Error in fread");
     2882                            int iColumn=indices[0];
     2883                            if (iColumn>=0) {
     2884                              nTightened++;
     2885                              nRead = fread(elements, sizeof(double), 4, fp);
     2886                              if (nRead != 4)
     2887                                throw("Error in fread");
     2888                              solver_->setColLower(iColumn,elements[0]);
     2889                              solver_->setColUpper(iColumn,elements[1]);
     2890                            }
     2891                          }
     2892                        }
     2893                      }
     2894                      if (fp)
     2895                        fclose(fp);
     2896                      char printBuffer[200];
     2897                      sprintf(printBuffer,"%d cuts read in of which %d were unique, %d bounds tightened",
     2898                             initialCuts,
     2899                             globalCuts_.sizeRowCuts()-initialGlobal,nTightened);
     2900                      messageHandler()->message(CBC_GENERAL, messages())
     2901                        << printBuffer << CoinMessageEol ;
     2902                      delete [] elements;
     2903                      delete [] indices;
     2904                    }
     2905                  }
    27192906                }
    27202907                feasible = solveWithCuts(cuts, maximumCutPassesAtRoot_,
    27212908                                         NULL);
     2909                if (multipleRootTries_&&
     2910                    (moreSpecialOptions_&134217728)!=0) {
     2911                  FILE * fp=NULL;
     2912                  size_t nRead;
     2913                  int numberColumns=solver_->getNumCols();
     2914                  int initialCuts=0;
     2915                  if ((moreSpecialOptions_&134217728)!=0) {
     2916                    // append so go down to end
     2917                    fp = fopen("global.cuts","r+b");
     2918                    if (fp) {
     2919                      int numCols;
     2920                      nRead = fread(&numCols, sizeof(int), 1, fp);
     2921                      if (nRead != 1)
     2922                        throw("Error in fread");
     2923                      if (numberColumns!=numCols) {
     2924                        printf("Mismatch on columns %d %d\n",numberColumns,numCols);
     2925                        fclose(fp);
     2926                        fp=NULL;
     2927                      }
     2928                    }
     2929                  }
     2930                  double * elements = new double [numberColumns+2];
     2931                  int * indices = new int [numberColumns];
     2932                  if (fp) {
     2933                    int numberEntries=1;
     2934                    while (numberEntries>0) {
     2935                      fpos_t position;
     2936                      fgetpos(fp, &position);
     2937                      nRead = fread(&numberEntries, sizeof(int), 1, fp);
     2938                      if (nRead != 1)
     2939                        throw("Error in fread");
     2940                      if (numberEntries>0) {
     2941                        initialCuts++;
     2942                        nRead = fread(elements, sizeof(double), numberEntries+2, fp);
     2943                        if (nRead != static_cast<size_t>(numberEntries+2))
     2944                          throw("Error in fread");
     2945                        nRead = fread(indices, sizeof(int), numberEntries, fp);
     2946                        if (nRead != static_cast<size_t>(numberEntries))
     2947                          throw("Error in fread");
     2948                      } else {
     2949                        // end
     2950                        fsetpos(fp, &position);
     2951                      }
     2952                    }
     2953                  } else {
     2954                    fp = fopen("global.cuts","wb");
     2955                    size_t nWrite;
     2956                    nWrite=fwrite(&numberColumns,sizeof(int),1,fp);
     2957                    if (nWrite != 1)
     2958                      throw("Error in fwrite");
     2959                  }
     2960                  size_t nWrite;
     2961                  // now append binding cuts
     2962                  int numberC=continuousSolver_->getNumRows();
     2963                  int numberRows=solver_->getNumRows();
     2964                  printf("Saving %d cuts (up from %d)\n",
     2965                         initialCuts+numberRows-numberC,initialCuts);
     2966                  const double * rowLower = solver_->getRowLower();
     2967                  const double * rowUpper = solver_->getRowUpper();
     2968                  // Row copy
     2969                  CoinPackedMatrix matrixByRow(*solver_->getMatrixByRow());
     2970                  const double * elementByRow = matrixByRow.getElements();
     2971                  const int * column = matrixByRow.getIndices();
     2972                  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
     2973                  const int * rowLength = matrixByRow.getVectorLengths();
     2974                  for (int iRow=numberC;iRow<numberRows;iRow++) {
     2975                    int n=rowLength[iRow];
     2976                    assert (n);
     2977                    CoinBigIndex start=rowStart[iRow];
     2978                    memcpy(elements,elementByRow+start,n*sizeof(double));
     2979                    memcpy(indices,column+start,n*sizeof(int));
     2980                    elements[n]=rowLower[iRow];
     2981                    elements[n+1]=rowUpper[iRow];
     2982                    nWrite=fwrite(&n,sizeof(int),1,fp);
     2983                    if (nWrite != 1)
     2984                      throw("Error in fwrite");
     2985                    nWrite=fwrite(elements,sizeof(double),n+2,fp);
     2986                    if (nWrite != static_cast<size_t>(n+2))
     2987                      throw("Error in fwrite");
     2988                    nWrite=fwrite(indices,sizeof(int),n,fp);
     2989                    if (nWrite != static_cast<size_t>(n))
     2990                      throw("Error in fwrite");
     2991                  }
     2992                  // eof marker
     2993                  int eofMarker=-1;
     2994                  nWrite=fwrite(&eofMarker,sizeof(int),1,fp);
     2995                  if (nWrite != 1)
     2996                    throw("Error in fwrite");
     2997                  fclose(fp);
     2998                  // do tighter bounds (? later extra to original columns)
     2999                  int nTightened=0;
     3000                  const double * lower = solver_->getColLower();
     3001                  const double * upper = solver_->getColUpper();
     3002                  const double * originalLower = continuousSolver_->getColLower();
     3003                  const double * originalUpper = continuousSolver_->getColUpper();
     3004                  double tolerance=1.0e-5;
     3005                  for (int i=0;i<numberColumns;i++) {
     3006                    if (lower[i]>originalLower[i]+tolerance) {
     3007                      nTightened++;
     3008                    }
     3009                    if (upper[i]<originalUpper[i]-tolerance) {
     3010                      nTightened++;
     3011                    }
     3012                  }
     3013                  if (nTightened) {
     3014                    fp = fopen("global.fix","wb");
     3015                    size_t nWrite;
     3016                    indices[0]=numberColumns;
     3017                    if (originalColumns_)
     3018                      indices[1]=COIN_INT_MAX;
     3019                    else
     3020                      indices[1]=-1;
     3021                    nWrite=fwrite(indices,sizeof(int),2,fp);
     3022                    if (nWrite != 2)
     3023                      throw("Error in fwrite");
     3024                    for (int i=0;i<numberColumns;i++) {
     3025                      int nTightened=0;
     3026                      if (lower[i]>originalLower[i]+tolerance) {
     3027                        nTightened++;
     3028                      }
     3029                      if (upper[i]<originalUpper[i]-tolerance) {
     3030                        nTightened++;
     3031                      }
     3032                      if (nTightened) {
     3033                        indices[0]=i;
     3034                        if (originalColumns_)
     3035                          indices[1]=originalColumns_[i];
     3036                        elements[0]=lower[i];
     3037                        elements[1]=upper[i];
     3038                        elements[2]=originalLower[i];
     3039                        elements[3]=originalUpper[i];
     3040                        nWrite=fwrite(indices,sizeof(int),2,fp);
     3041                        if (nWrite != 2)
     3042                          throw("Error in fwrite");
     3043                        nWrite=fwrite(elements,sizeof(double),4,fp);
     3044                        if (nWrite != 4)
     3045                          throw("Error in fwrite");
     3046                      }
     3047                    }
     3048                    // eof marker
     3049                    indices[0]=-1;
     3050                    nWrite=fwrite(indices,sizeof(int),2,fp);
     3051                    if (nWrite != 2)
     3052                      throw("Error in fwrite");
     3053                    fclose(fp);
     3054                  }
     3055                  delete [] elements;
     3056                  delete [] indices;
     3057                }
    27223058                if ((specialOptions_&524288) != 0 && !parentModel_
    27233059                        && storedRowCuts_) {
     
    29763312    }
    29773313#endif
     3314    if (!parentModel_&&(moreSpecialOptions_&268435456) != 0) {
     3315      // try idiotic idea
     3316      CbcObject * obj = new CbcIdiotBranch(this);
     3317      obj->setPriority(1); // temp
     3318      addObjects(1, &obj);
     3319      delete obj;
     3320    }
     3321   
    29783322    /*
    29793323      A hook to use clp to quickly explore some part of the tree.
     
    29853329        obj->setPriority(1);
    29863330        addObjects(1, &obj);
     3331        delete obj;
    29873332    }
    29883333    int saveNumberSolves = numberSolves_;
    29893334    int saveNumberIterations = numberIterations_;
    2990     if (fastNodeDepth_ >= 0 &&/*!parentModel_*/(specialOptions_&2048) == 0) {
     3335    if ((fastNodeDepth_ >= 0||(moreSpecialOptions_&33554432)!=0)
     3336        &&/*!parentModel_*/(specialOptions_&2048) == 0) {
    29913337        // add in a general depth object doClp
    29923338        int type = (fastNodeDepth_ <= 100) ? fastNodeDepth_ : -(fastNodeDepth_ - 100);
     3339        if ((moreSpecialOptions_&33554432)!=0)
     3340          type=12;
     3341        else
     3342          fastNodeDepth_ += 1000000;     // mark as done
    29933343        CbcObject * obj =
    29943344            new CbcGeneralDepth(this, type);
    29953345        addObjects(1, &obj);
    2996         // mark as done
    2997         fastNodeDepth_ += 1000000;
    29983346        delete obj;
    29993347        // fake number of objects
     
    40524400            if (!node || node->objectiveValue() > cutoff)
    40534401                continue;
    4054             // Do main work of solving node here
    4055             doOneNode(this, node, createdNode);
     4402            // Do main work of solving node here
     4403            doOneNode(this, node, createdNode);
    40564404#ifdef JJF_ZERO
    4057             if (node) {
    4058                 if (createdNode) {
    4059                     printf("Node %d depth %d, created %d depth %d\n",
    4060                            node->nodeNumber(), node->depth(),
    4061                            createdNode->nodeNumber(), createdNode->depth());
    4062                 } else {
    4063                     printf("Node %d depth %d,  no created node\n",
    4064                            node->nodeNumber(), node->depth());
    4065                 }
    4066             } else if (createdNode) {
    4067                 printf("Node exhausted, created %d depth %d\n",
    4068                        createdNode->nodeNumber(), createdNode->depth());
    4069             } else {
    4070                 printf("Node exhausted,  no created node\n");
    4071                 numberConsecutiveInfeasible = 2;
    4072             }
     4405            if (node) {
     4406              if (createdNode) {
     4407                printf("Node %d depth %d, created %d depth %d\n",
     4408                       node->nodeNumber(), node->depth(),
     4409                       createdNode->nodeNumber(), createdNode->depth());
     4410              } else {
     4411                printf("Node %d depth %d,  no created node\n",
     4412                       node->nodeNumber(), node->depth());
     4413              }
     4414            } else if (createdNode) {
     4415              printf("Node exhausted, created %d depth %d\n",
     4416                     createdNode->nodeNumber(), createdNode->depth());
     4417            } else {
     4418              printf("Node exhausted,  no created node\n");
     4419              numberConsecutiveInfeasible = 2;
     4420            }
    40734421#endif
    40744422            //if (createdNode)
     
    44534801        setCutoff(1.0e50) ; // As best solution should be worse than cutoff
    44544802        // change cutoff as constraint if wanted
    4455         if ((moreSpecialOptions_&16777216)!=0) {
    4456           if (solver_->getNumRows()>=numberRowsAtContinuous_)
    4457             solver_->setRowUpper(numberRowsAtContinuous_-1,1.0e50);
     4803        if (cutoffRowNumber_>=0) {
     4804          if (solver_->getNumRows()>cutoffRowNumber_)
     4805            solver_->setRowUpper(cutoffRowNumber_,1.0e50);
    44584806        }
    44594807        // also in continuousSolver_
     
    45184866    delete globalConflictCuts_;
    45194867    globalConflictCuts_=NULL;
    4520     if (!bestSolution_) {
     4868    if (!bestSolution_ && (specialOptions_&8388608)==0) {
    45214869        // make sure lp solver is infeasible
    45224870        int numberColumns = solver_->getNumCols();
     
    45454893    }
    45464894#endif
    4547     if (fastNodeDepth_ >= 1000000 && !parentModel_) {
    4548         // delete object off end
    4549         delete object_[numberObjects_];
     4895    if ((fastNodeDepth_ >= 1000000 || (moreSpecialOptions_&33554432)!=0)
     4896         && !parentModel_) {
     4897      // delete object off end
     4898      delete object_[numberObjects_];
     4899      if ((moreSpecialOptions_&33554432)==0)
    45504900        fastNodeDepth_ -= 1000000;
    45514901    }
     
    47285078        numberIntegers_(0),
    47295079        numberRowsAtContinuous_(0),
     5080        cutoffRowNumber_(-1),
    47305081        maximumNumberCuts_(0),
    47315082        phase_(0),
     
    47755126        ownObjects_(true),
    47765127        originalColumns_(NULL),
    4777         howOftenGlobalScan_(1),
     5128        howOftenGlobalScan_(3),
    47785129        numberGlobalViolations_(0),
    47795130        numberExtraIterations_(0),
     
    48945245        secondaryStatus_(-1),
    48955246        numberRowsAtContinuous_(0),
     5247        cutoffRowNumber_(-1),
    48965248        maximumNumberCuts_(0),
    48975249        phase_(0),
     
    49385290        ownObjects_(true),
    49395291        originalColumns_(NULL),
    4940         howOftenGlobalScan_(1),
     5292        howOftenGlobalScan_(3),
    49415293        numberGlobalViolations_(0),
    49425294        numberExtraIterations_(0),
     
    50655417}
    50665418
     5419static int * resizeInt(int * array,int oldLength, int newLength)
     5420{
     5421  if (!array)
     5422    return NULL;
     5423  assert (newLength>oldLength);
     5424  int * newArray = new int [newLength];
     5425  memcpy(newArray,array,oldLength*sizeof(int));
     5426  delete [] array;
     5427  memset(newArray+oldLength,0,(newLength-oldLength)*sizeof(int));
     5428  return newArray;
     5429}
     5430static double * resizeDouble(double * array,int oldLength, int newLength)
     5431{
     5432  if (!array)
     5433    return NULL;
     5434  assert (newLength>oldLength);
     5435  double * newArray = new double [newLength];
     5436  memcpy(newArray,array,oldLength*sizeof(double));
     5437  delete [] array;
     5438  memset(newArray+oldLength,0,(newLength-oldLength)*sizeof(double));
     5439  return newArray;
     5440}
    50675441/*
    50685442  Assign a solver to the model (model assumes ownership)
     
    50885462
    50895463{
    5090     // resize best solution if exists
    5091     if (bestSolution_ && solver && solver_) {
     5464    // resize stuff if exists
     5465    if (solver && solver_) {
    50925466        int nOld = solver_->getNumCols();
    50935467        int nNew = solver->getNumCols();
    50945468        if (nNew > nOld) {
    5095             double * temp = new double[nNew];
    5096             memcpy(temp, bestSolution_, nOld*sizeof(double));
    5097             memset(temp + nOld, 0, (nNew - nOld)*sizeof(double));
    5098             delete [] bestSolution_;
    5099             bestSolution_ = temp;
     5469          originalColumns_ = resizeInt(originalColumns_,nOld,nNew);
     5470          usedInSolution_ = resizeInt(usedInSolution_,nOld,nNew);
     5471          continuousSolution_ = resizeDouble(continuousSolution_,nOld,nNew);
     5472          hotstartSolution_ = resizeDouble(hotstartSolution_,nOld,nNew);
     5473          bestSolution_ = resizeDouble(bestSolution_,nOld,nNew);
     5474          currentSolution_ = resizeDouble(currentSolution_,nOld,nNew);
     5475          if (savedSolutions_) {
     5476            for (int i = 0; i < maximumSavedSolutions_; i++)
     5477              savedSolutions_[i] = resizeDouble(savedSolutions_[i],nOld,nNew);
     5478          }
    51005479        }
    51015480    }
     
    53885767    testSolution_ = currentSolution_;
    53895768    numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
     5769    cutoffRowNumber_ = rhs.cutoffRowNumber_;
    53905770    maximumNumberCuts_ = rhs.maximumNumberCuts_;
    53915771    phase_ = rhs.phase_;
     
    57156095        }
    57166096        numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
     6097        cutoffRowNumber_ = rhs.cutoffRowNumber_;
    57176098        maximumNumberCuts_ = rhs.maximumNumberCuts_;
    57186099        phase_ = rhs.phase_;
     
    60076388        int i;
    60086389        for (i = 0; i < numberCutGenerators_; i++) {
    6009             if (mode < 2)
     6390            if (mode < 2) {
    60106391                generator_[i] = new CbcCutGenerator(*rhs.generator_[i]);
    6011             else
     6392            } else {
    60126393                generator_[i] = new CbcCutGenerator(*rhs.virginGenerator_[i]);
     6394                // But copy across maximumTries
     6395                generator_[i]->setMaximumTries(rhs.generator_[i]->maximumTries());
     6396            }
    60136397            virginGenerator_[i] = new CbcCutGenerator(*rhs.virginGenerator_[i]);
    60146398        }
     
    60726456    strongInfo_[6] = rhs.strongInfo_[6];
    60736457    numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
     6458    cutoffRowNumber_ = rhs.cutoffRowNumber_;
    60746459    maximumDepth_ = rhs.maximumDepth_;
    60756460}
     
    62616646        heuristic_[where]->setHeuristicName(name) ;
    62626647    heuristic_[where]->setSeed(987654321 + where);
    6263     if (randomSeed_!=-1)
    6264       heuristic_[where]->setSeed(randomSeed_);
    62656648    numberHeuristics_++ ;
    62666649}
     
    62886671{
    62896672    int nNode = 0;
     6673    CbcNodeInfo * nodeInfo = node->nodeInfo();
    62906674    int numberColumns = getNumCols();
    6291     CbcNodeInfo * nodeInfo = node->nodeInfo();
    62926675
    62936676    /*
     
    65456928#         endif
    65466929                    addCuts[numberToAdd++] = addedCuts_[i];
     6930#if 1
     6931                    if ((specialOptions_&1) != 0) {
     6932                      const OsiRowCutDebugger * debugger = solver_->getRowCutDebugger() ;
     6933                      if (debugger)
     6934                        CoinAssert (!debugger->invalidCut(*addedCuts_[i]));
     6935                    }
     6936#endif
    65476937                } else {
    65486938#         ifdef CHECK_CUT_COUNTS
     
    66527042                    delete [] which;
    66537043                }
     7044                //#define CHECK_DEBUGGER
     7045#ifdef CHECK_DEBUGGER
     7046                if ((specialOptions_&1) != 0 ) {
     7047                  const OsiRowCutDebugger * debugger =
     7048                    solver_->getRowCutDebugger();
     7049                  if (debugger) {
     7050                    for (int j=0;j<numberToAdd;j++)
     7051                      CoinAssert (!debugger->invalidCut(*addCuts[j]));
     7052                    //addCuts[j]->print();
     7053                  }
     7054                }
     7055#endif
    66547056                //int n2=solver_->getNumRows();
    66557057                //for (int j=0;j<numberToAdd;j++)
     
    67097111CbcModel::synchronizeHandlers(int /*makeDefault*/)
    67107112{
     7113    bool defaultHandler = defaultHandler_;
    67117114    if (!defaultHandler_) {
    67127115        // Must have clone
     
    67157118    }
    67167119#ifdef COIN_HAS_CLP
    6717     OsiClpSolverInterface * solver;
    6718     solver = dynamic_cast<OsiClpSolverInterface *>(solver_) ;
    6719     if (solver) {
     7120    if (!defaultHandler) {
     7121      OsiClpSolverInterface * solver;
     7122      solver = dynamic_cast<OsiClpSolverInterface *>(solver_) ;
     7123      if (solver) {
    67207124        solver->passInMessageHandler(handler_);
    67217125        solver->getModelPtr()->passInMessageHandler(handler_);
    6722     }
    6723     solver = dynamic_cast<OsiClpSolverInterface *>(continuousSolver_) ;
    6724     if (solver) {
     7126      }
     7127      solver = dynamic_cast<OsiClpSolverInterface *>(continuousSolver_) ;
     7128      if (solver) {
    67257129        solver->passInMessageHandler(handler_);
    67267130        solver->getModelPtr()->passInMessageHandler(handler_);
     7131      }
    67277132    }
    67287133#endif
     
    69877392          if ((specialOptions_&1) != 0) {
    69887393            debugger = continuousSolver_->getRowCutDebugger() ;
    6989             if (debugger)
     7394            if (debugger) {
    69907395              if (debugger->invalidCut(*cut)) {
    69917396                continuousSolver_->applyRowCuts(1,cut);
     
    69937398              }
    69947399              CoinAssert (!debugger->invalidCut(*cut));
     7400            }
    69957401          }
    69967402        } else {
     
    72747680        allowZeroIterations = true;
    72757681    }
     7682    int saveNumberTries=numberTries;
    72767683    /*
    72777684      Is it time to scan the cuts in order to remove redundant cuts? If so, set
     
    72927699    // Really primalIntegerTolerance; relates to an illposed problem with various
    72937700    // integer solutions depending on integer tolerance.
    7294     double primalTolerance = 1.0e-7 ;
     7701    //double primalTolerance = 1.0e-7 ;
    72957702    // We may need to keep going on
    72967703    bool keepGoing = false;
     
    73427749        */
    73437750        int numberViolated = 0;
    7344         if (currentPassNumber_ == 1 && howOftenGlobalScan_ > 0 &&
     7751        if ((currentPassNumber_ == 1 ||!numberNodes_) && howOftenGlobalScan_ > 0 &&
    73457752                (numberNodes_ % howOftenGlobalScan_) == 0 &&
    73467753                (doCutsNow(1) || true)) {
    73477754            // global column cuts now done in node at top of tree
    7348             int numberCuts = globalCuts_.sizeRowCuts() ;
    7349             // possibly extend whichGenerator
    7350             resizeWhichGenerator(numberViolated, numberViolated + numberCuts);
    7351             for (int i = 0; i < numberCuts; i++) {
     7755            int numberCuts = numberCutGenerators_ ? globalCuts_.sizeRowCuts() : 0;
     7756            if (numberCuts) {
     7757              // possibly extend whichGenerator
     7758              resizeWhichGenerator(numberViolated, numberViolated + numberCuts);
     7759              // only add new cuts up to 10% of current elements
     7760              int numberElements = solver_->getNumElements();
     7761              int numberColumns = solver_->getNumCols();
     7762              int maximumAdd = CoinMax(numberElements/10,2*numberColumns)+100;
     7763              double * violations = new double[numberCuts];
     7764              int * which = new int[numberCuts];
     7765              int numberPossible=0;
     7766              for (int i = 0; i < numberCuts; i++) {
    73527767                OsiRowCut * thisCut = globalCuts_.rowCutPtr(i) ;
    7353                 if (thisCut->violated(cbcColSolution_) > primalTolerance ||
    7354                         thisCut->effectiveness() == COIN_DBL_MAX) {
     7768                double violation = thisCut->violated(cbcColSolution_);
     7769                if(thisCut->effectiveness() == COIN_DBL_MAX) {
     7770                  // see if already there
     7771                  int j;
     7772                  for (j = 0; j < currentNumberCuts_; j++) {
     7773                    if (addedCuts_[j]==thisCut)
     7774                      break;
     7775                  }
     7776                  if (j==currentNumberCuts_)
     7777                    violation = COIN_DBL_MAX;
     7778                  //else
     7779                  //printf("already done??\n");
     7780                }
     7781                if (violation > 0.005) {
     7782                  violations[numberPossible]=-violation;
     7783                  which[numberPossible++]=i;
     7784                }
     7785              }
     7786              CoinSort_2(violations,violations+numberPossible,which);
     7787              for (int i = 0; i < numberPossible; i++) {
     7788                int k=which[i];
     7789                OsiRowCut * thisCut = globalCuts_.rowCutPtr(k) ;
     7790                assert (thisCut->violated(cbcColSolution_) > 0.005/*primalTolerance*/ ||
     7791                        thisCut->effectiveness() == COIN_DBL_MAX);
     7792#define CHECK_DEBUGGER
     7793#ifdef CHECK_DEBUGGER
     7794                if ((specialOptions_&1) != 0 ) {
     7795                  const OsiRowCutDebugger * debugger =
     7796                    solver_->getRowCutDebuggerAlways();
     7797                  CoinAssert (!debugger->invalidCut(*thisCut));
     7798                }
     7799#endif
    73557800#if 0 //ndef NDEBUG
    7356                   printf("Global cut added - violation %g\n",
    7357                            thisCut->violated(cbcColSolution_)) ;
    7358 #endif
    7359                     whichGenerator_[numberViolated++] = -1;
     7801                printf("Global cut added - violation %g\n",
     7802                       thisCut->violated(cbcColSolution_)) ;
     7803#endif
     7804                whichGenerator_[numberViolated++] = -1;
    73607805#ifndef GLOBAL_CUTS_JUST_POINTERS
    7361                     theseCuts.insert(*thisCut) ;
     7806                theseCuts.insert(*thisCut) ;
    73627807#else
    7363                     theseCuts.insert(thisCut) ;
    7364 #endif
    7365                 }
    7366             }
    7367             numberGlobalViolations_ += numberViolated;
     7808                theseCuts.insert(thisCut) ;
     7809#endif
     7810                if (violations[i]!=-COIN_DBL_MAX)
     7811                  maximumAdd -= thisCut->row().getNumElements();
     7812                if (maximumAdd<0)
     7813                  break;
     7814              }
     7815              delete [] which;
     7816              delete [] violations;
     7817              numberGlobalViolations_ += numberViolated;
     7818            }
    73687819        }
    73697820        /*
     
    77218172            //solver_->setHintParam(OsiDoDualInResolve,true,OsiHintTry);
    77228173            if ( maximumSecondsReached() ) {
    7723                 numberTries = 0; // exit
     8174                numberTries = -1000; // exit
    77248175                feasible = false;
    77258176                break;
     
    78798330                    // maybe give it one more try
    78808331                    if (numberLastAttempts > 2 || currentDepth_ || experimentBreak < 2)
    7881                         break ;
     8332                        numberTries=0 ;
    78828333                    else
    78838334                        numberLastAttempts++;
     
    79138364            keepGoing=false;
    79148365        }
    7915         if (numberTries <=0 && feasible && !keepGoing && !parentModel_ && !numberNodes_) {
     8366        if (numberTries ==0 && feasible && !keepGoing && !parentModel_ && !numberNodes_) {
    79168367          for (int i = 0; i < numberCutGenerators_; i++) {
    79178368            if (generator_[i]->whetherCallAtEnd()
    79188369                &&!generator_[i]->whetherInMustCallAgainMode()) {
    7919               keepGoing= true;
    7920               break;
     8370              // give it some goes and switch off
     8371              numberTries=(saveNumberTries+4)/5;
     8372              generator_[i]->setWhetherCallAtEnd(false);
    79218373            }
    79228374          }
     
    79888440                }
    79898441                if (!feasible) { //node will be fathomed
     8442                    lockThread();
    79908443                    for (int i = 0; i < currentNumberCuts_; i++) {
    79918444                        // take off node
     
    79968449                        }
    79978450                    }
     8451                    unlockThread();
    79988452                }
    79998453            }
     
    81398593        */
    81408594        if (!numberNodes_) {
     8595          if (!parentModel_) {
     8596            //printf("%d global cuts\n",globalCuts_.sizeRowCuts()) ;
     8597            if ((specialOptions_&1) != 0) {
     8598              //specialOptions_ &= ~1;
     8599              int numberCuts = globalCuts_.sizeRowCuts();
     8600              const OsiRowCutDebugger *debugger =
     8601                continuousSolver_->getRowCutDebugger() ;
     8602              if (debugger) {
     8603                for (int i = 0; i < numberCuts; i++) {
     8604                  OsiRowCut * cut = globalCuts_.rowCutPtr(i) ;
     8605                  if (debugger->invalidCut(*cut)) {
     8606                    continuousSolver_->applyRowCuts(1,cut);
     8607                    continuousSolver_->writeMps("bad");
     8608                    printf("BAD cut\n");
     8609                  }
     8610                  //CoinAssert (!debugger->invalidCut(*cut));
     8611                }
     8612              }
     8613            }
     8614          }
    81418615          //solver_->writeMps("second");
    81428616            if (numberRowsAdded)
     
    82878761                willBeCutsInTree = 0;
    82888762        }
    8289         if (!numberNodes_)
     8763        if (!numberNodes_) {
    82908764            handler_->message(CBC_ROOT, messages_)
    82918765            << numberNewCuts_
     
    82938767            << currentPassNumber_
    82948768            << CoinMessageEol ;
     8769        }
    82958770        /*
    82968771          Count the number of cuts produced by each cut generator on this call. Not
     
    84038878                    howOften = -99; // switch off
    84048879            }
     8880            if (generator_[i]->maximumTries()!=-1)
     8881                howOften = -99; // switch off
    84058882            /*
    84068883              Below -99, this generator is switched off. There's no need to consider
    84078884              further. Then again, there was no point in persisting this far!
    84088885            */
    8409             if (howOften < -99)
    8410                 continue ;
     8886            if (howOften < -99) {
     8887              // may have been switched off - report
     8888              if (!numberNodes_) {
     8889                int n = generator_[i]->numberCutsInTotal();
     8890                if (n) {
     8891                  double average = 0.0;
     8892                  average = generator_[i]->numberElementsInTotal();
     8893                  average /= n;
     8894                  handler_->message(CBC_GENERATOR, messages_)
     8895                    << i
     8896                    << generator_[i]->cutGeneratorName()
     8897                    << n
     8898                    << average
     8899                    << generator_[i]->numberColumnCuts()
     8900                    << generator_[i]->numberCutsActive()
     8901                    + generator_[i]->numberColumnCuts();
     8902                  handler_->printing(generator_[i]->timing())
     8903                    << generator_[i]->timeInCutGenerator();
     8904                  handler_->message()
     8905                    << -100
     8906                    << CoinMessageEol ;
     8907                }
     8908              }
     8909              continue ;
     8910            }
    84118911            /*
    84128912              Adjust, if howOften is adjustable.
     
    87309230        specialOptions_ &= ~256; // mark as full scan done
    87319231    }
    8732 # ifdef COIN_HAS_CLP
    8733     if (!node && !parentModel_ && intParam_[CbcMaxNumNode] == -123456) {
    8734         OsiClpSolverInterface * clpSolver
    8735         = dynamic_cast<OsiClpSolverInterface *> (solver_);
    8736         if (clpSolver) {
    8737             clpSolver->lexSolve();
    8738         }
     9232# if 0 //def COIN_HAS_CLP
     9233    // check basis
     9234    OsiClpSolverInterface * clpSolver
     9235      = dynamic_cast<OsiClpSolverInterface *> (solver_);
     9236    if (clpSolver) {
     9237      ClpSimplex * simplex = clpSolver->getModelPtr();
     9238      int numberTotal=simplex->numberRows()+simplex->numberColumns();
     9239      int superbasic=0;
     9240      for (int i=0;i<numberTotal;i++) {
     9241        if (simplex->getStatus(i)==ClpSimplex::superBasic)
     9242          superbasic++;
     9243      }
     9244      if (superbasic) {
     9245        printf("%d superbasic!\n",superbasic);
     9246        clpSolver->resolve();
     9247        superbasic=0;
     9248        for (int i=0;i<numberTotal;i++) {
     9249          if (simplex->getStatus(i)==ClpSimplex::superBasic)
     9250            superbasic++;
     9251        }
     9252        assert (!superbasic);
     9253      }
    87399254    }
    87409255# endif
     
    87649279            }
    87659280        }
     9281        if (generator_[i]->whetherCallAtEnd())
     9282          generate=false;
    87669283        const OsiRowCutDebugger * debugger = NULL;
    87679284        bool onOptimalPath = false;
     
    88229339            }
    88239340#endif
    8824             if (mustResolve) {
     9341            if (mustResolve || (specialOptions_&1) != 0) {
    88259342                int returnCode = resolve(node ? node->nodeInfo() : NULL, 2);
    88269343                if (returnCode  == 0)
     
    88989415            if (thisCut->lb() > thisCut->ub())
    88999416                status = -1; // sub-problem is infeasible
    8900             if (thisCut->globallyValid()) {
     9417            if (thisCut->globallyValid()||!numberNodes_) {
    89019418                // add to global list
    89029419                OsiRowCut newCut(*thisCut);
     
    92239740        }
    92249741    }
    9225 #ifdef COIN_HAS_CLP
     9742  #ifdef COIN_HAS_CLP
    92269743    OsiClpSolverInterface * clpSolver
    92279744    = dynamic_cast<OsiClpSolverInterface *> (solver_);
     
    950010017#ifdef COIN_HAS_CLP
    950110018    if (clpSolver && !feasible) {
    9502         // make sure marked infeasible
     10019      // make sure marked infeasible
     10020      if (!clpSolver->isProvenDualInfeasible())
    950310021        clpSolver->getModelPtr()->setProblemStatus(1);
    950410022    }
     
    1051411032                int n = CoinMax(obj2->numberTimesDown(),
    1051511033                                obj2->numberTimesUp());
    10516                 if (n >= value)
    10517                     obj2->setNumberBeforeTrust(n + 1);
     11034                if (n >= value) {
     11035                  value = CoinMin(CoinMin(n+1,3*(value+1)/2),5*numberBeforeTrust_);
     11036                  obj2->setNumberBeforeTrust(value);
     11037                }
    1051811038            }
    1051911039        }
     
    1145211972            setCutoff(cutoff);
    1145311973            // change cutoff as constraint if wanted
    11454             if ((moreSpecialOptions_&16777216)!=0) {
    11455               if (solver_->getNumRows()>=numberRowsAtContinuous_) {
     11974            if (cutoffRowNumber_>=0) {
     11975              if (solver_->getNumRows()>cutoffRowNumber_) {
    1145611976                double offset;
    1145711977                solver_->getDblParam(OsiObjOffset, offset);
    11458                 solver_->setRowUpper(numberRowsAtContinuous_-1,cutoff+offset);
     11978                solver_->setRowUpper(cutoffRowNumber_,cutoff+offset);
    1145911979              }
    1146011980            }
     
    1163212152                setCutoff(cutoff);
    1163312153                // change cutoff as constraint if wanted
    11634                 if ((moreSpecialOptions_&16777216)!=0) {
    11635                   if (solver_->getNumRows()>=numberRowsAtContinuous_) {
     12154                if (cutoffRowNumber_>=0) {
     12155                  if (solver_->getNumRows()>cutoffRowNumber_) {
    1163612156                    double offset;
    1163712157                    solver_->getDblParam(OsiObjOffset, offset);
    11638                     solver_->setRowUpper(numberRowsAtContinuous_-1,cutoff+offset);
     12158                    solver_->setRowUpper(cutoffRowNumber_,cutoff+offset);
    1163912159                  }
    1164012160                }
     
    1235712877// Make partial cut into a global cut and save
    1235812878void
    12359 CbcModel::makePartialCut(const OsiRowCut * partialCut)
     12879CbcModel::makePartialCut(const OsiRowCut * partialCut,
     12880                         const OsiSolverInterface * solver)
    1236012881{
    1236112882  // get greedy cut
    1236212883  double bSum = partialCut->lb();
    1236312884  assert (bSum<0.0);
     12885  if (!solver)
     12886    solver=solver_;
    1236412887  int nConflict = partialCut->row().getNumElements();
    1236512888  const int * column = partialCut->row().getIndices();
    1236612889  const double * element = partialCut->row().getElements();
    1236712890  double * originalLower = topOfTree_->mutableLower();
    12368   const double * columnLower = solver_->getColLower();
     12891  const double * columnLower = solver->getColLower();
    1236912892  double * originalUpper = topOfTree_->mutableUpper();
    12370   const double * columnUpper = solver_->getColUpper();
     12893  const double * columnUpper = solver->getColUpper();
    1237112894  int nC=nConflict;
    1237212895  while (nConflict) {
     
    1240412927  printf("CUTa has %d (started at %d) - final bSum %g\n",nConflict,nC,bSum);
    1240512928  if (nConflict>1) {
     12929    if ((specialOptions_&1) != 0) {
     12930      const OsiRowCutDebugger *debugger = continuousSolver_->getRowCutDebugger() ;
     12931      if (debugger) {
     12932        if (debugger->invalidCut(newCut)) {
     12933          continuousSolver_->applyRowCuts(1,&newCut);
     12934          continuousSolver_->writeMps("bad");
     12935        }
     12936        CoinAssert (!debugger->invalidCut(newCut));
     12937      }
     12938    }
    1240612939    newCut.setGloballyValidAsInteger(2);
    1240712940    newCut.mutableRow().setTestForDuplicateIndex(false);
     
    1314613679            int totalNodes = numberNodes_ + numberExtraNodes_;
    1314713680            int totalIterations = numberIterations_ + numberExtraIterations_;
     13681            bool diving=false;
     13682            if ((moreSpecialOptions_&33554432)!=0) {
     13683              testDepth=COIN_INT_MAX;
     13684              if (oldNode&&(oldNode->depth()==-2||oldNode->depth()==4))
     13685                diving=true;
     13686            }
    1314813687            if (totalNodes*40 < totalIterations || numberNodes_ < 1000) {
    1314913688                doClp = false;
     
    1315313692                //     totalNodes,totalIterations,doClp ? 'Y' : 'N');
    1315413693            }
    13155             if (oldNode && fastNodeDepth_ >= 0 && oldNode->depth() >= testDepth &&/*!parentModel_*/(specialOptions_&2048) == 0
    13156                     && doClp && !cuts.sizeRowCuts()) {
     13694            if (oldNode &&
     13695                ((fastNodeDepth_ >= 0 && oldNode->depth() >= testDepth && doClp)||diving) &&/*!parentModel_*/(specialOptions_&2048) == 0
     13696                && !cuts.sizeRowCuts()) {
    1315713697                OsiClpSolverInterface * clpSolver
    1315813698                = dynamic_cast<OsiClpSolverInterface *> (solver_);
     
    1316513705#endif
    1316613706#ifdef COIN_HAS_CLP
    13167             int save;
     13707            int save=0;
    1316813708            OsiClpSolverInterface * clpSolver
    1316913709              = dynamic_cast<OsiClpSolverInterface *> (solver_);
     
    1334513885                assert (numberProblems);
    1334613886                int nProbMinus1 = numberProblems - 1;
     13887                lockThread();
    1334713888                for (int i = 0; i < currentNumberCuts_; i++) {
    1334813889                    if (addedCuts_[i])
    1334913890                        addedCuts_[i]->increment(nProbMinus1) ;
    1335013891                }
     13892                unlockThread();
    1335113893                for (int i = 0; i < numberProblems; i++) {
    1335213894                    double objectiveValue;
     
    1356714109        setCutoff(cutoff) ;
    1356814110        // change cutoff as constraint if wanted
    13569         if ((moreSpecialOptions_&16777216)!=0) {
    13570           if (solver_->getNumRows()>=numberRowsAtContinuous_) {
     14111        if (cutoffRowNumber_>=0) {
     14112          if (solver_->getNumRows()>cutoffRowNumber_) {
    1357114113            double offset;
    1357214114            solver_->getDblParam(OsiObjOffset, offset);
    13573             solver_->setRowUpper(numberRowsAtContinuous_-1,cutoff+offset);
     14115            solver_->setRowUpper(cutoffRowNumber_,cutoff+offset);
    1357414116          }
    1357514117        }
     
    1410314645{
    1410414646    int foundSolution = 0;
     14647    int saveNumberCutGenerators=numberCutGenerators_;
     14648    if ((moreSpecialOptions_&33554432)!=0 && (specialOptions_&2048)==0) {
     14649      if (node&&(node->depth()==-2||node->depth()==4))
     14650        numberCutGenerators_=0; // so can dive and branch
     14651    }
    1410514652    int currentNumberCuts = 0 ;
    1410614653    currentNode_ = node; // so can be accessed elsewhere
     
    1420014747            CbcBranchingObject * branch = static_cast <CbcBranchingObject *>(branch2) ;
    1420114748#endif
     14749#if 1
    1420214750            branch->setModel(this);
    1420314751            branchesLeft = node->branch(NULL); // old way
     14752#else
     14753            branchesLeft = node->branch(solver_);
     14754#endif
    1420414755            if (parallelMode() >= 0)
    1420514756                branch->setModel(baseModel);
     
    1460615157                  solver_->writeMpsNative("infeas2.mps", NULL, NULL, 2);
    1460715158                  solver_->getRowCutDebuggerAlways()->printOptimalSolution(*solver_);
    14608                   assert (solver_->getRowCutDebugger()) ;
     15159                  const OsiRowCutDebugger * debugger = solver_->getRowCutDebugger() ;
     15160                  assert (debugger) ;
     15161                  int numberRows0=continuousSolver_->getNumRows();
     15162                  int numberRows=solver_->getNumRows();
     15163                  const CoinPackedMatrix * rowCopy = solver_->getMatrixByRow();
     15164                  const int * rowLength = rowCopy->getVectorLengths();
     15165                  const double * elements = rowCopy->getElements();
     15166                  const int * column = rowCopy->getIndices();
     15167                  const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
     15168                  const double * rowLower = solver_->getRowLower();
     15169                  const double * rowUpper = solver_->getRowUpper();
     15170                  for (int iRow=numberRows0;iRow<numberRows;iRow++) {
     15171                    OsiRowCut rc;
     15172                    rc.setLb(rowLower[iRow]);
     15173                    rc.setUb(rowUpper[iRow]);
     15174                    CoinBigIndex start = rowStart[iRow];
     15175                    rc.setRow(rowLength[iRow],column+start,elements+start,false);
     15176                    CoinAssert (!debugger->invalidCut(rc));
     15177                  }
    1460915178                  assert (feasible);
    1461015179                }
     
    1487715446            if (parallelMode() >= 0)
    1487815447                newNode = new CbcNode() ;
     15448#if 0
     15449            // Try diving
     15450            if (parallelMode() >= 0 && (specialOptions_&2048) == 0) {
     15451              // See if any diving heuristics set to do dive+save
     15452              CbcHeuristicDive * dive=NULL;
     15453              for (int i = 0; i < numberHeuristics_; i++) {
     15454                CbcHeuristicDive * possible = dynamic_cast<CbcHeuristicDive *>(heuristic_[i]);
     15455                if (possible&&possible->maxSimplexIterations()==COIN_INT_MAX) {
     15456                  // if more than one then rotate later?
     15457                  //if (possible->canHeuristicRun()) {
     15458                  if (node->depth()==0||node->depth()==5) {
     15459                    dive=possible;
     15460                    break;
     15461                  }
     15462                }
     15463              }
     15464              if (dive) {
     15465                int numberNodes;
     15466                CbcSubProblem ** nodes=NULL;
     15467                int branchState=dive->fathom(this,numberNodes,nodes);
     15468                if (branchState) {
     15469                  printf("new solution\n");
     15470                }
     15471                if (0) {
     15472                  for (int iNode=0;iNode<numberNodes;iNode++) {
     15473                    //tree_->push(nodes[iNode]) ;
     15474                  }
     15475                  assert (node->nodeInfo());
     15476                  if (node->nodeInfo()->numberBranchesLeft()) {
     15477                    tree_->push(node) ;
     15478                  } else {
     15479                    node->setActive(false);
     15480                  }
     15481                }
     15482                delete [] nodes;
     15483              }
     15484            }
     15485            // end try diving
     15486#endif
    1487915487            // Set objective value (not so obvious if NLP etc)
    1488015488            setObjectiveValue(newNode, node);
     
    1526815876        unlockThread();
    1526915877    }
     15878    numberCutGenerators_=saveNumberCutGenerators;
    1527015879    return foundSolution;
    1527115880}
     
    1552516134    for (int i = 0; i < numberHeuristics_; i++) {
    1552616135        CbcHeuristicDive * heuristic = dynamic_cast<CbcHeuristicDive *> (heuristic_[i]);
    15527         if (heuristic) {
     16136        if (heuristic && heuristic->maxSimplexIterations()!=COIN_INT_MAX) {
    1552816137            heuristic->setMaxSimplexIterations(nTree);
    1552916138            heuristic->setMaxSimplexIterationsAtRoot(nRoot);
     
    1689417503    char general[200];
    1689517504    if (clpSolver) {
    16896       clpSolver->getModelPtr()->dual(); // to get more randomness
     17505      clpSolver->getModelPtr()->dual(); // probably not needed
    1689717506      clpSolver->setWarmStart(NULL);
    16898       sprintf(general,"Starting multiple root solver - %d iterations in initialSolve",
    16899               clpSolver->getIterationCount());
     17507      sprintf(general,"Starting multiple root solver");
    1690017508    } else {
    1690117509#endif
     
    1691617524    return NULL;
    1691717525}
    16918 
    16919 
     17526OsiRowCut *
     17527CbcModel::conflictCut(const OsiSolverInterface * solver, bool & localCuts)
     17528{
     17529  OsiRowCut * cut=NULL;
     17530  localCuts=false;
     17531# ifdef COIN_HAS_CLP
     17532  const OsiClpSolverInterface * clpSolver
     17533    = dynamic_cast<const OsiClpSolverInterface *> (solver);
     17534  if (clpSolver&&topOfTree_) {
     17535    int debugMode=0;
     17536    const double * originalLower = topOfTree_->lower();
     17537    const double * originalUpper = topOfTree_->upper();
     17538    int typeCut = 1;
     17539    ClpSimplex * simplex = clpSolver->getModelPtr();
     17540    assert(simplex->status()==1);
     17541    if(simplex->ray()) {
     17542      {
     17543        int numberRows=simplex->numberRows();
     17544        double * saveRay=CoinCopyOfArray(simplex->ray(),numberRows);
     17545#define SAFE_RAY
     17546#ifdef SAFE_RAY
     17547        ClpSimplex & tempSimplex=*simplex;
     17548#else
     17549        ClpSimplex tempSimplex=*simplex;
     17550#endif
     17551        int logLevel=simplex->logLevel();
     17552        tempSimplex.setLogLevel(63);
     17553        tempSimplex.scaling(0);
     17554        tempSimplex.dual();
     17555        tempSimplex.setLogLevel(logLevel);
     17556        if (!tempSimplex.numberIterations()) {
     17557          double * ray = tempSimplex.ray();
     17558          int nBad=0;
     17559          for (int i=0;i<numberRows;i++) {
     17560            if (fabs(ray[i]-saveRay[i])>1.0e-3) {
     17561              if (debugMode)
     17562                printf("row %d true %g bad %g - diff %g\n",
     17563                       i,ray[i],saveRay[i],ray[i]-saveRay[i]);
     17564              nBad++;
     17565            }
     17566          }
     17567          if (nBad)
     17568            printf("%d mismatch crunch ray values\n",nBad);
     17569        }
     17570        delete [] saveRay;
     17571      }
     17572     // make sure we use non-scaled versions
     17573      ClpPackedMatrix * saveMatrix = simplex->swapScaledMatrix(NULL);
     17574      double * saveScale = simplex->swapRowScale(NULL);
     17575      //printf("Could do normal cut\n");
     17576      // could use existing arrays
     17577      int numberRows=simplex->numberRows();
     17578      int numberColumns=simplex->numberColumns();
     17579      double * farkas = new double [2*numberColumns+numberRows];
     17580      double * bound = farkas + numberColumns;
     17581      double * effectiveRhs = bound + numberColumns;
     17582      // sign as internally for dual - so swap if primal
     17583      /*const*/ double * ray = simplex->ray();
     17584      // have to get rid of local cut rows
     17585      if (whichGenerator_) {
     17586        const int * whichGenerator = whichGenerator_ - numberRowsAtContinuous_;
     17587        int badRows=0;
     17588        for (int iRow=numberRowsAtContinuous_;iRow<numberRows;iRow++) {
     17589          int iType=whichGenerator[iRow];
     17590          if ((iType>=0&&iType<10000)||iType<-1) {
     17591            if (fabs(ray[iRow])>1.0e-10) {
     17592              badRows++;
     17593            } else {
     17594              ray[iRow]=0.0;
     17595            }
     17596          }
     17597        }
     17598        if (badRows) {
     17599          if ((debugMode&1)!=0)
     17600            printf("%d rows from local cuts\n",badRows);
     17601          localCuts=true;
     17602        }
     17603      }
     17604      // get farkas row
     17605      memset(farkas,0,(2*numberColumns+numberRows)*sizeof(double));
     17606      simplex->transposeTimes(-1.0,ray,farkas);
     17607      //const char * integerInformation = simplex->integerType_;
     17608      //assert (integerInformation);
     17609
     17610      int sequenceOut = simplex->sequenceOut();
     17611      // Put nonzero bounds in bound
     17612      const double * columnLower = simplex->columnLower();
     17613      const double * columnUpper = simplex->columnUpper();
     17614      int numberBad=0;
     17615      for (int i=0;i<numberColumns;i++) {
     17616        double value = farkas[i];
     17617        double boundValue=0.0;
     17618        if (simplex->getStatus(i)==ClpSimplex::basic) {
     17619          // treat as zero if small
     17620          if (fabs(value)<1.0e-8) {
     17621            value=0.0;
     17622            farkas[i]=0.0;
     17623          }
     17624          if (value) {
     17625            //printf("basic %d direction %d farkas %g\n",
     17626            //     i,simplex->directionOut(),value);
     17627            if (value<0.0)
     17628              boundValue=columnLower[i];
     17629            else
     17630              boundValue=columnUpper[i];
     17631          }
     17632        } else if (fabs(value)>1.0e-10) {
     17633          if (value<0.0)
     17634            boundValue=columnLower[i];
     17635          else
     17636            boundValue=columnUpper[i];
     17637        }
     17638        bound[i]=boundValue;
     17639        if (fabs(boundValue)>1.0e10)
     17640          numberBad++;
     17641      }
     17642      const double * rowLower = simplex->rowLower();
     17643      const double * rowUpper = simplex->rowUpper();
     17644      //int pivotRow = simplex->spareIntArray_[3];
     17645      //bool badPivot=pivotRow<0;
     17646      for (int i=0;i<numberRows;i++) {
     17647        double value = ray[i];
     17648        double rhsValue=0.0;
     17649        if (simplex->getRowStatus(i)==ClpSimplex::basic) {
     17650          // treat as zero if small
     17651          if (fabs(value)<1.0e-8) {
     17652            value=0.0;
     17653            ray[i]=0.0;
     17654          }
     17655          if (value) {
     17656            //printf("row basic %d direction %d ray %g\n",
     17657            //     i,simplex->directionOut(),value);
     17658            if (value<0.0)
     17659              rhsValue=rowLower[i];
     17660            else
     17661              rhsValue=rowUpper[i];
     17662          }
     17663        } else if (fabs(value)>1.0e-10) {
     17664          if (value<0.0)
     17665            rhsValue=rowLower[i];
     17666          else
     17667              rhsValue=rowUpper[i];
     17668        }
     17669        effectiveRhs[i]=rhsValue;
     17670      }
     17671      simplex->times(-1.0,bound,effectiveRhs);
     17672      simplex->swapRowScale(saveScale);
     17673      simplex->swapScaledMatrix(saveMatrix);
     17674      double bSum=0.0;
     17675      for (int i=0;i<numberRows;i++) {
     17676        bSum += effectiveRhs[i]*ray[i];
     17677      }
     17678      if (numberBad||bSum>-1.0e-4) {
     17679#ifndef NDEBUG
     17680        printf("bad BOUND bSum %g  - %d bad\n",
     17681               bSum,numberBad);
     17682#endif
     17683      } else {
     17684        const char * integerInformation = simplex->integerInformation();
     17685        assert (integerInformation);
     17686        int * conflict = new int[numberColumns];
     17687        double * sort = new double [numberColumns];
     17688        double relax=0.0;
     17689        int nConflict=0;
     17690        int nOriginal=0;
     17691        int nFixed=0;
     17692        for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     17693          if (integerInformation[iColumn]) {
     17694            if ((debugMode&1)!=0)
     17695            printf("%d status %d %g <= %g <=%g (orig %g, %g) farkas %g\n",
     17696                   iColumn,simplex->getStatus(iColumn),columnLower[iColumn],
     17697                   simplex->primalColumnSolution()[iColumn],columnUpper[iColumn],
     17698                   originalLower[iColumn],originalUpper[iColumn],
     17699                   farkas[iColumn]);
     17700            double gap = originalUpper[iColumn]-originalLower[iColumn];
     17701            if (!gap)
     17702              continue;
     17703            if (gap==columnUpper[iColumn]-columnLower[iColumn])
     17704              nOriginal++;
     17705            if (columnUpper[iColumn]==columnLower[iColumn])
     17706              nFixed++;
     17707            if (fabs(farkas[iColumn])<1.0e-15) {
     17708              farkas[iColumn]=0.0;
     17709              continue;
     17710            }
     17711            // temp
     17712            if (gap>=20000.0&&false) {
     17713              // can't use
     17714              if (farkas[iColumn]<0.0) {
     17715                assert(originalLower[iColumn]-columnLower[iColumn]<=0.0);
     17716                // farkas is negative - relax lower bound all way
     17717                relax +=
     17718                  farkas[iColumn]*(originalLower[iColumn]-columnLower[iColumn]);
     17719              } else {
     17720                assert(originalUpper[iColumn]-columnUpper[iColumn]>=0.0);
     17721                // farkas is positive - relax upper bound all way
     17722                relax +=
     17723                  farkas[iColumn]*(originalUpper[iColumn]-columnUpper[iColumn]);
     17724              }
     17725              continue;
     17726            }
     17727            if (originalLower[iColumn]==columnLower[iColumn]) {
     17728              if (farkas[iColumn]>0.0&&(simplex->getStatus(iColumn)==ClpSimplex::atUpperBound
     17729                                        ||simplex->getStatus(iColumn)==ClpSimplex::isFixed
     17730                                        ||iColumn==sequenceOut)) {
     17731                // farkas is positive - add to list
     17732                gap=originalUpper[iColumn]-columnUpper[iColumn];
     17733                if (gap) {
     17734                  sort[nConflict]=-farkas[iColumn]*gap;
     17735                  conflict[nConflict++]=iColumn;
     17736                }
     17737                //assert (gap>columnUpper[iColumn]-columnLower[iColumn]);
     17738              }
     17739            } else if (originalUpper[iColumn]==columnUpper[iColumn]) {
     17740              if (farkas[iColumn]<0.0&&(simplex->getStatus(iColumn)==ClpSimplex::atLowerBound
     17741                                        ||simplex->getStatus(iColumn)==ClpSimplex::isFixed
     17742                                        ||iColumn==sequenceOut)) {
     17743                // farkas is negative - add to list
     17744                gap=columnLower[iColumn]-originalLower[iColumn];
     17745                if (gap) {
     17746                  sort[nConflict]=farkas[iColumn]*gap;
     17747                  conflict[nConflict++]=iColumn;
     17748                }
     17749                //assert (gap>columnUpper[iColumn]-columnLower[iColumn]);
     17750              }
     17751            } else {
     17752              // can't use
     17753              if (farkas[iColumn]<0.0) {
     17754                assert(originalLower[iColumn]-columnLower[iColumn]<=0.0);
     17755                // farkas is negative - relax lower bound all way
     17756                relax +=
     17757                  farkas[iColumn]*(originalLower[iColumn]-columnLower[iColumn]);
     17758              } else {
     17759                assert(originalUpper[iColumn]-columnUpper[iColumn]>=0.0);
     17760                // farkas is positive - relax upper bound all way
     17761                relax +=
     17762                  farkas[iColumn]*(originalUpper[iColumn]-columnUpper[iColumn]);
     17763              }
     17764            }
     17765            assert(relax>=0.0);
     17766          } else {
     17767            // not integer - but may have been got at
     17768            double gap = originalUpper[iColumn]-originalLower[iColumn];
     17769            if (gap>columnUpper[iColumn]-columnLower[iColumn]) {
     17770              // can't use
     17771              if (farkas[iColumn]<0.0) {
     17772                assert(originalLower[iColumn]-columnLower[iColumn]<=0.0);
     17773                // farkas is negative - relax lower bound all way
     17774                relax +=
     17775                  farkas[iColumn]*(originalLower[iColumn]-columnLower[iColumn]);
     17776              } else {
     17777                assert(originalUpper[iColumn]-columnUpper[iColumn]>=0.0);
     17778                // farkas is positive - relax upper bound all way
     17779                relax +=
     17780                  farkas[iColumn]*(originalUpper[iColumn]-columnUpper[iColumn]);
     17781              }
     17782            }
     17783          }
     17784        }
     17785        if (relax+bSum>-1.0e-4||!nConflict) {
     17786          if (relax+bSum>-1.0e-4) {
     17787#ifndef NDEBUG
     17788            printf("General integers relax bSum to %g\n",relax+bSum);
     17789#endif
     17790          } else {
     17791            printf("All variables relaxed and still infeasible - what does this mean?\n");
     17792            int nR=0;
     17793            for (int i=0;i<numberRows;i++) {
     17794              if (fabs(ray[i])>1.0e-10)
     17795                nR++;
     17796              else
     17797                ray[i]=0.0;
     17798            }
     17799            int nC=0;
     17800            for (int i=0;i<numberColumns;i++) {
     17801              if (fabs(farkas[i])>1.0e-10)
     17802                nC++;
     17803              else
     17804                farkas[i]=0.0;
     17805            }
     17806            if (nR<3&&nC<5) {
     17807              printf("BAD %d nonzero rows, %d nonzero columns\n",nR,nC);
     17808            }
     17809          }
     17810        } else {
     17811          printf("BOUNDS violation bSum %g (relaxed %g) - %d at original bounds, %d fixed - %d in conflict\n",bSum,
     17812                 relax+bSum,nOriginal,nFixed,nConflict);
     17813          CoinSort_2(sort,sort+nConflict,conflict);
     17814          int nC=nConflict;
     17815          bSum+=relax;
     17816          double saveBsum = bSum;
     17817          while (nConflict) {
     17818            //int iColumn=conflict[nConflict-1];
     17819            double change=-sort[nConflict-1];
     17820            if (bSum+change>-1.0e-4)
     17821              break;
     17822            nConflict--;
     17823            bSum += change;
     17824          }
     17825          if (!nConflict) {
     17826            int nR=0;
     17827            for (int i=0;i<numberRows;i++) {
     17828              if (fabs(ray[i])>1.0e-10)
     17829                nR++;
     17830              else
     17831                ray[i]=0.0;
     17832            }
     17833            int nC=0;
     17834            for (int i=0;i<numberColumns;i++) {
     17835              if (fabs(farkas[i])>1.0e-10)
     17836                nC++;
     17837              else
     17838                farkas[i]=0.0;
     17839            }
     17840            if (nR<3&&nC<5) {
     17841              printf("BAD2 %d nonzero rows, %d nonzero columns\n",nR,nC);
     17842            }
     17843          }
     17844          // no point doing if no reduction (or big?) ?
     17845          if (nConflict<nC+1&&nConflict<500) {
     17846            cut=new OsiRowCut();
     17847            cut->setUb(COIN_DBL_MAX);
     17848            if (!typeCut) {
     17849              double lo=1.0;
     17850              for (int i=0;i<nConflict;i++) {
     17851                int iColumn = conflict[i];
     17852                if (originalLower[iColumn]==columnLower[iColumn]) {
     17853                  // must be at least one higher
     17854                  sort[i]=1.0;
     17855                  lo += originalLower[iColumn];
     17856                } else {
     17857                  // must be at least one lower
     17858                  sort[i]=-1.0;
     17859                  lo -= originalUpper[iColumn];
     17860                }
     17861              }
     17862              cut->setLb(lo);
     17863              cut->setRow(nConflict,conflict,sort);
     17864              printf("CUT has %d (started at %d) - final bSum %g\n",nConflict,nC,bSum);
     17865            } else {
     17866              // just save for use later
     17867              // first take off small
     17868              int nC2=nC;
     17869              while (nC2) {
     17870                //int iColumn=conflict[nConflict-1];
     17871                double change=-sort[nC2-1];
     17872                if (saveBsum+change>-1.0e-4||change>1.0e-4)
     17873                  break;
     17874                nC2--;
     17875                saveBsum += change;
     17876              }
     17877              cut->setLb(saveBsum);
     17878              for (int i=0;i<nC2;i++) {
     17879                int iColumn = conflict[i];
     17880                sort[i]=farkas[iColumn];
     17881              }
     17882              cut->setRow(nC2,conflict,sort);
     17883              printf("Stem CUT has %d (greedy %d - with small %d) - saved bSum %g final greedy bSum %g\n",
     17884                     nC2,nConflict,nC,saveBsum,bSum);
     17885            }
     17886          }
     17887        }
     17888        delete [] conflict;
     17889        delete [] sort;
     17890      }
     17891      delete [] farkas;
     17892    } else {
     17893      printf("No dual ray\n");
     17894    }
     17895  }
     17896#endif
     17897  return cut;
     17898}
     17899
  • trunk/Cbc/src/CbcModel.hpp

    r1876 r1880  
    359359    void makeGlobalCut(const OsiColCut & cut);
    360360    /// Make partial cut into a global cut and save
    361     void makePartialCut(const OsiRowCut * cut);
     361  void makePartialCut(const OsiRowCut * cut, const OsiSolverInterface * solver=NULL);
    362362    /// Make partial cuts into global cuts
    363363    void makeGlobalCuts();
     
    443443    */
    444444    void AddIntegers();
    445 
    446445    /**
    447446      Save copy of the model.
     
    868867    void setOriginalColumns(const int * originalColumns,
    869868                            int numberGood=COIN_INT_MAX) ;
     869    /// Create conflict cut (well - most of)
     870    OsiRowCut * conflictCut(const OsiSolverInterface * solver, bool & localCuts);
    870871
    871872    /** Set the print frequency.
     
    17351736    inline void setDefaultHandler(bool yesNo) {
    17361737        defaultHandler_ = yesNo;
     1738    }
     1739    /// Check default handler
     1740    inline bool defaultHandler() const {
     1741        return defaultHandler_;
    17371742    }
    17381743    //@}
     
    18191824        22 bit (4194304) - Conflict analysis
    18201825        23 bit (8388608) - Conflict analysis - temporary bit
    1821         24 bit (16777216) - Add cutoff as LP constraint
     1826        24 bit (16777216) - Add cutoff as LP constraint (out)
     1827        25 bit (33554432) - diving/reordering
     1828        26 bit (67108864) - load global cuts from file
     1829        27 bit (134217728) - append binding global cuts to file
     1830        28 bit (268435456) - idiot branching
     1831        29 bit (536870912) - don't make fake objective
    18221832    */
    18231833    inline void setMoreSpecialOptions(int value) {
     
    18281838        return moreSpecialOptions_;
    18291839    }
    1830   /// Set time method
     1840    /// Set cutoff as constraint
     1841    inline void setCutoffAsConstraint(bool yesNo) {
     1842      cutoffRowNumber_ = (yesNo) ? -2 : -1;
     1843    }
     1844    /// Set time method
    18311845    inline void setUseElapsedTime(bool yesNo) {
    18321846        if (yesNo)
     
    18391853        return (moreSpecialOptions_&131072)!=0;
    18401854    }
     1855    /// Get useful temporary pointer
     1856    inline void * temporaryPointer() const
     1857    { return temporaryPointer_;}
     1858    /// Set useful temporary pointer
     1859    inline void setTemporaryPointer(void * pointer)
     1860    { temporaryPointer_=pointer;}
    18411861    /// Go to dantzig pivot selection if easy problem (clp only)
    18421862#ifdef COIN_HAS_CLP
     
    24552475    /// Number of rows at continuous
    24562476    int numberRowsAtContinuous_;
     2477    /**
     2478       -1 - cutoff as constraint not activated
     2479       -2 - waiting to activate
     2480       >=0 - activated
     2481     */
     2482    int cutoffRowNumber_;
    24572483    /// Maximum number of cuts
    24582484    int maximumNumberCuts_;
     
    26222648    /// Arrays with analysis results
    26232649    double * analyzeResults_;
     2650    /// Useful temporary pointer
     2651    void * temporaryPointer_;
    26242652    /// Number of nodes infeasible by normal branching (before cuts)
    26252653    int numberInfeasibleNodes_;
  • trunk/Cbc/src/CbcNode.cpp

    r1839 r1880  
    16331633                    break;
    16341634                }
    1635             }
     1635            } else {
     1636              obj->initializeForBranching(model);
     1637            }
    16361638        }
    16371639    }
     
    21322134                    if (dynamicObject) {
    21332135                        // Use this object's numberBeforeTrust
    2134                         int numberBeforeTrust = dynamicObject->numberBeforeTrust();
     2136                        int numberBeforeTrustThis = dynamicObject->numberBeforeTrust();
    21352137                        iColumn = dynamicObject->columnNumber();
    21362138                        gotDown = false;
    21372139                        numberThisDown = dynamicObject->numberTimesDown();
    2138                         if (numberThisDown >= numberBeforeTrust)
     2140                        if (numberThisDown >= numberBeforeTrustThis)
    21392141                            gotDown = true;
    21402142                        gotUp = false;
    21412143                        numberThisUp = dynamicObject->numberTimesUp();
    2142                         if (numberThisUp >= numberBeforeTrust)
     2144                        if (numberThisUp >= numberBeforeTrustThis)
    21432145                            gotUp = true;
    21442146                        if (!depth_ && false) {
     
    25592561                            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
    25602562                        // Use this object's numberBeforeTrust
    2561                         int numberBeforeTrust = dynamicObject->numberBeforeTrust();
     2563                        int numberBeforeTrustThis = dynamicObject->numberBeforeTrust();
    25622564                        int iSequence = dynamicObject->columnNumber();
    25632565                        double value = saveSolution[iSequence];
     
    25672569                        int numberThisDown = dynamicObject->numberTimesDown();
    25682570                        int numberThisUp = dynamicObject->numberTimesUp();
    2569                         if (!numberBeforeTrust) {
     2571                        if (!numberBeforeTrustThis) {
    25702572                            // override
    25712573                            downEstimate[iObject] = downPenalty;
     
    25772579                            min1 = 0.8 * min1 + 0.2 * max1;
    25782580                            sort[j] = - min1;
    2579                         } else if (numberThisDown < numberBeforeTrust ||
    2580                                    numberThisUp < numberBeforeTrust) {
    2581                             double invTrust = 1.0 / static_cast<double> (numberBeforeTrust);
    2582                             if (numberThisDown < numberBeforeTrust) {
     2581                        } else if (numberThisDown < numberBeforeTrustThis ||
     2582                                   numberThisUp < numberBeforeTrustThis) {
     2583                            double invTrust = 1.0 / static_cast<double> (numberBeforeTrustThis);
     2584                            if (numberThisDown < numberBeforeTrustThis) {
    25832585                                double fraction = numberThisDown * invTrust;
    25842586                                downEstimate[iObject] = fraction * downEstimate[iObject] + (1.0 - fraction) * downPenalty;
    25852587                            }
    2586                             if (numberThisUp < numberBeforeTrust) {
     2588                            if (numberThisUp < numberBeforeTrustThis) {
    25872589                                double fraction = numberThisUp * invTrust;
    25882590                                upEstimate[iObject] = fraction * upEstimate[iObject] + (1.0 - fraction) * upPenalty;
     
    26172619                                delete branch;
    26182620                            }
    2619                             if (number >= numberBeforeTrust)
    2620                                 dynamicObject->setNumberBeforeTrust(number + 1);
     2621                            if (number >= numberBeforeTrustThis)
     2622                              dynamicObject->setNumberBeforeTrust(CoinMin(number + 1,5*numberBeforeTrust));
    26212623                            numberFixed++;
    26222624                        }
     
    28142816                numberTest = numberToDo;
    28152817                numberStrong=numberToDo;
    2816                 skipAll=false;
     2818                skipAll=0;
    28172819                searchStrategy=0;
    28182820                solver->setIntParam(OsiMaxNumIterationHotStart, 100000);
     2821                //printf("Strong branching type %d\n",strongType);
    28192822              }
    28202823            }
     
    28812884                if (model->messageHandler()->logLevel() > 3 && numberBeforeTrust && dynamicObject)
    28822885                    dynamicObject->print(1, choice.possibleBranch->value());
     2886                if (strongType)
     2887                  canSkip=0;
    28832888                if (skipAll < 0)
    2884                     canSkip = true;
    2885                 if (strongType)
    2886                   canSkip=false;
     2889                    canSkip = 1;
    28872890                if (!canSkip) {
    28882891                    if (!doneHotStart) {
     
    28962899                        if (!solver->isProvenOptimal()) {
    28972900                          skipAll=-2;
    2898                           canSkip = true;
     2901                          canSkip = 1;
    28992902                        }
    29002903                        xMark++;
     
    30893092                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
    30903093#endif
     3094                        if (!solver->isProvenOptimal()) {
     3095                          skipAll=-2;
     3096                          canSkip = 1;
     3097                        }
     3098                        xMark++;
    30913099                    }
    30923100#if 0 //def DO_ALL_AT_ROOT
     
    32773285                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
    32783286#endif
     3287                        if (!solver->isProvenOptimal()) {
     3288                          skipAll=-2;
     3289                          canSkip = 1;
     3290                        }
     3291                        xMark++;
    32793292                    }
    32803293
     
    34263439                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
    34273440#endif
     3441                        if (!solver->isProvenOptimal()) {
     3442                          skipAll=-2;
     3443                          canSkip = 1;
     3444                        }
     3445                        xMark++;
    34283446                        // may be infeasible (if other way stopped on iterations)
    34293447                        if (goneInfeasible) {
     
    34733491                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
    34743492#endif
     3493                        if (!solver->isProvenOptimal()) {
     3494                          skipAll=-2;
     3495                          canSkip = 1;
     3496                        }
     3497                        xMark++;
    34753498                        // may be infeasible (if other way stopped on iterations)
    34763499                        if (goneInfeasible) {
     
    36503673    }
    36513674    if (numberUnfinished*10 < numberStrongDone &&
    3652             numberStrongIterations*20 < model->getIterationCount()&&
     3675        model->numberStrongIterations()*20 < model->getIterationCount()&&
    36533676                                !auxiliaryInfo->solutionAddsCuts()) {
    36543677        //printf("increasing trust\n");
     
    45134536    double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
    45144537    if (thisOne->whichSolution() >= 0) {
    4515         ClpNode * nodeInfo = thisOne->nodeInfo(thisOne->whichSolution());
    4516         nodeInfo->applyNode(simplex, 2);
     4538        ClpNode * nodeInfo=NULL;
     4539        if ((model->moreSpecialOptions()&33554432)==0) {
     4540          nodeInfo = thisOne->nodeInfo(thisOne->whichSolution());
     4541          nodeInfo->applyNode(simplex, 2);
     4542        } else {
     4543          // from diving
     4544          CbcSubProblem ** nodes = reinterpret_cast<CbcSubProblem **>
     4545            (model->temporaryPointer());
     4546          assert (nodes);
     4547          int numberDo=thisOne->numberNodes()-1;
     4548          for (int iNode=0;iNode<numberDo;iNode++)
     4549            nodes[iNode]->apply(solver,1);
     4550          nodes[numberDo]->apply(solver,9+16);
     4551        }
    45174552        int saveLogLevel = simplex->logLevel();
    45184553        simplex->setLogLevel(0);
     
    45234558        if (simplex->status()) {
    45244559            //simplex->writeMps("bad7.mps",2);
    4525             if (nodeInfo->objectiveValue() > cutoff - 1.0e-2)
     4560            if (nodeInfo) {
     4561              if (nodeInfo->objectiveValue() > cutoff - 1.0e-2)
    45264562                goodSolution = false;
    4527             else
     4563              else
    45284564                assert (!simplex->status());
     4565            } else {
     4566              // debug diving
     4567              assert (!simplex->status());
     4568            }
    45294569        }
    45304570        if (goodSolution) {
     
    45844624        } else {
    45854625            branch_ = thisOne->createCbcBranch(solver, &usefulInfo, preferredWay);
    4586             // Set to firat one (and change when re-pushing)
    4587             CbcGeneralBranchingObject * branch =
     4626            if (branch_) {
     4627              // Set to first one (and change when re-pushing)
     4628              CbcGeneralBranchingObject * branch =
    45884629                dynamic_cast <CbcGeneralBranchingObject *> (branch_);
    4589             branch->state(objectiveValue_, sumInfeasibilities_,
    4590                           numberUnsatisfied_, 0);
    4591             branch->setNode(this);
    4592             anyAction = 0;
     4630              branch->state(objectiveValue_, sumInfeasibilities_,
     4631                            numberUnsatisfied_, 0);
     4632              branch->setNode(this);
     4633              anyAction = 0;
     4634            } else {
     4635              anyAction = -2; // mark as infeasible
     4636            }
    45934637        }
    45944638    } else {
  • trunk/Cbc/src/CbcObject.hpp

    r1854 r1880  
    252252    /// Redoes data when sequence numbers change
    253253    virtual void redoSequenceEtc(CbcModel * , int , const int * ) {}
     254    /// Initialize for branching
     255    virtual void initializeForBranching(CbcModel * ) {}
    254256
    255257protected:
  • trunk/Cbc/src/CbcSolver.cpp

    r1876 r1880  
    651651    parameters_[whichParam(CBC_PARAM_STR_MIXEDCUTS, numberParameters_, parameters_)].setCurrentOption("ifmove");
    652652    parameters_[whichParam(CBC_PARAM_STR_FLOWCUTS, numberParameters_, parameters_)].setCurrentOption("ifmove");
    653     parameters_[whichParam(CBC_PARAM_STR_TWOMIRCUTS, numberParameters_, parameters_)].setCurrentOption("root");
     653    parameters_[whichParam(CBC_PARAM_STR_TWOMIRCUTS, numberParameters_, parameters_)].setCurrentOption("ifmove");
    654654    parameters_[whichParam(CBC_PARAM_STR_LANDPCUTS, numberParameters_, parameters_)].setCurrentOption("off");
    655655    parameters_[whichParam(CBC_PARAM_STR_RESIDCUTS, numberParameters_, parameters_)].setCurrentOption("off");
     
    11101110            }
    11111111        }
    1112         input[i] = '\0';
     1112        input[i++] = '\0';
    11131113        argv[j+1] = CoinStrdup(input + saveI);
    11141114        while (input[i] == ' ')
     
    12741274        int * prioritiesIn = NULL;
    12751275        std::vector< std::pair< std::string, double > > mipStart;
     1276        std::vector< std::pair< std::string, double > > mipStartBefore;
    12761277        int numberSOS = 0;
    12771278        int * sosStart = NULL;
     
    15941595        //GMIGen.setLimit(100);
    15951596        // set default action (0=off,1=on,2=root)
    1596         // Off
     1597        // Off 
    15971598        int GMIAction = 0;
    15981599
     
    16181619        twomirGen.setMaxElements(250);
    16191620        // set default action (0=off,1=on,2=root)
    1620         int twomirAction = 2;
     1621        int twomirAction = 3;
    16211622#ifndef DEBUG_MALLOC
    16221623        CglLandP landpGen;
     1624        landpGen.validator().setMinViolation(1.0e-4);
    16231625#endif
    16241626        // set default action (0=off,1=on,2=root)
     
    16901692            iParam = whichParam(CBC_PARAM_STR_PROBINGCUTS, numberParameters_, parameters_);
    16911693            parameters_[iParam].setCurrentOption("on");
    1692             probingAction = 1;
    1693             parameters_[iParam].setCurrentOption("forceOnStrong");
    1694             probingAction = 8;
     1694            probingAction = 3;
     1695            //parameters_[iParam].setCurrentOption("forceOnStrong");
     1696            //probingAction = 8;
    16951697        }
    16961698        std::string field;
     
    33143316                                            assert (solver3);
    33153317                                            solution = solver3->bestSolution();
    3316                                             double bestObjectiveValue = solver3->bestObjectiveValue();
    3317                                             linkSolver->setBestObjectiveValue(bestObjectiveValue);
    3318                                             linkSolver->setBestSolution(solution, solver3->getNumCols());
     3318                                            double bestObjectiveValue = solver3->bestObjectiveValue();
     3319                                            linkSolver->setBestObjectiveValue(bestObjectiveValue);
     3320                                            if (solution) {
     3321                                              linkSolver->setBestSolution(solution, solver3->getNumCols());
     3322                                            }
    33193323                                            CbcHeuristicDynamic3 dynamic(model_);
    33203324                                            dynamic.setHeuristicName("dynamic pass thru");
    33213325                                            model_.addHeuristic(&dynamic);
    33223326                                            // if convex
    3323                                             if ((linkSolver->specialOptions2()&4) != 0) {
     3327                                            if ((linkSolver->specialOptions2()&4) != 0 && solution) {
    33243328                                                int numberColumns = coinModel->numberColumns();
    33253329                                                assert (linkSolver->objectiveVariable() == numberColumns);
     
    36683672                                    preProcess = 0;
    36693673                            }
     3674                            if (mipStartBefore.size())
     3675                              {
     3676                                CbcModel tempModel=*babModel_;
     3677                                std::vector< std::string > colNames;
     3678                                for ( int i=0 ; (i<babModel_->solver()->getNumCols()) ; ++i )
     3679                                  colNames.push_back( model_.solver()->getColName(i) );
     3680                                std::vector< double > x( babModel_->getNumCols(), 0.0 );
     3681                                double obj;
     3682                                int status = computeCompleteSolution( &tempModel, colNames, mipStartBefore, &x[0], obj );
     3683                                // set cutoff
     3684                                if (!status)
     3685                                  babModel_->setCutoff(CoinMin(babModel_->getCutoff(),obj+1.0e-4));
     3686                              }
    36703687                            if (preProcess && type == CBC_PARAM_ACTION_BAB) {
    36713688#ifndef CBC_OTHER_SOLVER
     
    37553772                                        process.passInProhibited(prohibited, numberColumns);
    37563773                                        delete [] prohibited;
     3774                                    }
     3775                                    if (0) {
     3776                                     
     3777                                      // Special integers
     3778                                      int numberColumns = saveSolver->getNumCols();
     3779                                      char * prohibited = new char[numberColumns];
     3780                                      memset(prohibited, 0, numberColumns);
     3781                                      const CoinPackedMatrix * matrix = saveSolver->getMatrixByCol();
     3782                                      const int * columnLength = matrix->getVectorLengths();
     3783                                      int numberProhibited=0;
     3784                                      for (int iColumn = numberColumns-1; iColumn >=0; iColumn--) {
     3785                                        if (!saveSolver->isInteger(iColumn)||
     3786                                            columnLength[iColumn]>1)
     3787                                          break;
     3788                                        numberProhibited++;
     3789                                        prohibited[iColumn] = 1;
     3790                                      }
     3791                                      if (numberProhibited) {
     3792                                        process.passInProhibited(prohibited, numberColumns);
     3793                                        printf("**** Treating last %d integers as special - give high priority?\n",numberProhibited);
     3794                                      }
     3795                                      delete [] prohibited;
    37573796                                    }
    37583797                                    if (!model_.numberObjects() && true) {
     
    41524191                                    delete [] rowAdd;
    41534192                                    delete [] objectiveNew;
    4154                           }
     4193                                  }
    41554194                                  delete [] which;
    41564195                                  delete [] obj;
     
    42834322                                        else if (useCosts == 5)
    42844323                                            dsort[n++] = -columnLength[iColumn];
     4324                                        else if (useCosts == 6)
     4325                                            dsort[n++] = (columnLength[iColumn]==1) ? -1.0 : 0.0;
     4326                                        else if (useCosts == 7)
     4327                                            dsort[n++] = (objective[iColumn]) ? -1.0 : 0.0;
    42854328                                    }
    42864329                                }
     
    43324375                            int numberGenerators = 0;
    43334376                            int translate[] = { -100, -1, -99, -98, 1, -1098, -999, 1, 1, 1, -1};
     4377                            int maximumSlowPasses =
     4378                              parameters_[whichParam(CBC_PARAM_INT_MAX_SLOW_CUTS,
     4379                                                     numberParameters_, parameters_)].intValue();
    43344380                            if (probingAction) {
    43354381                                int numberColumns = babModel_->solver()->getNumCols();
     
    44164462                                  int when = laGomory/3;
    44174463                                  char atEnd = (when<2) ? 1 : 0;
    4418                                   int gomoryTypeMajor = 0;
     4464                                  int gomoryTypeMajor = 10;
    44194465                                  if (when<3) {
    44204466                                    // normal as well
     
    44234469                                    switches[numberGenerators++] = 0;
    44244470                                    if (when==2)
    4425                                       gomoryTypeMajor=10;
     4471                                      gomoryTypeMajor=20;
    44264472                                  } else {
    44274473                                    when--; // so on
     
    44374483                                    accuracyFlag[numberGenerators] = 3;
    44384484                                    doAtEnd[numberGenerators]=atEnd;
     4485                                    if (atEnd) {
     4486                                      babModel_->cutGenerator(numberGenerators)->setMaximumTries(99999999);
     4487                                      babModel_->cutGenerator(numberGenerators)->setHowOften(1);
     4488                                    }
    44394489                                    switches[numberGenerators++] = 0;
    44404490                                  }
     
    44454495                                    accuracyFlag[numberGenerators] = 3;
    44464496                                    doAtEnd[numberGenerators]=atEnd;
     4497                                    if (atEnd) {
     4498                                      babModel_->cutGenerator(numberGenerators)->setMaximumTries(99999999);
     4499                                      babModel_->cutGenerator(numberGenerators)->setHowOften(1);
     4500                                    }
    44474501                                    switches[numberGenerators++] = 0;
    44484502                                  }
     
    44634517                                babModel_->addCutGenerator(&redsplitGen, translate[redsplitAction], "Reduce-and-split");
    44644518                                accuracyFlag[numberGenerators] = 5;
    4465                                 switches[numberGenerators++] = 1;
     4519                                // slow ? - just do a few times
     4520                                if (redsplitAction!=1) {
     4521                                  babModel_->cutGenerator(numberGenerators)->setMaximumTries(maximumSlowPasses);
     4522                                  babModel_->cutGenerator(numberGenerators)->setHowOften(10);
     4523                                }
     4524                                switches[numberGenerators++] = 1;
    44664525                            }
    44674526                            if (redsplit2Action && !complicatedInteger) {
     
    44754534                                babModel_->addCutGenerator(&redsplit2Gen, translate[redsplit2Action], "Reduce-and-split(2)");
    44764535                                accuracyFlag[numberGenerators] = 5;
     4536                                // slow ? - just do a few times
     4537                                if (redsplit2Action!=1) {
     4538                                  babModel_->cutGenerator(numberGenerators)->setHowOften(maximumSlowPasses);
     4539                                  babModel_->cutGenerator(numberGenerators)->setMaximumTries(maximumSlowPasses);
     4540                                  babModel_->cutGenerator(numberGenerators)->setHowOften(5);
     4541                                }
     4542                               
    44774543                                switches[numberGenerators++] = 1;
    44784544                            }
     
    44894555                                  GMIAction=2;
    44904556                                  doAtEnd[numberGenerators]=1;
     4557                                  babModel_->cutGenerator(numberGenerators)->setMaximumTries(99999999);
     4558                                  babModel_->cutGenerator(numberGenerators)->setHowOften(1);
    44914559                                }
    44924560                                accuracyFlag[numberGenerators] = 5;
     
    45184586                                    twomirGen.setMaxElements(2000);
    45194587                                }
    4520                                 babModel_->addCutGenerator(&twomirGen, translate[twomirAction], "TwoMirCuts");
    4521                                 accuracyFlag[numberGenerators] = 4;
    4522                                 switches[numberGenerators++] = 1;
     4588                                int laTwomir = parameters_[whichParam(CBC_PARAM_STR_LATWOMIRCUTS, numberParameters_, parameters_)].currentOptionAsInteger();
     4589                                int twomirType = translate[twomirAction];
     4590                                if (!laTwomir) {
     4591                                  // Normal
     4592                                  babModel_->addCutGenerator(&twomirGen, translate[twomirAction], "TwoMirCuts");
     4593                                  accuracyFlag[numberGenerators] = 4;
     4594                                  switches[numberGenerators++] = 1;
     4595                                } else {
     4596                                  laTwomir--;
     4597                                  int type = (laTwomir % 3)+1;
     4598                                  int when = laTwomir/3;
     4599                                  char atEnd = (when<2) ? 1 : 0;
     4600                                  int twomirTypeMajor = 10;
     4601                                  if (when<3) {
     4602                                    // normal as well
     4603                                    babModel_->addCutGenerator(&twomirGen, translate[twomirAction], "TwoMirCuts");
     4604                                    accuracyFlag[numberGenerators] = 4;
     4605                                    switches[numberGenerators++] = 1;
     4606                                    if (when==2)
     4607                                      twomirTypeMajor=10;
     4608                                  } else {
     4609                                    when--; // so on
     4610                                    twomirTypeMajor=20;
     4611                                  }
     4612                                  if (!when)
     4613                                    twomirType=-99; // root
     4614                                  twomirGen.passInOriginalSolver(babModel_->solver());
     4615                                  if ((type&1) !=0) {
     4616                                    // clean
     4617                                    twomirGen.setTwomirType(twomirTypeMajor+1);
     4618                                    babModel_->addCutGenerator(&twomirGen, twomirType, "TwoMirCutsL1");
     4619                                    accuracyFlag[numberGenerators] = 4;
     4620                                    doAtEnd[numberGenerators]=atEnd;
     4621                                    switches[numberGenerators++] = atEnd ? 0 : 1;
     4622                                  }
     4623                                  if ((type&2) !=0) {
     4624                                    // simple
     4625                                    twomirGen.setTwomirType(twomirTypeMajor+2);
     4626                                    babModel_->addCutGenerator(&twomirGen, twomirType, "TwoMirCutsL2");
     4627                                    accuracyFlag[numberGenerators] = 4;
     4628                                    doAtEnd[numberGenerators]=atEnd;
     4629                                    switches[numberGenerators++] = atEnd ? 0 : 1;
     4630                                  }
     4631                                }
    45234632                            }
    45244633#ifndef DEBUG_MALLOC
     
    45264635                                babModel_->addCutGenerator(&landpGen, translate[landpAction], "LiftAndProject");
    45274636                                accuracyFlag[numberGenerators] = 5;
     4637                                // slow ? - just do a few times
     4638                                if (landpAction!=1) {
     4639                                  babModel_->cutGenerator(numberGenerators)->setMaximumTries(maximumSlowPasses);
     4640                                  babModel_->cutGenerator(numberGenerators)->setHowOften(10);
     4641                                }
    45284642                                switches[numberGenerators++] = 1;
    45294643                            }
     
    45534667                                CbcCutGenerator * generator = babModel_->cutGenerator(iGenerator);
    45544668                                int howOften = generator->howOften();
    4555                                 if (howOften == -98 || howOften == -99)
     4669                                if (howOften == -98 || howOften == -99 || generator->maximumTries()>0)
    45564670                                    generator->setSwitchOffIfLessThan(switches[iGenerator]);
    45574671                                // Use if any at root as more likely later and fairly cheap
     
    45614675                                if (doAtEnd[iGenerator]) {
    45624676                                    generator->setWhetherCallAtEnd(true);
    4563                                     generator->setMustCallAgain(true);
     4677                                    //generator->setMustCallAgain(true);
    45644678                                }
    45654679                                generator->setTiming(true);
     
    48124926#endif
    48134927                                const int * originalColumns = preProcess ? process.originalColumns() : NULL;
    4814 
    48154928                                if (mipStart.size())
    48164929                                {
     4930                                   std::vector< std::string > colNames;
    48174931                                   if (preProcess)
    48184932                                   {
    4819                                       std::vector< std::string > colNames;
     4933                                     for ( int i=0 ; (i<babModel_->solver()->getNumCols()) ; ++i ) {
     4934                                       int iColumn = babModel_->originalColumns()[i];
     4935                                       if (iColumn>=0) {
     4936                                         colNames.push_back( model_.solver()->getColName( iColumn ) );
     4937                                       } else {
     4938                                         // created variable
     4939                                         char newName[15];
     4940                                         sprintf(newName,"C%7.7d",i);
     4941                                         colNames.push_back( newName );
     4942                                       }
     4943                                     }
     4944                                   } else {
    48204945                                      for ( int i=0 ; (i<babModel_->solver()->getNumCols()) ; ++i )
    4821                                          colNames.push_back( model_.solver()->getColName( babModel_->originalColumns()[i] ) );
    4822                                       //printf("--- %s %d\n", babModel_->solver()->getColName(0).c_str(), babModel_->solver()->getColNames().size() );
    4823                                       //printf("-- SIZES of models %d %d %d\n", model_.getNumCols(),  babModel_->solver()->getNumCols(), babModel_->solver()->getColNames().size() );
    4824                                       std::vector< double > x( babModel_->getNumCols(), 0.0 );
    4825                                       double obj;
    4826                                       int status = computeCompleteSolution( babModel_, colNames, mipStart, &x[0], obj );
    4827                                       if (!status)
    4828                                          babModel_->setBestSolution( &x[0], x.size(), obj, false );
    4829                                    }
    4830                                 }
     4946                                         colNames.push_back( model_.solver()->getColName(i) );
     4947                                   }
     4948                                   //printf("--- %s %d\n", babModel_->solver()->getColName(0).c_str(), babModel_->solver()->getColNames().size() );
     4949                                   //printf("-- SIZES of models %d %d %d\n", model_.getNumCols(),  babModel_->solver()->getNumCols(), babModel_->solver()->getColNames().size() );
     4950                                   std::vector< double > x( babModel_->getNumCols(), 0.0 );
     4951                                   double obj;
     4952                                   int status = computeCompleteSolution( babModel_, colNames, mipStart, &x[0], obj );
     4953                                   if (!status)
     4954                                     babModel_->setBestSolution( &x[0], x.size(), obj, false );
     4955                                }
    48314956
    48324957                                if (solutionIn && useSolution >= 0) {
     
    49135038                                        // extend arrays in case SOS
    49145039                                        assert (originalColumns);
    4915                                         int n = originalColumns[numberColumns-1] + 1;
     5040                                        int n = CoinMin(truncateColumns,numberColumns);
     5041                                        n = originalColumns[n-1] + 1;
    49165042                                        n = CoinMax(n, CoinMax(numberColumns, numberOriginalColumns));
    49175043                                        int * newColumn = new int[n];
     
    59486074                                }
    59496075#endif
    5950                                 int multipleRoot = parameters_[whichParam(CBC_PARAM_INT_MULTIPLEROOTS, numberParameters_, parameters_)].intValue();
    5951                                 babModel_->setMultipleRootTries(multipleRoot);
    59526076                                int specialOptions = parameters_[whichParam(CBC_PARAM_INT_STRONG_STRATEGY, numberParameters_, parameters_)].intValue();
    59536077                                if (specialOptions>=0)
    59546078                                  babModel_->setStrongStrategy(specialOptions);
     6079                                int jParam = whichParam(CBC_PARAM_STR_CUTOFF_CONSTRAINT,
     6080                                                        numberParameters_, parameters_);
     6081                                if(parameters_[jParam].currentOptionAsInteger())
     6082                                  babModel_->setCutoffAsConstraint(true);
     6083                                int multipleRoot = parameters_[whichParam(CBC_PARAM_INT_MULTIPLEROOTS, numberParameters_, parameters_)].intValue();
     6084                                if (multipleRoot<10000) {
     6085                                  babModel_->setMultipleRootTries(multipleRoot);
     6086                                } else {
     6087                                  // will be doing repeated solves and saves
     6088                                  int numberGoes=multipleRoot/10000;
     6089                                  multipleRoot-=10000*numberGoes;
     6090                                  int moreOptions=babModel_->moreSpecialOptions();
     6091                                  if (numberGoes<100) {
     6092                                    remove("global.cuts");
     6093                                    remove("global.fix");
     6094                                    moreOptions |= (67108864|134217728);
     6095                                  } else {
     6096                                    moreOptions |= 67108864*(numberGoes/100);
     6097                                    numberGoes=numberGoes%100;
     6098                                  }
     6099                                  babModel_->setMultipleRootTries(multipleRoot);
     6100                                  babModel_->setMoreSpecialOptions(moreOptions);
     6101                                  int numberColumns=babModel_->getNumCols();
     6102                                  double * bestValues=new double [numberGoes];
     6103                                  double ** bestSolutions=new double * [numberGoes];
     6104                                  int * which=new int[numberGoes];
     6105                                  int numberSolutions=0;
     6106                                  sprintf(generalPrint,"Starting %d passes each with %d solvers",
     6107                                          numberGoes, multipleRoot%10);
     6108                                  generalMessageHandler->message(CLP_GENERAL, generalMessages)
     6109                                    << generalPrint
     6110                                    << CoinMessageEol;
     6111                                  for (int iGo=0;iGo<numberGoes;iGo++) {
     6112                                    sprintf(generalPrint,"Starting pass %d",iGo+1);
     6113                                    generalMessageHandler->message(CLP_GENERAL, generalMessages)
     6114                                      << generalPrint
     6115                                      << CoinMessageEol;
     6116                                    CbcModel tempModel=*babModel_;
     6117                                    tempModel.setMaximumNodes(0);
     6118                                    // switch off cuts if none generated
     6119                                    int numberGenerators = tempModel.numberCutGenerators();
     6120                                    for (int iGenerator = 0; iGenerator < numberGenerators; iGenerator++) {
     6121                                      CbcCutGenerator * generator = tempModel.cutGenerator(iGenerator);
     6122                                      generator->setSwitchOffIfLessThan(1);
     6123                                    }
     6124                                    // random
     6125                                    tempModel.setRandomSeed(tempModel.getRandomSeed()+100000000*(iGo+1+5*numberGoes));
     6126                                    for (int i=0;i<tempModel.numberHeuristics();i++)
     6127                                      tempModel.heuristic(i)->setSeed(tempModel.heuristic(i)->getSeed()+100000000*iGo);
     6128#ifndef CBC_OTHER_SOLVER
     6129                                    OsiClpSolverInterface * solver = dynamic_cast<OsiClpSolverInterface *> (tempModel.solver());
     6130                                    ClpSimplex * simplex = solver->getModelPtr();
     6131                                    int solverSeed=simplex->randomNumberGenerator()->getSeed();
     6132                                    simplex->setRandomSeed(solverSeed+100000000*(iGo+1));
     6133#endif
     6134                                    tempModel.branchAndBound();
     6135                                    if (tempModel.bestSolution()) {
     6136                                      bestSolutions[numberSolutions]=
     6137                                        CoinCopyOfArray(tempModel.bestSolution(),
     6138                                                        numberColumns);
     6139                                      bestValues[numberSolutions]=-tempModel.getMinimizationObjValue();
     6140                                      which[numberSolutions]=numberSolutions;
     6141                                      numberSolutions++;
     6142                                    }
     6143                                  }
     6144                                  // allow solutions
     6145                                  double sense = babModel_->solver()->getObjSense();;
     6146                                  CoinSort_2(bestValues,bestValues+numberSolutions,which);
     6147                                  babModel_->setMoreSpecialOptions(moreOptions&(~16777216));
     6148                                  for (int i=0;i<numberSolutions;i++) {
     6149                                    int k=which[i];
     6150                                    if (bestValues[i]<babModel_->getCutoff()) {
     6151                                      babModel_->setBestSolution(bestSolutions[k],numberColumns,
     6152                                                                 -bestValues[i]*sense,true);
     6153                                      babModel_->incrementUsed(bestSolutions[k]);
     6154                                    }
     6155                                    delete [] bestSolutions[k];
     6156                                  }
     6157                                  babModel_->setMoreSpecialOptions(moreOptions);
     6158                                  if (numberSolutions)
     6159                                    sprintf(generalPrint,"Ending major passes - best solution %g",-bestValues[numberSolutions-1]);
     6160                                  else
     6161                                    sprintf(generalPrint,"Ending major passes - no solution found");
     6162                                    generalMessageHandler->message(CLP_GENERAL, generalMessages)
     6163                                      << generalPrint
     6164                                      << CoinMessageEol;
     6165                                  delete [] which;
     6166                                  delete [] bestValues;
     6167                                  delete [] bestSolutions;
     6168                                }
    59556169                                if (biLinearProblem)
    59566170                                  babModel_->setSpecialOptions(babModel_->specialOptions() &(~(512|32768)));
     
    75197733                            double msObj;
    75207734                            readMIPStart( &model_, fileName.c_str(), mipStart, msObj );
     7735                            // copy to before preprocess if has .before.
     7736                            if (strstr(fileName.c_str(),".before.")) {
     7737                              mipStartBefore = mipStart;
     7738                              sprintf(generalPrint,"file %s will be used before preprocessing.",fileName.c_str() );
     7739                              generalMessageHandler->message(CLP_GENERAL, generalMessages)
     7740                                << generalPrint
     7741                                << CoinMessageEol;
     7742                            }
    75217743                        } else {
    75227744#ifndef DISALLOW_PRINTING
  • trunk/Cbc/src/CbcSolverAnalyze.cpp

    r1854 r1880  
    307307        int logLevel = generalMessageHandler->logLevel();
    308308        CbcModel model(*solver);
    309         model.passInMessageHandler(generalMessageHandler);
     309        if (!model.defaultHandler())
     310          model.passInMessageHandler(generalMessageHandler);
    310311        if (noPrinting_)
    311312            model.setLogLevel(0);
  • trunk/Cbc/src/CbcSolverHeuristics.cpp

    r1854 r1880  
    14561456
    14571457    if (useDIVING > 0) {
     1458        int majorIterations=64;
    14581459        int diveOptions2 = parameters_[whichParam(CBC_PARAM_INT_DIVEOPT, numberParameters_, parameters_)].intValue();
    14591460        int diveOptions;
    14601461        if (diveOptions2 > 99) {
    14611462            // switch on various active set stuff
    1462             diveOptions = diveOptions2 % 100;
    1463             diveOptions2 -= diveOptions;
     1463            diveOptions = diveOptions2%100;
     1464            diveOptions2 /= 100;
    14641465        } else {
    14651466            diveOptions = diveOptions2;
     
    14721473            heuristicDV.setHeuristicName("DiveVectorLength");
    14731474            heuristicDV.setWhen(diveOptions);
     1475            if (diveOptions2) {
     1476              heuristicDV.setMaxIterations(majorIterations);
     1477              heuristicDV.setPercentageToFix(0.0);
     1478              heuristicDV.setMaxSimplexIterations(COIN_INT_MAX);
     1479              heuristicDV.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
     1480            }
    14741481            model->addHeuristic(&heuristicDV) ;
    14751482        }
     
    14781485            heuristicDG.setHeuristicName("DiveGuided");
    14791486            heuristicDG.setWhen(diveOptions);
     1487            if (diveOptions2) {
     1488              heuristicDG.setMaxIterations(majorIterations);
     1489              heuristicDG.setPercentageToFix(0.0);
     1490              heuristicDG.setMaxSimplexIterations(COIN_INT_MAX);
     1491              heuristicDG.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
     1492            }
    14801493            model->addHeuristic(&heuristicDG) ;
    14811494        }
     
    14841497            heuristicDF.setHeuristicName("DiveFractional");
    14851498            heuristicDF.setWhen(diveOptions);
     1499            if (diveOptions2) {
     1500              heuristicDF.setMaxIterations(majorIterations);
     1501              heuristicDF.setPercentageToFix(0.0);
     1502              heuristicDF.setMaxSimplexIterations(COIN_INT_MAX);
     1503              heuristicDF.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
     1504            }
    14861505            model->addHeuristic(&heuristicDF) ;
    14871506        }
     
    14901509            heuristicDC.setHeuristicName("DiveCoefficient");
    14911510            heuristicDC.setWhen(diveOptions);
     1511            if (diveOptions2) {
     1512              heuristicDC.setMaxIterations(majorIterations);
     1513              heuristicDC.setPercentageToFix(0.0);
     1514              heuristicDC.setMaxSimplexIterations(COIN_INT_MAX);
     1515              heuristicDC.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
     1516            }
    14921517            model->addHeuristic(&heuristicDC) ;
    14931518        }
     
    14961521            heuristicDL.setHeuristicName("DiveLineSearch");
    14971522            heuristicDL.setWhen(diveOptions);
     1523            if (diveOptions2) {
     1524              heuristicDL.setMaxIterations(majorIterations);
     1525              heuristicDL.setPercentageToFix(0.0);
     1526              heuristicDL.setMaxSimplexIterations(COIN_INT_MAX);
     1527              heuristicDL.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
     1528            }
    14981529            model->addHeuristic(&heuristicDL) ;
    14991530        }
     
    15011532            CbcHeuristicDivePseudoCost heuristicDP(*model);
    15021533            heuristicDP.setHeuristicName("DivePseudoCost");
    1503             heuristicDP.setWhen(diveOptions + diveOptions2);
     1534            heuristicDP.setWhen(diveOptions /*+ diveOptions2*/);
     1535            if (diveOptions2) {
     1536              heuristicDP.setMaxIterations(majorIterations);
     1537              heuristicDP.setPercentageToFix(0.0);
     1538              heuristicDP.setMaxSimplexIterations(COIN_INT_MAX);
     1539              heuristicDP.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
     1540            }
    15041541            model->addHeuristic(&heuristicDP) ;
    15051542        }
  • trunk/Cbc/src/CbcSubProblem.cpp

    r1854 r1880  
    2929// Default Constructor
    3030CbcSubProblem::CbcSubProblem()
    31         : objectiveValue_(0.0),
    32         sumInfeasibilities_(0.0),
    33         variables_(NULL),
    34         newBounds_(NULL),
    35         status_(NULL),
    36         depth_(0),
    37         numberChangedBounds_(0),
    38         numberInfeasibilities_(0)
     31  : objectiveValue_(0.0),
     32    sumInfeasibilities_(0.0),
     33    branchValue_(0.0),
     34    djValue_(0.0),
     35    variables_(NULL),
     36    newBounds_(NULL),
     37    status_(NULL),
     38    depth_(0),
     39    numberChangedBounds_(0),
     40    numberInfeasibilities_(0),
     41    problemStatus_(0),
     42    branchVariable_(0)
    3943{
    4044}
     
    4650                              const unsigned char * status,
    4751                              int depth)
    48         : objectiveValue_(0.0),
    49         sumInfeasibilities_(0.0),
    50         variables_(NULL),
    51         newBounds_(NULL),
    52         status_(NULL),
    53         depth_(depth),
    54         numberChangedBounds_(0),
    55         numberInfeasibilities_(0)
     52  : objectiveValue_(0.0),
     53    sumInfeasibilities_(0.0),
     54    branchValue_(0.0),
     55    djValue_(0.0),
     56    variables_(NULL),
     57    newBounds_(NULL),
     58    status_(NULL),
     59    depth_(depth),
     60    numberChangedBounds_(0),
     61    numberInfeasibilities_(0),
     62    problemStatus_(0),
     63    branchVariable_(0)
    5664{
    5765    const double * lower = solver->getColLower();
     
    110118CbcSubProblem::CbcSubProblem ( const CbcSubProblem & rhs)
    111119        : objectiveValue_(rhs.objectiveValue_),
    112         sumInfeasibilities_(rhs.sumInfeasibilities_),
    113         variables_(NULL),
    114         newBounds_(NULL),
    115         status_(NULL),
    116         depth_(rhs.depth_),
    117         numberChangedBounds_(rhs.numberChangedBounds_),
    118         numberInfeasibilities_(rhs.numberInfeasibilities_)
     120          sumInfeasibilities_(rhs.sumInfeasibilities_),
     121          branchValue_(rhs.branchValue_),
     122          djValue_(rhs.djValue_),
     123          variables_(NULL),
     124          newBounds_(NULL),
     125          status_(NULL),
     126          depth_(rhs.depth_),
     127          numberChangedBounds_(rhs.numberChangedBounds_),
     128          numberInfeasibilities_(rhs.numberInfeasibilities_),
     129          problemStatus_(rhs.problemStatus_),
     130          branchVariable_(rhs.branchVariable_)
    119131{
    120132    if (numberChangedBounds_) {
     
    137149        objectiveValue_ = rhs.objectiveValue_;
    138150        sumInfeasibilities_ = rhs.sumInfeasibilities_;
     151        branchValue_ = rhs.branchValue_;
     152        djValue_ = rhs.djValue_;
    139153        depth_ = rhs.depth_;
    140154        numberChangedBounds_ = rhs.numberChangedBounds_;
    141155        numberInfeasibilities_ = rhs.numberInfeasibilities_;
     156        problemStatus_ = rhs.problemStatus_;
     157        branchVariable_ = rhs.branchVariable_;
    142158        if (numberChangedBounds_) {
    143159            variables_ = CoinCopyOfArray(rhs.variables_, numberChangedBounds_);
     
    155171    return *this;
    156172}
     173// Take over
     174void
     175CbcSubProblem::takeOver( CbcSubProblem & rhs, bool cleanUp)
     176{
     177    if (this != &rhs) {
     178        delete [] variables_;
     179        delete [] newBounds_;
     180        delete status_;
     181        objectiveValue_ = rhs.objectiveValue_;
     182        sumInfeasibilities_ = rhs.sumInfeasibilities_;
     183        branchValue_ = rhs.branchValue_;
     184        djValue_ = rhs.djValue_;
     185        depth_ = rhs.depth_;
     186        numberChangedBounds_ = rhs.numberChangedBounds_;
     187        numberInfeasibilities_ = rhs.numberInfeasibilities_;
     188        problemStatus_ = rhs.problemStatus_;
     189        branchVariable_ = rhs.branchVariable_;
     190        variables_ = rhs.variables_;
     191        newBounds_ = rhs.newBounds_;
     192        rhs.variables_ = NULL;
     193        rhs.newBounds_ = NULL;
     194        status_ = rhs.status_;
     195        rhs.status_ = NULL;
     196        if (cleanUp) {
     197          delete [] variables_;
     198          delete [] newBounds_;
     199          variables_ = new int [1];
     200          newBounds_ = new double [1];
     201          // swap way and make only fix
     202          numberChangedBounds_=1;
     203          if ((problemStatus_&1)==0) {
     204            // last way was down
     205            newBounds_[0] = ceil(branchValue_);
     206            variables_[0] = branchVariable_;
     207          } else {
     208            // last way was up
     209            newBounds_[0] = floor(branchValue_);
     210            variables_[0] = branchVariable_ | 0x80000000;
     211          }
     212        }
     213    }
     214}
    157215
    158216// Destructor
     
    169227    int i;
    170228    if ((what&1) != 0) {
     229    printf("CbcSubapply depth %d column %d way %d bvalue %g obj %g\n",
     230                 this->depth_,this->branchVariable_,this->problemStatus_,
     231                 this->branchValue_,this->objectiveValue_);
     232    printf("current bounds %g <= %g <= %g\n",solver->getColLower()[branchVariable_],branchValue_,solver->getColUpper()[branchVariable_]);
    171233#ifndef NDEBUG
    172234        int nSame = 0;
     
    228290#endif
    229291#endif
     292    printf("new bounds %g <= %g <= %g\n",solver->getColLower()[branchVariable_],branchValue_,solver->getColUpper()[branchVariable_]);
    230293    }
    231294#ifdef JJF_ZERO
     
    258321        assert (clpSolver);
    259322        clpSolver->setBasis(*status_);
    260         delete status_;
    261         status_ = NULL;
    262     }
    263 }
    264 
     323        if ((what&16)==0) {
     324          delete status_;
     325          status_ = NULL;
     326        }
     327    }
     328}
     329
  • trunk/Cbc/src/CbcSubProblem.hpp

    r1854 r1880  
    4040    virtual ~CbcSubProblem ();
    4141
     42    /// Take over
     43  void takeOver ( CbcSubProblem &, bool cleanup);
    4244    /// Apply subproblem (1=bounds, 2=basis, 3=both)
    4345    void apply(OsiSolverInterface * model, int what = 3) const;
     
    4850    /// Sum of infeasibilities
    4951    double sumInfeasibilities_;
     52    /// Branch value
     53    double branchValue_;
     54    /// Dj on branching variable at end
     55    double djValue_;
    5056    /** Which variable (top bit if upper bound changing)
    5157        next bit if changing on down branch only */
     
    6167    /// Number of infeasibilities
    6268    int numberInfeasibilities_;
     69    /** Status 1 bit going up on first, 2 bit set first branch infeasible on second, 4 bit redundant branch,
     70        bits after 256 give reason for stopping (just last node)
     71        0 - solution
     72        1 - infeasible
     73        2 - maximum depth
     74        >2 - error or max time or something
     75    */
     76    int problemStatus_;
     77    /// Variable branched on
     78    int branchVariable_;
    6379};
    6480
  • trunk/Cbc/src/CbcThread.cpp

    r1877 r1880  
    10131013        saveObjects_[iObject]->updateBefore(object[iObject]);
    10141014    }
     1015    //#define FAKE_PARALLEL
     1016#ifndef FAKE_PARALLEL
    10151017    for (iThread = 0; iThread < numberThreads_; iThread++) {
    10161018        children_[iThread].setReturnCode( 0);
     
    10311033    for (iThread = 0; iThread < numberThreads_; iThread++)
    10321034        children_[iThread].setReturnCode(-1);
     1035#else
     1036    // wait
     1037    bool finished = false;
     1038    double time = getTime();
     1039    for (iThread = 0; iThread < numberThreads_; iThread++) {
     1040        children_[iThread].setReturnCode( 0);
     1041        children_[iThread].signal();
     1042        while (!finished) {
     1043          children_[numberThreads_].waitNano( 1000000); // millisecond
     1044          finished = (children_[iThread].returnCode() >0);
     1045        }
     1046        children_[iThread].setReturnCode(-1);
     1047        finished=false;
     1048    }
     1049#endif
    10331050    children_[numberThreads_].incrementTimeInThread(getTime() - time);
    10341051    // Unmark marked
  • trunk/Cbc/src/unitTestClp.cpp

    r1668 r1880  
    677677      << " -- (" << model->getNodeCount() << " n / "
    678678      << model->getIterationCount() << " i / "
    679       << timeOfSolution << " s) (subtotal " << timeTaken << " s)"
     679      << timeOfSolution << " s) (subtotal " << timeTaken << " seconds)"
    680680      << std::endl;
    681681    delete model;
Note: See TracChangeset for help on using the changeset viewer.