Changeset 833 for trunk


Ignore:
Timestamp:
Nov 1, 2007 4:12:08 PM (12 years ago)
Author:
forrest
Message:

changes to heuristics and allow collection of more statistics

Location:
trunk/Cbc/src
Files:
12 edited

Legend:

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

    r687 r833  
    336336    status_(0),
    337337    intDecrease_(0),
    338     branchingValue_(0.0)
     338    branchingValue_(0.0),
     339    originalObjective_(COIN_DBL_MAX),
     340    cutoff_(COIN_DBL_MAX)
    339341{
    340342}
     
    353355    status_(status),
    354356    intDecrease_(intDecrease),
    355     branchingValue_(branchingValue)
     357    branchingValue_(branchingValue),
     358    originalObjective_(COIN_DBL_MAX),
     359    cutoff_(COIN_DBL_MAX)
    356360{
    357361}
     
    370374    status_(rhs.status_),
    371375    intDecrease_(rhs.intDecrease_),
    372     branchingValue_(rhs.branchingValue_)
     376    branchingValue_(rhs.branchingValue_),
     377    originalObjective_(rhs.originalObjective_),
     378    cutoff_(rhs.cutoff_)
    373379{
    374380}
     
    386392    intDecrease_ = rhs.intDecrease_;
    387393    branchingValue_ = rhs.branchingValue_;
     394    originalObjective_ = rhs.originalObjective_;
     395    cutoff_ = rhs.cutoff_;
    388396  }
    389397  return *this;
  • trunk/Cbc/src/CbcBranchBase.hpp

    r765 r833  
    548548  /// Branching value
    549549  double branchingValue_;
     550  /// Objective value before branching
     551  double originalObjective_;
     552  /// Current cutoff
     553  double cutoff_;
    550554
    551555};
  • trunk/Cbc/src/CbcBranchDynamic.cpp

    r702 r833  
    1717#include "CoinSort.hpp"
    1818#include "CoinError.hpp"
    19 
     19#ifdef COIN_DEVELOP
     20typedef struct {
     21  char where_;
     22  char status_;
     23  unsigned short sequence_;
     24  int numberUp_;
     25  int numberUpInf_;
     26  float sumUp_;
     27  float upEst_; // or change in obj in update
     28  int numberDown_;
     29  int numberDownInf_;
     30  float sumDown_;
     31  float downEst_; // or movement in value in update
     32} History;
     33History * history=NULL;
     34int numberHistory=0;
     35int maxHistory=0;
     36bool getHistoryStatistics_=true;
     37static void increaseHistory() {
     38  if (numberHistory==maxHistory) {
     39    maxHistory = 100+(3*maxHistory)/2;
     40    History * temp = new History [maxHistory];
     41    memcpy(temp,history,numberHistory*sizeof(History));
     42    delete [] history;
     43    history=temp;
     44  }
     45}
     46static bool addRecord(History newOne) {
     47  //if (!getHistoryStatistics_)
     48    return false;
     49  bool fromCompare=false;
     50  int i;
     51  for ( i=numberHistory-1;i>=0;i--) {
     52    if (newOne.sequence_!=history[i].sequence_)
     53      continue;
     54    if (newOne.where_!=history[i].where_)
     55      continue;
     56    if (newOne.numberUp_!=history[i].numberUp_)
     57      continue;
     58    if (newOne.sumUp_!=history[i].sumUp_)
     59      continue;
     60    if (newOne.numberUpInf_!=history[i].numberUpInf_)
     61      continue;
     62    if (newOne.upEst_!=history[i].upEst_)
     63      continue;
     64    if (newOne.numberDown_!=history[i].numberDown_)
     65      continue;
     66    if (newOne.sumDown_!=history[i].sumDown_)
     67      continue;
     68    if (newOne.numberDownInf_!=history[i].numberDownInf_)
     69      continue;
     70    if (newOne.downEst_!=history[i].downEst_)
     71      continue;
     72    // If B knock out previous B
     73    if (newOne.where_=='C') {
     74      fromCompare=true;
     75      if (newOne.status_=='B') {
     76        int j;
     77        for (j=i-1;j>=0;j--) {
     78          if (history[j].where_=='C') {
     79            if (history[j].status_=='I') {
     80              break;
     81            } else if (history[j].status_=='B') {
     82              history[j].status_=' ';
     83              break;
     84            }
     85          }
     86        }
     87      }
     88      break;
     89    }
     90  }
     91  if (i==-1||fromCompare) {
     92    //add
     93    increaseHistory();
     94    history[numberHistory++]=newOne;
     95    return true;
     96  } else {
     97    return false;
     98  }
     99}
     100#endif
    20101/** Default Constructor
    21102
     
    98179  sumDownChange_ = 2.0;
    99180  numberTimesDown_ = 2;
    100 #if 0
     181#define TYPE2 0
     182#if TYPE2==0
    101183  // No
    102184  sumUpCost_ = 0.0;
     
    159241  sumDownChange_ = 2.0;
    160242  numberTimesDown_ = 2;
    161 #if 0
     243#if TYPE2==0
    162244  // No
    163245  sumUpCost_ = 0.0;
     
    389471    below = above -1;
    390472  }
     473  // was 1 - but that looks flakey
    391474#define INFEAS 1
    392475#if INFEAS==1
     
    400483#endif
    401484  double sum;
    402   int number;
     485  double number;
    403486  double downCost = CoinMax(value-below,0.0);
     487#if TYPE2==0
    404488  sum = sumDownCost_;
    405489  number = numberTimesDown_;
    406490#if INFEAS==1
     491  sum += numberTimesDownInfeasible_*CoinMax(distanceToCutoff/(downCost+1.0e-12),sumDownCost_);
     492  number += numberTimesDownInfeasible_;
     493#endif
     494#elif TYPE2==1
     495  sum = sumDownCost_;
     496  number = sumDownChange_;
     497#if INFEAS==1
     498  sum += numberTimesDownInfeasible_*CoinMax(distanceToCutoff/(downCost+1.0e-12),sumDownCost_);
     499  number += numberTimesDownInfeasible_;
     500#endif
     501#elif TYPE2==2
     502  abort();
     503#if INFEAS==1
    407504  sum += numberTimesDownInfeasible_*(distanceToCutoff/(downCost+1.0e-12));
    408505  number += numberTimesDownInfeasible_;
    409506#endif
    410   if (number>0)
    411     downCost *= sum / (double) number;
     507#endif
     508  if (number>0.0)
     509    downCost *= sum / number;
    412510  else
    413511    downCost  *=  downDynamicPseudoCost_;
    414512  double upCost = CoinMax((above-value),0.0);
     513#if TYPE2==0
    415514  sum = sumUpCost_;
    416515  number = numberTimesUp_;
    417516#if INFEAS==1
     517  sum += numberTimesUpInfeasible_*CoinMax(distanceToCutoff/(upCost+1.0e-12),sumUpCost_);
     518  number += numberTimesUpInfeasible_;
     519#endif
     520#elif TYPE2==1
     521  sum = sumUpCost_;
     522  number = sumUpChange_;
     523#if INFEAS==1
     524  sum += numberTimesUpInfeasible_*CoinMax(distanceToCutoff/(upCost+1.0e-12),sumUpCost_);
     525  number += numberTimesUpInfeasible_;
     526#endif
     527#elif TYPE2==1
     528  abort();
     529#if INFEAS==1
    418530  sum += numberTimesUpInfeasible_*(distanceToCutoff/(upCost+1.0e-12));
    419531  number += numberTimesUpInfeasible_;
    420532#endif
    421   if (number>0)
    422     upCost *= sum / (double) number;
     533#endif
     534  if (number>0.0)
     535    upCost *= sum / number;
    423536  else
    424537    upCost  *=  upDynamicPseudoCost_;
     
    439552    return 0.0;
    440553  } else {
    441     int stateOfSearch = model_->stateOfSearch();
     554    int stateOfSearch = model_->stateOfSearch()%10;
    442555    double returnValue=0.0;
    443556    double minValue = CoinMin(downCost,upCost);
    444557    double maxValue = CoinMax(downCost,upCost);
    445     if (stateOfSearch<=1||model_->currentNode()->depth()<=10/* was ||maxValue>0.2*distanceToCutoff*/) {
    446       // no solution
     558    char where;
     559    // was <= 10
     560    //if (stateOfSearch<=1||model_->currentNode()->depth()<=-10 /* was ||maxValue>0.2*distanceToCutoff*/) {
     561    if (stateOfSearch<=2) {
     562      // no branching solution
     563      where='i';
    447564      returnValue = WEIGHT_BEFORE*minValue + (1.0-WEIGHT_BEFORE)*maxValue;
    448565      if (0) {
     
    471588    } else {
    472589      // some solution
     590      where='I';
    473591      returnValue = WEIGHT_AFTER*minValue + (1.0-WEIGHT_AFTER)*maxValue;
    474592    }
     
    500618      returnValue *= 1.0e-3;
    501619    }
     620#ifdef COIN_DEVELOP
     621    History hist;
     622    hist.where_=where;
     623    hist.status_=' ';
     624    hist.sequence_=columnNumber_;
     625    hist.numberUp_=numberTimesUp_;
     626    hist.numberUpInf_=numberTimesUpInfeasible_;
     627    hist.sumUp_=sumUpCost_;
     628    hist.upEst_=upCost;
     629    hist.numberDown_=numberTimesDown_;
     630    hist.numberDownInf_=numberTimesDownInfeasible_;
     631    hist.sumDown_=sumDownCost_;
     632    hist.downEst_=downCost;
     633    if (stateOfSearch)
     634      addRecord(hist);
     635#endif
    502636    return CoinMax(returnValue,1.0e-15);
    503637  }
     
    741875  way = - way; // because after branch so moved on
    742876  double value = branchingObject->value();
    743   return CbcObjectUpdateData (this, way,
    744                                   change, iStatus,
    745                                   originalUnsatisfied-unsatisfied,value);
     877  CbcObjectUpdateData newData (this, way,
     878                               change, iStatus,
     879                               originalUnsatisfied-unsatisfied,value);
     880  newData.originalObjective_ = originalValue;
     881  // Solvers know about direction
     882  double direction = solver->getObjSense();
     883  solver->getDblParam(OsiDualObjectiveLimit,newData.cutoff_);
     884  newData.cutoff_ *= direction;
     885  return newData;
    746886}
    747887// Update object by CbcObjectUpdateData
     
    753893  double value = data.branchingValue_;
    754894  double change = data.change_;
    755 #define TYPE2 1
    756895#define TYPERATIO 0.9
     896#define MINIMUM_MOVEMENT 0.0
     897#ifdef COIN_DEVELOP
     898  History hist;
     899  hist.where_='U'; // need to tell if hot
     900#endif
     901  double movement=0.0;
    757902  if (way<0) {
    758903    // down
     904    movement = value-floor(value);
    759905    if (feasible) {
    760       //printf("(down change %g value down %g ",change,value-floor(value));
     906#ifdef COIN_DEVELOP
     907      hist.status_='D';
     908#endif
     909      movement = CoinMax(movement,MINIMUM_MOVEMENT);
     910      //printf("(down change %g value down %g ",change,movement);
    761911      incrementNumberTimesDown();
    762       addToSumDownChange(1.0e-30+value-floor(value));
     912      addToSumDownChange(1.0e-30+movement);
    763913      addToSumDownDecrease(data.intDecrease_);
    764914#if TYPE2==0
    765       addToSumDownCost(change/(1.0e-30+(value-floor(value))));
     915      addToSumDownCost(change/(1.0e-30+movement));
    766916      setDownDynamicPseudoCost(sumDownCost()/(double) numberTimesDown());
    767917#elif TYPE2==1
     
    769919      setDownDynamicPseudoCost(sumDownCost()/sumDownChange());
    770920#elif TYPE2==2
    771       addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(value-floor(value))));
     921      addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    772922      setDownDynamicPseudoCost(sumDownCost()*(TYPERATIO/sumDownChange()+(1.0-TYPERATIO)/(double) numberTimesDown()));
    773923#endif
    774924    } else {
    775       //printf("(down infeasible value down %g ",change,value-floor(value));
     925#ifdef COIN_DEVELOP
     926      hist.status_='d';
     927#endif
     928      //printf("(down infeasible value down %g ",change,movement);
    776929      incrementNumberTimesDownInfeasible();
    777930#if INFEAS==2
     
    782935        change = distanceToCutoff*2.0;
    783936      else
    784         change = downDynamicPseudoCost()*(value-floor(value))*10.0;
     937        change = downDynamicPseudoCost()*movement*10.0;
    785938      change = CoinMax(1.0e-12*(1.0+fabs(originalValue)),change);
    786939      incrementNumberTimesDown();
    787       addToSumDownChange(1.0e-30+value-floor(value));
     940      addToSumDownChange(1.0e-30+movement);
    788941      addToSumDownDecrease(data.intDecrease_);
    789942#if TYPE2==0
    790       addToSumDownCost(change/(1.0e-30+(value-floor(value))));
     943      addToSumDownCost(change/(1.0e-30+movement));
    791944      setDownDynamicPseudoCost(sumDownCost()/(double) numberTimesDown());
    792945#elif TYPE2==1
     
    794947      setDownDynamicPseudoCost(sumDownCost()/sumDownChange());
    795948#elif TYPE2==2
    796       addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(value-floor(value))));
     949      addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    797950      setDownDynamicPseudoCost(sumDownCost()*(TYPERATIO/sumDownChange()+(1.0-TYPERATIO)/(double) numberTimesDown()));
    798951#endif
    799952#endif
    800953    }
     954#if INFEAS==1
     955    double sum = sumDownCost_;
     956    int number = numberTimesDown_;
     957    double originalValue = data.originalObjective_;
     958    assert (originalValue!=COIN_DBL_MAX);
     959    double distanceToCutoff =  data.cutoff_  - originalValue;
     960    if (distanceToCutoff>1.0e20)
     961      distanceToCutoff=10.0+fabs(originalValue);
     962    sum += numberTimesDownInfeasible_*distanceToCutoff;
     963    number += numberTimesDownInfeasible_;
     964    setDownDynamicPseudoCost(sum/(double) number);
     965#endif
    801966  } else {
    802967    // up
     968    movement = ceil(value)-value;
    803969    if (feasible) {
    804       //printf("(up change %g value down %g ",change,ceil(value)-value);
     970#ifdef COIN_DEVELOP
     971      hist.status_='U';
     972#endif
     973      movement = CoinMax(movement,MINIMUM_MOVEMENT);
     974      //printf("(up change %g value down %g ",change,movement);
    805975      incrementNumberTimesUp();
    806       addToSumUpChange(1.0e-30+ceil(value)-value);
     976      addToSumUpChange(1.0e-30+movement);
    807977      addToSumUpDecrease(data.intDecrease_);
    808978#if TYPE2==0
    809       addToSumUpCost(change/(1.0e-30+(ceil(value)-value)));
     979      addToSumUpCost(change/(1.0e-30+movement));
    810980      setUpDynamicPseudoCost(sumUpCost()/(double) numberTimesUp());
    811981#elif TYPE2==1
     
    813983      setUpDynamicPseudoCost(sumUpCost()/sumUpChange());
    814984#elif TYPE2==2
    815       addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(ceil(value)-value)));
     985      addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    816986      setUpDynamicPseudoCost(sumUpCost()*(TYPERATIO/sumUpChange()+(1.0-TYPERATIO)/(double) numberTimesUp()));
    817987#endif
    818988    } else {
    819       //printf("(up infeasible value down %g ",change,ceil(value)-value);
     989#ifdef COIN_DEVELOP
     990      hist.status_='u';
     991#endif
     992      //printf("(up infeasible value down %g ",change,movement);
    820993      incrementNumberTimesUpInfeasible();
    821994#if INFEAS==2
     
    826999        change = distanceToCutoff*2.0;
    8271000      else
    828         change = upDynamicPseudoCost()*(ceil(value)-value)*10.0;
     1001        change = upDynamicPseudoCost()*movement*10.0;
    8291002      change = CoinMax(1.0e-12*(1.0+fabs(originalValue)),change);
    8301003      incrementNumberTimesUp();
    831       addToSumUpChange(1.0e-30+ceil(value)-value);
     1004      addToSumUpChange(1.0e-30+movement);
    8321005      addToSumUpDecrease(data.intDecrease_);
    8331006#if TYPE2==0
    834       addToSumUpCost(change/(1.0e-30+(ceil(value)-value)));
     1007      addToSumUpCost(change/(1.0e-30+movement));
    8351008      setUpDynamicPseudoCost(sumUpCost()/(double) numberTimesUp());
    8361009#elif TYPE2==1
     
    8381011      setUpDynamicPseudoCost(sumUpCost()/sumUpChange());
    8391012#elif TYPE2==2
    840       addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(ceil(value)-value)));
     1013      addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    8411014      setUpDynamicPseudoCost(sumUpCost()*(TYPERATIO/sumUpChange()+(1.0-TYPERATIO)/(double) numberTimesUp()));
    8421015#endif
    8431016#endif
    8441017    }
    845   }
     1018#if INFEAS==1
     1019    double sum = sumUpCost_;
     1020    int number = numberTimesUp_;
     1021    double originalValue = data.originalObjective_;
     1022    assert (originalValue!=COIN_DBL_MAX);
     1023    double distanceToCutoff =  data.cutoff_  - originalValue;
     1024    if (distanceToCutoff>1.0e20)
     1025      distanceToCutoff=10.0+fabs(originalValue);
     1026    sum += numberTimesUpInfeasible_*distanceToCutoff;
     1027    number += numberTimesUpInfeasible_;
     1028    setUpDynamicPseudoCost(sum/(double) number);
     1029#endif
     1030  }
     1031#ifdef COIN_DEVELOP
     1032  hist.sequence_=columnNumber_;
     1033  hist.numberUp_=numberTimesUp_;
     1034  hist.numberUpInf_=numberTimesUpInfeasible_;
     1035  hist.sumUp_=sumUpCost_;
     1036  hist.upEst_=change;
     1037  hist.numberDown_=numberTimesDown_;
     1038  hist.numberDownInf_=numberTimesDownInfeasible_;
     1039  hist.sumDown_=sumDownCost_;
     1040  hist.downEst_=movement;
     1041  addRecord(hist);
     1042#endif
    8461043  //print(1,0.5);
    8471044}
     
    10201217  assert (object_);
    10211218  assert (info.possibleBranch==this);
     1219    info.upMovement = object_->upDynamicPseudoCost()*(ceil(value_)-value_);
     1220    info.downMovement = object_->downDynamicPseudoCost()*(value_-floor(value_));
     1221    info.numIntInfeasUp  -= (int) (object_->sumUpDecrease()/
     1222                                   (1.0e-12+(double) object_->numberTimesUp()));
     1223    info.numObjInfeasUp = 0;
     1224    info.finishedUp = false;
     1225    info.numItersUp = 0;
     1226    info.numIntInfeasDown  -= (int) (object_->sumDownDecrease()/
     1227                                   (1.0e-12+(double) object_->numberTimesDown()));
     1228    info.numObjInfeasDown = 0;
     1229    info.finishedDown = false;
     1230    info.numItersDown = 0;
     1231    info.fix =0;
    10221232  if (object_->numberTimesUp()<object_->numberBeforeTrust()||
    10231233      object_->numberTimesDown()<object_->numberBeforeTrust()) {
    10241234    return 0;
    10251235  } else {
    1026     info.upMovement = object_->upDynamicPseudoCost()*(ceil(value_)-value_);
    1027     info.downMovement = object_->downDynamicPseudoCost()*(value_-floor(value_));
    1028     info.numIntInfeasUp  -= (int) (object_->sumUpDecrease()/
    1029                                    ((double) object_->numberTimesUp()));
    1030     info.numObjInfeasUp = 0;
    1031     info.finishedUp = false;
    1032     info.numItersUp = 0;
    1033     info.numIntInfeasDown  -= (int) (object_->sumDownDecrease()/
    1034                                    ((double) object_->numberTimesDown()));
    1035     info.numObjInfeasDown = 0;
    1036     info.finishedDown = false;
    1037     info.numItersDown = 0;
    1038     info.fix =0;
    10391236    return 1;
    10401237  }
     
    10871284  bestNumberDown_ = 0;
    10881285  bestObject_ = NULL;
     1286#ifdef COIN_DEVELOP
     1287  History hist;
     1288  hist.where_='C';
     1289  hist.status_='I';
     1290  hist.sequence_=55555;
     1291  hist.numberUp_=0;
     1292  hist.numberUpInf_=0;
     1293  hist.sumUp_=0.0;
     1294  hist.upEst_=0.0;
     1295  hist.numberDown_=0;
     1296  hist.numberDownInf_=0;
     1297  hist.sumDown_=0.0;
     1298  hist.downEst_=0.0;
     1299  addRecord(hist);
     1300#endif
    10891301}
    10901302
     
    11561368    // down
    11571369    if (feasible) {
    1158       //printf("(down change %g value down %g ",change,value-floor(value));
     1370      double movement = value-floor(value);
     1371      movement = CoinMax(movement,MINIMUM_MOVEMENT);
     1372      //printf("(down change %g value down %g ",change,movement);
    11591373      object->incrementNumberTimesDown();
    1160       object->addToSumDownChange(1.0e-30+value-floor(value));
     1374      object->addToSumDownChange(1.0e-30+movement);
    11611375      object->addToSumDownDecrease(originalUnsatisfied-unsatisfied);
    11621376#if TYPE2==0
    1163       object->addToSumDownCost(change/(1.0e-30+(value-floor(value))));
     1377      object->addToSumDownCost(change/(1.0e-30+movement));
    11641378      object->setDownDynamicPseudoCost(object->sumDownCost()/(double) object->numberTimesDown());
    11651379#elif TYPE2==1
     
    11671381      object->setDownDynamicPseudoCost(object->sumDownCost()/object->sumDownChange());
    11681382#elif TYPE2==2
    1169       object->addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(value-floor(value))));
     1383      object->addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    11701384      object->setDownDynamicPseudoCost(object->sumDownCost()*(TYPERATIO/object->sumDownChange()+(1.0-TYPERATIO)/(double) object->numberTimesDown()));
    11711385#endif
    11721386    } else {
    1173       //printf("(down infeasible value down %g ",change,value-floor(value));
     1387      //printf("(down infeasible value down %g ",change,movement);
    11741388      object->incrementNumberTimesDownInfeasible();
    11751389#if INFEAS==2
     
    11801394        change = distanceToCutoff*2.0;
    11811395      else
    1182         change = object->downDynamicPseudoCost()*(value-floor(value))*10.0;
     1396        change = object->downDynamicPseudoCost()*movement*10.0;
    11831397      change = CoinMax(1.0e-12*(1.0+fabs(originalValue)),change);
    11841398      object->incrementNumberTimesDown();
    1185       object->addToSumDownChange(1.0e-30+value-floor(value));
     1399      object->addToSumDownChange(1.0e-30+movement);
    11861400      object->addToSumDownDecrease(originalUnsatisfied-unsatisfied);
    11871401#if TYPE2==0
    1188       object->addToSumDownCost(change/(1.0e-30+(value-floor(value))));
     1402      object->addToSumDownCost(change/(1.0e-30+movement));
    11891403      object->setDownDynamicPseudoCost(object->sumDownCost()/(double) object->numberTimesDown());
    11901404#elif TYPE2==1
     
    11921406      object->setDownDynamicPseudoCost(object->sumDownCost()/object->sumDownChange());
    11931407#elif TYPE2==2
    1194       object->addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(value-floor(value))));
     1408      object->addToSumDownCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    11951409      object->setDownDynamicPseudoCost(object->sumDownCost()*(TYPERATIO/object->sumDownChange()+(1.0-TYPERATIO)/(double) object->numberTimesDown()));
    11961410#endif
     
    12001414    // up
    12011415    if (feasible) {
    1202       //printf("(up change %g value down %g ",change,ceil(value)-value);
     1416      double movement = ceil(value)-value;
     1417      movement = CoinMax(movement,MINIMUM_MOVEMENT);
     1418      //printf("(up change %g value down %g ",change,movement);
    12031419      object->incrementNumberTimesUp();
    1204       object->addToSumUpChange(1.0e-30+ceil(value)-value);
     1420      object->addToSumUpChange(1.0e-30+movement);
    12051421      object->addToSumUpDecrease(unsatisfied-originalUnsatisfied);
    12061422#if TYPE2==0
    1207       object->addToSumUpCost(change/(1.0e-30+(ceil(value)-value)));
     1423      object->addToSumUpCost(change/(1.0e-30+movement));
    12081424      object->setUpDynamicPseudoCost(object->sumUpCost()/(double) object->numberTimesUp());
    12091425#elif TYPE2==1
     
    12111427      object->setUpDynamicPseudoCost(object->sumUpCost()/object->sumUpChange());
    12121428#elif TYPE2==2
    1213       object->addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(ceil(value)-value)));
     1429      object->addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    12141430      object->setUpDynamicPseudoCost(object->sumUpCost()*(TYPERATIO/object->sumUpChange()+(1.0-TYPERATIO)/(double) object->numberTimesUp()));
    12151431#endif
    12161432    } else {
    1217       //printf("(up infeasible value down %g ",change,ceil(value)-value);
     1433      //printf("(up infeasible value down %g ",change,movement);
    12181434      object->incrementNumberTimesUpInfeasible();
    12191435#if INFEAS==2
     
    12241440        change = distanceToCutoff*2.0;
    12251441      else
    1226         change = object->upDynamicPseudoCost()*(ceil(value)-value)*10.0;
     1442        change = object->upDynamicPseudoCost()*movement*10.0;
    12271443      change = CoinMax(1.0e-12*(1.0+fabs(originalValue)),change);
    12281444      object->incrementNumberTimesUp();
    1229       object->addToSumUpChange(1.0e-30+ceil(value)-value);
     1445      object->addToSumUpChange(1.0e-30+movement);
    12301446      object->addToSumUpDecrease(unsatisfied-originalUnsatisfied);
    12311447#if TYPE2==0
    1232       object->addToSumUpCost(change/(1.0e-30+(ceil(value)-value)));
     1448      object->addToSumUpCost(change/(1.0e-30+movement));
    12331449      object->setUpDynamicPseudoCost(object->sumUpCost()/(double) object->numberTimesUp());
    12341450#elif TYPE2==1
     
    12361452      object->setUpDynamicPseudoCost(object->sumUpCost()/object->sumUpChange());
    12371453#elif TYPE2==2
    1238       object->addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+(ceil(value)-value)));
     1454      object->addToSumUpCost(change*TYPERATIO+(1.0-TYPERATIO)*change/(1.0e-30+movement));
    12391455      object->setUpDynamicPseudoCost(object->sumUpCost()*(TYPERATIO/object->sumUpChange()+(1.0-TYPERATIO)/(double) object->numberTimesUp()));
    12401456#endif
     
    12611477{
    12621478  CbcModel * model = thisOne->model();
    1263   int stateOfSearch = model->stateOfSearch();
     1479  int stateOfSearch = model->stateOfSearch()%10;
    12641480  int betterWay=0;
    12651481  double value=0.0;
     
    12691485    bestNumberDown_=COIN_INT_MAX;
    12701486  }
    1271   if (stateOfSearch<=1&&thisOne->model()->currentNode()->depth()>=8) {
     1487  if (stateOfSearch<=2) {
    12721488#define TRY_STUFF 1
    12731489#ifdef TRY_STUFF
     
    13651581    }
    13661582  }
     1583#ifdef COIN_DEVELOP
     1584  History hist;
     1585  {
     1586    CbcDynamicPseudoCostBranchingObject * branchingObject =
     1587      dynamic_cast<CbcDynamicPseudoCostBranchingObject *>(thisOne);
     1588    assert (branchingObject);
     1589    CbcSimpleIntegerDynamicPseudoCost *  object = branchingObject->object();
     1590    assert (object);
     1591    hist.where_='C';
     1592    hist.status_=' ';
     1593    hist.sequence_=object->columnNumber();
     1594    hist.numberUp_=object->numberTimesUp();
     1595    hist.numberUpInf_=numInfUp;
     1596    hist.sumUp_=object->sumUpCost();
     1597    hist.upEst_=changeUp;
     1598    hist.numberDown_=object->numberTimesDown();
     1599    hist.numberDownInf_=numInfDown;
     1600    hist.sumDown_=object->sumDownCost();
     1601    hist.downEst_=changeDown;
     1602  }
     1603#endif
    13671604  if (betterWay) {
     1605#ifdef COIN_DEVELOP
     1606    hist.status_='B';
     1607#endif
    13681608    // maybe change better way
    13691609    CbcDynamicPseudoCostBranchingObject * branchingObject =
     
    13881628      betterWay = thisOne->object()->preferredWay();
    13891629  }
     1630#ifdef COIN_DEVELOP
     1631  addRecord(hist);
     1632#endif
    13901633  return betterWay;
    13911634}
     
    14011644  return bestCriterion_;
    14021645}
     1646#ifdef COIN_DEVELOP
     1647void printHistory(const char * file)
     1648{
     1649  if (!numberHistory)
     1650    return;
     1651  FILE * fp = fopen(file,"w");
     1652  assert(fp);
     1653  unsigned short numberIntegers=0;
     1654  int i;
     1655  for (i=0;i<numberHistory;i++) {
     1656    if (history[i].where_!='C'||history[i].status_!='I')
     1657      numberIntegers = CoinMax(numberIntegers,history[i].sequence_);
     1658  }
     1659  numberIntegers++;
     1660  for (int iC=0;iC<numberIntegers;iC++) {
     1661    int n=0;
     1662    for (i=0;i<numberHistory;i++) {
     1663      if (history[i].sequence_==iC) {
     1664        if (!n)
     1665          fprintf(fp,"XXX %d\n",iC);
     1666        n++;
     1667        fprintf(fp,"%c%c up %8d %8d %12.5f %12.5f down  %8d %8d %12.5f %12.5f\n",
     1668               history[i].where_,
     1669               history[i].status_,
     1670               history[i].numberUp_,
     1671               history[i].numberUpInf_,
     1672               history[i].sumUp_,
     1673               history[i].upEst_,
     1674               history[i].numberDown_,
     1675               history[i].numberDownInf_,
     1676               history[i].sumDown_,
     1677               history[i].downEst_);
     1678      }
     1679    }
     1680  }
     1681  fclose(fp);
     1682}
     1683#endif
  • trunk/Cbc/src/CbcHeuristic.cpp

    r776 r833  
    1818#include "CbcMessage.hpp"
    1919#include "CbcHeuristic.hpp"
     20#include "CbcHeuristicFPump.hpp"
    2021#include "CbcStrategy.hpp"
    2122#include "CglPreProcess.hpp"
    2223#include "OsiAuxInfo.hpp"
     24#include "OsiPresolve.hpp"
    2325
    2426// Default Constructor
     
    2729   when_(2),
    2830   numberNodes_(200),
     31   feasibilityPumpOptions_(-1),
    2932   fractionSmall_(1.0),
    3033   heuristicName_("Unknown")
     
    3942  when_(2),
    4043  numberNodes_(200),
     44  feasibilityPumpOptions_(-1),
    4145  fractionSmall_(1.0),
    4246  heuristicName_("Unknown")
     
    5054  when_(rhs.when_),
    5155  numberNodes_(rhs.numberNodes_),
     56  feasibilityPumpOptions_(rhs.feasibilityPumpOptions_),
    5257  fractionSmall_(rhs.fractionSmall_),
    5358  heuristicName_(rhs.heuristicName_)
     
    6267    when_ = rhs.when_;
    6368    numberNodes_ = rhs.numberNodes_;
     69    feasibilityPumpOptions_ = rhs.feasibilityPumpOptions_;
    6470    fractionSmall_ = rhs.fractionSmall_;
    6571    heuristicName_ = rhs.heuristicName_ ;
     
    109115  model_ = model;
    110116}
     117#ifdef COIN_DEVELOP
     118extern bool getHistoryStatistics_;
     119#endif
    111120// Do mini branch and bound (return 1 if solution)
    112121int
     
    131140  }
    132141#endif
     142#ifdef COIN_DEVELOP
     143  getHistoryStatistics_=false;
     144#endif
     145  int logLevel = model_->logLevel();
     146  char generalPrint[100];
     147  // Do presolve to see if possible
     148  int numberColumns = solver->getNumCols();
     149  char * reset = NULL;
     150  int returnCode=1;
     151  {
     152    int saveLogLevel = solver->messageHandler()->logLevel();
     153    if (saveLogLevel==1)
     154      solver->messageHandler()->setLogLevel(0);
     155    OsiPresolve * pinfo = new OsiPresolve();
     156    int presolveActions=0;
     157    // Allow dual stuff on integers
     158    presolveActions=1;
     159    // Do not allow all +1 to be tampered with
     160    //if (allPlusOnes)
     161    //presolveActions |= 2;
     162    // allow transfer of costs
     163    // presolveActions |= 4;
     164    pinfo->setPresolveActions(presolveActions);
     165    OsiSolverInterface * presolvedModel = pinfo->presolvedModel(*solver,1.0e-8,true,2);
     166    delete pinfo;
     167    // see if too big
     168    double before = 2*solver->getNumRows()+solver->getNumCols();
     169    if (presolvedModel) {
     170      int afterRows = presolvedModel->getNumRows();
     171      int afterCols = presolvedModel->getNumCols();
     172      delete presolvedModel;
     173      double after = 2*afterRows+afterCols;
     174      if (after>fractionSmall_*before&&after>300) {
     175        // Need code to try again to compress further using used
     176        const int * used =  model_->usedInSolution();
     177        int maxUsed=0;
     178        int iColumn;
     179        const double * lower = solver->getColLower();
     180        const double * upper = solver->getColUpper();
     181        for (iColumn=0;iColumn<numberColumns;iColumn++) {
     182          if (upper[iColumn]>lower[iColumn]) {
     183            if (solver->isBinary(iColumn))
     184              maxUsed = CoinMax(maxUsed,used[iColumn]);
     185          }
     186        }
     187        if (maxUsed) {
     188          reset = new char [numberColumns];
     189          int nFix=0;
     190          for (iColumn=0;iColumn<numberColumns;iColumn++) {
     191            reset[iColumn]=0;
     192            if (upper[iColumn]>lower[iColumn]) {
     193              if (solver->isBinary(iColumn)&&used[iColumn]==maxUsed) {
     194                bool setValue=true;
     195                if (maxUsed==1) {
     196                  double randomNumber = CoinDrand48();
     197                  if (randomNumber>0.3)
     198                    setValue=false;
     199                }
     200                if (setValue) {
     201                  reset[iColumn]=1;
     202                  solver->setColLower(iColumn,1.0);
     203                  nFix++;
     204                }
     205              }
     206            }
     207          }
     208          pinfo = new OsiPresolve();
     209          presolveActions=0;
     210          // Allow dual stuff on integers
     211          presolveActions=1;
     212          // Do not allow all +1 to be tampered with
     213          //if (allPlusOnes)
     214          //presolveActions |= 2;
     215          // allow transfer of costs
     216          // presolveActions |= 4;
     217          pinfo->setPresolveActions(presolveActions);
     218          presolvedModel = pinfo->presolvedModel(*solver,1.0e-8,true,2);
     219          delete pinfo;
     220          if(presolvedModel) {
     221            // see if too big
     222            int afterRows2 = presolvedModel->getNumRows();
     223            int afterCols2 = presolvedModel->getNumCols();
     224            delete presolvedModel;
     225            double after = 2*afterRows2+afterCols2;
     226            if (after>fractionSmall_*before&&after>300) {
     227              sprintf(generalPrint,"Full problem %d rows %d columns, reduced to %d rows %d columns - %d fixed gives %d, %d - still too large",
     228                      solver->getNumRows(),solver->getNumCols(),
     229                      afterRows,afterCols,nFix,afterRows2,afterCols2);
     230            } else {
     231              sprintf(generalPrint,"Full problem %d rows %d columns, reduced to %d rows %d columns - %d fixed gives %d, %d - ok now",
     232                      solver->getNumRows(),solver->getNumCols(),
     233                      afterRows,afterCols,nFix,afterRows2,afterCols2);
     234            }
     235            model_->messageHandler()->message(CBC_FPUMP1,model_->messages())
     236              << generalPrint
     237              <<CoinMessageEol;
     238          } else {
     239            returnCode=-1; // infeasible
     240          }
     241        }
     242      }
     243    } else {
     244      returnCode=-1; // infeasible
     245    }
     246    solver->messageHandler()->setLogLevel(saveLogLevel);
     247  }
     248  if (returnCode==-1) {
     249    delete [] reset;
     250#ifdef COIN_DEVELOP
     251    getHistoryStatistics_=true;
     252#endif
     253    return returnCode;
     254  }
    133255  // Reduce printout
    134256  solver->setHintParam(OsiDoReducePrint,true,OsiHintTry);
     
    136258  solver->setDblParam(OsiDualObjectiveLimit,cutoff*solver->getObjSense());
    137259  solver->initialSolve();
    138   int returnCode=1;
    139   int logLevel = model_->logLevel();
    140260  if (solver->isProvenOptimal()) {
    141261    CglPreProcess process;
     
    151271    } else {
    152272      // see if too big
    153       double before = solver->getNumRows()+solver->getNumCols();
    154       double after = solver2->getNumRows()+solver2->getNumCols();
    155       char generalPrint[100];
     273      double before = 2*solver->getNumRows()+solver->getNumCols();
     274      double after = 2*solver2->getNumRows()+solver2->getNumCols();
    156275      if (after>fractionSmall_*before&&after>300) {
    157276        sprintf(generalPrint,"Full problem %d rows %d columns, reduced to %d rows %d columns - too large",
     
    161280          << generalPrint
    162281          <<CoinMessageEol;
    163         return -1;
     282        returnCode = -1;
    164283      } else {
    165284        sprintf(generalPrint,"Full problem %d rows %d columns, reduced to %d rows %d columns",
     
    170289          <<CoinMessageEol;
    171290      }
    172       solver2->resolve();
    173       CbcModel model(*solver2);
    174       if (logLevel<=1)
    175         model.setLogLevel(0);
    176       else
    177         model.setLogLevel(logLevel);
    178       model.setCutoff(cutoff);
    179       model.setMaximumNodes(numberNodes);
    180       model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
    181       // Lightweight
    182       CbcStrategyDefaultSubTree strategy(model_,true,5,1,0);
    183       model.setStrategy(strategy);
    184       model.solver()->setIntParam(OsiMaxNumIterationHotStart,10);
    185       // Do search
    186       if (logLevel>1)
    187         model_->messageHandler()->message(CBC_START_SUB,model_->messages())
    188           << name
    189           << model.getMaximumNodes()
    190           <<CoinMessageEol;
    191       // probably faster to use a basis to get integer solutions
    192       model.setSpecialOptions(2);
     291      if (returnCode==1) {
     292        solver2->resolve();
     293        CbcModel model(*solver2);
     294        if (logLevel<=1)
     295          model.setLogLevel(0);
     296        else
     297          model.setLogLevel(logLevel);
     298        if (feasibilityPumpOptions_>=0) {
     299          CbcHeuristicFPump heuristic4;
     300          int pumpTune=feasibilityPumpOptions_;
     301          if (pumpTune>0) {
     302            /*
     303              >=10000000 for using obj
     304              >=1000000 use as accumulate switch
     305              >=1000 use index+1 as number of large loops
     306              >=100 use 0.05 objvalue as increment
     307              >=10 use +0.1 objvalue for cutoff (add)
     308              1 == fix ints at bounds, 2 fix all integral ints, 3 and continuous at bounds
     309              4 and static continuous, 5 as 3 but no internal integers
     310              6 as 3 but all slack basis!
     311            */
     312            double value = solver2->getObjSense()*solver2->getObjValue();
     313            int w = pumpTune/10;
     314            int c = w % 10;
     315            w /= 10;
     316            int i = w % 10;
     317            w /= 10;
     318            int r = w;
     319            int accumulate = r/1000;
     320            r -= 1000*accumulate;
     321            if (accumulate>=10) {
     322              int which = accumulate/10;
     323              accumulate -= 10*which;
     324              which--;
     325              // weights and factors
     326              double weight[]={0.1,0.1,0.5,0.5,1.0,1.0,5.0,5.0};
     327              double factor[] = {0.1,0.5,0.1,0.5,0.1,0.5,0.1,0.5};
     328              heuristic4.setInitialWeight(weight[which]);
     329              heuristic4.setWeightFactor(factor[which]);
     330            }
     331            // fake cutoff
     332            if (c) {
     333              double cutoff;
     334              solver2->getDblParam(OsiDualObjectiveLimit,cutoff);
     335              cutoff = CoinMin(cutoff,value + 0.1*fabs(value)*c);
     336              heuristic4.setFakeCutoff(cutoff);
     337            }
     338            if (i||r) {
     339              // also set increment
     340              //double increment = (0.01*i+0.005)*(fabs(value)+1.0e-12);
     341              double increment = 0.0;
     342              heuristic4.setAbsoluteIncrement(increment);
     343              heuristic4.setAccumulate(accumulate);
     344              heuristic4.setMaximumRetries(r+1);
     345            }
     346            pumpTune = pumpTune%100;
     347            if (pumpTune==6)
     348              pumpTune =13;
     349            heuristic4.setWhen(pumpTune+10);
     350          }
     351          heuristic4.setHeuristicName("feasibility pump");
     352          model.addHeuristic(&heuristic4);
     353        }
     354        model.setCutoff(cutoff);
     355        model.setMaximumNodes(numberNodes);
     356        model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
     357        // Lightweight
     358        CbcStrategyDefaultSubTree strategy(model_,true,5,1,0);
     359        model.setStrategy(strategy);
     360        model.solver()->setIntParam(OsiMaxNumIterationHotStart,10);
     361        // Do search
     362        if (logLevel>1)
     363          model_->messageHandler()->message(CBC_START_SUB,model_->messages())
     364            << name
     365            << model.getMaximumNodes()
     366            <<CoinMessageEol;
     367        // probably faster to use a basis to get integer solutions
     368        model.setSpecialOptions(2);
    193369#ifdef CBC_THREAD
    194       if (model_->getNumberThreads()>0&&(model_->getThreadMode()&1)!=0) {
    195         // See if at root node
    196         bool atRoot = model_->getNodeCount()==0;
    197         int passNumber = model_->getCurrentPassNumber();
    198         if (atRoot&&passNumber==1)
    199           model.setNumberThreads(model_->getNumberThreads());
    200       }
    201 #endif
    202       model.setMaximumCutPassesAtRoot(CoinMin(20,model_->getMaximumCutPassesAtRoot()));
    203       model.setParentModel(*model_);
    204       model.setOriginalColumns(process.originalColumns());
    205       model.branchAndBound();
    206       if (logLevel>1)
    207         model_->messageHandler()->message(CBC_END_SUB,model_->messages())
    208           << name
    209           <<CoinMessageEol;
    210       if (model.getMinimizationObjValue()<CoinMin(cutoff,1.0e30)) {
    211         // solution
    212         returnCode=model.isProvenOptimal() ? 3 : 1;
    213         // post process
    214 #ifdef COIN_HAS_CLP
    215         OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model.solver());
    216         if (clpSolver) {
    217           ClpSimplex * lpSolver = clpSolver->getModelPtr();
    218           lpSolver->setSpecialOptions(lpSolver->specialOptions()|0x01000000); // say is Cbc (and in branch and bound)
     370        if (model_->getNumberThreads()>0&&(model_->getThreadMode()&1)!=0) {
     371          // See if at root node
     372          bool atRoot = model_->getNodeCount()==0;
     373          int passNumber = model_->getCurrentPassNumber();
     374          if (atRoot&&passNumber==1)
     375            model.setNumberThreads(model_->getNumberThreads());
    219376        }
    220377#endif
    221         process.postProcess(*model.solver());
    222         if (solver->isProvenOptimal()) {
    223           // Solution now back in solver
    224           int numberColumns = solver->getNumCols();
    225           memcpy(newSolution,solver->getColSolution(),
    226                  numberColumns*sizeof(double));
    227           newSolutionValue = model.getMinimizationObjValue();
    228         } else {
    229           // odd - but no good
    230           returnCode=0; // so will be infeasible
    231         }
    232       } else {
     378        model.setMaximumCutPassesAtRoot(CoinMin(20,model_->getMaximumCutPassesAtRoot()));
     379        model.setParentModel(*model_);
     380        model.setOriginalColumns(process.originalColumns());
     381        model.branchAndBound();
     382        if (logLevel>1)
     383          model_->messageHandler()->message(CBC_END_SUB,model_->messages())
     384            << name
     385            <<CoinMessageEol;
     386        if (model.getMinimizationObjValue()<CoinMin(cutoff,1.0e30)) {
     387          // solution
     388          returnCode=model.isProvenOptimal() ? 3 : 1;
     389          // post process
     390#ifdef COIN_HAS_CLP
     391          OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model.solver());
     392          if (clpSolver) {
     393            ClpSimplex * lpSolver = clpSolver->getModelPtr();
     394            lpSolver->setSpecialOptions(lpSolver->specialOptions()|0x01000000); // say is Cbc (and in branch and bound)
     395          }
     396#endif
     397          process.postProcess(*model.solver());
     398          if (solver->isProvenOptimal()) {
     399            // Solution now back in solver
     400            int numberColumns = solver->getNumCols();
     401            memcpy(newSolution,solver->getColSolution(),
     402                   numberColumns*sizeof(double));
     403            newSolutionValue = model.getMinimizationObjValue();
     404          } else {
     405            // odd - but no good
     406            returnCode=0; // so will be infeasible
     407          }
     408        } else {
    233409        // no good
    234         returnCode=model.isProvenInfeasible() ? 2 : 0; // so will be infeasible
    235       }
    236       if (model.status()==5)
    237         returnCode=-2; // stop
     410          returnCode=model.isProvenInfeasible() ? 2 : 0; // so will be infeasible
     411        }
     412        if (model.status()==5)
     413          returnCode=-2; // stop
     414      }
    238415    }
    239416  } else {
     
    241418  }
    242419  model_->setLogLevel(logLevel);
     420  if (reset) {
     421    for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     422      if (reset[iColumn])
     423        solver->setColLower(iColumn,0.0);
     424    }
     425    delete [] reset;
     426  }
     427#ifdef COIN_DEVELOP
     428  getHistoryStatistics_=true;
     429#endif
    243430  return returnCode;
    244431}
  • trunk/Cbc/src/CbcHeuristic.hpp

    r765 r833  
    7979  inline int numberNodes() const
    8080  { return numberNodes_;}
     81  /// Sets feasibility pump options (-1 is off)
     82  inline void setFeasibilityPumpOptions(int value)
     83  { feasibilityPumpOptions_=value;}
     84  /// Gets feasibility pump options (-1 is off)
     85  inline int feasibilityPumpOptions() const
     86  { return feasibilityPumpOptions_;}
    8187  /// Just set model - do not do anything else
    8288  inline void setModelOnly(CbcModel * model)
     
    123129  /// Number of nodes in any sub tree
    124130  int numberNodes_;
     131  /// Feasibility pump options (-1 is off)
     132  int feasibilityPumpOptions_;
    125133  /// Fraction of new(rows+columns)/old(rows+columns) before doing small branch and bound
    126134  double fractionSmall_;
  • trunk/Cbc/src/CbcHeuristicFPump.cpp

    r827 r833  
    3131   initialWeight_(0.0),
    3232   weightFactor_(0.1),
     33   artificialCost_(COIN_DBL_MAX),
    3334   maximumPasses_(100),
    3435   maximumRetries_(1),
     
    5253   initialWeight_(0.0),
    5354   weightFactor_(0.1),
     55   artificialCost_(COIN_DBL_MAX),
    5456   maximumPasses_(100),
    5557   maximumRetries_(1),
     
    139141  initialWeight_(rhs.initialWeight_),
    140142  weightFactor_(rhs.weightFactor_),
     143  artificialCost_(rhs.artificialCost_),
    141144  maximumPasses_(rhs.maximumPasses_),
    142145  maximumRetries_(rhs.maximumRetries_),
     
    161164    initialWeight_ = rhs.initialWeight_;
    162165    weightFactor_ = rhs.weightFactor_;
     166    artificialCost_ = rhs.artificialCost_;
    163167    maximumPasses_ = rhs.maximumPasses_;
    164168    maximumRetries_ = rhs.maximumRetries_;
     
    311315  int numberSolutions=0;
    312316  OsiSolverInterface * solver = NULL;
     317  double artificialFactor = 0.00001;
    313318  while (!exitAll) {
    314319    int numberPasses=0;
     320    artificialFactor *= 10.0;
    315321    numberTries++;
    316322    // Clone solver - otherwise annoys root node computations
     
    418424    // Get amount for original objective
    419425    double scaleFactor = 0.0;
    420     for (i=0;i<numberColumns;i++)
    421       scaleFactor += saveObjective[i]*saveObjective[i];
     426#ifdef COIN_DEVELOP
     427    double largestCost=0.0;
     428    int nArtificial=0;
     429#endif
     430    for (i=0;i<numberColumns;i++) {
     431      double value = saveObjective[i];
     432      scaleFactor += value*value;
     433#ifdef COIN_DEVELOP
     434      largestCost=CoinMax(largestCost,fabs(value));
     435      if (value*direction>=artificialCost_)
     436        nArtificial++;
     437#endif
     438    }
    422439    if (scaleFactor)
    423440      scaleFactor = (initialWeight_*sqrt((double) numberIntegers))/sqrt(scaleFactor);
     441#ifdef COIN_DEVELOP
     442    if (scaleFactor)
     443      printf("Using %g fraction of original objective - largest %g - %d artificials\n",scaleFactor,
     444             largestCost,nArtificial);
     445#endif
    424446    // 5. MAIN WHILE LOOP
    425447    bool newLineNeeded=false;
     
    664686          // below so we can keep original code and allow for objective
    665687          int iColumn = i;
     688          // Special code for "artificials"
     689          if (direction*saveObjective[iColumn]>=artificialCost_) {
     690            //solver->setObjCoeff(iColumn,scaleFactor*saveObjective[iColumn]);
     691            solver->setObjCoeff(iColumn,(artificialFactor*saveObjective[iColumn])/artificialCost_);
     692          }
    666693          if(!solver->isBinary(iColumn)&&!doGeneral)
    667694            continue;
    668695          // deal with fixed variables (i.e., upper=lower)
    669           if (fabs(lower[iColumn]-upper[iColumn]) < primalTolerance) {
     696          if (fabs(lower[iColumn]-upper[iColumn]) < primalTolerance||!solver->isInteger(iColumn)) {
    670697            //if (lower[iColumn] > 1. - primalTolerance) solver->setObjCoeff(iColumn,-costValue);
    671698            //else                                       solver->setObjCoeff(iColumn,costValue);
  • trunk/Cbc/src/CbcHeuristicFPump.hpp

    r765 r833  
    9898  inline double weightFactor() const
    9999  { return weightFactor_;}
     100  /// Set threshold cost for using original cost - even on continuous (default infinity)
     101  inline void setArtificialCost(double value)
     102  { artificialCost_ = value;}
     103  /// Get threshold cost for using original cost - even on continuous (default infinity)
     104  inline double artificialCost() const
     105  { return artificialCost_;}
    100106  /// Set maximum passes (default 100)
    101107  inline void setMaximumPasses(int value)
     
    150156  /// Initial weight for true objective
    151157  double initialWeight_;
    152   /// factor for decreasing weight
     158  /// Factor for decreasing weight
    153159  double weightFactor_;
     160  /// Threshold cost for using original cost - even on continuous
     161  double artificialCost_;
    154162  /// Maximum number of passes
    155163  int maximumPasses_;
  • trunk/Cbc/src/CbcHeuristicRINS.cpp

    r759 r833  
    217217      numberTries_++;
    218218      if ((numberTries_%10)==0&&numberSuccesses_*3<numberTries_)
    219         howOften_ += howOften_/10;
     219        howOften_ += howOften_/2;
    220220    }
    221221
     
    235235  memset(used_,0,numberColumns);
    236236}
    237 
    238  
     237// Default Constructor
     238CbcHeuristicRENS::CbcHeuristicRENS()
     239  :CbcHeuristic()
     240{
     241  numberTries_=0;
     242}
     243
     244// Constructor with model - assumed before cuts
     245
     246CbcHeuristicRENS::CbcHeuristicRENS(CbcModel & model)
     247  :CbcHeuristic(model)
     248{
     249  numberTries_=0;
     250}
     251
     252// Destructor
     253CbcHeuristicRENS::~CbcHeuristicRENS ()
     254{
     255}
     256
     257// Clone
     258CbcHeuristic *
     259CbcHeuristicRENS::clone() const
     260{
     261  return new CbcHeuristicRENS(*this);
     262}
     263
     264// Assignment operator
     265CbcHeuristicRENS &
     266CbcHeuristicRENS::operator=( const CbcHeuristicRENS& rhs)
     267{
     268  if (this!=&rhs) {
     269    CbcHeuristic::operator=(rhs);
     270    numberTries_ = rhs.numberTries_;
     271  }
     272  return *this;
     273}
     274
     275// Copy constructor
     276CbcHeuristicRENS::CbcHeuristicRENS(const CbcHeuristicRENS & rhs)
     277:
     278  CbcHeuristic(rhs),
     279  numberTries_(rhs.numberTries_)
     280{
     281}
     282// Resets stuff if model changes
     283void
     284CbcHeuristicRENS::resetModel(CbcModel * model)
     285{
     286}
     287int
     288CbcHeuristicRENS::solution(double & solutionValue,
     289                         double * betterSolution)
     290{
     291  int returnCode=0;
     292  if (numberTries_)
     293    return 0;
     294  numberTries_++;
     295  OsiSolverInterface * solver = model_->solver();
     296 
     297  int numberIntegers = model_->numberIntegers();
     298  const int * integerVariable = model_->integerVariable();
     299 
     300  const double * currentSolution = solver->getColSolution();
     301  OsiSolverInterface * newSolver = model_->continuousSolver()->clone();
     302  const double * colLower = newSolver->getColLower();
     303  const double * colUpper = newSolver->getColUpper();
     304
     305  double primalTolerance;
     306  solver->getDblParam(OsiPrimalTolerance,primalTolerance);
     307   
     308  int i;
     309  int numberFixed=0;
     310  int numberTightened=0;
     311
     312  for (i=0;i<numberIntegers;i++) {
     313    int iColumn=integerVariable[i];
     314    double value = currentSolution[iColumn];
     315    double lower = colLower[iColumn];
     316    double upper = colUpper[iColumn];
     317    value = CoinMax(value,lower);
     318    value = CoinMin(value,upper);
     319    if (fabs(value-floor(value+0.5))<1.0e-8) {
     320      value = floor(value+0.5);
     321      newSolver->setColLower(iColumn,value);
     322      newSolver->setColUpper(iColumn,value);
     323      numberFixed++;
     324    } else if (colUpper[iColumn]-colLower[iColumn]>=2.0) {
     325      numberTightened++;
     326      newSolver->setColLower(iColumn,floor(value));
     327      newSolver->setColUpper(iColumn,ceil(value));
     328    }
     329  }
     330  if (numberFixed>numberIntegers/5) {
     331#ifdef COIN_DEVELOP
     332    printf("%d integers fixed and %d tightened\n",numberFixed,numberTightened);
     333#endif
     334    returnCode = smallBranchAndBound(newSolver,numberNodes_,betterSolution,solutionValue,
     335                                     model_->getCutoff(),"CbcHeuristicRENS");
     336    if (returnCode<0)
     337      returnCode=0; // returned on size
     338    //printf("return code %d",returnCode);
     339    if ((returnCode&2)!=0) {
     340      // could add cut
     341      returnCode &= ~2;
     342      //printf("could add cut with %d elements (if all 0-1)\n",nFix);
     343    } else {
     344      //printf("\n");
     345    }
     346  }
     347 
     348  delete newSolver;
     349  return returnCode;
     350}
     351// update model
     352void CbcHeuristicRENS::setModel(CbcModel * model)
     353{
     354  model_ = model;
     355}
     356
     357 
  • trunk/Cbc/src/CbcHeuristicRINS.hpp

    r765 r833  
    7575};
    7676
     77/** LocalSearch class
     78 */
     79
     80class CbcHeuristicRENS : public CbcHeuristic {
     81public:
     82
     83  // Default Constructor
     84  CbcHeuristicRENS ();
     85
     86  /* Constructor with model - assumed before cuts
     87     Initial version does not do Lps
     88  */
     89  CbcHeuristicRENS (CbcModel & model);
     90 
     91  // Copy constructor
     92  CbcHeuristicRENS ( const CbcHeuristicRENS &);
     93   
     94  // Destructor
     95  ~CbcHeuristicRENS ();
     96 
     97  /// Clone
     98  virtual CbcHeuristic * clone() const;
     99
     100
     101  /// Assignment operator
     102  CbcHeuristicRENS & operator=(const CbcHeuristicRENS& rhs);
     103
     104  /// Resets stuff if model changes
     105  virtual void resetModel(CbcModel * model);
     106
     107  /// update model (This is needed if cliques update matrix etc)
     108  virtual void setModel(CbcModel * model);
     109 
     110  using CbcHeuristic::solution ;
     111  /** returns 0 if no solution, 1 if valid solution.
     112      Sets solution values if good, sets objective value (only if good)
     113      This does Relaxation Extension Neighborhood Search
     114  */
     115  virtual int solution(double & objectiveValue,
     116                       double * newSolution);
     117
     118protected:
     119  // Data
     120  /// Number of tries
     121  int numberTries_;
     122};
     123
    77124
    78125#endif
  • trunk/Cbc/src/CbcModel.cpp

    r822 r833  
    6161#include "CglDuplicateRow.hpp"
    6262#include "CglStored.hpp"
     63#include "CglClique.hpp"
    6364
    6465#include "CoinTime.hpp"
     
    603604        printf("Problem type is %d\n",problemType_);
    604605#endif
     606      }
     607    }
     608    // But try again
     609    if (continuousMultiplier<1.0) {
     610      memset(rhs,0,numberRows*sizeof(double));
     611      int * count = new int [numberRows];
     612      memset(count,0,numberRows*sizeof(int));
     613      for (iColumn=0;iColumn<numberColumns;iColumn++) {
     614        CoinBigIndex start = columnStart[iColumn];
     615        CoinBigIndex end = start + columnLength[iColumn];
     616        if (upper[iColumn]==lower[iColumn]) {
     617          for (CoinBigIndex j=start;j<end;j++) {
     618            int iRow = row[j];
     619            rhs[iRow] += lower[iColumn]*element[j];
     620          }
     621        } else if (solver_->isInteger(iColumn)) {
     622          for (CoinBigIndex j=start;j<end;j++) {
     623            int iRow = row[j];
     624            if (fabs(element[j]-floor(element[j]+0.5))>1.0e-10)
     625              rhs[iRow]  = COIN_DBL_MAX;
     626          }
     627        } else {
     628          for (CoinBigIndex j=start;j<end;j++) {
     629            int iRow = row[j];
     630            count[iRow]++;
     631            if (fabs(element[j])!=1.0)
     632              rhs[iRow]  = COIN_DBL_MAX;
     633          }
     634        }
     635      }
     636      // now look at continuous
     637      bool allGood=true;
     638      double direction = solver_->getObjSense() ;
     639      int numberObj=0;
     640      for (iColumn=0;iColumn<numberColumns;iColumn++) {
     641        if (upper[iColumn]>lower[iColumn]) {
     642          double objValue = objective[iColumn]*direction;
     643          if (objValue&&!solver_->isInteger(iColumn)) {
     644            numberObj++;
     645            CoinBigIndex start = columnStart[iColumn];
     646            CoinBigIndex end = start + columnLength[iColumn];
     647            if (objValue>0.0) {
     648              // wants to be as low as possible
     649              if (lower[iColumn]<-1.0e10||fabs(lower[iColumn]-floor(lower[iColumn]+0.5))>1.0e-10) {
     650                allGood=false;
     651                break;
     652              } else if (upper[iColumn]<1.0e10&&fabs(upper[iColumn]-floor(upper[iColumn]+0.5))>1.0e-10) {
     653                allGood=false;
     654                break;
     655              }
     656              bool singletonRow=true;
     657              bool equality=false;
     658              for (CoinBigIndex j=start;j<end;j++) {
     659                int iRow = row[j];
     660                if (count[iRow]>1)
     661                  singletonRow=false;
     662                else if (rowLower[iRow]==rowUpper[iRow])
     663                  equality=true;
     664                if (fabs(rhs[iRow])>1.0e20||fabs(rhs[iRow]-floor(rhs[iRow]+0.5))>1.0e-10
     665                    ||fabs(element[j])!=1.0) {
     666                  // no good
     667                  allGood=false;
     668                  break;
     669                }
     670                if (element[j]>0.0) {
     671                  if (rowLower[iRow]>-1.0e20&&fabs(rowLower[iRow]-floor(rowLower[iRow]+0.5))>1.0e-10) {
     672                    // no good
     673                    allGood=false;
     674                    break;
     675                  }
     676                } else {
     677                  if (rowUpper[iRow]<1.0e20&&fabs(rowUpper[iRow]-floor(rowUpper[iRow]+0.5))>1.0e-10) {
     678                    // no good
     679                    allGood=false;
     680                    break;
     681                  }
     682                }
     683              }
     684              if (!singletonRow&&end>start+1&&!equality)
     685                allGood=false;
     686              if (!allGood)
     687                break;
     688            } else {
     689              // wants to be as high as possible
     690              if (upper[iColumn]>1.0e10||fabs(upper[iColumn]-floor(upper[iColumn]+0.5))>1.0e-10) {
     691                allGood=false;
     692                break;
     693              } else if (lower[iColumn]>-1.0e10&&fabs(lower[iColumn]-floor(lower[iColumn]+0.5))>1.0e-10) {
     694                allGood=false;
     695                break;
     696              }
     697              bool singletonRow=true;
     698              bool equality=false;
     699              for (CoinBigIndex j=start;j<end;j++) {
     700                int iRow = row[j];
     701                if (count[iRow]>1)
     702                  singletonRow=false;
     703                else if (rowLower[iRow]==rowUpper[iRow])
     704                  equality=true;
     705                if (fabs(rhs[iRow])>1.0e20||fabs(rhs[iRow]-floor(rhs[iRow]+0.5))>1.0e-10
     706                    ||fabs(element[j])!=1.0) {
     707                  // no good
     708                  allGood=false;
     709                  break;
     710                }
     711                if (element[j]<0.0) {
     712                  if (rowLower[iRow]>-1.0e20&&fabs(rowLower[iRow]-floor(rowLower[iRow]+0.5))>1.0e-10) {
     713                    // no good
     714                    allGood=false;
     715                    break;
     716                  }
     717                } else {
     718                  if (rowUpper[iRow]<1.0e20&&fabs(rowUpper[iRow]-floor(rowUpper[iRow]+0.5))>1.0e-10) {
     719                    // no good
     720                    allGood=false;
     721                    break;
     722                  }
     723                }
     724              }
     725              if (!singletonRow&&end>start+1&&!equality)
     726                allGood=false;
     727              if (!allGood)
     728                break;
     729            }
     730          }
     731        }
     732      }
     733      if (allGood) {
     734        if (numberObj)
     735          printf("YYYY analysis says all continuous with costs will be integer\n");
     736        continuousMultiplier=1.0;
    605737      }
    606738    }
     
    9211053  if (eventHandler)
    9221054    eventHandler->setModel(this);
     1055#ifdef CLIQUE_ANALYSIS
    9231056  // set up for probing
    924   //probingInfo_ = new CglTreeProbingInfo(solver_);
     1057  probingInfo_ = new CglTreeProbingInfo(solver_);
     1058#else
    9251059  probingInfo_=NULL;
     1060#endif
    9261061
    9271062  // Try for dominated columns
     
    14021537  //solverCharacteristics_->setSolver(solver_);
    14031538  if (feasible) {
     1539    if (probingInfo_) {
     1540      int number01 = probingInfo_->numberIntegers();
     1541      //const fixEntry * entry = probingInfo_->fixEntries();
     1542      const int * toZero = probingInfo_->toZero();
     1543      //const int * toOne = probingInfo_->toOne();
     1544      //const int * integerVariable = probingInfo_->integerVariable();
     1545      if (toZero[number01]) {
     1546        for (int i = 0;i<numberCutGenerators_;i++) {
     1547          CglFakeClique * clique = dynamic_cast<CglFakeClique*>(generator_[i]->generator());
     1548          if (clique) {
     1549            OsiSolverInterface * fakeSolver = probingInfo_->analyze(*solver_,1);
     1550            if (fakeSolver) {
     1551              printf("Probing fake solver has %d rows\n",fakeSolver->getNumRows());
     1552              //if (fakeSolver)
     1553              //fakeSolver->writeMps("bad");
     1554              if (generator_[i]->numberCutsInTotal())
     1555                generator_[i]->setHowOften(1);
     1556            }
     1557            clique->assignSolver(fakeSolver);
     1558            //stored->setProbingInfo(probingInfo_);
     1559            break;
     1560          }
     1561        }
     1562      }
     1563      delete probingInfo_;
     1564      probingInfo_=NULL;
     1565    }
    14041566    newNode = new CbcNode ;
    14051567    // Set objective value (not so obvious if NLP etc)
     
    19002062      }
    19012063    }
    1902     // If no solution but many nodes - signal change in strategy
    1903     if (numberNodes_>2*numberObjects_+1000&&stateOfSearch_!=2)
    1904       stateOfSearch_=3;
    19052064    // See if can stop on gap
    19062065    double testGap = CoinMax(dblParam_[CbcAllowableGap],
     
    21992358            OsiCuts theseCuts;
    22002359            // reset probing info
    2201             if (probingInfo_)
    2202               probingInfo_->initializeFixing();
     2360            //if (probingInfo_)
     2361            //probingInfo_->initializeFixing();
    22032362            for (int i = 0;i<numberCutGenerators_;i++) {
    22042363              bool generate = generator_[i]->normal();
     
    60726231      generator_[i]->incrementNumberCutsInTotal(countRowCuts[i]);
    60736232      generator_[i]->incrementNumberCutsActive(count[i]);
     6233      CglStored * stored = dynamic_cast<CglStored*>(generator_[i]->generator());
     6234      if (stored&&!countRowCuts[i])
     6235        continue;
    60746236      if (handler_->logLevel()>1||!numberNodes_) {
    60756237        handler_->message(CBC_GENERATOR,messages_)
     
    77597921        int lastNumberCuts=0;
    77607922        // reset probing info
    7761         if (probingInfo_)
    7762           probingInfo_->initializeFixing();
     7923        //if (probingInfo_)
     7924        //probingInfo_->initializeFixing();
    77637925        for (i=0;i<numberCutGenerators_;i++)
    77647926          {
     
    79098071        numberHeuristicSolutions_++;
    79108072      numberSolutions_++;
    7911       if (numberHeuristicSolutions_==numberSolutions_)
    7912         stateOfSearch_ = 1;
    7913       else
    7914         stateOfSearch_ = 2;
    79158073
    79168074      if (how!=CBC_ROUNDING) {
     
    79398097      int lastNumberCuts=0;
    79408098      // reset probing info
    7941       if (probingInfo_)
    7942         probingInfo_->initializeFixing();
     8099      //if (probingInfo_)
     8100      //probingInfo_->initializeFixing();
    79438101      for (i=0;i<numberCutGenerators_;i++) {
    79448102        bool generate = generator_[i]->atSolution();
     
    80798237     
    80808238        numberSolutions_++;
    8081         if (numberNodes_== 0 || numberHeuristicSolutions_==numberSolutions_)
    8082           stateOfSearch_ = 1;
    8083         else
    8084           stateOfSearch_ = 2;
    80858239       
    80868240        if (how!=CBC_ROUNDING) {
     
    82958449  int iGen;
    82968450  // reset probing info
    8297   if (probingInfo_)
    8298     probingInfo_->initializeFixing();
     8451  //if (probingInfo_)
     8452  //probingInfo_->initializeFixing();
    82998453  for (iGen=0;iGen<numberCutGenerators_;iGen++) {
    83008454    generator = dynamic_cast<CglProbing*>(generator_[iGen]->generator());
     
    97539907                       OsiSolverBranch * & branches)
    97549908{
     9909  // Set state of search
     9910  /*
     9911    0 - outside CbcNode
     9912    1 - no solutions
     9913    2 - all heuristic solutions
     9914    3 - a solution reached by branching (could be strong)
     9915    4 - no solution but many nodes
     9916       add 10 if depth >= K
     9917  */
     9918  stateOfSearch_=1;
     9919  if (numberSolutions_>0) {
     9920    if (numberHeuristicSolutions_==numberSolutions_)
     9921      stateOfSearch_ = 3;
     9922    else
     9923      stateOfSearch_ = 3;
     9924  } if (numberNodes_>2*numberObjects_+1000) {
     9925    stateOfSearch_=4;
     9926  }
     9927  //stateOfSearch_=3;
     9928  if (currentNode_&&currentNode_->depth()>=8)
     9929    stateOfSearch_ +=10;
    97559930  int anyAction =-1 ;
    97569931  resolved = false ;
     
    990110076    }
    990210077  }
     10078  stateOfSearch_ =0; // outside chooseBranch
    990310079  return anyAction;
    990410080}
     
    1076410940        OsiCuts theseCuts;
    1076510941        // reset probing info
    10766         if (probingInfo_)
    10767           probingInfo_->initializeFixing();
     10942        //if (probingInfo_)
     10943        //probingInfo_->initializeFixing(solver_);
    1076810944        for (int i = 0;i<numberCutGenerators_;i++) {
    1076910945          bool generate = generator_[i]->normal();
  • trunk/Cbc/src/CbcNode.cpp

    r822 r833  
    22022202  if (model->hotstartSolution())
    22032203    return -3;
    2204 #define RANGING
     2204  //#define RANGING
    22052205#ifdef RANGING
    22062206  // Pass number
     
    24232423        toOne = probingInfo->toOne();
    24242424        backward = probingInfo->backward();
    2425         if (!toZero[number01]||number01<numberObjects) {
     2425        if (!toZero[number01]||number01<numberObjects||true) {
    24262426          // no info
    24272427          probingInfo=NULL;
     
    29282928          skipAll=false;
    29292929        }
    2930         model->setStateOfSearch(2); // use min min
     2930        //model->setStateOfSearch(2); // use min min
    29312931      }
    29322932      // could adjust using average iterations per branch
     
    30823082        // see if can skip strong branching
    30833083        int canSkip = choice.possibleBranch->fillStrongInfo(choice);
     3084#if 0
    30843085        if (!newWay) {
    30853086          if ((maxChange>distanceToCutoff2)&&(!doQuickly||(numberTest>0&&searchStrategy!=2)))
     
    31043105          }
    31053106        }
     3107#else
     3108        if (((numberTest2<=0&&numberTest<=0)||skipAll)&&sort[iDo]>distanceToCutoff) {
     3109          //canSkip=1; // always skip
     3110          if (iDo>20) {
     3111            delete choice.possibleBranch;
     3112            choice.possibleBranch=NULL;
     3113            break; // give up anyway
     3114          }
     3115        }
     3116#endif
    31063117        if (model->messageHandler()->logLevel()>3&&numberBeforeTrust)
    31073118          dynamicObject->print(1,choice.possibleBranch->value());
     
    41264137    int jInt = back[iColumn];
    41274138    newLower[jInt]=upperValue;
    4128     if (choice.finishedDown||!fastIterations)
     4139    if (choice.finishedDown)
    41294140      objLower[jInt]=choice.downMovement+objectiveValue_;
    41304141    else
    41314142      objLower[jInt]=objectiveValue_;
    41324143    newUpper[jInt]=lowerValue;
    4133     if (choice.finishedUp||!fastIterations)
     4144    if (choice.finishedUp)
    41344145      objUpper[jInt]=choice.upMovement+objectiveValue_;
    41354146    else
  • trunk/Cbc/src/CbcSolver.cpp

    r824 r833  
    472472  parameters_[whichParam(COMBINE,numberParameters_,parameters_)].setCurrentOption("on");
    473473  parameters_[whichParam(RINS,numberParameters_,parameters_)].setCurrentOption("off");
     474  parameters_[whichParam(RENS,numberParameters_,parameters_)].setCurrentOption("off");
    474475  parameters_[whichParam(LOCALTREE,numberParameters_,parameters_)].setCurrentOption("off");
    475476  parameters_[whichParam(COSTSTRATEGY,numberParameters_,parameters_)].setCurrentOption("off");
     
    730731  return new CbcStopNow(*this);
    731732}
    732 //#undef NEW_STYLE_SOLVER
     733//#undef NEW_STYLE_SOLVER 
    733734//#define NEW_STYLE_SOLVER
    734735//#undef COIN_HAS_ASL
     
    11481149             1 return fixed non-presolved solver
    11491150             2 as one but use presolve Inside this
    1150              3 do heuristics and set best solution
    1151              4 do BAB and just set best solution
     1151             3 use presolve and fix ones with large cost
     1152             ? do heuristics and set best solution
     1153             ? do BAB and just set best solution
    11521154             10+ then use lastSolution and relax a few
    11531155             -2 cleanup afterwards if using 2
     
    11601162  if (doAction==11&&!lastSolution)
    11611163    lastSolution = model.bestSolution();
    1162   assert (((doAction==1||doAction==2)&&!lastSolution)||(doAction==11&&lastSolution));
     1164  assert (((doAction>=1&&doAction<=3)&&!lastSolution)||(doAction==11&&lastSolution));
    11631165  double fractionIntFixed = dextra[3];
    11641166  double fractionFixed = dextra[4];
     1167  double fixAbove = dextra[2];
    11651168  double time1 = CoinCpuTime();
    11661169  OsiSolverInterface * originalSolver = model.solver();
     
    11711174  ClpSimplex * lpSolver;
    11721175  ClpPresolve pinfo;
    1173   if (doAction==2) {
    1174     doAction=1;
    1175     lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10);
    1176     assert (lpSolver);
     1176  assert(originalSolver->getObjSense()>0);
     1177  if (doAction==2||doAction==3) {
     1178    double * saveLB=NULL;
     1179    double * saveUB=NULL;
     1180    int numberColumns = originalLpSolver->numberColumns();
     1181    if (fixAbove>0.0) {
     1182      double time1 = CoinCpuTime();
     1183      originalClpSolver->initialSolve();
     1184      printf("first solve took %g seconds\n",CoinCpuTime()-time1);
     1185      double * columnLower = originalLpSolver->columnLower() ;
     1186      double * columnUpper = originalLpSolver->columnUpper() ;
     1187      const double * solution = originalLpSolver->primalColumnSolution();
     1188      saveLB = CoinCopyOfArray(columnLower,numberColumns);
     1189      saveUB = CoinCopyOfArray(columnUpper,numberColumns);
     1190      const double * objective = originalLpSolver->getObjCoefficients() ;
     1191      int iColumn;
     1192      int nFix=0;
     1193      int nArt=0;
     1194      for (iColumn=0;iColumn<numberColumns;iColumn++) {
     1195        if (objective[iColumn]>fixAbove) {
     1196          if (solution[iColumn]<columnLower[iColumn]+1.0e-8) {
     1197            columnUpper[iColumn]=columnLower[iColumn];
     1198            nFix++;
     1199          } else {
     1200            nArt++;
     1201          }
     1202        } else if (objective[iColumn]<-fixAbove) {
     1203          if (solution[iColumn]>columnUpper[iColumn]-1.0e-8) {
     1204            columnLower[iColumn]=columnUpper[iColumn];
     1205            nFix++;
     1206          } else {
     1207            nArt++;
     1208          }
     1209        }
     1210      }
     1211      printf("%d artificials fixed, %d left as in solution\n",nFix,nArt);
     1212      lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10);
     1213      if (!lpSolver||doAction==2) {
     1214        // take off fixing in original
     1215        memcpy(columnLower,saveLB,numberColumns*sizeof(double));
     1216        memcpy(columnUpper,saveUB,numberColumns*sizeof(double));
     1217      }
     1218      delete [] saveLB;
     1219      delete [] saveUB;
     1220      if (!lpSolver) {
     1221        // try again
     1222        pinfo.destroyPresolve();
     1223        lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10);
     1224        assert (lpSolver);
     1225      }
     1226    } else {
     1227      lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10);
     1228      assert (lpSolver);
     1229    }
    11771230    clpSolver = new OsiClpSolverInterface(lpSolver,true);
    11781231    assert(lpSolver == clpSolver->getModelPtr());
    1179     //int numberColumns = lpSolver->numberColumns();
    1180     int numberColumns = lpSolver->numberColumns();
     1232    numberColumns = lpSolver->numberColumns();
    11811233    originalColumns = CoinCopyOfArray(pinfo.originalColumns(),numberColumns);
     1234    doAction=1;
    11821235  } else {
    11831236    OsiSolverInterface * solver = originalSolver->clone();
     
    31013154  parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption("on");
    31023155  parameters[whichParam(RINS,numberParameters,parameters)].setCurrentOption("off");
     3156  parameters[whichParam(RENS,numberParameters,parameters)].setCurrentOption("off");
    31033157  parameters[whichParam(LOCALTREE,numberParameters,parameters)].setCurrentOption("off");
    31043158  parameters[whichParam(COSTSTRATEGY,numberParameters,parameters)].setCurrentOption("off");
     
    31273181  int useCombine = parameters_[whichParam(COMBINE,numberParameters_,parameters_)].currentOptionAsInteger();
    31283182  int useRINS = parameters_[whichParam(RINS,numberParameters_,parameters_)].currentOptionAsInteger();
     3183  int useRENS = parameters_[whichParam(RENS,numberParameters_,parameters_)].currentOptionAsInteger();
    31293184  // FPump done first as it only works if no solution
    31303185  int kType = (type<3) ? type : 1;
     
    31363191    if (dextra3)
    31373192      heuristic4.setFractionSmall(dextra3);
     3193    double dextra1 = parameters_[whichParam(DEXTRA1,numberParameters_,parameters_)].doubleValue();
     3194    if (dextra1)
     3195      heuristic4.setArtificialCost(dextra1);
    31383196    heuristic4.setMaximumPasses(parameters_[whichParam(FPUMPITS,numberParameters_,parameters_)].intValue());
    31393197    int pumpTune=parameters_[whichParam(FPUMPTUNE,numberParameters_,parameters_)].intValue();
     
    32533311    anyToDo=true;
    32543312  }
     3313  if (useRENS>=kType) {
     3314    CbcHeuristicRENS heuristic6(*model);
     3315    heuristic6.setHeuristicName("RENS");
     3316    heuristic6.setFractionSmall(0.5);
     3317    heuristic6.setFeasibilityPumpOptions(1008003);
     3318    int nodes []={-2,-1,200,1000};
     3319    heuristic6.setNumberNodes(nodes[useRENS]);
     3320    model->addHeuristic(&heuristic6) ;
     3321    anyToDo=true;
     3322  }
    32553323  if (type==2&&anyToDo) {
    32563324    // Do heuristics
     
    34543522    int cutPass=-1234567;
    34553523    int cutPassInTree=-1234567;
    3456     int tunePreProcess=5;
     3524    int tunePreProcess=0;
    34573525    int testOsiParameters=-1;
    34583526    // 0 normal, 1 from ampl or MIQP etc (2 allows cuts)
     
    38053873    int redsplitAction=0;
    38063874
    3807     CglClique cliqueGen(false,true);
     3875    CglFakeClique cliqueGen(NULL,false);
     3876    //CglClique cliqueGen(false,true);
    38083877    cliqueGen.setStarCliqueReport(false);
    38093878    cliqueGen.setRowCliqueReport(false);
     
    43784447            case RINS:
    43794448              break;
     4449            case RENS:
     4450              break;
    43804451            case CUTSSTRATEGY:
    43814452              gomoryAction = action;
     
    49535024                // Just ones which affect >= extra3
    49545025                int extra3 = parameters_[whichParam(EXTRA3,numberParameters_,parameters_)].intValue();
    4955                 /* 3 is fraction of integer variables fixed (0.97)
     5026                /* 2 is cost above which to fix if feasible
     5027                   3 is fraction of integer variables fixed (0.97)
    49565028                   4 is fraction of all variables fixed (0.0)
    49575029                */
     
    53185390                    <<CoinMessageEol;
    53195391                }
    5320                 if (!complicatedInteger&&clpSolver->tightenPrimalBounds()!=0) {
     5392                if (!complicatedInteger&&preProcess==0&&clpSolver->tightenPrimalBounds()!=0) {
    53215393#ifndef DISALLOW_PRINTING
    53225394                  std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
     
    55145586                  generator1.setMaxLookRoot(50);
    55155587                  generator1.setRowCuts(3);
     5588                  if ((tunePreProcess&1)!=0) {
     5589                    // heavy probing
     5590                    generator1.setMaxPassRoot(2);
     5591                    generator1.setMaxElements(300);
     5592                  }
    55165593                  // Add in generators
    55175594                  process.addCutGenerator(&generator1);
    5518                   int translate[]={9999,0,0,-1,2,3,-2,9999,4};
     5595                  int translate[]={9999,0,0,-3,2,3,-2,9999,4,4};
    55195596                  process.passInMessageHandler(babModel_->messageHandler());
    55205597                  //process.messageHandler()->setLogLevel(babModel_->logLevel());
     
    58495926              int switches[20];
    58505927              int numberGenerators=0;
    5851               int translate[]={-100,-1,-99,-98,1,1,1,1};
     5928              int translate[]={-100,-1,-99,-98,1,1,1,1,-1};
    58525929              if (probingAction) {
    58535930                if (probingAction==5||probingAction==7)
     
    58595936                  // How far to follow the consequences
    58605937                  probingGen.setMaxLook(50);
     5938                  probingGen.setMaxLookRoot(50);
     5939                }
     5940                if (probingAction==8) {
     5941                  probingGen.setMaxPassRoot(2);
     5942                  probingGen.setMaxProbeRoot(babModel_->solver()->getNumCols());
    58615943                  probingGen.setMaxLookRoot(50);
    58625944                }
     
    58765958                switches[numberGenerators++]=0;
    58775959              }
     5960#ifdef CLIQUE_ANALYSIS
     5961              if (miplib&&!storedAmpl.sizeRowCuts()) {
     5962                printf("looking at probing\n");
     5963                babModel_->addCutGenerator(&storedAmpl,1,"Stored");
     5964              }
     5965#endif
    58785966              if (knapsackAction) {
    58795967                babModel_->addCutGenerator(&knapsackGen,translate[knapsackAction],"Knapsack");
     
    69707058                  return returnCode;
    69717059                }
     7060                {
     7061                  int extra1 = parameters_[whichParam(EXTRA1,numberParameters_,parameters_)].intValue();
     7062                  if (extra1!=-1) {
     7063                    if (extra1<10000) {
     7064                      babModel_->setSearchStrategy(extra1);
     7065                      printf("XXXXX searchStrategy %d\n",extra1);
     7066                    } else {
     7067                      babModel_->setNumberAnalyzeIterations(extra1);
     7068                      printf("XXXXX analyze %d\n",extra1);
     7069                    }
     7070                  }
     7071                }
     7072#ifdef CLIQUE_ANALYSIS
     7073                if (!storedAmpl.sizeRowCuts()) {
     7074                  printf("looking at probing\n");
     7075                  babModel_->addCutGenerator(&storedAmpl,1,"Stored");
     7076                }
     7077#endif
    69727078                babModel_->branchAndBound(statistics);
     7079#ifdef COIN_DEVELOP
     7080                void printHistory(const char * file);
     7081                printHistory("branch.log");
     7082#endif
    69737083#ifdef NEW_STYLE_SOLVER
    69747084                returnCode = callBack_->callBack(babModel_,4);
     
    70557165                for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
    70567166                  CbcCutGenerator * generator = babModel_->cutGenerator(iGenerator);
     7167                  CglStored * stored = dynamic_cast<CglStored*>(generator->generator());
     7168                  if (stored&&!generator->numberCutsInTotal())
     7169                    continue;
    70577170                  sprintf(generalPrint,"%s was tried %d times and created %d cuts of which %d were active after adding rounds of cuts",
    70587171                          generator->cutGeneratorName(),
     
    72217334                  clpSolver = dynamic_cast< OsiClpSolverInterface*> (babModel_->solver());
    72227335                  lpSolver = clpSolver->getModelPtr();
    7223                   int numberColumns = originalCoinModel_->numberColumns();
    72247336                  int numberColumns2 = lpSolver->numberColumns();
    7225                   assert (numberColumns2>numberColumns);
     7337                  assert (numberColumns2>originalCoinModel_->numberColumns());
    72267338                  double * primalSolution = new double [numberColumns2];
    72277339                  memset(primalSolution,0,numberColumns2*sizeof(double));
Note: See TracChangeset for help on using the changeset viewer.