Changeset 1293 for branches/sandbox/Cbc/src/CbcBranchActual.cpp
 Timestamp:
 Nov 10, 2009 3:03:28 PM (10 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

branches/sandbox/Cbc/src/CbcBranchActual.cpp
r1286 r1293 1 /* $Id$ */2 // Copyright (C) 2002, International Business Machines3 // Corporation and others. All Rights Reserved.4 #if defined(_MSC_VER)5 // Turn off compiler warning about long names6 # pragma warning(disable:4786)7 #endif8 #include <cassert>9 #include <cstdlib>10 #include <cmath>11 #include <cfloat>12 //#define CBC_DEBUG13 1 14 #include "CoinTypes.hpp"15 #include "OsiSolverInterface.hpp"16 #include "OsiSolverBranch.hpp"17 #include "CbcModel.hpp"18 #include "CbcMessage.hpp"19 #include "CbcBranchActual.hpp"20 #include "CoinSort.hpp"21 #include "CoinError.hpp"22 23 //##############################################################################24 25 // Default Constructor26 CbcClique::CbcClique ()27 : CbcObject(),28 numberMembers_(0),29 numberNonSOSMembers_(0),30 members_(NULL),31 type_(NULL),32 cliqueType_(1),33 slack_(1)34 {35 }36 37 // Useful constructor (which are integer indices)38 CbcClique::CbcClique (CbcModel * model, int cliqueType, int numberMembers,39 const int * which, const char * type, int identifier, int slack)40 : CbcObject(model)41 {42 id_ = identifier;43 numberMembers_ = numberMembers;44 if (numberMembers_) {45 members_ = new int[numberMembers_];46 memcpy(members_, which, numberMembers_*sizeof(int));47 type_ = new char[numberMembers_];48 if (type) {49 memcpy(type_, type, numberMembers_*sizeof(char));50 } else {51 for (int i = 0; i < numberMembers_; i++)52 type_[i] = 1;53 }54 } else {55 members_ = NULL;56 type_ = NULL;57 }58 // Find out how many non sos59 int i;60 numberNonSOSMembers_ = 0;61 for (i = 0; i < numberMembers_; i++)62 if (!type_[i])63 numberNonSOSMembers_++;64 cliqueType_ = cliqueType;65 slack_ = slack;66 }67 68 // Copy constructor69 CbcClique::CbcClique ( const CbcClique & rhs)70 : CbcObject(rhs)71 {72 numberMembers_ = rhs.numberMembers_;73 numberNonSOSMembers_ = rhs.numberNonSOSMembers_;74 if (numberMembers_) {75 members_ = new int[numberMembers_];76 memcpy(members_, rhs.members_, numberMembers_*sizeof(int));77 type_ = new char[numberMembers_];78 memcpy(type_, rhs.type_, numberMembers_*sizeof(char));79 } else {80 members_ = NULL;81 type_ = NULL;82 }83 cliqueType_ = rhs.cliqueType_;84 slack_ = rhs.slack_;85 }86 87 // Clone88 CbcObject *89 CbcClique::clone() const90 {91 return new CbcClique(*this);92 }93 94 // Assignment operator95 CbcClique &96 CbcClique::operator=( const CbcClique & rhs)97 {98 if (this != &rhs) {99 CbcObject::operator=(rhs);100 delete [] members_;101 delete [] type_;102 numberMembers_ = rhs.numberMembers_;103 numberNonSOSMembers_ = rhs.numberNonSOSMembers_;104 if (numberMembers_) {105 members_ = new int[numberMembers_];106 memcpy(members_, rhs.members_, numberMembers_*sizeof(int));107 type_ = new char[numberMembers_];108 memcpy(type_, rhs.type_, numberMembers_*sizeof(char));109 } else {110 members_ = NULL;111 type_ = NULL;112 }113 cliqueType_ = rhs.cliqueType_;114 slack_ = rhs.slack_;115 }116 return *this;117 }118 119 // Destructor120 CbcClique::~CbcClique ()121 {122 delete [] members_;123 delete [] type_;124 }125 double126 CbcClique::infeasibility(const OsiBranchingInformation * /*info*/,127 int &preferredWay) const128 {129 int numberUnsatis = 0, numberFree = 0;130 int j;131 const int * integer = model_>integerVariable();132 OsiSolverInterface * solver = model_>solver();133 const double * solution = model_>testSolution();134 const double * lower = solver>getColLower();135 const double * upper = solver>getColUpper();136 double largestValue = 0.0;137 double integerTolerance =138 model_>getDblParam(CbcModel::CbcIntegerTolerance);139 double * sort = new double[numberMembers_];140 141 double slackValue = 0.0;142 for (j = 0; j < numberMembers_; j++) {143 int sequence = members_[j];144 int iColumn = integer[sequence];145 double value = solution[iColumn];146 value = CoinMax(value, lower[iColumn]);147 value = CoinMin(value, upper[iColumn]);148 double nearest = floor(value + 0.5);149 double distance = fabs(value  nearest);150 if (distance > integerTolerance) {151 if (!type_[j])152 value = 1.0  value; // non SOS153 // if slack then choose that154 if (j == slack_ && value > 0.05)155 slackValue = value;156 largestValue = CoinMax(value, largestValue);157 sort[numberUnsatis++] = value;158 } else if (upper[iColumn] > lower[iColumn]) {159 numberFree++;160 }161 }162 preferredWay = 1;163 double otherWay = 0.0;164 if (numberUnsatis) {165 // sort166 std::sort(sort, sort + numberUnsatis);167 for (j = 0; j < numberUnsatis; j++) {168 if ((j&1) != 0)169 otherWay += sort[j];170 }171 // Need to think more172 double value = 0.2 * numberUnsatis + 0.01 * (numberMembers_  numberFree);173 if (fabs(largestValue  0.5) < 0.1) {174 // close to half175 value += 0.1;176 }177 if (slackValue) {178 // branching on slack179 value += slackValue;180 }181 // scale other way182 otherWay *= value / (1.0  otherWay);183 delete [] sort;184 return value;185 } else {186 delete [] sort;187 return 0.0; // satisfied188 }189 }190 191 // This looks at solution and sets bounds to contain solution192 void193 CbcClique::feasibleRegion()194 {195 int j;196 const int * integer = model_>integerVariable();197 OsiSolverInterface * solver = model_>solver();198 const double * solution = model_>testSolution();199 const double * lower = solver>getColLower();200 const double * upper = solver>getColUpper();201 #ifndef NDEBUG202 double integerTolerance =203 model_>getDblParam(CbcModel::CbcIntegerTolerance);204 #endif205 for (j = 0; j < numberMembers_; j++) {206 int sequence = members_[j];207 int iColumn = integer[sequence];208 double value = solution[iColumn];209 value = CoinMax(value, lower[iColumn]);210 value = CoinMin(value, upper[iColumn]);211 double nearest = floor(value + 0.5);212 #ifndef NDEBUG213 double distance = fabs(value  nearest);214 assert(distance <= integerTolerance);215 #endif216 solver>setColLower(iColumn, nearest);217 solver>setColUpper(iColumn, nearest);218 }219 }220 // Redoes data when sequence numbers change221 void222 CbcClique::redoSequenceEtc(CbcModel * model, int numberColumns, const int * originalColumns)223 {224 model_ = model;225 int n2 = 0;226 for (int j = 0; j < numberMembers_; j++) {227 int iColumn = members_[j];228 int i;229 for (i = 0; i < numberColumns; i++) {230 if (originalColumns[i] == iColumn)231 break;232 }233 if (i < numberColumns) {234 members_[n2] = i;235 type_[n2++] = type_[j];236 }237 }238 if (n2 < numberMembers_) {239 //printf("** SOS number of members reduced from %d to %d!\n",numberMembers_,n2);240 numberMembers_ = n2;241 }242 // Find out how many non sos243 int i;244 numberNonSOSMembers_ = 0;245 for (i = 0; i < numberMembers_; i++)246 if (!type_[i])247 numberNonSOSMembers_++;248 }249 CbcBranchingObject *250 CbcClique::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int way)251 {252 int numberUnsatis = 0;253 int j;254 int nUp = 0;255 int nDown = 0;256 int numberFree = numberMembers_;257 const int * integer = model_>integerVariable();258 //OsiSolverInterface * solver = model_>solver();259 const double * solution = model_>testSolution();260 const double * lower = solver>getColLower();261 const double * upper = solver>getColUpper();262 int * upList = new int[numberMembers_];263 int * downList = new int[numberMembers_];264 double * sort = new double[numberMembers_];265 double integerTolerance =266 model_>getDblParam(CbcModel::CbcIntegerTolerance);267 268 double slackValue = 0.0;269 for (j = 0; j < numberMembers_; j++) {270 int sequence = members_[j];271 int iColumn = integer[sequence];272 double value = solution[iColumn];273 value = CoinMax(value, lower[iColumn]);274 value = CoinMin(value, upper[iColumn]);275 double nearest = floor(value + 0.5);276 double distance = fabs(value  nearest);277 if (distance > integerTolerance) {278 if (!type_[j])279 value = 1.0  value; // non SOS280 // if slack then choose that281 if (j == slack_ && value > 0.05)282 slackValue = value;283 value = value; // for sort284 upList[numberUnsatis] = j;285 sort[numberUnsatis++] = value;286 } else if (upper[iColumn] > lower[iColumn]) {287 upList[numberFree] = j;288 }289 }290 assert (numberUnsatis);291 if (!slackValue) {292 // sort293 CoinSort_2(sort, sort + numberUnsatis, upList);294 // put first in up etc295 int kWay = 1;296 for (j = 0; j < numberUnsatis; j++) {297 if (kWay > 0)298 upList[nUp++] = upList[j];299 else300 downList[nDown++] = upList[j];301 kWay = kWay;302 }303 for (j = numberFree; j < numberMembers_; j++) {304 if (kWay > 0)305 upList[nUp++] = upList[j];306 else307 downList[nDown++] = upList[j];308 kWay = kWay;309 }310 } else {311 // put slack to 0 in first way312 nUp = 1;313 upList[0] = slack_;314 for (j = 0; j < numberUnsatis; j++) {315 downList[nDown++] = upList[j];316 }317 for (j = numberFree; j < numberMembers_; j++) {318 downList[nDown++] = upList[j];319 }320 }321 // create object322 CbcBranchingObject * branch;323 if (numberMembers_ <= 64)324 branch = new CbcCliqueBranchingObject(model_, this, way,325 nDown, downList, nUp, upList);326 else327 branch = new CbcLongCliqueBranchingObject(model_, this, way,328 nDown, downList, nUp, upList);329 delete [] upList;330 delete [] downList;331 delete [] sort;332 return branch;333 }334 335 //##############################################################################336 337 // Default Constructor338 CbcSOS::CbcSOS ()339 : CbcObject(),340 members_(NULL),341 weights_(NULL),342 shadowEstimateDown_(1.0),343 shadowEstimateUp_(1.0),344 downDynamicPseudoRatio_(0.0),345 upDynamicPseudoRatio_(0.0),346 numberTimesDown_(0),347 numberTimesUp_(0),348 numberMembers_(0),349 sosType_(1),350 integerValued_(false)351 {352 }353 354 // Useful constructor (which are indices)355 CbcSOS::CbcSOS (CbcModel * model, int numberMembers,356 const int * which, const double * weights, int identifier, int type)357 : CbcObject(model),358 shadowEstimateDown_(1.0),359 shadowEstimateUp_(1.0),360 downDynamicPseudoRatio_(0.0),361 upDynamicPseudoRatio_(0.0),362 numberTimesDown_(0),363 numberTimesUp_(0),364 numberMembers_(numberMembers),365 sosType_(type)366 {367 id_ = identifier;368 integerValued_ = type == 1;369 if (integerValued_) {370 // check all members integer371 OsiSolverInterface * solver = model>solver();372 if (solver) {373 for (int i = 0; i < numberMembers_; i++) {374 if (!solver>isInteger(which[i]))375 integerValued_ = false;376 }377 } else {378 // can't tell379 integerValued_ = false;380 }381 }382 if (numberMembers_) {383 members_ = new int[numberMembers_];384 weights_ = new double[numberMembers_];385 memcpy(members_, which, numberMembers_*sizeof(int));386 if (weights) {387 memcpy(weights_, weights, numberMembers_*sizeof(double));388 } else {389 for (int i = 0; i < numberMembers_; i++)390 weights_[i] = i;391 }392 // sort so weights increasing393 CoinSort_2(weights_, weights_ + numberMembers_, members_);394 double last = COIN_DBL_MAX;395 int i;396 for (i = 0; i < numberMembers_; i++) {397 double possible = CoinMax(last + 1.0e10, weights_[i]);398 weights_[i] = possible;399 last = possible;400 }401 } else {402 members_ = NULL;403 weights_ = NULL;404 }405 assert (sosType_ > 0 && sosType_ < 3);406 }407 408 // Copy constructor409 CbcSOS::CbcSOS ( const CbcSOS & rhs)410 : CbcObject(rhs)411 {412 shadowEstimateDown_ = rhs.shadowEstimateDown_;413 shadowEstimateUp_ = rhs.shadowEstimateUp_;414 downDynamicPseudoRatio_ = rhs.downDynamicPseudoRatio_;415 upDynamicPseudoRatio_ = rhs.upDynamicPseudoRatio_;416 numberTimesDown_ = rhs.numberTimesDown_;417 numberTimesUp_ = rhs.numberTimesUp_;418 numberMembers_ = rhs.numberMembers_;419 sosType_ = rhs.sosType_;420 integerValued_ = rhs.integerValued_;421 if (numberMembers_) {422 members_ = new int[numberMembers_];423 weights_ = new double[numberMembers_];424 memcpy(members_, rhs.members_, numberMembers_*sizeof(int));425 memcpy(weights_, rhs.weights_, numberMembers_*sizeof(double));426 } else {427 members_ = NULL;428 weights_ = NULL;429 }430 }431 432 // Clone433 CbcObject *434 CbcSOS::clone() const435 {436 return new CbcSOS(*this);437 }438 439 // Assignment operator440 CbcSOS &441 CbcSOS::operator=( const CbcSOS & rhs)442 {443 if (this != &rhs) {444 CbcObject::operator=(rhs);445 delete [] members_;446 delete [] weights_;447 shadowEstimateDown_ = rhs.shadowEstimateDown_;448 shadowEstimateUp_ = rhs.shadowEstimateUp_;449 downDynamicPseudoRatio_ = rhs.downDynamicPseudoRatio_;450 upDynamicPseudoRatio_ = rhs.upDynamicPseudoRatio_;451 numberTimesDown_ = rhs.numberTimesDown_;452 numberTimesUp_ = rhs.numberTimesUp_;453 numberMembers_ = rhs.numberMembers_;454 sosType_ = rhs.sosType_;455 integerValued_ = rhs.integerValued_;456 if (numberMembers_) {457 members_ = new int[numberMembers_];458 weights_ = new double[numberMembers_];459 memcpy(members_, rhs.members_, numberMembers_*sizeof(int));460 memcpy(weights_, rhs.weights_, numberMembers_*sizeof(double));461 } else {462 members_ = NULL;463 weights_ = NULL;464 }465 }466 return *this;467 }468 469 // Destructor470 CbcSOS::~CbcSOS ()471 {472 delete [] members_;473 delete [] weights_;474 }475 double476 CbcSOS::infeasibility(const OsiBranchingInformation * info,477 int &preferredWay) const478 {479 int j;480 int firstNonZero = 1;481 int lastNonZero = 1;482 OsiSolverInterface * solver = model_>solver();483 const double * solution = model_>testSolution();484 //const double * lower = solver>getColLower();485 const double * upper = solver>getColUpper();486 //double largestValue=0.0;487 double integerTolerance =488 model_>getDblParam(CbcModel::CbcIntegerTolerance);489 double weight = 0.0;490 double sum = 0.0;491 492 // check bounds etc493 double lastWeight = 1.0e100;494 for (j = 0; j < numberMembers_; j++) {495 int iColumn = members_[j];496 if (lastWeight >= weights_[j]  1.0e7)497 throw CoinError("Weights too close together in SOS", "infeasibility", "CbcSOS");498 double value = CoinMax(0.0, solution[iColumn]);499 sum += value;500 if (value > integerTolerance && upper[iColumn]) {501 // Possibly due to scaling a fixed variable might slip through502 if (value > upper[iColumn]) {503 value = upper[iColumn];504 // Could change to #ifdef CBC_DEBUG505 #ifndef NDEBUG506 if (model_>messageHandler()>logLevel() > 2)507 printf("** Variable %d (%d) has value %g and upper bound of %g\n",508 iColumn, j, value, upper[iColumn]);509 #endif510 }511 weight += weights_[j] * value;512 if (firstNonZero < 0)513 firstNonZero = j;514 lastNonZero = j;515 }516 }517 preferredWay = 1;518 if (lastNonZero  firstNonZero >= sosType_) {519 // find where to branch520 assert (sum > 0.0);521 weight /= sum;522 if (info>defaultDual_ >= 0.0 && info>usefulRegion_ && info>columnStart_) {523 assert (sosType_ == 1);524 int iWhere;525 for (iWhere = firstNonZero; iWhere < lastNonZero  1; iWhere++) {526 if (weight < weights_[iWhere+1]) {527 break;528 }529 }530 int jColumnDown = members_[iWhere];531 int jColumnUp = members_[iWhere+1];532 int n = 0;533 CoinBigIndex j;534 double objMove = info>objective_[jColumnDown];535 for (j = info>columnStart_[jColumnDown];536 j < info>columnStart_[jColumnDown] + info>columnLength_[jColumnDown]; j++) {537 double value = info>elementByColumn_[j];538 int iRow = info>row_[j];539 info>indexRegion_[n++] = iRow;540 info>usefulRegion_[iRow] = value;541 }542 for (iWhere = firstNonZero; iWhere < lastNonZero; iWhere++) {543 int jColumn = members_[iWhere];544 double solValue = info>solution_[jColumn];545 if (!solValue)546 continue;547 objMove = info>objective_[jColumn] * solValue;548 for (j = info>columnStart_[jColumn];549 j < info>columnStart_[jColumn] + info>columnLength_[jColumn]; j++) {550 double value = info>elementByColumn_[j] * solValue;551 int iRow = info>row_[j];552 double oldValue = info>usefulRegion_[iRow];553 if (!oldValue) {554 info>indexRegion_[n++] = iRow;555 } else {556 value += oldValue;557 if (!value)558 value = 1.0e100;559 }560 info>usefulRegion_[iRow] = value;561 }562 }563 const double * pi = info>pi_;564 const double * activity = info>rowActivity_;565 const double * lower = info>rowLower_;566 const double * upper = info>rowUpper_;567 double tolerance = info>primalTolerance_;568 double direction = info>direction_;569 shadowEstimateDown_ = objMove * direction;570 bool infeasible = false;571 for (int k = 0; k < n; k++) {572 int iRow = info>indexRegion_[k];573 double movement = info>usefulRegion_[iRow];574 // not this time info>usefulRegion_[iRow]=0.0;575 if (lower[iRow] < 1.0e20)576 assert (pi[iRow] <= 1.0e3);577 if (upper[iRow] > 1.0e20)578 assert (pi[iRow] >= 1.0e3);579 double valueP = pi[iRow] * direction;580 // if move makes infeasible then make at least default581 double newValue = activity[iRow] + movement;582 if (newValue > upper[iRow] + tolerance  newValue < lower[iRow]  tolerance) {583 shadowEstimateDown_ += fabs(movement) * CoinMax(fabs(valueP), info>defaultDual_);584 infeasible = true;585 }586 }587 if (shadowEstimateDown_ < info>integerTolerance_) {588 if (!infeasible) {589 shadowEstimateDown_ = 1.0e10;590 #ifdef COIN_DEVELOP591 printf("zero pseudoShadowPrice\n");592 #endif593 } else594 shadowEstimateDown_ = info>integerTolerance_;595 }596 // And other way597 // take off598 objMove = info>objective_[jColumnDown];599 for (j = info>columnStart_[jColumnDown];600 j < info>columnStart_[jColumnDown] + info>columnLength_[jColumnDown]; j++) {601 double value = info>elementByColumn_[j];602 int iRow = info>row_[j];603 double oldValue = info>usefulRegion_[iRow];604 if (!oldValue) {605 info>indexRegion_[n++] = iRow;606 } else {607 value += oldValue;608 if (!value)609 value = 1.0e100;610 }611 info>usefulRegion_[iRow] = value;612 }613 // add on614 objMove += info>objective_[jColumnUp];615 for (j = info>columnStart_[jColumnUp];616 j < info>columnStart_[jColumnUp] + info>columnLength_[jColumnUp]; j++) {617 double value = info>elementByColumn_[j];618 int iRow = info>row_[j];619 double oldValue = info>usefulRegion_[iRow];620 if (!oldValue) {621 info>indexRegion_[n++] = iRow;622 } else {623 value += oldValue;624 if (!value)625 value = 1.0e100;626 }627 info>usefulRegion_[iRow] = value;628 }629 shadowEstimateUp_ = objMove * direction;630 infeasible = false;631 for (int k = 0; k < n; k++) {632 int iRow = info>indexRegion_[k];633 double movement = info>usefulRegion_[iRow];634 info>usefulRegion_[iRow] = 0.0;635 if (lower[iRow] < 1.0e20)636 assert (pi[iRow] <= 1.0e3);637 if (upper[iRow] > 1.0e20)638 assert (pi[iRow] >= 1.0e3);639 double valueP = pi[iRow] * direction;640 // if move makes infeasible then make at least default641 double newValue = activity[iRow] + movement;642 if (newValue > upper[iRow] + tolerance  newValue < lower[iRow]  tolerance) {643 shadowEstimateUp_ += fabs(movement) * CoinMax(fabs(valueP), info>defaultDual_);644 infeasible = true;645 }646 }647 if (shadowEstimateUp_ < info>integerTolerance_) {648 if (!infeasible) {649 shadowEstimateUp_ = 1.0e10;650 #ifdef COIN_DEVELOP651 printf("zero pseudoShadowPrice\n");652 #endif653 } else654 shadowEstimateUp_ = info>integerTolerance_;655 }656 // adjust657 double downCost = shadowEstimateDown_;658 double upCost = shadowEstimateUp_;659 if (numberTimesDown_)660 downCost *= downDynamicPseudoRatio_ /661 static_cast<double> (numberTimesDown_);662 if (numberTimesUp_)663 upCost *= upDynamicPseudoRatio_ /664 static_cast<double> (numberTimesUp_);665 #define WEIGHT_AFTER 0.7666 #define WEIGHT_BEFORE 0.1667 int stateOfSearch = model_>stateOfSearch() % 10;668 double returnValue = 0.0;669 double minValue = CoinMin(downCost, upCost);670 double maxValue = CoinMax(downCost, upCost);671 if (stateOfSearch <= 2) {672 // no branching solution673 returnValue = WEIGHT_BEFORE * minValue + (1.0  WEIGHT_BEFORE) * maxValue;674 } else {675 returnValue = WEIGHT_AFTER * minValue + (1.0  WEIGHT_AFTER) * maxValue;676 }677 #ifdef PRINT_SHADOW678 printf("%d id  down %d %g up %d %g shadow %g, %g returned %g\n",679 id_, numberTimesDown_, downDynamicPseudoRatio_,680 numberTimesUp_, upDynamicPseudoRatio_, shadowEstimateDown_,681 shadowEstimateUp_, returnValue);682 #endif683 return returnValue;684 } else {685 double value = lastNonZero  firstNonZero + 1;686 value *= 0.5 / static_cast<double> (numberMembers_);687 return value;688 }689 } else {690 return 0.0; // satisfied691 }692 }693 694 // This looks at solution and sets bounds to contain solution695 void696 CbcSOS::feasibleRegion()697 {698 int j;699 int firstNonZero = 1;700 int lastNonZero = 1;701 OsiSolverInterface * solver = model_>solver();702 const double * solution = model_>testSolution();703 //const double * lower = solver>getColLower();704 const double * upper = solver>getColUpper();705 double integerTolerance =706 model_>getDblParam(CbcModel::CbcIntegerTolerance);707 double weight = 0.0;708 double sum = 0.0;709 710 for (j = 0; j < numberMembers_; j++) {711 int iColumn = members_[j];712 double value = CoinMax(0.0, solution[iColumn]);713 sum += value;714 if (value > integerTolerance && upper[iColumn]) {715 weight += weights_[j] * value;716 if (firstNonZero < 0)717 firstNonZero = j;718 lastNonZero = j;719 }720 }721 assert (lastNonZero  firstNonZero < sosType_) ;722 for (j = 0; j < firstNonZero; j++) {723 int iColumn = members_[j];724 solver>setColUpper(iColumn, 0.0);725 }726 for (j = lastNonZero + 1; j < numberMembers_; j++) {727 int iColumn = members_[j];728 solver>setColUpper(iColumn, 0.0);729 }730 }731 // Redoes data when sequence numbers change732 void733 CbcSOS::redoSequenceEtc(CbcModel * model, int numberColumns, const int * originalColumns)734 {735 model_ = model;736 int n2 = 0;737 for (int j = 0; j < numberMembers_; j++) {738 int iColumn = members_[j];739 int i;740 for (i = 0; i < numberColumns; i++) {741 if (originalColumns[i] == iColumn)742 break;743 }744 if (i < numberColumns) {745 members_[n2] = i;746 weights_[n2++] = weights_[j];747 }748 }749 if (n2 < numberMembers_) {750 //printf("** SOS number of members reduced from %d to %d!\n",numberMembers_,n2);751 numberMembers_ = n2;752 }753 }754 CbcBranchingObject *755 CbcSOS::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int way)756 {757 int j;758 const double * solution = model_>testSolution();759 double integerTolerance =760 model_>getDblParam(CbcModel::CbcIntegerTolerance);761 //OsiSolverInterface * solver = model_>solver();762 const double * upper = solver>getColUpper();763 int firstNonFixed = 1;764 int lastNonFixed = 1;765 int firstNonZero = 1;766 int lastNonZero = 1;767 double weight = 0.0;768 double sum = 0.0;769 for (j = 0; j < numberMembers_; j++) {770 int iColumn = members_[j];771 if (upper[iColumn]) {772 double value = CoinMax(0.0, solution[iColumn]);773 sum += value;774 if (firstNonFixed < 0)775 firstNonFixed = j;776 lastNonFixed = j;777 if (value > integerTolerance) {778 weight += weights_[j] * value;779 if (firstNonZero < 0)780 firstNonZero = j;781 lastNonZero = j;782 }783 }784 }785 assert (lastNonZero  firstNonZero >= sosType_) ;786 // find where to branch787 assert (sum > 0.0);788 weight /= sum;789 int iWhere;790 double separator = 0.0;791 for (iWhere = firstNonZero; iWhere < lastNonZero; iWhere++)792 if (weight < weights_[iWhere+1])793 break;794 if (sosType_ == 1) {795 // SOS 1796 separator = 0.5 * (weights_[iWhere] + weights_[iWhere+1]);797 } else {798 // SOS 2799 if (iWhere == firstNonFixed)800 iWhere++;;801 if (iWhere == lastNonFixed  1)802 iWhere = lastNonFixed  2;803 separator = weights_[iWhere+1];804 }805 // create object806 CbcBranchingObject * branch;807 branch = new CbcSOSBranchingObject(model_, this, way, separator);808 branch>setOriginalObject(this);809 return branch;810 }811 /* Pass in information on branch just done and create CbcObjectUpdateData instance.812 If object does not need data then backward pointer will be NULL.813 Assumes can get information from solver */814 CbcObjectUpdateData815 CbcSOS::createUpdateInformation(const OsiSolverInterface * solver,816 const CbcNode * node,817 const CbcBranchingObject * branchingObject)818 {819 double originalValue = node>objectiveValue();820 int originalUnsatisfied = node>numberUnsatisfied();821 double objectiveValue = solver>getObjValue() * solver>getObjSense();822 int unsatisfied = 0;823 int i;824 //might be base model  doesn't matter825 int numberIntegers = model_>numberIntegers();;826 const double * solution = solver>getColSolution();827 //const double * lower = solver>getColLower();828 //const double * upper = solver>getColUpper();829 double change = CoinMax(0.0, objectiveValue  originalValue);830 int iStatus;831 if (solver>isProvenOptimal())832 iStatus = 0; // optimal833 else if (solver>isIterationLimitReached()834 && !solver>isDualObjectiveLimitReached())835 iStatus = 2; // unknown836 else837 iStatus = 1; // infeasible838 839 bool feasible = iStatus != 1;840 if (feasible) {841 double integerTolerance =842 model_>getDblParam(CbcModel::CbcIntegerTolerance);843 const int * integerVariable = model_>integerVariable();844 for (i = 0; i < numberIntegers; i++) {845 int j = integerVariable[i];846 double value = solution[j];847 double nearest = floor(value + 0.5);848 if (fabs(value  nearest) > integerTolerance)849 unsatisfied++;850 }851 }852 int way = branchingObject>way();853 way =  way; // because after branch so moved on854 double value = branchingObject>value();855 CbcObjectUpdateData newData (this, way,856 change, iStatus,857 originalUnsatisfied  unsatisfied, value);858 newData.originalObjective_ = originalValue;859 // Solvers know about direction860 double direction = solver>getObjSense();861 solver>getDblParam(OsiDualObjectiveLimit, newData.cutoff_);862 newData.cutoff_ *= direction;863 return newData;864 }865 // Update object by CbcObjectUpdateData866 void867 CbcSOS::updateInformation(const CbcObjectUpdateData & data)868 {869 bool feasible = data.status_ != 1;870 int way = data.way_;871 //double value = data.branchingValue_;872 double originalValue = data.originalObjective_;873 double change = data.change_;874 if (way < 0) {875 // down876 if (!feasible) {877 double distanceToCutoff = 0.0;878 //double objectiveValue = model_>getCurrentMinimizationObjValue();879 distanceToCutoff = model_>getCutoff()  originalValue;880 if (distanceToCutoff < 1.0e20)881 change = distanceToCutoff * 2.0;882 else883 change = (downDynamicPseudoRatio_ * shadowEstimateDown_ + 1.0e3) * 10.0;884 }885 change = CoinMax(1.0e12 * (1.0 + fabs(originalValue)), change);886 #ifdef PRINT_SHADOW887 if (numberTimesDown_)888 printf("Updating id %d  down change %g (true %g)  ndown %d estimated change %g  raw shadow estimate %g\n",889 id_, change, data.change_, numberTimesDown_, shadowEstimateDown_*890 (downDynamicPseudoRatio_ / ((double) numberTimesDown_)),891 shadowEstimateDown_);892 else893 printf("Updating id %d  down change %g (true %g)  shadow estimate %g\n",894 id_, change, data.change_, shadowEstimateDown_);895 #endif896 numberTimesDown_++;897 downDynamicPseudoRatio_ += change / shadowEstimateDown_;898 } else {899 // up900 if (!feasible) {901 double distanceToCutoff = 0.0;902 //double objectiveValue = model_>getCurrentMinimizationObjValue();903 distanceToCutoff = model_>getCutoff()  originalValue;904 if (distanceToCutoff < 1.0e20)905 change = distanceToCutoff * 2.0;906 else907 change = (upDynamicPseudoRatio_ * shadowEstimateUp_ + 1.0e3) * 10.0;908 }909 change = CoinMax(1.0e12 * (1.0 + fabs(originalValue)), change);910 #ifdef PRINT_SHADOW911 if (numberTimesUp_)912 printf("Updating id %d  up change %g (true %g)  nup %d estimated change %g  raw shadow estimate %g\n",913 id_, change, data.change_, numberTimesUp_, shadowEstimateUp_*914 (upDynamicPseudoRatio_ / ((double) numberTimesUp_)),915 shadowEstimateUp_);916 else917 printf("Updating id %d  up change %g (true %g)  shadow estimate %g\n",918 id_, change, data.change_, shadowEstimateUp_);919 #endif920 numberTimesUp_++;921 upDynamicPseudoRatio_ += change / shadowEstimateUp_;922 }923 }924 925 /* Create an OsiSolverBranch object926 927 This returns NULL if branch not represented by bound changes928 */929 OsiSolverBranch *930 CbcSOS::solverBranch() const931 {932 int j;933 const double * solution = model_>testSolution();934 double integerTolerance =935 model_>getDblParam(CbcModel::CbcIntegerTolerance);936 OsiSolverInterface * solver = model_>solver();937 const double * upper = solver>getColUpper();938 int firstNonFixed = 1;939 int lastNonFixed = 1;940 int firstNonZero = 1;941 int lastNonZero = 1;942 double weight = 0.0;943 double sum = 0.0;944 double * fix = new double[numberMembers_];945 int * which = new int[numberMembers_];946 for (j = 0; j < numberMembers_; j++) {947 int iColumn = members_[j];948 // fix all on one side or other (even if fixed)949 fix[j] = 0.0;950 which[j] = iColumn;951 if (upper[iColumn]) {952 double value = CoinMax(0.0, solution[iColumn]);953 sum += value;954 if (firstNonFixed < 0)955 firstNonFixed = j;956 lastNonFixed = j;957 if (value > integerTolerance) {958 weight += weights_[j] * value;959 if (firstNonZero < 0)960 firstNonZero = j;961 lastNonZero = j;962 }963 }964 }965 assert (lastNonZero  firstNonZero >= sosType_) ;966 // find where to branch967 assert (sum > 0.0);968 weight /= sum;969 // down branch fixes ones above weight to 0970 int iWhere;971 int iDownStart = 0;972 int iUpEnd = 0;973 for (iWhere = firstNonZero; iWhere < lastNonZero; iWhere++)974 if (weight < weights_[iWhere+1])975 break;976 if (sosType_ == 1) {977 // SOS 1978 iUpEnd = iWhere + 1;979 iDownStart = iUpEnd;980 } else {981 // SOS 2982 if (iWhere == firstNonFixed)983 iWhere++;;984 if (iWhere == lastNonFixed  1)985 iWhere = lastNonFixed  2;986 iUpEnd = iWhere + 1;987 iDownStart = iUpEnd + 1;988 }989 //990 OsiSolverBranch * branch = new OsiSolverBranch();991 branch>addBranch(1, 0, NULL, NULL, numberMembers_  iDownStart, which + iDownStart, fix);992 branch>addBranch(1, 0, NULL, NULL, iUpEnd, which, fix);993 delete [] fix;994 delete [] which;995 return branch;996 }997 // Construct an OsiSOS object998 OsiSOS *999 CbcSOS::osiObject(const OsiSolverInterface * solver) const1000 {1001 OsiSOS * obj = new OsiSOS(solver, numberMembers_, members_, weights_, sosType_);1002 obj>setPriority(priority());1003 return obj;1004 }1005 1006 //##############################################################################1007 1008 /** Default Constructor1009 1010 Equivalent to an unspecified binary variable.1011 */1012 CbcSimpleInteger::CbcSimpleInteger ()1013 : CbcObject(),1014 originalLower_(0.0),1015 originalUpper_(1.0),1016 breakEven_(0.5),1017 columnNumber_(1),1018 preferredWay_(0)1019 {1020 }1021 1022 /** Useful constructor1023 1024 Loads actual upper & lower bounds for the specified variable.1025 */1026 CbcSimpleInteger::CbcSimpleInteger ( CbcModel * model, int iColumn, double breakEven)1027 : CbcObject(model)1028 {1029 columnNumber_ = iColumn ;1030 originalLower_ = model>solver()>getColLower()[columnNumber_] ;1031 originalUpper_ = model>solver()>getColUpper()[columnNumber_] ;1032 breakEven_ = breakEven;1033 assert (breakEven_ > 0.0 && breakEven_ < 1.0);1034 preferredWay_ = 0;1035 }1036 1037 1038 // Copy constructor1039 CbcSimpleInteger::CbcSimpleInteger ( const CbcSimpleInteger & rhs)1040 : CbcObject(rhs)1041 1042 {1043 columnNumber_ = rhs.columnNumber_;1044 originalLower_ = rhs.originalLower_;1045 originalUpper_ = rhs.originalUpper_;1046 breakEven_ = rhs.breakEven_;1047 preferredWay_ = rhs.preferredWay_;1048 }1049 1050 // Clone1051 CbcObject *1052 CbcSimpleInteger::clone() const1053 {1054 return new CbcSimpleInteger(*this);1055 }1056 1057 // Assignment operator1058 CbcSimpleInteger &1059 CbcSimpleInteger::operator=( const CbcSimpleInteger & rhs)1060 {1061 if (this != &rhs) {1062 CbcObject::operator=(rhs);1063 columnNumber_ = rhs.columnNumber_;1064 originalLower_ = rhs.originalLower_;1065 originalUpper_ = rhs.originalUpper_;1066 breakEven_ = rhs.breakEven_;1067 preferredWay_ = rhs.preferredWay_;1068 }1069 return *this;1070 }1071 1072 // Destructor1073 CbcSimpleInteger::~CbcSimpleInteger ()1074 {1075 }1076 // Construct an OsiSimpleInteger object1077 OsiSimpleInteger *1078 CbcSimpleInteger::osiObject() const1079 {1080 OsiSimpleInteger * obj = new OsiSimpleInteger(columnNumber_,1081 originalLower_, originalUpper_);1082 obj>setPriority(priority());1083 return obj;1084 }1085 double1086 CbcSimpleInteger::infeasibility(const OsiBranchingInformation * info,1087 int &preferredWay) const1088 {1089 double value = info>solution_[columnNumber_];1090 value = CoinMax(value, info>lower_[columnNumber_]);1091 value = CoinMin(value, info>upper_[columnNumber_]);1092 double nearest = floor(value + (1.0  breakEven_));1093 assert (breakEven_ > 0.0 && breakEven_ < 1.0);1094 if (nearest > value)1095 preferredWay = 1;1096 else1097 preferredWay = 1;1098 if (preferredWay_)1099 preferredWay = preferredWay_;1100 double weight = fabs(value  nearest);1101 // normalize so weight is 0.5 at break even1102 if (nearest < value)1103 weight = (0.5 / breakEven_) * weight;1104 else1105 weight = (0.5 / (1.0  breakEven_)) * weight;1106 if (fabs(value  nearest) <= info>integerTolerance_)1107 return 0.0;1108 else1109 return weight;1110 }1111 double1112 CbcSimpleInteger::feasibleRegion(OsiSolverInterface * solver, const OsiBranchingInformation * info) const1113 {1114 double value = info>solution_[columnNumber_];1115 #ifdef COIN_DEVELOP1116 if (fabs(value  floor(value + 0.5)) > 1.0e5)1117 printf("value for %d away from integer %g\n", columnNumber_, value);1118 #endif1119 double newValue = CoinMax(value, info>lower_[columnNumber_]);1120 newValue = CoinMin(newValue, info>upper_[columnNumber_]);1121 newValue = floor(newValue + 0.5);1122 solver>setColLower(columnNumber_, newValue);1123 solver>setColUpper(columnNumber_, newValue);1124 return fabs(value  newValue);1125 }1126 1127 /* Create an OsiSolverBranch object1128 1129 This returns NULL if branch not represented by bound changes1130 */1131 OsiSolverBranch *1132 CbcSimpleInteger::solverBranch(OsiSolverInterface * /*solver*/,1133 const OsiBranchingInformation * info) const1134 {1135 double value = info>solution_[columnNumber_];1136 value = CoinMax(value, info>lower_[columnNumber_]);1137 value = CoinMin(value, info>upper_[columnNumber_]);1138 assert (info>upper_[columnNumber_] > info>lower_[columnNumber_]);1139 #ifndef NDEBUG1140 double nearest = floor(value + 0.5);1141 assert (fabs(value  nearest) > info>integerTolerance_);1142 #endif1143 OsiSolverBranch * branch = new OsiSolverBranch();1144 branch>addBranch(columnNumber_, value);1145 return branch;1146 }1147 // Creates a branching object1148 CbcBranchingObject *1149 CbcSimpleInteger::createCbcBranch(OsiSolverInterface * /*solver*/,1150 const OsiBranchingInformation * info, int way)1151 {1152 CbcIntegerBranchingObject * branch = new CbcIntegerBranchingObject(model_, 0, 1, 0.5);1153 fillCreateBranch(branch, info, way);1154 return branch;1155 }1156 // Fills in a created branching object1157 void1158 CbcSimpleInteger::fillCreateBranch(CbcIntegerBranchingObject * branch, const OsiBranchingInformation * info, int way)1159 {1160 branch>setOriginalObject(this);1161 double value = info>solution_[columnNumber_];1162 value = CoinMax(value, info>lower_[columnNumber_]);1163 value = CoinMin(value, info>upper_[columnNumber_]);1164 assert (info>upper_[columnNumber_] > info>lower_[columnNumber_]);1165 if (!info>hotstartSolution_ && priority_ != 999) {1166 #ifndef NDEBUG1167 double nearest = floor(value + 0.5);1168 assert (fabs(value  nearest) > info>integerTolerance_);1169 #endif1170 } else if (info>hotstartSolution_) {1171 double targetValue = info>hotstartSolution_[columnNumber_];1172 if (way > 0)1173 value = targetValue  0.1;1174 else1175 value = targetValue + 0.1;1176 } else {1177 if (value <= info>lower_[columnNumber_])1178 value += 0.1;1179 else if (value >= info>upper_[columnNumber_])1180 value = 0.1;1181 }1182 assert (value >= info>lower_[columnNumber_] &&1183 value <= info>upper_[columnNumber_]);1184 branch>fillPart(columnNumber_, way, value);1185 }1186 /* Column number if single column object 1 otherwise,1187 so returns >= 01188 Used by heuristics1189 */1190 int1191 CbcSimpleInteger::columnNumber() const1192 {1193 return columnNumber_;1194 }1195 /* Reset variable bounds to their original values.1196 1197 Bounds may be tightened, so it may be good to be able to set this info in object.1198 */1199 void1200 CbcSimpleInteger::resetBounds(const OsiSolverInterface * solver)1201 {1202 originalLower_ = solver>getColLower()[columnNumber_] ;1203 originalUpper_ = solver>getColUpper()[columnNumber_] ;1204 }1205 1206 /* Change column numbers after preprocessing1207 */1208 void1209 CbcSimpleInteger::resetSequenceEtc(int /*numberColumns*/,1210 const int * originalColumns)1211 {1212 //assert (numberColumns>0);1213 int iColumn;1214 #if 01215 for (iColumn = 0; iColumn < numberColumns; iColumn++) {1216 if (columnNumber_ == originalColumns[iColumn])1217 break;1218 }1219 assert (iColumn < numberColumns);1220 #else1221 iColumn = originalColumns[columnNumber_];1222 assert (iColumn >= 0);1223 #endif1224 columnNumber_ = iColumn;1225 }1226 // This looks at solution and sets bounds to contain solution1227 /** More precisely: it first forces the variable within the existing1228 bounds, and then tightens the bounds to fix the variable at the1229 nearest integer value.1230 */1231 void1232 CbcSimpleInteger::feasibleRegion()1233 {1234 abort();1235 }1236 1237 //##############################################################################1238 1239 // Default Constructor1240 CbcIntegerBranchingObject::CbcIntegerBranchingObject()1241 : CbcBranchingObject()1242 {1243 down_[0] = 0.0;1244 down_[1] = 0.0;1245 up_[0] = 0.0;1246 up_[1] = 0.0;1247 #ifdef FUNNY_BRANCHING1248 variables_ = NULL;1249 newBounds_ = NULL;1250 numberExtraChangedBounds_ = 0;1251 #endif1252 }1253 // Useful constructor1254 CbcIntegerBranchingObject::CbcIntegerBranchingObject (CbcModel * model,1255 int variable, int way , double value)1256 : CbcBranchingObject(model, variable, way, value)1257 {1258 int iColumn = variable;1259 assert (model_>solver()>getNumCols() > 0);1260 down_[0] = model_>solver()>getColLower()[iColumn];1261 down_[1] = floor(value_);1262 up_[0] = ceil(value_);1263 up_[1] = model>getColUpper()[iColumn];1264 #ifdef FUNNY_BRANCHING1265 variables_ = NULL;1266 newBounds_ = NULL;1267 numberExtraChangedBounds_ = 0;1268 #endif1269 }1270 // Does part of constructor1271 void1272 CbcIntegerBranchingObject::fillPart (int variable,1273 int way , double value)1274 {1275 //originalObject_=NULL;1276 branchIndex_ = 0;1277 value_ = value;1278 numberBranches_ = 2;1279 //model_= model;1280 //originalCbcObject_=NULL;1281 variable_ = variable;1282 way_ = way;1283 int iColumn = variable;1284 down_[0] = model_>solver()>getColLower()[iColumn];1285 down_[1] = floor(value_);1286 up_[0] = ceil(value_);1287 up_[1] = model_>getColUpper()[iColumn];1288 }1289 // Useful constructor for fixing1290 CbcIntegerBranchingObject::CbcIntegerBranchingObject (CbcModel * model,1291 int variable, int way,1292 double lowerValue,1293 double upperValue)1294 : CbcBranchingObject(model, variable, way, lowerValue)1295 {1296 setNumberBranchesLeft(1);1297 down_[0] = lowerValue;1298 down_[1] = upperValue;1299 up_[0] = lowerValue;1300 up_[1] = upperValue;1301 #ifdef FUNNY_BRANCHING1302 variables_ = NULL;1303 newBounds_ = NULL;1304 numberExtraChangedBounds_ = 0;1305 #endif1306 }1307 1308 1309 // Copy constructor1310 CbcIntegerBranchingObject::CbcIntegerBranchingObject ( const CbcIntegerBranchingObject & rhs) : CbcBranchingObject(rhs)1311 {1312 down_[0] = rhs.down_[0];1313 down_[1] = rhs.down_[1];1314 up_[0] = rhs.up_[0];1315 up_[1] = rhs.up_[1];1316 #ifdef FUNNY_BRANCHING1317 numberExtraChangedBounds_ = rhs.numberExtraChangedBounds_;1318 int size = numberExtraChangedBounds_ * (sizeof(double) + sizeof(int));1319 char * temp = new char [size];1320 newBounds_ = (double *) temp;1321 variables_ = (int *) (newBounds_ + numberExtraChangedBounds_);1322 1323 int i ;1324 for (i = 0; i < numberExtraChangedBounds_; i++) {1325 variables_[i] = rhs.variables_[i];1326 newBounds_[i] = rhs.newBounds_[i];1327 }1328 #endif1329 }1330 1331 // Assignment operator1332 CbcIntegerBranchingObject &1333 CbcIntegerBranchingObject::operator=( const CbcIntegerBranchingObject & rhs)1334 {1335 if (this != &rhs) {1336 CbcBranchingObject::operator=(rhs);1337 down_[0] = rhs.down_[0];1338 down_[1] = rhs.down_[1];1339 up_[0] = rhs.up_[0];1340 up_[1] = rhs.up_[1];1341 #ifdef FUNNY_BRANCHING1342 delete [] newBounds_;1343 numberExtraChangedBounds_ = rhs.numberExtraChangedBounds_;1344 int size = numberExtraChangedBounds_ * (sizeof(double) + sizeof(int));1345 char * temp = new char [size];1346 newBounds_ = (double *) temp;1347 variables_ = (int *) (newBounds_ + numberExtraChangedBounds_);1348 1349 int i ;1350 for (i = 0; i < numberExtraChangedBounds_; i++) {1351 variables_[i] = rhs.variables_[i];1352 newBounds_[i] = rhs.newBounds_[i];1353 }1354 #endif1355 }1356 return *this;1357 }1358 CbcBranchingObject *1359 CbcIntegerBranchingObject::clone() const1360 {1361 return (new CbcIntegerBranchingObject(*this));1362 }1363 1364 1365 // Destructor1366 CbcIntegerBranchingObject::~CbcIntegerBranchingObject ()1367 {1368 // for debugging threads1369 way_ = 23456789;1370 #ifdef FUNNY_BRANCHING1371 delete [] newBounds_;1372 #endif1373 }1374 1375 /*1376 Perform a branch by adjusting the bounds of the specified variable. Note1377 that each arm of the branch advances the object to the next arm by1378 advancing the value of way_.1379 1380 Providing new values for the variable's lower and upper bounds for each1381 branching direction gives a little bit of additional flexibility and will1382 be easily extensible to multiway branching.1383 Returns change in guessed objective on next branch1384 */1385 double1386 CbcIntegerBranchingObject::branch()1387 {1388 // for debugging threads1389 if (way_ < 1  way_ > 100000) {1390 printf("way %d, left %d, iCol %d, variable %d\n",1391 way_, numberBranchesLeft(),1392 originalCbcObject_>columnNumber(), variable_);1393 assert (way_ != 23456789);1394 }1395 decrementNumberBranchesLeft();1396 if (down_[1] == COIN_DBL_MAX)1397 return 0.0;1398 int iColumn = originalCbcObject_>columnNumber();1399 assert (variable_ == iColumn);1400 double olb, oub ;1401 olb = model_>solver()>getColLower()[iColumn] ;1402 oub = model_>solver()>getColUpper()[iColumn] ;1403 //#define TIGHTEN_BOUNDS1404 #ifndef TIGHTEN_BOUNDS1405 #ifdef COIN_DEVELOP1406 if (olb != down_[0]  oub != up_[1]) {1407 if (way_ > 0)1408 printf("branching up on var %d: [%g,%g] => [%g,%g]  other [%g,%g]\n",1409 iColumn, olb, oub, up_[0], up_[1], down_[0], down_[1]) ;1410 else1411 printf("branching down on var %d: [%g,%g] => [%g,%g]  other [%g,%g]\n",1412 iColumn, olb, oub, down_[0], down_[1], up_[0], up_[1]) ;1413 }1414 #endif1415 #endif1416 if (way_ < 0) {1417 #ifdef CBC_DEBUG1418 { double olb, oub ;1419 olb = model_>solver()>getColLower()[iColumn] ;1420 oub = model_>solver()>getColUpper()[iColumn] ;1421 printf("branching down on var %d: [%g,%g] => [%g,%g]\n",1422 iColumn, olb, oub, down_[0], down_[1]) ;1423 }1424 #endif1425 #ifndef TIGHTEN_BOUNDS1426 model_>solver()>setColLower(iColumn, down_[0]);1427 #else1428 model_>solver()>setColLower(iColumn, CoinMax(down_[0], olb));1429 #endif1430 model_>solver()>setColUpper(iColumn, down_[1]);1431 //#define CBC_PRINT21432 #ifdef CBC_PRINT21433 printf("%d branching down has bounds %g %g", iColumn, down_[0], down_[1]);1434 #endif1435 #ifdef FUNNY_BRANCHING1436 // branch  do extra bounds1437 for (int i = 0; i < numberExtraChangedBounds_; i++) {1438 int variable = variables_[i];1439 if ((variable&0x40000000) != 0) {1440 // for going down1441 int k = variable & 0x3fffffff;1442 assert (k != iColumn);1443 if ((variable&0x80000000) == 0) {1444 // lower bound changing1445 #ifdef CBC_PRINT21446 printf(" extra for %d changes lower from %g to %g",1447 k, model_>solver()>getColLower()[k], newBounds_[i]);1448 #endif1449 model_>solver()>setColLower(k, newBounds_[i]);1450 } else {1451 // upper bound changing1452 #ifdef CBC_PRINT21453 printf(" extra for %d changes upper from %g to %g",1454 k, model_>solver()>getColUpper()[k], newBounds_[i]);1455 #endif1456 model_>solver()>setColUpper(k, newBounds_[i]);1457 }1458 }1459 }1460 #endif1461 #ifdef CBC_PRINT21462 printf("\n");1463 #endif1464 way_ = 1;1465 } else {1466 #ifdef CBC_DEBUG1467 { double olb, oub ;1468 olb = model_>solver()>getColLower()[iColumn] ;1469 oub = model_>solver()>getColUpper()[iColumn] ;1470 printf("branching up on var %d: [%g,%g] => [%g,%g]\n",1471 iColumn, olb, oub, up_[0], up_[1]) ;1472 }1473 #endif1474 model_>solver()>setColLower(iColumn, up_[0]);1475 #ifndef TIGHTEN_BOUNDS1476 model_>solver()>setColUpper(iColumn, up_[1]);1477 #else1478 model_>solver()>setColUpper(iColumn, CoinMin(up_[1], oub));1479 #endif1480 #ifdef CBC_PRINT21481 printf("%d branching up has bounds %g %g", iColumn, up_[0], up_[1]);1482 #endif1483 #ifdef FUNNY_BRANCHING1484 // branch  do extra bounds1485 for (int i = 0; i < numberExtraChangedBounds_; i++) {1486 int variable = variables_[i];1487 if ((variable&0x40000000) == 0) {1488 // for going up1489 int k = variable & 0x3fffffff;1490 assert (k != iColumn);1491 if ((variable&0x80000000) == 0) {1492 // lower bound changing1493 #ifdef CBC_PRINT21494 printf(" extra for %d changes lower from %g to %g",1495 k, model_>solver()>getColLower()[k], newBounds_[i]);1496 #endif1497 model_>solver()>setColLower(k, newBounds_[i]);1498 } else {1499 // upper bound changing1500 #ifdef CBC_PRINT21501 printf(" extra for %d changes upper from %g to %g",1502 k, model_>solver()>getColUpper()[k], newBounds_[i]);1503 #endif1504 model_>solver()>setColUpper(k, newBounds_[i]);1505 }1506 }1507 }1508 #endif1509 #ifdef CBC_PRINT21510 printf("\n");1511 #endif1512 way_ = 1; // Swap direction1513 }1514 double nlb = model_>solver()>getColLower()[iColumn];1515 double nub = model_>solver()>getColUpper()[iColumn];1516 if (nlb < olb) {1517 #ifndef NDEBUG1518 printf("bad lb change for column %d from %g to %g\n", iColumn, olb, nlb);1519 #endif1520 model_>solver()>setColLower(iColumn, CoinMin(olb, nub));1521 nlb = olb;1522 }1523 if (nub > oub) {1524 #ifndef NDEBUG1525 printf("bad ub change for column %d from %g to %g\n", iColumn, oub, nub);1526 #endif1527 model_>solver()>setColUpper(iColumn, CoinMax(oub, nlb));1528 }1529 #ifndef NDEBUG1530 if (nlb < olb + 1.0e8 && nub > oub  1.0e8 && false)1531 printf("bad null change for column %d  bounds %g,%g\n", iColumn, olb, oub);1532 #endif1533 return 0.0;1534 }1535 /* Update bounds in solver as in 'branch' and update given bounds.1536 branchState is 1 for 'down' +1 for 'up' */1537 void1538 CbcIntegerBranchingObject::fix(OsiSolverInterface * /*solver*/,1539 double * lower, double * upper,1540 int branchState) const1541 {1542 int iColumn = originalCbcObject_>columnNumber();1543 assert (variable_ == iColumn);1544 if (branchState < 0) {1545 model_>solver()>setColLower(iColumn, down_[0]);1546 lower[iColumn] = down_[0];1547 model_>solver()>setColUpper(iColumn, down_[1]);1548 upper[iColumn] = down_[1];1549 } else {1550 model_>solver()>setColLower(iColumn, up_[0]);1551 lower[iColumn] = up_[0];1552 model_>solver()>setColUpper(iColumn, up_[1]);1553 upper[iColumn] = up_[1];1554 }1555 }1556 #ifdef FUNNY_BRANCHING1557 // Deactivate bounds for branching1558 void1559 CbcIntegerBranchingObject::deactivate()1560 {1561 down_[1] = COIN_DBL_MAX;1562 }1563 int1564 CbcIntegerBranchingObject::applyExtraBounds(int iColumn, double lower, double upper, int way)1565 {1566 // branch  do bounds1567 1568 int i;1569 int found = 0;1570 if (variable_ == iColumn) {1571 printf("odd applyExtra %d\n", iColumn);1572 if (way < 0) {1573 down_[0] = CoinMax(lower, down_[0]);1574 down_[1] = CoinMin(upper, down_[1]);1575 assert (down_[0] <= down_[1]);1576 } else {1577 up_[0] = CoinMax(lower, up_[0]);1578 up_[1] = CoinMin(upper, up_[1]);1579 assert (up_[0] <= up_[1]);1580 }1581 return 0;1582 }1583 int check = (way < 0) ? 0x40000000 : 0;1584 double newLower = lower;1585 double newUpper = upper;1586 for (i = 0; i < numberExtraChangedBounds_; i++) {1587 int variable = variables_[i];1588 if ((variable&0x40000000) == check) {1589 int k = variable & 0x3fffffff;1590 if (k == iColumn) {1591 if ((variable&0x80000000) == 0) {1592 // lower bound changing1593 found = 1;1594 newBounds_[i] = CoinMax(lower, newBounds_[i]);1595 newLower = newBounds_[i];1596 } else {1597 // upper bound changing1598 found = 2;1599 newBounds_[i] = CoinMin(upper, newBounds_[i]);1600 newUpper = newBounds_[i];1601 }1602 }1603 }1604 }1605 int nAdd = 0;1606 if ((found&2) == 0) {1607 // need to add new upper1608 nAdd++;1609 }1610 if ((found&1) == 0) {1611 // need to add new lower1612 nAdd++;1613 }1614 if (nAdd) {1615 int size = (numberExtraChangedBounds_ + nAdd) * (sizeof(double) + sizeof(int));1616 char * temp = new char [size];1617 double * newBounds = (double *) temp;1618 int * variables = (int *) (newBounds + numberExtraChangedBounds_ + nAdd);1619 1620 int i ;1621 for (i = 0; i < numberExtraChangedBounds_; i++) {1622 variables[i] = variables_[i];1623 newBounds[i] = newBounds_[i];1624 }1625 delete [] newBounds_;1626 newBounds_ = newBounds;1627 variables_ = variables;1628 if ((found&2) == 0) {1629 // need to add new upper1630 int variable = iColumn  0x80000000;1631 variables_[numberExtraChangedBounds_] = variable;1632 newBounds_[numberExtraChangedBounds_++] = newUpper;1633 }1634 if ((found&1) == 0) {1635 // need to add new lower1636 int variable = iColumn;1637 variables_[numberExtraChangedBounds_] = variable;1638 newBounds_[numberExtraChangedBounds_++] = newLower;1639 }1640 }1641 1642 return (newUpper >= newLower) ? 0 : 1;1643 }1644 #endif1645 // Print what would happen1646 void1647 CbcIntegerBranchingObject::print()1648 {1649 int iColumn = originalCbcObject_>columnNumber();1650 assert (variable_ == iColumn);1651 if (way_ < 0) {1652 {1653 double olb, oub ;1654 olb = model_>solver()>getColLower()[iColumn] ;1655 oub = model_>solver()>getColUpper()[iColumn] ;1656 printf("CbcInteger would branch down on var %d (int var %d): [%g,%g] => [%g,%g]\n",1657 iColumn, variable_, olb, oub, down_[0], down_[1]) ;1658 }1659 } else {1660 {1661 double olb, oub ;1662 olb = model_>solver()>getColLower()[iColumn] ;1663 oub = model_>solver()>getColUpper()[iColumn] ;1664 printf("CbcInteger would branch up on var %d (int var %d): [%g,%g] => [%g,%g]\n",1665 iColumn, variable_, olb, oub, up_[0], up_[1]) ;1666 }1667 }1668 }1669 1670 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the1671 same type and must have the same original object, but they may have1672 different feasible regions.1673 Return the appropriate CbcRangeCompare value (first argument being the1674 sub/superset if that's the case). In case of overlap (and if \c1675 replaceIfOverlap is true) replace the current branching object with one1676 whose feasible region is the overlap.1677 */1678 CbcRangeCompare1679 CbcIntegerBranchingObject::compareBranchingObject1680 (const CbcBranchingObject* brObj, const bool replaceIfOverlap)1681 {1682 const CbcIntegerBranchingObject* br =1683 dynamic_cast<const CbcIntegerBranchingObject*>(brObj);1684 assert(br);1685 double* thisBd = way_ < 0 ? down_ : up_;1686 const double* otherBd = br>way_ < 0 ? br>down_ : br>up_;1687 return CbcCompareRanges(thisBd, otherBd, replaceIfOverlap);1688 }1689 1690 //##############################################################################1691 1692 /** Default Constructor1693 1694 Equivalent to an unspecified binary variable.1695 */1696 CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost ()1697 : CbcSimpleInteger(),1698 downPseudoCost_(1.0e5),1699 upPseudoCost_(1.0e5),1700 upDownSeparator_(1.0),1701 method_(0)1702 {1703 }1704 1705 /** Useful constructor1706 1707 Loads actual upper & lower bounds for the specified variable.1708 */1709 CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost (CbcModel * model,1710 int iColumn, double breakEven)1711 : CbcSimpleInteger(model, iColumn, breakEven)1712 {1713 const double * cost = model>getObjCoefficients();1714 double costValue = CoinMax(1.0e5, fabs(cost[iColumn]));1715 // treat as if will cost what it says up1716 upPseudoCost_ = costValue;1717 // and balance at breakeven1718 downPseudoCost_ = ((1.0  breakEven_) * upPseudoCost_) / breakEven_;1719 upDownSeparator_ = 1.0;1720 method_ = 0;1721 }1722 1723 /** Useful constructor1724 1725 Loads actual upper & lower bounds for the specified variable.1726 */1727 CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost (CbcModel * model,1728 int iColumn, double downPseudoCost,1729 double upPseudoCost)1730 : CbcSimpleInteger(model, iColumn)1731 {1732 downPseudoCost_ = CoinMax(1.0e10, downPseudoCost);1733 upPseudoCost_ = CoinMax(1.0e10, upPseudoCost);1734 breakEven_ = upPseudoCost_ / (upPseudoCost_ + downPseudoCost_);1735 upDownSeparator_ = 1.0;1736 method_ = 0;1737 }1738 // Useful constructor  passed and model index and pseudo costs1739 CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost (CbcModel * model,1740 int /*dummy*/,1741 int iColumn,1742 double downPseudoCost, double upPseudoCost)1743 {1744 *this = CbcSimpleIntegerPseudoCost(model, iColumn, downPseudoCost, upPseudoCost);1745 columnNumber_ = iColumn;1746 }1747 1748 // Copy constructor1749 CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost ( const CbcSimpleIntegerPseudoCost & rhs)1750 : CbcSimpleInteger(rhs),1751 downPseudoCost_(rhs.downPseudoCost_),1752 upPseudoCost_(rhs.upPseudoCost_),1753 upDownSeparator_(rhs.upDownSeparator_),1754 method_(rhs.method_)1755 1756 {1757 }1758 1759 // Clone1760 CbcObject *1761 CbcSimpleIntegerPseudoCost::clone() const1762 {1763 return new CbcSimpleIntegerPseudoCost(*this);1764 }1765 1766 // Assignment operator1767 CbcSimpleIntegerPseudoCost &1768 CbcSimpleIntegerPseudoCost::operator=( const CbcSimpleIntegerPseudoCost & rhs)1769 {1770 if (this != &rhs) {1771 CbcSimpleInteger::operator=(rhs);1772 downPseudoCost_ = rhs.downPseudoCost_;1773 upPseudoCost_ = rhs.upPseudoCost_;1774 upDownSeparator_ = rhs.upDownSeparator_;1775 method_ = rhs.method_;1776 }1777 return *this;1778 }1779 1780 // Destructor1781 CbcSimpleIntegerPseudoCost::~CbcSimpleIntegerPseudoCost ()1782 {1783 }1784 CbcBranchingObject *1785 CbcSimpleIntegerPseudoCost::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int way)1786 {1787 //OsiSolverInterface * solver = model_>solver();1788 const double * solution = model_>testSolution();1789 const double * lower = solver>getColLower();1790 const double * upper = solver>getColUpper();1791 double value = solution[columnNumber_];1792 value = CoinMax(value, lower[columnNumber_]);1793 value = CoinMin(value, upper[columnNumber_]);1794 #ifndef NDEBUG1795 double nearest = floor(value + 0.5);1796 double integerTolerance =1797 model_>getDblParam(CbcModel::CbcIntegerTolerance);1798 assert (upper[columnNumber_] > lower[columnNumber_]);1799 #endif1800 if (!model_>hotstartSolution()) {1801 assert (fabs(value  nearest) > integerTolerance);1802 } else {1803 const double * hotstartSolution = model_>hotstartSolution();1804 double targetValue = hotstartSolution[columnNumber_];1805 if (way > 0)1806 value = targetValue  0.1;1807 else1808 value = targetValue + 0.1;1809 }1810 CbcIntegerPseudoCostBranchingObject * newObject =1811 new CbcIntegerPseudoCostBranchingObject(model_, columnNumber_, way,1812 value);1813 double up = upPseudoCost_ * (ceil(value)  value);1814 double down = downPseudoCost_ * (value  floor(value));1815 double changeInGuessed = up  down;1816 if (way > 0)1817 changeInGuessed =  changeInGuessed;1818 changeInGuessed = CoinMax(0.0, changeInGuessed);1819 //if (way>0)1820 //changeInGuessed += 1.0e8; // bias to stay up1821 newObject>setChangeInGuessed(changeInGuessed);1822 newObject>setOriginalObject(this);1823 return newObject;1824 }1825 double1826 CbcSimpleIntegerPseudoCost::infeasibility(const OsiBranchingInformation * /*info*/,1827 int &preferredWay) const1828 {1829 OsiSolverInterface * solver = model_>solver();1830 const double * solution = model_>testSolution();1831 const double * lower = solver>getColLower();1832 const double * upper = solver>getColUpper();1833 if (upper[columnNumber_] == lower[columnNumber_]) {1834 // fixed1835 preferredWay = 1;1836 return 0.0;1837 }1838 double value = solution[columnNumber_];1839 value = CoinMax(value, lower[columnNumber_]);1840 value = CoinMin(value, upper[columnNumber_]);1841 /*printf("%d %g %g %g %g\n",columnNumber_,value,lower[columnNumber_],1842 solution[columnNumber_],upper[columnNumber_]);*/1843 double nearest = floor(value + 0.5);1844 double integerTolerance =1845 model_>getDblParam(CbcModel::CbcIntegerTolerance);1846 double below = floor(value + integerTolerance);1847 double above = below + 1.0;1848 if (above > upper[columnNumber_]) {1849 above = below;1850 below = above  1;1851 }1852 double downCost = CoinMax((value  below) * downPseudoCost_, 0.0);1853 double upCost = CoinMax((above  value) * upPseudoCost_, 0.0);1854 if (downCost >= upCost)1855 preferredWay = 1;1856 else1857 preferredWay = 1;1858 // See if up down choice set1859 if (upDownSeparator_ > 0.0) {1860 preferredWay = (value  below >= upDownSeparator_) ? 1 : 1;1861 }1862 if (preferredWay_)1863 preferredWay = preferredWay_;1864 if (fabs(value  nearest) <= integerTolerance) {1865 return 0.0;1866 } else {1867 // can't get at model so 1,2 don't make sense1868 assert(method_ < 1  method_ > 2);1869 if (!method_)1870 return CoinMin(downCost, upCost);1871 else1872 return CoinMax(downCost, upCost);1873 }1874 }1875 1876 // Return "up" estimate1877 double1878 CbcSimpleIntegerPseudoCost::upEstimate() const1879 {1880 OsiSolverInterface * solver = model_>solver();1881 const double * solution = model_>testSolution();1882 const double * lower = solver>getColLower();1883 const double * upper = solver>getColUpper();1884 double value = solution[columnNumber_];1885 value = CoinMax(value, lower[columnNumber_]);1886 value = CoinMin(value, upper[columnNumber_]);1887 if (upper[columnNumber_] == lower[columnNumber_]) {1888 // fixed1889 return 0.0;1890 }1891 double integerTolerance =1892 model_>getDblParam(CbcModel::CbcIntegerTolerance);1893 double below = floor(value + integerTolerance);1894 double above = below + 1.0;1895 if (above > upper[columnNumber_]) {1896 above = below;1897 below = above  1;1898 }1899 double upCost = CoinMax((above  value) * upPseudoCost_, 0.0);1900 return upCost;1901 }1902 // Return "down" estimate1903 double1904 CbcSimpleIntegerPseudoCost::downEstimate() const1905 {1906 OsiSolverInterface * solver = model_>solver();1907 const double * solution = model_>testSolution();1908 const double * lower = solver>getColLower();1909 const double * upper = solver>getColUpper();1910 double value = solution[columnNumber_];1911 value = CoinMax(value, lower[columnNumber_]);1912 value = CoinMin(value, upper[columnNumber_]);1913 if (upper[columnNumber_] == lower[columnNumber_]) {1914 // fixed1915 return 0.0;1916 }1917 double integerTolerance =1918 model_>getDblParam(CbcModel::CbcIntegerTolerance);1919 double below = floor(value + integerTolerance);1920 double above = below + 1.0;1921 if (above > upper[columnNumber_]) {1922 above = below;1923 below = above  1;1924 }1925 double downCost = CoinMax((value  below) * downPseudoCost_, 0.0);1926 return downCost;1927 }1928 1929 //##############################################################################1930 1931 // Default Constructor1932 CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject()1933 : CbcIntegerBranchingObject()1934 {1935 changeInGuessed_ = 1.0e5;1936 }1937 1938 // Useful constructor1939 CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject (CbcModel * model,1940 int variable, int way , double value)1941 : CbcIntegerBranchingObject(model, variable, way, value)1942 {1943 }1944 // Useful constructor for fixing1945 CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject (CbcModel * model,1946 int variable, int way,1947 double lowerValue,1948 double /*upperValue*/)1949 : CbcIntegerBranchingObject(model, variable, way, lowerValue)1950 {1951 changeInGuessed_ = 1.0e100;1952 }1953 1954 1955 // Copy constructor1956 CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject (1957 const CbcIntegerPseudoCostBranchingObject & rhs)1958 : CbcIntegerBranchingObject(rhs)1959 {1960 changeInGuessed_ = rhs.changeInGuessed_;1961 }1962 1963 // Assignment operator1964 CbcIntegerPseudoCostBranchingObject &1965 CbcIntegerPseudoCostBranchingObject::operator=( const CbcIntegerPseudoCostBranchingObject & rhs)1966 {1967 if (this != &rhs) {1968 CbcIntegerBranchingObject::operator=(rhs);1969 changeInGuessed_ = rhs.changeInGuessed_;1970 }1971 return *this;1972 }1973 CbcBranchingObject *1974 CbcIntegerPseudoCostBranchingObject::clone() const1975 {1976 return (new CbcIntegerPseudoCostBranchingObject(*this));1977 }1978 1979 1980 // Destructor1981 CbcIntegerPseudoCostBranchingObject::~CbcIntegerPseudoCostBranchingObject ()1982 {1983 }1984 1985 /*1986 Perform a branch by adjusting the bounds of the specified variable. Note1987 that each arm of the branch advances the object to the next arm by1988 advancing the value of way_.1989 1990 Providing new values for the variable's lower and upper bounds for each1991 branching direction gives a little bit of additional flexibility and will1992 be easily extensible to multiway branching.1993 Returns change in guessed objective on next branch1994 */1995 double1996 CbcIntegerPseudoCostBranchingObject::branch()1997 {1998 CbcIntegerBranchingObject::branch();1999 return changeInGuessed_;2000 }2001 2002 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the2003 same type and must have the same original object, but they may have2004 different feasible regions.2005 Return the appropriate CbcRangeCompare value (first argument being the2006 sub/superset if that's the case). In case of overlap (and if \c2007 replaceIfOverlap is true) replace the current branching object with one2008 whose feasible region is the overlap.2009 */2010 CbcRangeCompare2011 CbcIntegerPseudoCostBranchingObject::compareBranchingObject2012 (const CbcBranchingObject* brObj, const bool replaceIfOverlap)2013 {2014 const CbcIntegerPseudoCostBranchingObject* br =2015 dynamic_cast<const CbcIntegerPseudoCostBranchingObject*>(brObj);2016 assert(br);2017 double* thisBd = way_ < 0 ? down_ : up_;2018 const double* otherBd = br>way_ < 0 ? br>down_ : br>up_;2019 return CbcCompareRanges(thisBd, otherBd, replaceIfOverlap);2020 }2021 2022 2023 //##############################################################################2024 2025 // Default Constructor2026 CbcCliqueBranchingObject::CbcCliqueBranchingObject()2027 : CbcBranchingObject()2028 {2029 clique_ = NULL;2030 downMask_[0] = 0;2031 downMask_[1] = 0;2032 upMask_[0] = 0;2033 upMask_[1] = 0;2034 }2035 2036 // Useful constructor2037 CbcCliqueBranchingObject::CbcCliqueBranchingObject (CbcModel * model,2038 const CbcClique * clique,2039 int way ,2040 int numberOnDownSide, const int * down,2041 int numberOnUpSide, const int * up)2042 : CbcBranchingObject(model, clique>id(), way, 0.5)2043 {2044 clique_ = clique;2045 downMask_[0] = 0;2046 downMask_[1] = 0;2047 upMask_[0] = 0;2048 upMask_[1] = 0;2049 int i;2050 for (i = 0; i < numberOnDownSide; i++) {2051 int sequence = down[i];2052 int iWord = sequence >> 5;2053 int iBit = sequence  32 * iWord;2054 unsigned int k = 1 << iBit;2055 downMask_[iWord] = k;2056 }2057 for (i = 0; i < numberOnUpSide; i++) {2058 int sequence = up[i];2059 int iWord = sequence >> 5;2060 int iBit = sequence  32 * iWord;2061 unsigned int k = 1 << iBit;2062 upMask_[iWord] = k;2063 }2064 }2065 2066 // Copy constructor2067 CbcCliqueBranchingObject::CbcCliqueBranchingObject ( const CbcCliqueBranchingObject & rhs) : CbcBranchingObject(rhs)2068 {2069 clique_ = rhs.clique_;2070 downMask_[0] = rhs.downMask_[0];2071 downMask_[1] = rhs.downMask_[1];2072 upMask_[0] = rhs.upMask_[0];2073 upMask_[1] = rhs.upMask_[1];2074 }2075 2076 // Assignment operator2077 CbcCliqueBranchingObject &2078 CbcCliqueBranchingObject::operator=( const CbcCliqueBranchingObject & rhs)2079 {2080 if (this != &rhs) {2081 CbcBranchingObject::operator=(rhs);2082 clique_ = rhs.clique_;2083 downMask_[0] = rhs.downMask_[0];2084 downMask_[1] = rhs.downMask_[1];2085 upMask_[0] = rhs.upMask_[0];2086 upMask_[1] = rhs.upMask_[1];2087 }2088 return *this;2089 }2090 CbcBranchingObject *2091 CbcCliqueBranchingObject::clone() const2092 {2093 return (new CbcCliqueBranchingObject(*this));2094 }2095 2096 2097 // Destructor2098 CbcCliqueBranchingObject::~CbcCliqueBranchingObject ()2099 {2100 }2101 double2102 CbcCliqueBranchingObject::branch()2103 {2104 decrementNumberBranchesLeft();2105 int iWord;2106 int numberMembers = clique_>numberMembers();2107 const int * which = clique_>members();2108 const int * integerVariables = model_>integerVariable();2109 int numberWords = (numberMembers + 31) >> 5;2110 // *** for way  up means fix all those in down section2111 if (way_ < 0) {2112 #ifdef FULL_PRINT2113 printf("Down Fix ");2114 #endif2115 for (iWord = 0; iWord < numberWords; iWord++) {2116 int i;2117 for (i = 0; i < 32; i++) {2118 unsigned int k = 1 << i;2119 if ((upMask_[iWord]&k) != 0) {2120 int iColumn = which[i+32*iWord];2121 #ifdef FULL_PRINT2122 printf("%d ", i + 32*iWord);2123 #endif2124 // fix weak way2125 if (clique_>type(i + 32*iWord))2126 model_>solver()>setColUpper(integerVariables[iColumn], 0.0);2127 else2128 model_>solver()>setColLower(integerVariables[iColumn], 1.0);2129 }2130 }2131 }2132 way_ = 1; // Swap direction2133 } else {2134 #ifdef FULL_PRINT2135 printf("Up Fix ");2136 #endif2137 for (iWord = 0; iWord < numberWords; iWord++) {2138 int i;2139 for (i = 0; i < 32; i++) {2140 unsigned int k = 1 << i;2141 if ((downMask_[iWord]&k) != 0) {2142 int iColumn = which[i+32*iWord];2143 #ifdef FULL_PRINT2144 printf("%d ", i + 32*iWord);2145 #endif2146 // fix weak way2147 if (clique_>type(i + 32*iWord))2148 model_>solver()>setColUpper(integerVariables[iColumn], 0.0);2149 else2150 model_>solver()>setColLower(integerVariables[iColumn], 1.0);2151 }2152 }2153 }2154 way_ = 1; // Swap direction2155 }2156 #ifdef FULL_PRINT2157 printf("\n");2158 #endif2159 return 0.0;2160 }2161 // Print what would happen2162 void2163 CbcCliqueBranchingObject::print()2164 {2165 int iWord;2166 int numberMembers = clique_>numberMembers();2167 const int * which = clique_>members();2168 const int * integerVariables = model_>integerVariable();2169 int numberWords = (numberMembers + 31) >> 5;2170 // *** for way  up means fix all those in down section2171 if (way_ < 0) {2172 printf("Clique  Down Fix ");2173 for (iWord = 0; iWord < numberWords; iWord++) {2174 int i;2175 for (i = 0; i < 32; i++) {2176 unsigned int k = 1 << i;2177 if ((upMask_[iWord]&k) != 0) {2178 int iColumn = which[i+32*iWord];2179 printf("%d ", integerVariables[iColumn]);2180 }2181 }2182 }2183 } else {2184 printf("Clique  Up Fix ");2185 for (iWord = 0; iWord < numberWords; iWord++) {2186 int i;2187 for (i = 0; i < 32; i++) {2188 unsigned int k = 1 << i;2189 if ((downMask_[iWord]&k) != 0) {2190 int iColumn = which[i+32*iWord];2191 printf("%d ", integerVariables[iColumn]);2192 }2193 }2194 }2195 }2196 printf("\n");2197 }2198 2199 static inline int2200 CbcCompareCliques(const CbcClique* cl0, const CbcClique* cl1)2201 {2202 if (cl0>cliqueType() < cl1>cliqueType()) {2203 return 1;2204 }2205 if (cl0>cliqueType() > cl1>cliqueType()) {2206 return 1;2207 }2208 if (cl0>numberMembers() != cl1>numberMembers()) {2209 return cl0>numberMembers()  cl1>numberMembers();2210 }2211 if (cl0>numberNonSOSMembers() != cl1>numberNonSOSMembers()) {2212 return cl0>numberNonSOSMembers()  cl1>numberNonSOSMembers();2213 }2214 return memcmp(cl0>members(), cl1>members(),2215 cl0>numberMembers() * sizeof(int));2216 }2217 2218 /** Compare the original object of \c this with the original object of \c2219 brObj. Assumes that there is an ordering of the original objects.2220 This method should be invoked only if \c this and brObj are of the same2221 type.2222 Return negative/0/positive depending on whether \c this is2223 smaller/same/larger than the argument.2224 */2225 int2226 CbcCliqueBranchingObject::compareOriginalObject2227 (const CbcBranchingObject* brObj) const2228 {2229 const CbcCliqueBranchingObject* br =2230 dynamic_cast<const CbcCliqueBranchingObject*>(brObj);2231 assert(br);2232 return CbcCompareCliques(clique_, br>clique_);2233 }2234 2235 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the2236 same type and must have the same original object, but they may have2237 different feasible regions.2238 Return the appropriate CbcRangeCompare value (first argument being the2239 sub/superset if that's the case). In case of overlap (and if \c2240 replaceIfOverlap is true) replace the current branching object with one2241 whose feasible region is the overlap.2242 */2243 CbcRangeCompare2244 CbcCliqueBranchingObject::compareBranchingObject2245 (const CbcBranchingObject* brObj, const bool /*replaceIfOverlap*/)2246 {2247 const CbcCliqueBranchingObject* br =2248 dynamic_cast<const CbcCliqueBranchingObject*>(brObj);2249 assert(br);2250 unsigned int* thisMask = way_ < 0 ? upMask_ : downMask_;2251 const unsigned int* otherMask = br>way_ < 0 ? br>upMask_ : br>downMask_;2252 const CoinUInt64 cl0 =2253 (static_cast<CoinUInt64>(thisMask[0]) << 32)  thisMask[1];2254 const CoinUInt64 cl1 =2255 (static_cast<CoinUInt64>(otherMask[0]) << 32)  otherMask[1];2256 if (cl0 == cl1) {2257 return CbcRangeSame;2258 }2259 const CoinUInt64 cl_intersection = (cl0 & cl1);2260 if (cl_intersection == cl0) {2261 return CbcRangeSuperset;2262 }2263 if (cl_intersection == cl1) {2264 return CbcRangeSubset;2265 }2266 const CoinUInt64 cl_xor = (cl0 ^ cl1);2267 if (cl_intersection == 0 && cl_xor == 0) {2268 return CbcRangeDisjoint;2269 }2270 const CoinUInt64 cl_union = (cl0  cl1);2271 thisMask[0] = static_cast<unsigned int>(cl_union >> 32);2272 thisMask[1] = static_cast<unsigned int>(cl_union & 0xffffffff);2273 return CbcRangeOverlap;2274 }2275 2276 //##############################################################################2277 2278 // Default Constructor2279 CbcLongCliqueBranchingObject::CbcLongCliqueBranchingObject()2280 : CbcBranchingObject()2281 {2282 clique_ = NULL;2283 downMask_ = NULL;2284 upMask_ = NULL;2285 }2286 2287 // Useful constructor2288 CbcLongCliqueBranchingObject::CbcLongCliqueBranchingObject (CbcModel * model,2289 const CbcClique * clique,2290 int way ,2291 int numberOnDownSide, const int * down,2292 int numberOnUpSide, const int * up)2293 : CbcBranchingObject(model, clique>id(), way, 0.5)2294 {2295 clique_ = clique;2296 int numberMembers = clique_>numberMembers();2297 int numberWords = (numberMembers + 31) >> 5;2298 downMask_ = new unsigned int [numberWords];2299 upMask_ = new unsigned int [numberWords];2300 memset(downMask_, 0, numberWords*sizeof(unsigned int));2301 memset(upMask_, 0, numberWords*sizeof(unsigned int));2302 int i;2303 for (i = 0; i < numberOnDownSide; i++) {2304 int sequence = down[i];2305 int iWord = sequence >> 5;2306 int iBit = sequence  32 * iWord;2307 unsigned int k = 1 << iBit;2308 downMask_[iWord] = k;2309 }2310 for (i = 0; i < numberOnUpSide; i++) {2311 int sequence = up[i];2312 int iWord = sequence >> 5;2313 int iBit = sequence  32 * iWord;2314 unsigned int k = 1 << iBit;2315 upMask_[iWord] = k;2316 }2317 }2318 2319 // Copy constructor2320 CbcLongCliqueBranchingObject::CbcLongCliqueBranchingObject ( const CbcLongCliqueBranchingObject & rhs) : CbcBranchingObject(rhs)2321 {2322 clique_ = rhs.clique_;2323 if (rhs.downMask_) {2324 int numberMembers = clique_>numberMembers();2325 int numberWords = (numberMembers + 31) >> 5;2326 downMask_ = new unsigned int [numberWords];2327 memcpy(downMask_, rhs.downMask_, numberWords*sizeof(unsigned int));2328 upMask_ = new unsigned int [numberWords];2329 memcpy(upMask_, rhs.upMask_, numberWords*sizeof(unsigned int));2330 } else {2331 downMask_ = NULL;2332 upMask_ = NULL;2333 }2334 }2335 2336 // Assignment operator2337 CbcLongCliqueBranchingObject &2338 CbcLongCliqueBranchingObject::operator=( const CbcLongCliqueBranchingObject & rhs)2339 {2340 if (this != &rhs) {2341 CbcBranchingObject::operator=(rhs);2342 clique_ = rhs.clique_;2343 delete [] downMask_;2344 delete [] upMask_;2345 if (rhs.downMask_) {2346 int numberMembers = clique_>numberMembers();2347 int numberWords = (numberMembers + 31) >> 5;2348 downMask_ = new unsigned int [numberWords];2349 memcpy(downMask_, rhs.downMask_, numberWords*sizeof(unsigned int));2350 upMask_ = new unsigned int [numberWords];2351 memcpy(upMask_, rhs.upMask_, numberWords*sizeof(unsigned int));2352 } else {2353 downMask_ = NULL;2354 upMask_ = NULL;2355 }2356 }2357 return *this;2358 }2359 CbcBranchingObject *2360 CbcLongCliqueBranchingObject::clone() const2361 {2362 return (new CbcLongCliqueBranchingObject(*this));2363 }2364 2365 2366 // Destructor2367 CbcLongCliqueBranchingObject::~CbcLongCliqueBranchingObject ()2368 {2369 delete [] downMask_;2370 delete [] upMask_;2371 }2372 double2373 CbcLongCliqueBranchingObject::branch()2374 {2375 decrementNumberBranchesLeft();2376 int iWord;2377 int numberMembers = clique_>numberMembers();2378 const int * which = clique_>members();2379 const int * integerVariables = model_>integerVariable();2380 int numberWords = (numberMembers + 31) >> 5;2381 // *** for way  up means fix all those in down section2382 if (way_ < 0) {2383 #ifdef FULL_PRINT2384 printf("Down Fix ");2385 #endif2386 for (iWord = 0; iWord < numberWords; iWord++) {2387 int i;2388 for (i = 0; i < 32; i++) {2389 unsigned int k = 1 << i;2390 if ((upMask_[iWord]&k) != 0) {2391 int iColumn = which[i+32*iWord];2392 #ifdef FULL_PRINT2393 printf("%d ", i + 32*iWord);2394 #endif2395 // fix weak way2396 if (clique_>type(i + 32*iWord))2397 model_>solver()>setColUpper(integerVariables[iColumn], 0.0);2398 else2399 model_>solver()>setColLower(integerVariables[iColumn], 1.0);2400 }2401 }2402 }2403 way_ = 1; // Swap direction2404 } else {2405 #ifdef FULL_PRINT2406 printf("Up Fix ");2407 #endif2408 for (iWord = 0; iWord < numberWords; iWord++) {2409 int i;2410 for (i = 0; i < 32; i++) {2411 unsigned int k = 1 << i;2412 if ((downMask_[iWord]&k) != 0) {2413 int iColumn = which[i+32*iWord];2414 #ifdef FULL_PRINT2415 printf("%d ", i + 32*iWord);2416 #endif2417 // fix weak way2418 if (clique_>type(i + 32*iWord))2419 model_>solver()>setColUpper(integerVariables[iColumn], 0.0);2420 else2421 model_>solver()>setColLower(integerVariables[iColumn], 1.0);2422 }2423 }2424 }2425 way_ = 1; // Swap direction2426 }2427 #ifdef FULL_PRINT2428 printf("\n");2429 #endif2430 return 0.0;2431 }2432 void2433 CbcLongCliqueBranchingObject::print()2434 {2435 int iWord;2436 int numberMembers = clique_>numberMembers();2437 const int * which = clique_>members();2438 const int * integerVariables = model_>integerVariable();2439 int numberWords = (numberMembers + 31) >> 5;2440 // *** for way  up means fix all those in down section2441 if (way_ < 0) {2442 printf("Clique  Down Fix ");2443 for (iWord = 0; iWord < numberWords; iWord++) {2444 int i;2445 for (i = 0; i < 32; i++) {2446 unsigned int k = 1 << i;2447 if ((upMask_[iWord]&k) != 0) {2448 int iColumn = which[i+32*iWord];2449 printf("%d ", integerVariables[iColumn]);2450 }2451 }2452 }2453 } else {2454 printf("Clique  Up Fix ");2455 for (iWord = 0; iWord < numberWords; iWord++) {2456 int i;2457 for (i = 0; i < 32; i++) {2458 unsigned int k = 1 << i;2459 if ((downMask_[iWord]&k) != 0) {2460 int iColumn = which[i+32*iWord];2461 printf("%d ", integerVariables[iColumn]);2462 }2463 }2464 }2465 }2466 printf("\n");2467 }2468 2469 /** Compare the original object of \c this with the original object of \c2470 brObj. Assumes that there is an ordering of the original objects.2471 This method should be invoked only if \c this and brObj are of the same2472 type.2473 Return negative/0/positive depending on whether \c this is2474 smaller/same/larger than the argument.2475 */2476 int2477 CbcLongCliqueBranchingObject::compareOriginalObject2478 (const CbcBranchingObject* brObj) const2479 {2480 const CbcLongCliqueBranchingObject* br =2481 dynamic_cast<const CbcLongCliqueBranchingObject*>(brObj);2482 assert(br);2483 return CbcCompareCliques(clique_, br>clique_);2484 }2485 2486 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the2487 same type and must have the same original object, but they may have2488 different feasible regions.2489 Return the appropriate CbcRangeCompare value (first argument being the2490 sub/superset if that's the case). In case of overlap (and if \c2491 replaceIfOverlap is true) replace the current branching object with one2492 whose feasible region is the overlap.2493 */2494 CbcRangeCompare2495 CbcLongCliqueBranchingObject::compareBranchingObject2496 (const CbcBranchingObject* brObj, const bool /*replaceIfOverlap*/)2497 {2498 const CbcLongCliqueBranchingObject* br =2499 dynamic_cast<const CbcLongCliqueBranchingObject*>(brObj);2500 assert(br);2501 const int numberMembers = clique_>numberMembers();2502 const int numberWords = (numberMembers + 31) >> 5;2503 unsigned int* thisMask = way_ < 0 ? upMask_ : downMask_;2504 const unsigned int* otherMask = br>way_ < 0 ? br>upMask_ : br>downMask_;2505 2506 if (memcmp(thisMask, otherMask, numberWords * sizeof(unsigned int)) == 0) {2507 return CbcRangeSame;2508 }2509 bool canBeSuperset = true;2510 bool canBeSubset = true;2511 int i;2512 for (i = numberWords  1; i >= 0 && (canBeSuperset  canBeSubset); i) {2513 const unsigned int both = (thisMask[i] & otherMask[i]);2514 canBeSuperset &= (both == thisMask[i]);2515 canBeSubset &= (both == otherMask[i]);2516 }2517 if (canBeSuperset) {2518 return CbcRangeSuperset;2519 }2520 if (canBeSubset) {2521 return CbcRangeSubset;2522 }2523 2524 for (i = numberWords  1; i >= 0; i) {2525 if ((thisMask[i] ^ otherMask[i]) != 0) {2526 break;2527 }2528 }2529 if (i == 1) { // complement2530 return CbcRangeDisjoint;2531 }2532 // must be overlap2533 for (i = numberWords  1; i >= 0; i) {2534 thisMask[i] = otherMask[i];2535 }2536 return CbcRangeOverlap;2537 }2538 2539 //##############################################################################2540 2541 // Default Constructor2542 CbcSOSBranchingObject::CbcSOSBranchingObject()2543 : CbcBranchingObject(),2544 firstNonzero_(1),2545 lastNonzero_(1)2546 {2547 set_ = NULL;2548 separator_ = 0.0;2549 }2550 2551 // Useful constructor2552 CbcSOSBranchingObject::CbcSOSBranchingObject (CbcModel * model,2553 const CbcSOS * set,2554 int way ,2555 double separator)2556 : CbcBranchingObject(model, set>id(), way, 0.5)2557 {2558 set_ = set;2559 separator_ = separator;2560 computeNonzeroRange();2561 }2562 2563 // Copy constructor2564 CbcSOSBranchingObject::CbcSOSBranchingObject (const CbcSOSBranchingObject & rhs)2565 : CbcBranchingObject(rhs),2566 firstNonzero_(rhs.firstNonzero_),2567 lastNonzero_(rhs.lastNonzero_)2568 {2569 set_ = rhs.set_;2570 separator_ = rhs.separator_;2571 }2572 2573 // Assignment operator2574 CbcSOSBranchingObject &2575 CbcSOSBranchingObject::operator=( const CbcSOSBranchingObject & rhs)2576 {2577 if (this != &rhs) {2578 CbcBranchingObject::operator=(rhs);2579 set_ = rhs.set_;2580 separator_ = rhs.separator_;2581 firstNonzero_ = rhs.firstNonzero_;2582 lastNonzero_ = rhs.lastNonzero_;2583 }2584 return *this;2585 }2586 CbcBranchingObject *2587 CbcSOSBranchingObject::clone() const2588 {2589 return (new CbcSOSBranchingObject(*this));2590 }2591 2592 2593 // Destructor2594 CbcSOSBranchingObject::~CbcSOSBranchingObject ()2595 {2596 }2597 2598 void2599 CbcSOSBranchingObject::computeNonzeroRange()2600 {2601 const int numberMembers = set_>numberMembers();2602 const double * weights = set_>weights();2603 int i = 0;2604 if (way_ < 0) {2605 for ( i = 0; i < numberMembers; i++) {2606 if (weights[i] > separator_)2607 break;2608 }2609 assert (i < numberMembers);2610 firstNonzero_ = 0;2611 lastNonzero_ = i;2612 } else {2613 for ( i = 0; i < numberMembers; i++) {2614 if (weights[i] >= separator_)2615 break;2616 }2617 assert (i < numberMembers);2618 firstNonzero_ = i;2619 lastNonzero_ = numberMembers;2620 }2621 }2622 2623 double2624 CbcSOSBranchingObject::branch()2625 {2626 decrementNumberBranchesLeft();2627 int numberMembers = set_>numberMembers();2628 const int * which = set_>members();2629 const double * weights = set_>weights();2630 OsiSolverInterface * solver = model_>solver();2631 //const double * lower = solver>getColLower();2632 //const double * upper = solver>getColUpper();2633 // *** for way  up means fix all those in down section2634 if (way_ < 0) {2635 int i;2636 for ( i = 0; i < numberMembers; i++) {2637 if (weights[i] > separator_)2638 break;2639 }2640 assert (i < numberMembers);2641 for (; i < numberMembers; i++)2642 solver>setColUpper(which[i], 0.0);2643 way_ = 1; // Swap direction2644 } else {2645 int i;2646 for ( i = 0; i < numberMembers; i++) {2647 if (weights[i] >= separator_)2648 break;2649 else2650 solver>setColUpper(which[i], 0.0);2651 }2652 assert (i < numberMembers);2653 way_ = 1; // Swap direction2654 }2655 computeNonzeroRange();2656 return 0.0;2657 }2658 /* Update bounds in solver as in 'branch' and update given bounds.2659 branchState is 1 for 'down' +1 for 'up' */2660 void2661 CbcSOSBranchingObject::fix(OsiSolverInterface * solver,2662 double * /*lower*/, double * upper,2663 int branchState) const2664 {2665 int numberMembers = set_>numberMembers();2666 const int * which = set_>members();2667 const double * weights = set_>weights();2668 //const double * lower = solver>getColLower();2669 //const double * upper = solver>getColUpper();2670 // *** for way  up means fix all those in down section2671 if (branchState < 0) {2672 int i;2673 for ( i = 0; i < numberMembers; i++) {2674 if (weights[i] > separator_)2675 break;2676 }2677 assert (i < numberMembers);2678 for (; i < numberMembers; i++) {2679 solver>setColUpper(which[i], 0.0);2680 upper[which[i]] = 0.0;2681 }2682 } else {2683 int i;2684 for ( i = 0; i < numberMembers; i++) {2685 if (weights[i] >= separator_) {2686 break;2687 } else {2688 solver>setColUpper(which[i], 0.0);2689 upper[which[i]] = 0.0;2690 }2691 }2692 assert (i < numberMembers);2693 }2694 }2695 // Print what would happen2696 void2697 CbcSOSBranchingObject::print()2698 {2699 int numberMembers = set_>numberMembers();2700 const int * which = set_>members();2701 const double * weights = set_>weights();2702 OsiSolverInterface * solver = model_>solver();2703 //const double * lower = solver>getColLower();2704 const double * upper = solver>getColUpper();2705 int first = numberMembers;2706 int last = 1;2707 int numberFixed = 0;2708 int numberOther = 0;2709 int i;2710 for ( i = 0; i < numberMembers; i++) {2711 double bound = upper[which[i]];2712 if (bound) {2713 first = CoinMin(first, i);2714 last = CoinMax(last, i);2715 }2716 }2717 // *** for way  up means fix all those in down section2718 if (way_ < 0) {2719 printf("SOS Down");2720 for ( i = 0; i < numberMembers; i++) {2721 double bound = upper[which[i]];2722 if (weights[i] > separator_)2723 break;2724 else if (bound)2725 numberOther++;2726 }2727 assert (i < numberMembers);2728 for (; i < numberMembers; i++) {2729 double bound = upper[which[i]];2730 if (bound)2731 numberFixed++;2732 }2733 } else {2734 printf("SOS Up");2735 for ( i = 0; i < numberMembers; i++) {2736 double bound = upper[which[i]];2737 if (weights[i] >= separator_)2738 break;2739 else if (bound)2740 numberFixed++;2741 }2742 assert (i < numberMembers);2743 for (; i < numberMembers; i++) {2744 double bound = upper[which[i]];2745 if (bound)2746 numberOther++;2747 }2748 }2749 printf("  at %g, free range %d (%g) => %d (%g), %d would be fixed, %d other way\n",2750 separator_, which[first], weights[first], which[last], weights[last], numberFixed, numberOther);2751 }2752 2753 /** Compare the original object of \c this with the original object of \c2754 brObj. Assumes that there is an ordering of the original objects.2755 This method should be invoked only if \c this and brObj are of the same2756 type.2757 Return negative/0/positive depending on whether \c this is2758 smaller/same/larger than the argument.2759 */2760 int2761 CbcSOSBranchingObject::compareOriginalObject2762 (const CbcBranchingObject* brObj) const2763 {2764 const CbcSOSBranchingObject* br =2765 dynamic_cast<const CbcSOSBranchingObject*>(brObj);2766 assert(br);2767 const CbcSOS* s0 = set_;2768 const CbcSOS* s1 = br>set_;2769 if (s0>sosType() != s1>sosType()) {2770 return s0>sosType()  s1>sosType();2771 }2772 if (s0>numberMembers() != s1>numberMembers()) {2773 return s0>numberMembers()  s1>numberMembers();2774 }2775 const int memberCmp = memcmp(s0>members(), s1>members(),2776 s0>numberMembers() * sizeof(int));2777 if (memberCmp != 0) {2778 return memberCmp;2779 }2780 return memcmp(s0>weights(), s1>weights(),2781 s0>numberMembers() * sizeof(double));2782 }2783 2784 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the2785 same type and must have the same original object, but they may have2786 different feasible regions.2787 Return the appropriate CbcRangeCompare value (first argument being the2788 sub/superset if that's the case). In case of overlap (and if \c2789 replaceIfOverlap is true) replace the current branching object with one2790 whose feasible region is the overlap.2791 */2792 CbcRangeCompare2793 CbcSOSBranchingObject::compareBranchingObject2794 (const CbcBranchingObject* brObj, const bool replaceIfOverlap)2795 {2796 const CbcSOSBranchingObject* br =2797 dynamic_cast<const CbcSOSBranchingObject*>(brObj);2798 assert(br);2799 if (firstNonzero_ < br>firstNonzero_) {2800 if (lastNonzero_ >= br>lastNonzero_) {2801 return CbcRangeSuperset;2802 } else if (lastNonzero_ <= br>firstNonzero_) {2803 return CbcRangeDisjoint;2804 } else {2805 // overlap2806 if (replaceIfOverlap) {2807 firstNonzero_ = br>firstNonzero_;2808 }2809 return CbcRangeOverlap;2810 }2811 } else if (firstNonzero_ > br>firstNonzero_) {2812 if (lastNonzero_ <= br>lastNonzero_) {2813 return CbcRangeSubset;2814 } else if (firstNonzero_ >= br>lastNonzero_) {2815 return CbcRangeDisjoint;2816 } else {2817 // overlap2818 if (replaceIfOverlap) {2819 lastNonzero_ = br>lastNonzero_;2820 }2821 return CbcRangeOverlap;2822 }2823 } else {2824 if (lastNonzero_ == br>lastNonzero_) {2825 return CbcRangeSame;2826 }2827 return lastNonzero_ < br>lastNonzero_ ? CbcRangeSubset : CbcRangeSuperset;2828 }2829 return CbcRangeSame; // fake return2830 }2831 2832 //##############################################################################2833 2834 // Default Constructor2835 CbcBranchDefaultDecision::CbcBranchDefaultDecision()2836 : CbcBranchDecision()2837 {2838 bestCriterion_ = 0.0;2839 bestChangeUp_ = 0.0;2840 bestNumberUp_ = 0;2841 bestChangeDown_ = 0.0;2842 bestObject_ = NULL;2843 bestNumberDown_ = 0;2844 }2845 2846 // Copy constructor2847 CbcBranchDefaultDecision::CbcBranchDefaultDecision (2848 const CbcBranchDefaultDecision & rhs)2849 : CbcBranchDecision(rhs)2850 {2851 bestCriterion_ = rhs.bestCriterion_;2852 bestChangeUp_ = rhs.bestChangeUp_;2853 bestNumberUp_ = rhs.bestNumberUp_;2854 bestChangeDown_ = rhs.bestChangeDown_;2855 bestNumberDown_ = rhs.bestNumberDown_;2856 bestObject_ = rhs.bestObject_;2857 model_ = rhs.model_;2858 }2859 2860 CbcBranchDefaultDecision::~CbcBranchDefaultDecision()2861 {2862 }2863 2864 // Clone2865 CbcBranchDecision *2866 CbcBranchDefaultDecision::clone() const2867 {2868 return new CbcBranchDefaultDecision(*this);2869 }2870 2871 // Initialize i.e. before start of choosing at a node2872 void2873 CbcBranchDefaultDecision::initialize(CbcModel * model)2874 {2875 bestCriterion_ = 0.0;2876 bestChangeUp_ = 0.0;2877 bestNumberUp_ = 0;2878 bestChangeDown_ = 0.0;2879 bestNumberDown_ = 0;2880 bestObject_ = NULL;2881 model_ = model;2882 }2883 2884 2885 /*2886 Simple default decision algorithm. Compare based on infeasibility (numInfUp,2887 numInfDn) until a solution is found by search, then switch to change in2888 objective (changeUp, changeDn). Note that bestSoFar is remembered in2889 bestObject_, so the parameter bestSoFar is unused.2890 */2891 2892 int2893 CbcBranchDefaultDecision::betterBranch(CbcBranchingObject * thisOne,2894 CbcBranchingObject * /*bestSoFar*/,2895 double changeUp, int numInfUp,2896 double changeDn, int numInfDn)2897 {2898 bool beforeSolution = cbcModel()>getSolutionCount() ==2899 cbcModel()>getNumberHeuristicSolutions();;2900 int betterWay = 0;2901 if (beforeSolution) {2902 if (!bestObject_) {2903 bestNumberUp_ = COIN_INT_MAX;2904 bestNumberDown_ = COIN_INT_MAX;2905 }2906 // before solution  choose smallest number2907 // could add in depth as well2908 int bestNumber = CoinMin(bestNumberUp_, bestNumberDown_);2909 if (numInfUp < numInfDn) {2910 if (numInfUp < bestNumber) {2911 betterWay = 1;2912 } else if (numInfUp == bestNumber) {2913 if (changeUp < bestCriterion_)2914 betterWay = 1;2915 }2916 } else if (numInfUp > numInfDn) {2917 if (numInfDn < bestNumber) {2918 betterWay = 1;2919 } else if (numInfDn == bestNumber) {2920 if (changeDn < bestCriterion_)2921 betterWay = 1;2922 }2923 } else {2924 // up and down have same number2925 bool better = false;2926 if (numInfUp < bestNumber) {2927 better = true;2928 } else if (numInfUp == bestNumber) {2929 if (CoinMin(changeUp, changeDn) < bestCriterion_)2930 better = true;;2931 }2932 if (better) {2933 // see which way2934 if (changeUp <= changeDn)2935 betterWay = 1;2936 else2937 betterWay = 1;2938 }2939 }2940 } else {2941 if (!bestObject_) {2942 bestCriterion_ = 1.0;2943 }2944 // got a solution2945 if (changeUp <= changeDn) {2946 if (changeUp > bestCriterion_)2947 betterWay = 1;2948 } else {2949 if (changeDn > bestCriterion_)2950 betterWay = 1;2951 }2952 }2953 if (betterWay) {2954 bestCriterion_ = CoinMin(changeUp, changeDn);2955 bestChangeUp_ = changeUp;2956 bestNumberUp_ = numInfUp;2957 bestChangeDown_ = changeDn;2958 bestNumberDown_ = numInfDn;2959 bestObject_ = thisOne;2960 // See if user is overriding way2961 if (thisOne>object() && thisOne>object()>preferredWay())2962 betterWay = thisOne>object()>preferredWay();2963 }2964 return betterWay;2965 }2966 /* Sets or gets best criterion so far */2967 void2968 CbcBranchDefaultDecision::setBestCriterion(double value)2969 {2970 bestCriterion_ = value;2971 }2972 double2973 CbcBranchDefaultDecision::getBestCriterion() const2974 {2975 return bestCriterion_;2976 }2977 2978 /* Compare N branching objects. Return index of best2979 and sets way of branching in chosen object.2980 2981 This routine is used only after strong branching.2982 */2983 2984 int2985 CbcBranchDefaultDecision::bestBranch (CbcBranchingObject ** objects, int numberObjects,2986 int numberUnsatisfied,2987 double * changeUp, int * numberInfeasibilitiesUp,2988 double * changeDown, int * numberInfeasibilitiesDown,2989 double objectiveValue)2990 {2991 2992 int bestWay = 0;2993 int whichObject = 1;2994 if (numberObjects) {2995 CbcModel * model = cbcModel();2996 // at continuous2997 //double continuousObjective = model>getContinuousObjective();2998 //int continuousInfeasibilities = model>getContinuousInfeasibilities();2999 3000 // average cost to get rid of infeasibility3001 //double averageCostPerInfeasibility =3002 //(objectiveValuecontinuousObjective)/3003 //(double) (abs(continuousInfeasibilitiesnumberUnsatisfied)+1);3004 /* beforeSolution is :3005 0  before any solution3006 n  n heuristic solutions but no branched one3007 1  branched solution found3008 */3009 int numberSolutions = model>getSolutionCount();3010 double cutoff = model>getCutoff();3011 int method = 0;3012 int i;3013 if (numberSolutions) {3014 int numberHeuristic = model>getNumberHeuristicSolutions();3015 if (numberHeuristic < numberSolutions) {3016 method = 1;3017 } else {3018 method = 2;3019 // look further3020 for ( i = 0 ; i < numberObjects ; i++) {3021 int numberNext = numberInfeasibilitiesUp[i];3022 3023 if (numberNext < numberUnsatisfied) {3024 int numberUp = numberUnsatisfied  numberInfeasibilitiesUp[i];3025 double perUnsatisfied = changeUp[i] / static_cast<double> (numberUp);3026 double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;3027 if (estimatedObjective < cutoff)3028 method = 3;3029 }3030 numberNext = numberInfeasibilitiesDown[i];3031 if (numberNext < numberUnsatisfied) {3032 int numberDown = numberUnsatisfied  numberInfeasibilitiesDown[i];3033 double perUnsatisfied = changeDown[i] / static_cast<double> (numberDown);3034 double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;3035 if (estimatedObjective < cutoff)3036 method = 3;3037 }3038 }3039 }3040 method = 2;3041 } else {3042 method = 0;3043 }3044 // Uncomment next to force method 43045 //method=4;3046 /* Methods :3047 0  fewest infeasibilities3048 1  largest min change in objective3049 2  as 1 but use sum of changes if min close3050 3  predicted best solution3051 4  take cheapest up branch if infeasibilities same3052 */3053 int bestNumber = COIN_INT_MAX;3054 double bestCriterion = 1.0e50;3055 double alternativeCriterion = 1.0;3056 double bestEstimate = 1.0e100;3057 switch (method) {3058 case 0:3059 // could add in depth as well3060 for ( i = 0 ; i < numberObjects ; i++) {3061 int thisNumber = CoinMin(numberInfeasibilitiesUp[i], numberInfeasibilitiesDown[i]);3062 if (thisNumber <= bestNumber) {3063 int betterWay = 0;3064 if (numberInfeasibilitiesUp[i] < numberInfeasibilitiesDown[i]) {3065 if (numberInfeasibilitiesUp[i] < bestNumber) {3066 betterWay = 1;3067 } else {3068 if (changeUp[i] < bestCriterion)3069 betterWay = 1;3070 }3071 } else if (numberInfeasibilitiesUp[i] > numberInfeasibilitiesDown[i]) {3072 if (numberInfeasibilitiesDown[i] < bestNumber) {3073 betterWay = 1;3074 } else {3075 if (changeDown[i] < bestCriterion)3076 betterWay = 1;3077 }3078 } else {3079 // up and down have same number3080 bool better = false;3081 if (numberInfeasibilitiesUp[i] < bestNumber) {3082 better = true;3083 } else if (numberInfeasibilitiesUp[i] == bestNumber) {3084 if (CoinMin(changeUp[i], changeDown[i]) < bestCriterion)3085 better = true;;3086 }3087 if (better) {3088 // see which way3089 if (changeUp[i] <= changeDown[i])3090 betterWay = 1;3091 else3092 betterWay = 1;3093 }3094 }3095 if (betterWay) {3096 bestCriterion = CoinMin(changeUp[i], changeDown[i]);3097 bestNumber = thisNumber;3098 whichObject = i;3099 bestWay = betterWay;3100 }3101 }3102 }3103 break;3104 case 1:3105 for ( i = 0 ; i < numberObjects ; i++) {3106 int betterWay = 0;3107 if (changeUp[i] <= changeDown[i]) {3108 if (changeUp[i] > bestCriterion)3109 betterWay = 1;3110 } else {3111 if (changeDown[i] > bestCriterion)3112 betterWay = 1;3113 }3114 if (betterWay) {3115 bestCriterion = CoinMin(changeUp[i], changeDown[i]);3116 whichObject = i;3117 bestWay = betterWay;3118 }3119 }3120 break;3121 case 2:3122 for ( i = 0 ; i < numberObjects ; i++) {3123 double change = CoinMin(changeUp[i], changeDown[i]);3124 double sum = changeUp[i] + changeDown[i];3125 bool take = false;3126 if (change > 1.1*bestCriterion)3127 take = true;3128 else if (change > 0.9*bestCriterion && sum + change > bestCriterion + alternativeCriterion)3129 take = true;3130 if (take) {3131 if (changeUp[i] <= changeDown[i]) {3132 if (changeUp[i] > bestCriterion)3133 bestWay = 1;3134 } else {3135 if (changeDown[i] > bestCriterion)3136 bestWay = 1;3137 }3138 bestCriterion = change;3139 alternativeCriterion = sum;3140 whichObject = i;3141 }3142 }3143 break;3144 case 3:3145 for ( i = 0 ; i < numberObjects ; i++) {3146 int numberNext = numberInfeasibilitiesUp[i];3147 3148 if (numberNext < numberUnsatisfied) {3149 int numberUp = numberUnsatisfied  numberInfeasibilitiesUp[i];3150 double perUnsatisfied = changeUp[i] / static_cast<double> (numberUp);3151 double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;3152 if (estimatedObjective < bestEstimate) {3153 bestEstimate = estimatedObjective;3154 bestWay = 1;3155 whichObject = i;3156 }3157 }3158 numberNext = numberInfeasibilitiesDown[i];3159 if (numberNext < numberUnsatisfied) {3160 int numberDown = numberUnsatisfied  numberInfeasibilitiesDown[i];3161 double perUnsatisfied = changeDown[i] / static_cast<double> (numberDown);3162 double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;3163 if (estimatedObjective < bestEstimate) {3164 bestEstimate = estimatedObjective;3165 bestWay = 1;3166 whichObject = i;3167 }3168 }3169 }3170 break;3171 case 4:3172 // if number infeas same then cheapest up3173 // first get best number or when going down3174 // now choose smallest change up amongst equal number infeas3175 for ( i = 0 ; i < numberObjects ; i++) {3176 int thisNumber = CoinMin(numberInfeasibilitiesUp[i], numberInfeasibilitiesDown[i]);3177 if (thisNumber <= bestNumber) {3178 int betterWay = 0;3179 if (numberInfeasibilitiesUp[i] < numberInfeasibilitiesDown[i]) {3180 if (numberInfeasibilitiesUp[i] < bestNumber) {3181 betterWay = 1;3182 } else {3183 if (changeUp[i] < bestCriterion)3184 betterWay = 1;3185 }3186 } else if (numberInfeasibilitiesUp[i] > numberInfeasibilitiesDown[i]) {3187 if (numberInfeasibilitiesDown[i] < bestNumber) {3188 betterWay = 1;3189 } else {3190 if (changeDown[i] < bestCriterion)3191 betterWay = 1;3192 }3193 } else {3194 // up and down have same number3195 bool better = false;3196 if (numberInfeasibilitiesUp[i] < bestNumber) {3197 better = true;3198 } else if (numberInfeasibilitiesUp[i] == bestNumber) {3199 if (CoinMin(changeUp[i], changeDown[i]) < bestCriterion)3200 better = true;;3201 }3202 if (better) {3203 // see which way3204 if (changeUp[i] <= changeDown[i])3205 betterWay = 1;3206 else3207 betterWay = 1;3208 }3209 }3210 if (betterWay) {3211 bestCriterion = CoinMin(changeUp[i], changeDown[i]);3212 bestNumber = thisNumber;3213 whichObject = i;3214 bestWay = betterWay;3215 }3216 }3217 }3218 bestCriterion = 1.0e50;3219 for ( i = 0 ; i < numberObjects ; i++) {3220 int thisNumber = numberInfeasibilitiesUp[i];3221 if (thisNumber == bestNumber && changeUp) {3222 if (changeUp[i] < bestCriterion) {3223 bestCriterion = changeUp[i];3224 whichObject = i;3225 bestWay = 1;3226 }3227 }3228 }3229 break;3230 }3231 // set way in best3232 if (whichObject >= 0) {3233 CbcBranchingObject * bestObject = objects[whichObject];3234 if (bestObject>object() && bestObject>object()>preferredWay())3235 bestWay = bestObject>object()>preferredWay();3236 bestObject>way(bestWay);3237 } else {3238 printf("debug\n");3239 }3240 }3241 return whichObject;3242 }3243 3244 //##############################################################################3245 3246 // Default Constructor3247 CbcFollowOn::CbcFollowOn ()3248 : CbcObject(),3249 rhs_(NULL)3250 {3251 }3252 3253 // Useful constructor3254 CbcFollowOn::CbcFollowOn (CbcModel * model)3255 : CbcObject(model)3256 {3257 assert (model);3258 OsiSolverInterface * solver = model_>solver();3259 matrix_ = *solver>getMatrixByCol();3260 matrix_.removeGaps();3261 matrix_.setExtraGap(0.0);3262 matrixByRow_ = *solver>getMatrixByRow();3263 int numberRows = matrix_.getNumRows();3264 3265 rhs_ = new int[numberRows];3266 int i;3267 const double * rowLower = solver>getRowLower();3268 const double * rowUpper = solver>getRowUpper();3269 // Row copy3270 const double * elementByRow = matrixByRow_.getElements();3271 const int * column = matrixByRow_.getIndices();3272 const CoinBigIndex * rowStart = matrixByRow_.getVectorStarts();3273 const int * rowLength = matrixByRow_.getVectorLengths();3274 for (i = 0; i < numberRows; i++) {3275 rhs_[i] = 0;3276 double value = rowLower[i];3277 if (value == rowUpper[i]) {3278 if (floor(value) == value && value >= 1.0 && value < 10.0) {3279 // check elements3280 bool good = true;3281 for (int j = rowStart[i]; j < rowStart[i] + rowLength[i]; j++) {3282 int iColumn = column[j];3283 if (!solver>isBinary(iColumn))3284 good = false;3285 double elValue = elementByRow[j];3286 if (floor(elValue) != elValue  value < 1.0)3287 good = false;3288 }3289 if (good)3290 rhs_[i] = static_cast<int> (value);3291 }3292 }3293 }3294 }3295 3296 // Copy constructor3297 CbcFollowOn::CbcFollowOn ( const CbcFollowOn & rhs)3298 : CbcObject(rhs),3299 matrix_(rhs.matrix_),3300 matrixByRow_(rhs.matrixByRow_)3301 {3302 int numberRows = matrix_.getNumRows();3303 rhs_ = CoinCopyOfArray(rhs.rhs_, numberRows);3304 }3305 3306 // Clone3307 CbcObject *3308 CbcFollowOn::clone() const3309 {3310 return new CbcFollowOn(*this);3311 }3312 3313 // Assignment operator3314 CbcFollowOn &3315 CbcFollowOn::operator=( const CbcFollowOn & rhs)3316 {3317 if (this != &rhs) {3318 CbcObject::operator=(rhs);3319 delete [] rhs_;3320 matrix_ = rhs.matrix_;3321 matrixByRow_ = rhs.matrixByRow_;3322 int numberRows = matrix_.getNumRows();3323 rhs_ = CoinCopyOfArray(rhs.rhs_, numberRows);3324 }3325 return *this;3326 }3327 3328 // Destructor3329 CbcFollowOn::~CbcFollowOn ()3330 {3331 delete [] rhs_;3332 }3333 // As some computation is needed in more than one place  returns row3334 int3335 CbcFollowOn::gutsOfFollowOn(int & otherRow, int & preferredWay) const3336 {3337 int whichRow = 1;3338 otherRow = 1;3339 int numberRows = matrix_.getNumRows();3340 3341 int i;3342 // For sorting3343 int * sort = new int [numberRows];3344 int * isort = new int [numberRows];3345 // Column copy3346 //const double * element = matrix_.getElements();3347 const int * row = matrix_.getIndices();3348 const CoinBigIndex * columnStart = matrix_.getVectorStarts();3349 const int * columnLength = matrix_.getVectorLengths();3350 // Row copy3351 const double * elementByRow = matrixByRow_.getElements();3352 const int * column = matrixByRow_.getIndices();3353 const CoinBigIndex * rowStart = matrixByRow_.getVectorStarts();3354 const int * rowLength = matrixByRow_.getVectorLengths();3355 OsiSolverInterface * solver = model_>solver();3356 const double * columnLower = solver>getColLower();3357 const double * columnUpper = solver>getColUpper();3358 const double * solution = solver>getColSolution();3359 double integerTolerance = model_>getDblParam(CbcModel::CbcIntegerTolerance);3360 int nSort = 0;3361 for (i = 0; i < numberRows; i++) {3362 if (rhs_[i]) {3363 // check elements3364 double smallest = 1.0e10;3365 double largest = 0.0;3366 int rhsValue = rhs_[i];3367 int number1 = 0;3368 int numberUnsatisfied = 0;3369 for (int j = rowStart[i]; j < rowStart[i] + rowLength[i]; j++) {3370 int iColumn = column[j];3371 double value = elementByRow[j];3372 double solValue = solution[iColumn];3373 if (columnLower[iColumn] != columnUpper[iColumn]) {3374 smallest = CoinMin(smallest, value);3375 largest = CoinMax(largest, value);3376 if (value == 1.0)3377 number1++;3378 if (solValue < 1.0  integerTolerance && solValue > integerTolerance)3379 numberUnsatisfied++;3380 } else {3381 rhsValue = static_cast<int>(value * floor(solValue + 0.5));3382 }3383 }3384 if (numberUnsatisfied > 1) {3385 if (smallest < largest) {3386 // probably no good but check a few things3387 assert (largest <= rhsValue);3388 if (number1 == 1 && largest == rhsValue)3389 printf("could fix\n");3390 } else if (largest == rhsValue) {3391 sort[nSort] = i;3392 isort[nSort++] = numberUnsatisfied;3393 }3394 }3395 }3396 }3397 if (nSort > 1) {3398 CoinSort_2(isort, isort + nSort, sort);3399 CoinZeroN(isort, numberRows);3400 double * other = new double[numberRows];3401 CoinZeroN(other, numberRows);3402 int * which = new int[numberRows];3403 //#define COUNT3404 #ifndef COUNT3405 bool beforeSolution = model_>getSolutionCount() == 0;3406 #endif3407 for (int k = 0; k < nSort  1; k++) {3408 i = sort[k];3409 int numberUnsatisfied = 0;3410 int n = 0;3411 int j;3412 for (j = rowStart[i]; j < rowStart[i] + rowLength[i]; j++) {3413 int iColumn = column[j];3414 if (columnLower[iColumn] != columnUpper[iColumn]) {3415 double solValue = solution[iColumn]  columnLower[iColumn];3416 if (solValue < 1.0  integerTolerance && solValue > integerTolerance) {3417 numberUnsatisfied++;3418 for (int jj = columnStart[iColumn]; jj < columnStart[iColumn] + columnLength[iColumn]; jj++) {3419 int iRow = row[jj];3420 if (rhs_[iRow]) {3421 other[iRow] += solValue;3422 if (isort[iRow]) {3423 isort[iRow]++;3424 } else {3425 isort[iRow] = 1;3426 which[n++] = iRow;3427 }3428 }3429 }3430 }3431 }3432 }3433 double total = 0.0;3434 // Take out row3435 double sumThis = other[i];3436 other[i] = 0.0;3437 assert (numberUnsatisfied == isort[i]);3438 // find one nearest half if solution, one if before solution3439 int iBest = 1;3440 double dtarget = 0.5 * total;3441 #ifdef COUNT3442 int target = (numberUnsatisfied + 1) >> 1;3443 int best = numberUnsatisfied;3444 #else3445 double best;3446 if (beforeSolution)3447 best = dtarget;3448 else3449 best = 1.0e30;3450 #endif3451 for (j = 0; j < n; j++) {3452 int iRow = which[j];3453 double dvalue = other[iRow];3454 other[iRow] = 0.0;3455 #ifdef COUNT3456 int value = isort[iRow];3457 #endif3458 isort[iRow] = 0;3459 if (fabs(dvalue) < 1.0e8  fabs(sumThis  dvalue) < 1.0e8)3460 continue;3461 if (dvalue < integerTolerance  dvalue > 1.0  integerTolerance)3462 continue;3463 #ifdef COUNT3464 if (abs(value  target) < best && value != numberUnsatisfied) {3465 best = abs(value  target);3466 iBest = iRow;3467 if (dvalue < dtarget)3468 preferredWay = 1;3469 else3470 preferredWay = 1;3471 }3472 #else3473 if (beforeSolution) {3474 if (fabs(dvalue  dtarget) > best) {3475 best = fabs(dvalue  dtarget);3476 iBest = iRow;3477 if (dvalue < dtarget)3478 preferredWay = 1;3479 else3480 preferredWay = 1;3481 }3482 } else {3483 if (fabs(dvalue  dtarget) < best) {3484 best = fabs(dvalue  dtarget);3485 iBest = iRow;3486 if (dvalue < dtarget)3487 preferredWay = 1;3488 else3489 preferredWay = 1;3490 }3491 }3492 #endif3493 }3494 if (iBest >= 0) {3495 whichRow = i;3496 otherRow = iBest;3497 break;3498 }3499 }3500 delete [] which;3501 delete [] other;3502 }3503 delete [] sort;3504 delete [] isort;3505 return whichRow;3506 }3507 double3508 CbcFollowOn::infeasibility(const OsiBranchingInformation * /*info*/,3509 int &preferredWay) const3510 {3511 int otherRow = 0;3512 int whichRow = gutsOfFollowOn(otherRow, preferredWay);3513 if (whichRow < 0)3514 return 0.0;3515 else3516 return 2.0* model_>getDblParam(CbcModel::CbcIntegerTolerance);3517 }3518 3519 // This looks at solution and sets bounds to contain solution3520 void3521 CbcFollowOn::feasibleRegion()3522 {3523 }3524 3525 CbcBranchingObject *3526 CbcFollowOn::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int way)3527 {3528 int otherRow = 0;3529 int preferredWay;3530 int whichRow = gutsOfFollowOn(otherRow, preferredWay);3531 assert(way == preferredWay);3532 assert (whichRow >= 0);3533 int numberColumns = matrix_.getNumCols();3534 3535 // Column copy3536 //const double * element = matrix_.getElements();3537 const int * row = matrix_.getIndices();3538 const CoinBigIndex * columnStart = matrix_.getVectorStarts();3539 const int * columnLength = matrix_.getVectorLengths();3540 // Row copy3541 //const double * elementByRow = matrixByRow_.getElements();3542 const int * column = matrixByRow_.getIndices();3543 const CoinBigIndex * rowStart = matrixByRow_.getVectorStarts();3544 const int * rowLength = matrixByRow_.getVectorLengths();3545 //OsiSolverInterface * solver = model_>solver();3546 const double * columnLower = solver>getColLower();3547 const double * columnUpper = solver>getColUpper();3548 //const double * solution = solver>getColSolution();3549 int nUp = 0;3550 int nDown = 0;3551 int * upList = new int[numberColumns];3552 int * downList = new int[numberColumns];3553 int j;3554 for (j = rowStart[whichRow]; j < rowStart[whichRow] + rowLength[whichRow]; j++) {3555 int iColumn = column[j];3556 if (columnLower[iColumn] != columnUpper[iColumn]) {3557 bool up = true;3558 for (int jj = columnStart[iColumn]; jj < columnStart[iColumn] + columnLength[iColumn]; jj++) {3559 int iRow = row[jj];3560 if (iRow == otherRow) {3561 up = false;3562 break;3563 }3564 }3565 if (up)3566 upList[nUp++] = iColumn;3567 else3568 downList[nDown++] = iColumn;3569 }3570 }3571 //printf("way %d\n",way);3572 // create object3573 //printf("would fix %d down and %d up\n",nDown,nUp);3574 CbcBranchingObject * branch3575 = new CbcFixingBranchingObject(model_, way,3576 nDown, downList, nUp, upList);3577 delete [] upList;3578 delete [] downList;3579 return branch;3580 }3581 3582 //##############################################################################3583 3584 // Default Constructor3585 CbcFixingBranchingObject::CbcFixingBranchingObject()3586 : CbcBranchingObject()3587 {3588 numberDown_ = 0;3589 numberUp_ = 0;3590 downList_ = NULL;3591 upList_ = NULL;3592 }3593 3594 // Useful constructor3595 CbcFixingBranchingObject::CbcFixingBranchingObject (CbcModel * model,3596 int way ,3597 int numberOnDownSide, const int * down,3598 int numberOnUpSide, const int * up)3599 : CbcBranchingObject(model, 0, way, 0.5)3600 {3601 numberDown_ = numberOnDownSide;3602 numberUp_ = numberOnUpSide;3603 downList_ = CoinCopyOfArray(down, numberDown_);3604 upList_ = CoinCopyOfArray(up, numberUp_);3605 }3606 3607 // Copy constructor3608 CbcFixingBranchingObject::CbcFixingBranchingObject ( const CbcFixingBranchingObject & rhs) : CbcBranchingObject(rhs)3609 {3610 numberDown_ = rhs.numberDown_;3611 numberUp_ = rhs.numberUp_;3612 downList_ = CoinCopyOfArray(rhs.downList_, numberDown_);3613 upList_ = CoinCopyOfArray(rhs.upList_, numberUp_);3614 }3615 3616 // Assignment operator3617 CbcFixingBranchingObject &3618 CbcFixingBranchingObject::operator=( const CbcFixingBranchingObject & rhs)3619 {3620 if (this != &rhs) {3621 CbcBranchingObject::operator=(rhs);3622 delete [] downList_;3623 delete [] upList_;3624 numberDown_ = rhs.numberDown_;3625 numberUp_ = rhs.numberUp_;3626 downList_ = CoinCopyOfArray(rhs.downList_, numberDown_);3627 upList_ = CoinCopyOfArray(rhs.upList_, numberUp_);3628 }3629 return *this;3630 }3631 CbcBranchingObject *3632 CbcFixingBranchingObject::clone() const3633 {3634 return (new CbcFixingBranchingObject(*this));3635 }3636 3637 3638 // Destructor3639 CbcFixingBranchingObject::~CbcFixingBranchingObject ()3640 {3641 delete [] downList_;3642 delete [] upList_;3643 }3644 double3645 CbcFixingBranchingObject::branch()3646 {3647 decrementNumberBranchesLeft();3648 OsiSolverInterface * solver = model_>solver();3649 const double * columnLower = solver>getColLower();3650 int i;3651 // *** for way  up means fix all those in up section3652 if (way_ < 0) {3653 #ifdef FULL_PRINT3654 printf("Down Fix ");3655 #endif3656 //printf("Down Fix %d\n",numberDown_);3657 for (i = 0; i < numberDown_; i++) {3658 int iColumn = downList_[i];3659 model_>solver()>setColUpper(iColumn, columnLower[iColumn]);3660 #ifdef FULL_PRINT3661 printf("Setting bound on %d to lower bound\n", iColumn);3662 #endif3663 }3664 way_ = 1; // Swap direction3665 } else {3666 #ifdef FULL_PRINT3667 printf("Up Fix ");3668 #endif3669 //printf("Up Fix %d\n",numberUp_);3670 for (i = 0; i < numberUp_; i++) {3671 int iColumn = upList_[i];3672 model_>solver()>setColUpper(iColumn, columnLower[iColumn]);3673 #ifdef FULL_PRINT3674 printf("Setting bound on %d to lower bound\n", iColumn);3675 #endif3676 }3677 way_ = 1; // Swap direction3678 }3679 #ifdef FULL_PRINT3680 printf("\n");3681 #endif3682 return 0.0;3683 }3684 void3685 CbcFixingBranchingObject::print()3686 {3687 int i;3688 // *** for way  up means fix all those in up section3689 if (way_ < 0) {3690 printf("Down Fix ");3691 for (i = 0; i < numberDown_; i++) {3692 int iColumn = downList_[i];3693 printf("%d ", iColumn);3694 }3695 } else {3696 printf("Up Fix ");3697 for (i = 0; i < numberUp_; i++) {3698 int iColumn = upList_[i];3699 printf("%d ", iColumn);3700 }3701 }3702 printf("\n");3703 }3704 3705 /** Compare the original object of \c this with the original object of \c3706 brObj. Assumes that there is an ordering of the original objects.3707 This method should be invoked only if \c this and brObj are of the same3708 type.3709 Return negative/0/positive depending on whether \c this is3710 smaller/same/larger than the argument.3711 */3712 int3713 CbcFixingBranchingObject::compareOriginalObject3714 (const CbcBranchingObject* /*brObj*/) const3715 {3716 throw("must implement");3717 }3718 3719 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the3720 same type and must have the same original object, but they may have3721 different feasible regions.3722 Return the appropriate CbcRangeCompare value (first argument being the3723 sub/superset if that's the case). In case of overlap (and if \c3724 replaceIfOverlap is true) replace the current branching object with one3725 whose feasible region is the overlap.3726 */3727 CbcRangeCompare3728 CbcFixingBranchingObject::compareBranchingObject3729 (const CbcBranchingObject* /*brObj*/, const bool /*replaceIfOverlap*/)3730 {3731 #if 0 //ndef NDEBUG3732 const CbcFixingBranchingObject* br =3733 dynamic_cast<const CbcFixingBranchingObject*>(brObj);3734 assert(br);3735 #endif3736 // If two FixingBranchingObject's have the same base object then it's pretty3737 // much guaranteed3738 throw("must implement");3739 }3740 3741 //##############################################################################3742 3743 // Default Constructor3744 CbcNWay::CbcNWay ()3745 : CbcObject(),3746 numberMembers_(0),3747 members_(NULL),3748 consequence_(NULL)3749 {3750 }3751 3752 // Useful constructor (which are integer indices)3753 CbcNWay::CbcNWay (CbcModel * model, int numberMembers,3754 const int * which, int identifier)3755 : CbcObject(model)3756 {3757 id_ = identifier;3758 numberMembers_ = numberMembers;3759 if (numberMembers_) {3760 members_ = new int[numberMembers_];3761 memcpy(members_, which, numberMembers_*sizeof(int));3762 } else {3763 members_ = NULL;3764 }3765 consequence_ = NULL;3766 }3767 3768 // Copy constructor3769 CbcNWay::CbcNWay ( const CbcNWay & rhs)3770 : CbcObject(rhs)3771 {3772 numberMembers_ = rhs.numberMembers_;3773 consequence_ = NULL;3774 if (numberMembers_) {3775 members_ = new int[numberMembers_];3776 memcpy(members_, rhs.members_, numberMembers_*sizeof(int));3777 if (rhs.consequence_) {3778 consequence_ = new CbcConsequence * [numberMembers_];3779 for (int i = 0; i < numberMembers_; i++) {3780 if (rhs.consequence_[i])3781 consequence_[i] = rhs.consequence_[i]>clone();3782 else3783 consequence_[i] = NULL;3784 }3785 }3786 } else {3787 members_ = NULL;3788 }3789 }3790 3791 // Clone3792 CbcObject *3793 CbcNWay::clone() const3794 {3795 return new CbcNWay(*this);3796 }3797 3798 // Assignment operator3799 CbcNWay &3800 CbcNWay::operator=( const CbcNWay & rhs)3801 {3802 if (this != &rhs) {3803 CbcObject::operator=(rhs);3804 delete [] members_;3805 numberMembers_ = rhs.numberMembers_;3806 if (consequence_) {3807 for (int i = 0; i < numberMembers_; i++)3808 delete consequence_[i];3809 delete [] consequence_;3810 consequence_ = NULL;3811 }3812 if (numberMembers_) {3813 members_ = new int[numberMembers_];3814 memcpy(members_, rhs.members_, numberMembers_*sizeof(int));3815 } else {3816 members_ = NULL;3817 }3818 if (rhs.consequence_) {3819 consequence_ = new CbcConsequence * [numberMembers_];3820 for (int i = 0; i < numberMembers_; i++) {3821 if (rhs.consequence_[i])3822 consequence_[i] = rhs.consequence_[i]>clone();3823 else3824 consequence_[i] = NULL;3825 }3826 }3827 }3828 return *this;3829 }3830 3831 // Destructor3832 CbcNWay::~CbcNWay ()3833 {3834 delete [] members_;3835 if (consequence_) {3836 for (int i = 0; i < numberMembers_; i++)3837 delete consequence_[i];3838 delete [] consequence_;3839 }3840 }3841 // Set up a consequence for a single member3842 void3843 CbcNWay::setConsequence(int iColumn, const CbcConsequence & consequence)3844 {3845 if (!consequence_) {3846 consequence_ = new CbcConsequence * [numberMembers_];3847 for (int i = 0; i < numberMembers_; i++)3848 consequence_[i] = NULL;3849 }3850 for (int i = 0; i < numberMembers_; i++) {3851 if (members_[i] == iColumn) {3852 consequence_[i] = consequence.clone();3853 break;3854 }3855 }3856 }3857 3858 // Applies a consequence for a single member3859 void3860 CbcNWay::applyConsequence(int iSequence, int state) const3861 {3862 assert (state == 9999  state == 9999);3863 if (consequence_) {3864 CbcConsequence * consequence = consequence_[iSequence];3865 if (consequence)3866 consequence>applyToSolver(model_>solver(), state);3867 }3868 }3869 double3870 CbcNWay::infeasibility(const OsiBranchingInformation * /*info*/,3871 int &preferredWay) const3872 {3873 int numberUnsatis = 0;3874 int j;3875 OsiSolverInterface * solver = model_>solver();3876 const double * solution = model_>testSolution();3877 const double * lower = solver>getColLower();3878 const double * upper = solver>getColUpper();3879 double largestValue = 0.0;3880 3881 double integerTolerance =3882 model_>getDblParam(CbcModel::CbcIntegerTolerance);3883 3884 for (j = 0; j < numberMembers_; j++) {3885 int iColumn = members_[j];3886 double value = solution[iColumn];3887 value = CoinMax(value, lower[iColumn]);3888 value = CoinMin(value, upper[iColumn]);3889 double distance = CoinMin(value  lower[iColumn], upper[iColumn]  value);3890 if (distance > integerTolerance) {3891 numberUnsatis++;3892 largestValue = CoinMax(distance, largestValue);3893 }3894 }3895 preferredWay = 1;3896 if (numberUnsatis) {3897 return largestValue;3898 } else {3899 return 0.0; // satisfied3900 }3901 }3902 3903 // This looks at solution and sets bounds to contain solution3904 void3905 CbcNWay::feasibleRegion()3906 {3907 int j;3908 OsiSolverInterface * solver = model_>solver();3909 const double * solution = model_>testSolution();3910 const double * lower = solver>getColLower();3911 const double * upper = solver>getColUpper();3912 double integerTolerance =3913 model_>getDblParam(CbcModel::CbcIntegerTolerance);3914 for (j = 0; j < numberMembers_; j++) {3915 int iColumn = members_[j];3916 double value = solution[iColumn];3917 value = CoinMax(value, lower[iColumn]);3918 value = CoinMin(value, upper[iColumn]);3919 if (value >= upper[iColumn]  integerTolerance) {3920 solver>setColLower(iColumn, upper[iColumn]);3921 } else {3922 assert (value <= lower[iColumn] + integerTolerance);3923 solver>setColUpper(iColumn, lower[iColumn]);3924 }3925 }3926 }3927 // Redoes data when sequence numbers change3928 void3929 CbcNWay::redoSequenceEtc(CbcModel * model, int numberColumns, const int * originalColumns)3930 {3931 model_ = model;3932 int n2 = 0;3933 for (int j = 0; j < numberMembers_; j++) {3934 int iColumn = members_[j];3935 int i;3936 for (i = 0; i < numberColumns; i++) {3937 if (originalColumns[i] == iColumn)3938 break;3939 }3940 if (i < numberColumns) {3941 members_[n2] = i;3942 consequence_[n2++] = consequence_[j];3943 } else {3944 delete consequence_[j];3945 }3946 }3947 if (n2 < numberMembers_) {3948 printf("** NWay number of members reduced from %d to %d!\n", numberMembers_, n2);3949 numberMembers_ = n2;3950 }3951 }3952 CbcBranchingObject *3953 CbcNWay::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int /*way*/)3954 {3955 int numberFree = 0;3956 int j;3957 3958 //OsiSolverInterface * solver = model_>solver();3959 const double * solution = model_>testSolution();3960 const double * lower = solver>getColLower();3961 const double * upper = solver>getColUpper();3962 int * list = new int[numberMembers_];3963 double * sort = new double[numberMembers_];3964 3965 for (j = 0; j < numberMembers_; j++) {3966 int iColumn = members_[j];3967 double value = solution[iColumn];3968 value = CoinMax(value, lower[iColumn]);3969 value = CoinMin(value, upper[iColumn]);3970 if (upper[iColumn] > lower[iColumn]) {3971 double distance = upper[iColumn]  value;3972 list[numberFree] = j;3973 sort[numberFree++] = distance;3974 }3975 }3976 assert (numberFree);3977 // sort3978 CoinSort_2(sort, sort + numberFree, list);3979 // create object3980 CbcBranchingObject * branch;3981 branch = new CbcNWayBranchingObject(model_, this, numberFree, list);3982 branch>setOriginalObject(this);3983 delete [] list;3984 delete [] sort;3985 return branch;3986 }3987 3988 //##############################################################################3989 3990 // Default Constructor3991 CbcNWayBranchingObject::CbcNWayBranchingObject()3992 : CbcBranchingObject()3993 {3994 order_ = NULL;3995 object_ = NULL;3996 numberInSet_ = 0;3997 way_ = 0;3998 }3999 4000 // Useful constructor4001 CbcNWayBranchingObject::CbcNWayBranchingObject (CbcModel * model,4002 const CbcNWay * nway,4003 int number, const int * order)4004 : CbcBranchingObject(model, nway>id(), 1, 0.5)4005 {4006 numberBranches_ = number;4007 order_ = new int [number];4008 object_ = nway;4009 numberInSet_ = number;4010 memcpy(order_, order, number*sizeof(int));4011 }4012 4013 // Copy constructor4014 CbcNWayBranchingObject::CbcNWayBranchingObject ( const CbcNWayBranchingObject & rhs) : CbcBranchingObject(rhs)4015 {4016 numberInSet_ = rhs.numberInSet_;4017 object_ = rhs.object_;4018 if (numberInSet_) {4019 order_ = new int [numberInSet_];4020 memcpy(order_, rhs.order_, numberInSet_*sizeof(int));4021 } else {4022 order_ = NULL;4023 }4024 }4025 4026 // Assignment operator4027 CbcNWayBranchingObject &4028 CbcNWayBranchingObject::operator=( const CbcNWayBranchingObject & rhs)4029 {4030 if (this != &rhs) {4031 CbcBranchingObject::operator=(rhs);4032 object_ = rhs.object_;4033 delete [] order_;4034 numberInSet_ = rhs.numberInSet_;4035 if (numberInSet_) {4036 order_ = new int [numberInSet_];4037 memcpy(order_, rhs.order_, numberInSet_*sizeof(int));4038 } else {4039 order_ = NULL;4040 }4041 }4042 return *this;4043 }4044 CbcBranchingObject *4045 CbcNWayBranchingObject::clone() const4046 {4047 return (new CbcNWayBranchingObject(*this));4048 }4049 4050 4051 // Destructor4052 CbcNWayBranchingObject::~CbcNWayBranchingObject ()4053 {4054 delete [] order_;4055 }4056 double4057 CbcNWayBranchingObject::branch()4058 {4059 int which = branchIndex_;4060 branchIndex_++;4061 assert (numberBranchesLeft() >= 0);4062 if (which == 0) {4063 // first branch so way_ may mean something4064 assert (way_ == 1  way_ == 1);4065 if (way_ == 1)4066 which++;4067 } else if (which == 1) {4068 // second branch so way_ may mean something4069 assert (way_ == 1  way_ == 1);4070 if (way_ == 1)4071 which;4072 // switch way off4073 way_ = 0;4074 }4075 const double * lower = model_>solver()>getColLower();4076 const double * upper = model_>solver()>getColUpper();4077 const int * members = object_>members();4078 for (int j = 0; j < numberInSet_; j++) {4079 int iSequence = order_[j];4080 int iColumn = members[iSequence];4081 if (j != which) {4082 model_>solver()>setColUpper(iColumn, lower[iColumn]);4083 //model_>solver()>setColLower(iColumn,lower[iColumn]);4084 assert (lower[iColumn] > 1.0e20);4085 // apply any consequences4086 object_>applyConsequence(iSequence, 9999);4087 } else {4088 model_>solver()>setColLower(iColumn, upper[iColumn]);4089 //model_>solver()>setColUpper(iColumn,upper[iColumn]);4090 #ifdef FULL_PRINT4091 printf("Up Fix %d to %g\n", iColumn, upper[iColumn]);4092 #endif4093 assert (upper[iColumn] < 1.0e20);4094 // apply any consequences4095 object_>applyConsequence(iSequence, 9999);4096 }4097 }4098 return 0.0;4099 }4100 void4101 CbcNWayBranchingObject::print()4102 {4103 printf("NWay  Up Fix ");4104 const int * members = object_>members();4105 for (int j = 0; j < way_; j++) {4106 int iColumn = members[order_[j]];4107 printf("%d ", iColumn);4108 }4109 printf("\n");4110 }4111 4112 /** Compare the original object of \c this with the original object of \c4113 brObj. Assumes that there is an ordering of the original objects.4114 This method should be invoked only if \c this and brObj are of the same4115 type.4116 Return negative/0/positive depending on whether \c this is4117 smaller/same/larger than the argument.4118 */4119 int4120 CbcNWayBranchingObject::compareOriginalObject4121 (const CbcBranchingObject* /*brObj*/) const4122 {4123 throw("must implement");4124 }4125 4126 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the4127 same type and must have the same original object, but they may have4128 different feasible regions.4129 Return the appropriate CbcRangeCompare value (first argument being the4130 sub/superset if that's the case). In case of overlap (and if \c4131 replaceIfOverlap is true) replace the current branching object with one4132 whose feasible region is the overlap.4133 */4134 CbcRangeCompare4135 CbcNWayBranchingObject::compareBranchingObject4136 (const CbcBranchingObject* /*brObj*/, const bool /*replaceIfOverlap*/)4137 {4138 throw("must implement");4139 }4140 4141 //##############################################################################4142 4143 // Default Constructor4144 CbcFixVariable::CbcFixVariable ()4145 : CbcConsequence(),4146 numberStates_(0),4147 states_(NULL),4148 startLower_(NULL),4149 startUpper_(NULL),4150 newBound_(NULL),4151 variable_(NULL)4152 {4153 }4154 4155 // One useful Constructor4156 CbcFixVariable::CbcFixVariable (int numberStates, const int * states, const int * numberNewLower,4157 const int ** newLowerValue,4158 const int ** lowerColumn,4159 const int * numberNewUpper, const int ** newUpperValue,4160 const int ** upperColumn)4161 : CbcConsequence(),4162 states_(NULL),4163 startLower_(NULL),4164 startUpper_(NULL),4165 newBound_(NULL),4166 variable_(NULL)4167 {4168 // How much space4169 numberStates_ = numberStates;4170 if (numberStates_) {4171 states_ = new int[numberStates_];4172 memcpy(states_, states, numberStates_*sizeof(int));4173 int i;4174 int n = 0;4175 startLower_ = new int[numberStates_+1];4176 startUpper_ = new int[numberStates_+1];4177 startLower_[0] = 0;4178 //count4179 for (i = 0; i < numberStates_; i++) {4180 n += numberNewLower[i];4181 startUpper_[i] = n;4182 n += numberNewUpper[i];4183 startLower_[i+1] = n;4184 }4185 newBound_ = new double [n];4186 variable_ = new int [n];4187 n = 0;4188 for (i = 0; i < numberStates_; i++) {4189 int j;4190 int k;4191 const int * bound;4192 const int * variable;4193 k = numberNewLower[i];4194 bound = newLowerValue[i];4195 variable = lowerColumn[i];4196 for (j = 0; j < k; j++) {4197 newBound_[n] = bound[j];4198 variable_[n++] = variable[j];4199 }4200 k = numberNewUpper[i];4201 bound = newUpperValue[i];4202 variable = upperColumn[i];4203 for (j = 0; j < k; j++) {4204 newBound_[n] = bound[j];4205 variable_[n++] = variable[j];4206 }4207 }4208 }4209 }4210 4211 // Copy constructor4212 CbcFixVariable::CbcFixVariable ( const CbcFixVariable & rhs)4213 : CbcConsequence(rhs)4214 {4215 numberStates_ = rhs.numberStates_;4216 states_ = NULL;4217 startLower_ = NULL;4218 startUpper_ = NULL;4219 newBound_ = NULL;4220 variable_ = NULL;4221 if (numberStates_) {4222 states_ = CoinCopyOfArray(rhs.states_, numberStates_);4223 startLower_ = CoinCopyOfArray(rhs.startLower_, numberStates_ + 1);4224 startUpper_ = CoinCopyOfArray(rhs.startUpper_, numberStates_ + 1);4225 int n = startLower_[numberStates_];4226 newBound_ = CoinCopyOfArray(rhs.newBound_, n);4227 variable_ = CoinCopyOfArray(rhs.variable_, n);4228 }4229 }4230 4231 // Clone4232 CbcConsequence *4233 CbcFixVariable::clone() const4234 {4235 return new CbcFixVariable(*this);4236 }4237 4238 // Assignment operator4239 CbcFixVariable &4240 CbcFixVariable::operator=( const CbcFixVariable & rhs)4241 {4242 if (this != &rhs) {4243 CbcConsequence::operator=(rhs);4244 delete [] states_;4245 delete [] startLower_;4246 delete [] startUpper_;4247 delete [] newBound_;4248 delete [] variable_;4249 states_ = NULL;4250 startLower_ = NULL;4251 startUpper_ = NULL;4252 newBound_ = NULL;4253 variable_ = NULL;4254 numberStates_ = rhs.numberStates_;4255 if (numberStates_) {4256 states_ = CoinCopyOfArray(rhs.states_, numberStates_);4257 startLower_ = CoinCopyOfArray(rhs.startLower_, numberStates_ + 1);4258 startUpper_ = CoinCopyOfArray(rhs.startUpper_, numberStates_ + 1);4259 int n = startLower_[numberStates_];4260 newBound_ = CoinCopyOfArray(rhs.newBound_, n);4261 variable_ = CoinCopyOfArray(rhs.variable_, n);4262 }4263 }4264 return *this;4265 }4266 4267 // Destructor4268 CbcFixVariable::~CbcFixVariable ()4269 {4270 delete [] states_;4271 delete [] startLower_;4272 delete [] startUpper_;4273 delete [] newBound_;4274 delete [] variable_;4275 }4276 // Set up a startLower for a single member4277 void4278 CbcFixVariable::applyToSolver(OsiSolverInterface * solver, int state) const4279 {4280 assert (state == 9999  state == 9999);4281 // Find state4282 int find;4283 for (find = 0; find < numberStates_; find++)4284 if (states_[find] == state)4285 break;4286 if (find == numberStates_)4287 return;4288 int i;4289 // Set new lower bounds4290 for (i = startLower_[find]; i < startUpper_[find]; i++) {4291 int iColumn = variable_[i];4292 double value = newBound_[i];4293 double oldValue = solver>getColLower()[iColumn];4294 //printf("for %d old lower bound %g, new %g",iColumn,oldValue,value);4295 solver>setColLower(iColumn, CoinMax(value, oldValue));4296 //printf(" => %g\n",solver>getColLower()[iColumn]);4297 }4298 // Set new upper bounds4299 for (i = startUpper_[find]; i < startLower_[find+1]; i++) {4300 int iColumn = variable_[i];4301 double value = newBound_[i];4302 double oldValue = solver>getColUpper()[iColumn];4303 //printf("for %d old upper bound %g, new %g",iColumn,oldValue,value);4304 solver>setColUpper(iColumn, CoinMin(value, oldValue));4305 //printf(" => %g\n",solver>getColUpper()[iColumn]);4306 }4307 }4308 4309 //##############################################################################4310 4311 // Default Constructor4312 CbcDummyBranchingObject::CbcDummyBranchingObject(CbcModel * model)4313 : CbcBranchingObject(model, 0, 0, 0.5)4314 {4315 setNumberBranchesLeft(1);4316 }4317 4318 4319 // Copy constructor4320 CbcDummyBranchingObject::CbcDummyBranchingObject ( const CbcDummyBranchingObject & rhs) : CbcBranchingObject(rhs)4321 {4322 }4323 4324 // Assignment operator4325 CbcDummyBranchingObject &4326 CbcDummyBranchingObject::operator=( const CbcDummyBranchingObject & rhs)4327 {4328 if (this != &rhs) {4329 CbcBranchingObject::operator=(rhs);4330 }4331 return *this;4332 }4333 CbcBranchingObject *4334 CbcDummyBranchingObject::clone() const4335 {4336 return (new CbcDummyBranchingObject(*this));4337 }4338 4339 4340 // Destructor4341 CbcDummyBranchingObject::~CbcDummyBranchingObject ()4342 {4343 }4344 4345 /*4346 Perform a dummy branch4347 */4348 double4349 CbcDummyBranchingObject::branch()4350 {4351 decrementNumberBranchesLeft();4352 return 0.0;4353 }4354 // Print what would happen4355 void4356 CbcDummyBranchingObject::print()4357 {4358 printf("Dummy branch\n");4359 }4360 4361 // Default Constructor4362 CbcGeneral::CbcGeneral()4363 : CbcObject()4364 {4365 }4366 4367 // Constructor from model4368 CbcGeneral::CbcGeneral(CbcModel * model)4369 : CbcObject(model)4370 {4371 }4372 4373 4374 // Destructor4375 CbcGeneral::~CbcGeneral ()4376 {4377 }4378 4379 // Copy constructor4380 CbcGeneral::CbcGeneral ( const CbcGeneral & rhs)4381 : CbcObject(rhs)4382 {4383 }4384 #ifdef COIN_HAS_CLP4385 #include "OsiClpSolverInterface.hpp"4386 #include "CoinWarmStartBasis.hpp"4387 #include "ClpNode.hpp"4388 #include "CbcBranchDynamic.hpp"4389 // Assignment operator4390 CbcGeneral &4391 CbcGeneral::operator=( const CbcGeneral & rhs)4392 {4393 if (this != &rhs) {4394 CbcObject::operator=(rhs);4395 }4396 return *this;4397 }4398 // Infeasibility  large is 0.54399 double4400 CbcGeneral::infeasibility(const OsiBranchingInformation * /*info*/,4401 int &/*preferredWay*/) const4402 {4403 abort();4404 return 0.0;4405 }4406 CbcBranchingObject *4407 CbcGeneral::createCbcBranch(OsiSolverInterface * /*solver*/, const OsiBranchingInformation * /*info*/, int /*way*/)4408 {4409 abort();4410 return NULL;4411 }4412 4413 // Default Constructor4414 CbcGeneralDepth::CbcGeneralDepth ()4415 : CbcGeneral(),4416 maximumDepth_(0),4417 maximumNodes_(0),4418 whichSolution_(1),4419 numberNodes_(0),4420 nodeInfo_(NULL)4421 {4422 }4423 4424 // Useful constructor (which are integer indices)4425 CbcGeneralDepth::CbcGeneralDepth (CbcModel * model, int maximumDepth)4426 : CbcGeneral(model),4427 maximumDepth_(maximumDepth),4428 maximumNodes_(0),4429 whichSolution_(1),4430 numberNodes_(0),4431 nodeInfo_(NULL)4432 {4433 assert(maximumDepth_ < 1000000);4434 if (maximumDepth_ > 0)4435 maximumNodes_ = (1 << maximumDepth_) + 1 + maximumDepth_;4436 else if (maximumDepth_ < 0)4437 maximumNodes_ = 1 + 1  maximumDepth_;4438 else4439 maximumNodes_ = 0;4440 #define MAX_NODES 1004441 maximumNodes_ = CoinMin(maximumNodes_, 1 + maximumDepth_ + MAX_NODES);4442 if (maximumNodes_) {4443 nodeInfo_ = new ClpNodeStuff();4444 nodeInfo_>maximumNodes_ = maximumNodes_;4445 ClpNodeStuff * info = nodeInfo_;4446 // for reduced costs and duals4447 info>solverOptions_ = 7;4448 if (maximumDepth_ > 0) {4449 info>nDepth_ = maximumDepth_;4450 } else {4451 info>nDepth_ =  maximumDepth_;4452 info>solverOptions_ = 32;4453 }4454 ClpNode ** nodeInfo = new ClpNode * [maximumNodes_];4455 for (int i = 0; i < maximumNodes_; i++)4456 nodeInfo[i] = NULL;4457 info>nodeInfo_ = nodeInfo;4458 } else {4459 nodeInfo_ = NULL;4460 }4461 }4462 4463 // Copy constructor4464 CbcGeneralDepth::CbcGeneralDepth ( const CbcGeneralDepth & rhs)4465 : CbcGeneral(rhs)4466 {4467 maximumDepth_ = rhs.maximumDepth_;4468 maximumNodes_ = rhs.maximumNodes_;4469 whichSolution_ = 1;4470 numberNodes_ = 0;4471 if (maximumNodes_) {4472 assert (rhs.nodeInfo_);4473 nodeInfo_ = new ClpNodeStuff(*rhs.nodeInfo_);4474 nodeInfo_>maximumNodes_ = maximumNodes_;4475 ClpNodeStuff * info = nodeInfo_;4476 if (maximumDepth_ > 0) {4477 info>nDepth_ = maximumDepth_;4478 } else {4479 info>nDepth_ =  maximumDepth_;4480 info>solverOptions_ = 32;4481 }4482 if (!info>nodeInfo_) {4483 ClpNode ** nodeInfo = new ClpNode * [maximumNodes_];4484 for (int i = 0; i < maximumNodes_; i++)4485 nodeInfo[i] = NULL;4486 info>nodeInfo_ = nodeInfo;4487 }4488 } else {4489 nodeInfo_ = NULL;4490 }4491 }4492 4493 // Clone4494 CbcObject *4495 CbcGeneralDepth::clone() const4496 {4497 return new CbcGeneralDepth(*this);4498 }4499 4500 // Assignment operator4501 CbcGeneralDepth &4502 CbcGeneralDepth::operator=( const CbcGeneralDepth & rhs)4503 {4504 if (this != &rhs) {4505 CbcGeneral::operator=(rhs);4506 delete nodeInfo_;4507 maximumDepth_ = rhs.maximumDepth_;4508 maximumNodes_ = rhs.maximumNodes_;4509 whichSolution_ = 1;4510 numberNodes_ = 0;4511 if (maximumDepth_) {4512 assert (rhs.nodeInfo_);4513 nodeInfo_ = new ClpNodeStuff(*rhs.nodeInfo_);4514 nodeInfo_>maximumNodes_ = maximumNodes_;4515 } else {4516 nodeInfo_ = NULL;4517 }4518 }4519 return *this;4520 }4521 4522 // Destructor4523 CbcGeneralDepth::~CbcGeneralDepth ()4524 {4525 delete nodeInfo_;4526 }4527 // Infeasibility  large is 0.54528 double4529 CbcGeneralDepth::infeasibility(const OsiBranchingInformation * /*info*/,4530 int &/*preferredWay*/) const4531 {4532 whichSolution_ = 1;4533 // should use genuine OsiBranchingInformation usefulInfo = model_>usefulInformation();4534 // for now assume only called when correct4535 //if (usefulInfo.depth_>=4&&!model_>parentModel()4536 // &&(usefulInfo.depth_%2)==0) {4537 if (true) {4538 OsiSolverInterface * solver = model_>solver();4539 OsiClpSolverInterface * clpSolver4540 = dynamic_cast<OsiClpSolverInterface *> (solver);4541 if (clpSolver) {4542 ClpNodeStuff * info = nodeInfo_;4543 info>integerTolerance_ = model_>getIntegerTolerance();4544 info>integerIncrement_ = model_>getCutoffIncrement();4545 info>numberBeforeTrust_ = model_>numberBeforeTrust();4546 info>stateOfSearch_ = model_>stateOfSearch();4547 // Compute "small" change in branch4548 int nBranches = model_>getIntParam(CbcModel::CbcNumberBranches);4549 if (nBranches) {4550 double average = model_>getDblParam(CbcModel::CbcSumChange) / static_cast<double>(nBranches);4551 info>smallChange_ =4552 CoinMax(average * 1.0e5, model_>getDblParam(CbcModel::CbcSmallestChange));4553 info>smallChange_ = CoinMax(info>smallChange_, 1.0e8);4554 } else {4555 info>smallChange_ = 1.0e8;4556 }4557 int numberIntegers = model_>numberIntegers();4558 double * down = new double[numberIntegers];4559 double * up = new double[numberIntegers];4560 int * priority = new int[numberIntegers];4561 int * numberDown = new int[numberIntegers];4562 int * numberUp = new int[numberIntegers];4563 int * numberDownInfeasible = new int[numberIntegers];4564 int * numberUpInfeasible = new int[numberIntegers];4565 model_>fillPseudoCosts(down, up, priority, numberDown, numberUp,4566 numberDownInfeasible, numberUpInfeasible);4567 info>fillPseudoCosts(down, up, priority, numberDown, numberUp,4568 numberDownInfeasible,4569 numberUpInfeasible, numberIntegers);4570 info>presolveType_ = 1;4571 delete [] down;4572 delete [] up;4573 delete [] numberDown;4574 delete [] numberUp;4575 delete [] numberDownInfeasible;4576 delete [] numberUpInfeasible;4577 bool takeHint;4578 OsiHintStrength strength;4579 solver>getHintParam(OsiDoReducePrint, takeHint, strength);4580 ClpSimplex * simplex = clpSolver>getModelPtr();4581 int saveLevel = simplex>logLevel();4582 if (strength != OsiHintIgnore && takeHint && saveLevel == 1)4583 simplex>setLogLevel(0);4584 clpSolver>setBasis();4585 whichSolution_ = simplex>fathomMany(info);4586 //printf("FAT %d nodes, %d iterations\n",4587 //info>numberNodesExplored_,info>numberIterations_);4588 //printf("CbcBranch %d rows, %d columns\n",clpSolver>getNumRows(),4589 // clpSolver>getNumCols());4590 model_>incrementExtra(info>numberNodesExplored_,4591 info>numberIterations_);4592 // update pseudo costs4593 double smallest = 1.0e50;4594 double largest = 1.0;4595 OsiObject ** objects = model_>objects();4596 #ifndef NDEBUG4597 const int * integerVariable = model_>integerVariable();4598 #endif4599 for (int i = 0; i < numberIntegers; i++) {4600 #ifndef NDEBUG4601 CbcSimpleIntegerDynamicPseudoCost * obj =4602 dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(objects[i]) ;4603 assert (obj && obj>columnNumber() == integerVariable[i]);4604 #else4605 CbcSimpleIntegerDynamicPseudoCost * obj =4606 static_cast <CbcSimpleIntegerDynamicPseudoCost *>(objects[i]) ;4607 #endif4608 if (info>numberUp_[i] > 0) {4609 if (info>downPseudo_[i] > largest)4610 largest = info>downPseudo_[i];4611 if (info>downPseudo_[i] < smallest)4612 smallest = info>downPseudo_[i];4613 if (info>upPseudo_[i] > largest)4614 largest = info>upPseudo_[i];4615 if (info>upPseudo_[i] < smallest)4616 smallest = info>upPseudo_[i];4617 obj>updateAfterMini(info>numberDown_[i],4618 info>numberDownInfeasible_[i],4619 info>downPseudo_[i],4620 info>numberUp_[i],4621 info>numberUpInfeasible_[i],4622 info>upPseudo_[i]);4623 }4624 }4625 //printf("range of costs %g to %g\n",smallest,largest);4626 simplex>setLogLevel(saveLevel);4627 numberNodes_ = info>nNodes_;4628 int numberDo = numberNodes_;4629 if (whichSolution_ >= 0)4630 numberDo;4631 if (numberDo > 0) {4632 return 0.5;4633 } else {4634 // no solution4635 return COIN_DBL_MAX; // say infeasible4636 }4637 } else {4638 return 1.0;4639 }4640 } else {4641 return 1.0;4642 }4643 }4644 4645 // This looks at solution and sets bounds to contain solution4646 void4647 CbcGeneralDepth::feasibleRegion()4648 {4649 // Other stuff should have done this4650 }4651 // Redoes data when sequence numbers change4652 void4653 CbcGeneralDepth::redoSequenceEtc(CbcModel * /*model*/,4654 int /*numberColumns*/,4655 const int * /*originalColumns*/)4656 {4657 }4658 4659 //#define CHECK_PATH4660 #ifdef CHECK_PATH4661 extern const double * debuggerSolution_Z;4662 extern int numberColumns_Z;4663 extern int gotGoodNode_Z;4664 #endif4665 CbcBranchingObject *4666 CbcGeneralDepth::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int /*way*/)4667 {4668 int numberDo = numberNodes_;4669 if (whichSolution_ >= 0)4670 numberDo;4671 assert (numberDo > 0);4672 // create object4673 CbcGeneralBranchingObject * branch = new CbcGeneralBranchingObject(model_);4674 // skip solution4675 branch>numberSubProblems_ = numberDo;4676 // If parentBranch_ back in then will have to be 2*4677 branch>numberSubLeft_ = numberDo;4678 branch>setNumberBranches(numberDo);4679 CbcSubProblem * sub = new CbcSubProblem[numberDo];4680 int iProb = 0;4681 branch>subProblems_ = sub;4682 branch>numberRows_ = model_>solver()>getNumRows();4683 int iNode;4684 //OsiSolverInterface * solver = model_>solver();4685 OsiClpSolverInterface * clpSolver4686 = dynamic_cast<OsiClpSolverInterface *> (solver);4687 assert (clpSolver);4688 ClpSimplex * simplex = clpSolver>getModelPtr();4689 int numberColumns = simplex>numberColumns();4690 double * lowerBefore = CoinCopyOfArray(simplex>getColLower(),4691 numberColumns);4692 double * upperBefore = CoinCopyOfArray(simplex>getColUpper(),4693 numberColumns);4694 ClpNodeStuff * info = nodeInfo_;4695 double * weight = new double[numberNodes_];4696 int * whichNode = new int [numberNodes_];4697 // Sort4698 for (iNode = 0; iNode < numberNodes_; iNode++) {4699 if (iNode != whichSolution_) {4700 double objectiveValue = info>nodeInfo_[iNode]>objectiveValue();4701 double sumInfeasibilities = info>nodeInfo_[iNode]>sumInfeasibilities();4702 int numberInfeasibilities = info>nodeInfo_[iNode]>numberInfeasibilities();4703 double thisWeight = 0.0;4704 #if 14705 // just closest4706 thisWeight = 1.0e9 * numberInfeasibilities;4707 thisWeight += sumInfeasibilities;4708 thisWeight += 1.0e7 * objectiveValue;4709 // Try estimate4710 thisWeight = info>nodeInfo_[iNode]>estimatedSolution();4711 #else4712 thisWeight = 1.0e3 * numberInfeasibilities;4713 thisWeight += 1.0e5 * sumInfeasibilities;4714 thisWeight += objectiveValue;4715 #endif4716 whichNode[iProb] = iNode;4717 weight[iProb++] = thisWeight;4718 }4719 }4720 assert (iProb == numberDo);4721 CoinSort_2(weight, weight + numberDo, whichNode);4722 for (iProb = 0; iProb < numberDo; iProb++) {4723 iNode = whichNode[iProb];4724 ClpNode * node = info>nodeInfo_[iNode];4725 // move bounds4726 node>applyNode(simplex, 3);4727 // create subproblem4728 sub[iProb] = CbcSubProblem(clpSolver, lowerBefore, upperBefore,4729 node>statusArray(), node>depth());4730 sub[iProb].objectiveValue_ = node>objectiveValue();4731 sub[iProb].sumInfeasibilities_ = node>sumInfeasibilities();4732 sub[iProb].numberInfeasibilities_ = node>numberInfeasibilities();4733 #ifdef CHECK_PATH4734 if (simplex>numberColumns() == numberColumns_Z) {4735 bool onOptimal = true;4736 const double * columnLower = simplex>columnLower();4737 const double * columnUpper = simplex>columnUpper();4738 for (int i = 0; i < numberColumns_Z; i++) {4739 if (iNode == gotGoodNode_Z)4740 printf("good %d %d %g %g\n", iNode, i, columnLower[i], columnUpper[i]);4741 if (columnUpper[i] < debuggerSolution_Z[i]  columnLower[i] > debuggerSolution_Z[i] && simplex>isInteger(i)) {4742 onOptimal = false;4743 break;4744 }4745 }4746 if (onOptimal) {4747 printf("adding to node %x as %d  objs\n", this, iProb);4748 for (int j = 0; j <= iProb; j++)4749 printf("%d %g\n", j, sub[j].objectiveValue_);4750 }4751 }4752 #endif4753 }4754 delete [] weight;4755 delete [] whichNode;4756 const double * lower = solver>getColLower();4757 const double * upper = solver>getColUpper();4758 // restore bounds4759 for ( int j = 0; j < numberColumns; j++) {4760 if (lowerBefore[j] != lower[j])4761 solver>setColLower(j, lowerBefore[j]);4762 if (upperBefore[j] != upper[j])4763 solver>setColUpper(j, upperBefore[j]);4764 }4765 delete [] upperBefore;4766 delete [] lowerBefore;4767 return branch;4768 }4769 4770 // Default Constructor4771 CbcGeneralBranchingObject::CbcGeneralBranchingObject()4772 : CbcBranchingObject(),4773 subProblems_(NULL),4774 node_(NULL),4775 numberSubProblems_(0),4776 numberSubLeft_(0),4777 whichNode_(1),4778 numberRows_(0)4779 {4780 // printf("CbcGeneral %x default constructor\n",this);4781 }4782 4783 // Useful constructor4784 CbcGeneralBranchingObject::CbcGeneralBranchingObject (CbcModel * model)4785 : CbcBranchingObject(model, 1, 1, 0.5),4786 subProblems_(NULL),4787 node_(NULL),4788 numberSubProblems_(0),4789 numberSubLeft_(0),4790 whichNode_(1),4791 numberRows_(0)4792 {4793 //printf("CbcGeneral %x useful constructor\n",this);4794 }4795 4796 // Copy constructor4797 CbcGeneralBranchingObject::CbcGeneralBranchingObject ( const CbcGeneralBranchingObject & rhs)4798 : CbcBranchingObject(rhs),4799 subProblems_(NULL),4800 node_(rhs.node_),4801 numberSubProblems_(rhs.numberSubProblems_),4802 numberSubLeft_(rhs.numberSubLeft_),4803 whichNode_(rhs.whichNode_),4804 numberRows_(rhs.numberRows_)4805 {4806 abort();4807 if (numberSubProblems_) {4808 subProblems_ = new CbcSubProblem[numberSubProblems_];4809 for (int i = 0; i < numberSubProblems_; i++)4810 subProblems_[i] = rhs.subProblems_[i];4811 }4812 }4813 4814 // Assignment operator4815 CbcGeneralBranchingObject &4816 CbcGeneralBranchingObject::operator=( const CbcGeneralBranchingObject & rhs)4817 {4818 if (this != &rhs) {4819 abort();4820 CbcBranchingObject::operator=(rhs);4821 delete [] subProblems_;4822 numberSubProblems_ = rhs.numberSubProblems_;4823 numberSubLeft_ = rhs.numberSubLeft_;4824 whichNode_ = rhs.whichNode_;4825 numberRows_ = rhs.numberRows_;4826 if (numberSubProblems_) {4827 subProblems_ = new CbcSubProblem[numberSubProblems_];4828 for (int i = 0; i < numberSubProblems_; i++)4829 subProblems_[i] = rhs.subProblems_[i];4830 } else {4831 subProblems_ = NULL;4832 }4833 node_ = rhs.node_;4834 }4835 return *this;4836 }4837 CbcBranchingObject *4838 CbcGeneralBranchingObject::clone() const4839 {4840 return (new CbcGeneralBranchingObject(*this));4841 }4842 4843 4844 // Destructor4845 CbcGeneralBranchingObject::~CbcGeneralBranchingObject ()4846 {4847 //printf("CbcGeneral %x destructor\n",this);4848 delete [] subProblems_;4849 }4850 bool doingDoneBranch = false;4851 double4852 CbcGeneralBranchingObject::branch()4853 {4854 double cutoff = model_>getCutoff();4855 //printf("GenB %x whichNode %d numberLeft %d which %d\n",4856 // this,whichNode_,numberBranchesLeft(),branchIndex());4857 if (whichNode_ < 0) {4858 assert (node_);4859 bool applied = false;4860 while (numberBranchesLeft()) {4861 int which = branchIndex();4862 decrementNumberBranchesLeft();4863 CbcSubProblem * thisProb = subProblems_ + which;4864 if (thisProb>objectiveValue_ < cutoff) {4865 //printf("branch %x (sub %x) which now %d\n",this,4866 // subProblems_,which);4867 OsiSolverInterface * solver = model_>solver();4868 thisProb>apply(solver);4869 OsiClpSolverInterface * clpSolver4870 = dynamic_cast<OsiClpSolverInterface *> (solver);4871 assert (clpSolver);4872 // Move status to basis4873 clpSolver>setWarmStart(NULL);4874 //ClpSimplex * simplex = clpSolver>getModelPtr();4875 node_>setObjectiveValue(thisProb>objectiveValue_);4876 node_>setSumInfeasibilities(thisProb>sumInfeasibilities_);4877 node_>setNumberUnsatisfied(thisProb>numberInfeasibilities_);4878 applied = true;4879 doingDoneBranch = true;4880 break;4881 } else if (numberBranchesLeft()) {4882 node_>nodeInfo()>branchedOn() ;4883 }4884 }4885 if (!applied) {4886 // no good one4887 node_>setObjectiveValue(cutoff + 1.0e20);4888 node_>setSumInfeasibilities(1.0);4889 node_>setNumberUnsatisfied(1);4890 assert (whichNode_ < 0);4891 }4892 } else {4893 decrementNumberBranchesLeft();4894 CbcSubProblem * thisProb = subProblems_ + whichNode_;4895 assert (thisProb>objectiveValue_ < cutoff);4896 OsiSolverInterface * solver = model_>solver();4897 thisProb>apply(solver);4898 //OsiClpSolverInterface * clpSolver4899 //= dynamic_cast<OsiClpSolverInterface *> (solver);4900 //assert (clpSolver);4901 // Move status to basis4902 //clpSolver>setWarmStart(NULL);4903 }4904 return 0.0;4905 }4906 /* Double checks in case node can change its mind!4907 Can change objective etc */4908 void4909 CbcGeneralBranchingObject::checkIsCutoff(double cutoff)4910 {4911 assert (node_);4912 int first = branchIndex();4913 int last = first + numberBranchesLeft();4914 for (int which = first; which < last; which++) {4915 CbcSubProblem * thisProb = subProblems_ + which;4916 if (thisProb>objectiveValue_ < cutoff) {4917 node_>setObjectiveValue(thisProb>objectiveValue_);4918 node_>setSumInfeasibilities(thisProb>sumInfeasibilities_);4919 node_>setNumberUnsatisfied(thisProb>numberInfeasibilities_);4920 break;4921 }4922 }4923 }4924 // Print what would happen4925 void4926 CbcGeneralBranchingObject::print()4927 {4928 //printf("CbcGeneralObject has %d subproblems\n",numberSubProblems_);4929 }4930 // Fill in current objective etc4931 void4932 CbcGeneralBranchingObject::state(double & objectiveValue,4933 double & sumInfeasibilities,4934 int & numberUnsatisfied, int which) const4935 {4936 assert (which >= 0 && which < numberSubProblems_);4937 const CbcSubProblem * thisProb = subProblems_ + which;4938 objectiveValue = thisProb>objectiveValue_;4939 sumInfeasibilities = thisProb>sumInfeasibilities_;4940 numberUnsatisfied = thisProb>numberInfeasibilities_;4941 }4942 /** Compare the original object of \c this with the original object of \c4943 brObj. Assumes that there is an ordering of the original objects.4944 This method should be invoked only if \c this and brObj are of the same4945 type.4946 Return negative/0/positive depending on whether \c this is4947 smaller/same/larger than the argument.4948 */4949 int4950 CbcGeneralBranchingObject::compareOriginalObject4951 (const CbcBranchingObject* /*brObj*/) const4952 {4953 throw("must implement");4954 }4955 4956 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the4957 same type and must have the same original object, but they may have4958 different feasible regions.4959 Return the appropriate CbcRangeCompare value (first argument being the4960 sub/superset if that's the case). In case of overlap (and if \c4961 replaceIfOverlap is true) replace the current branching object with one4962 whose feasible region is the overlap.4963 */4964 CbcRangeCompare4965 CbcGeneralBranchingObject::compareBranchingObject4966 (const CbcBranchingObject* /*brObj*/, const bool /*replaceIfOverlap*/)4967 {4968 throw("must implement");4969 }4970 4971 // Default Constructor4972 CbcOneGeneralBranchingObject::CbcOneGeneralBranchingObject()4973 : CbcBranchingObject(),4974 object_(NULL),4975 whichOne_(1)4976 {4977 //printf("CbcOneGeneral %x default constructor\n",this);4978 }4979 4980 // Useful constructor4981 CbcOneGeneralBranchingObject::CbcOneGeneralBranchingObject (CbcModel * model,4982 CbcGeneralBranchingObject * object,4983 int whichOne)4984 : CbcBranchingObject(model, 1, 1, 0.5),4985 object_(object),4986 whichOne_(whichOne)4987 {4988 //printf("CbcOneGeneral %x useful constructor object %x %d left\n",this,4989 // object_,object_>numberSubLeft_);4990 numberBranches_ = 1;4991 }4992 4993 // Copy constructor4994 CbcOneGeneralBranchingObject::CbcOneGeneralBranchingObject ( const CbcOneGeneralBranchingObject & rhs)4995 : CbcBranchingObject(rhs),4996 object_(rhs.object_),4997 whichOne_(rhs.whichOne_)4998 {4999 }5000 5001 // Assignment operator5002 CbcOneGeneralBranchingObject &5003 CbcOneGeneralBranchingObject::operator=( const CbcOneGeneralBranchingObject & rhs)5004 {5005 if (this != &rhs) {5006 CbcBranchingObject::operator=(rhs);5007 object_ = rhs.object_;5008 whichOne_ = rhs.whichOne_;5009 }5010 return *this;5011 }5012 CbcBranchingObject *5013 CbcOneGeneralBranchingObject::clone() const5014 {5015 return (new CbcOneGeneralBranchingObject(*this));5016 }5017 5018 5019 // Destructor5020 CbcOneGeneralBranchingObject::~CbcOneGeneralBranchingObject ()5021 {5022 //printf("CbcOneGeneral %x destructor object %x %d left\n",this,5023 // object_,object_>numberSubLeft_);5024 assert (object_>numberSubLeft_ > 0 &&5025 object_>numberSubLeft_ < 1000000);5026 if (!object_>decrementNumberLeft()) {5027 // printf("CbcGeneral %x yy destructor\n",object_);5028 delete object_;5029 }5030 }5031 double5032 CbcOneGeneralBranchingObject::branch()5033 {5034 assert (numberBranchesLeft());5035 decrementNumberBranchesLeft();5036 assert (!numberBranchesLeft());5037 object_>setWhichNode(whichOne_);5038 object_>branch();5039 return 0.0;5040 }5041 /* Double checks in case node can change its mind!5042 Can change objective etc */5043 void5044 CbcOneGeneralBranchingObject::checkIsCutoff(double /*cutoff*/)5045 {5046 assert (numberBranchesLeft());5047 }5048 // Print what would happen5049 void5050 CbcOneGeneralBranchingObject::print()5051 {5052 //printf("CbcOneGeneralObject has 1 subproblem\n");5053 }5054 /** Compare the original object of \c this with the original object of \c5055 brObj. Assumes that there is an ordering of the original objects.5056 This method should be invoked only if \c this and brObj are of the same5057 type.5058 Return negative/0/positive depending on whether \c this is5059 smaller/same/larger than the argument.5060 */5061 int5062 CbcOneGeneralBranchingObject::compareOriginalObject5063 (const CbcBranchingObject* /*brObj*/) const5064 {5065 throw("must implement");5066 }5067 5068 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the5069 same type and must have the same original object, but they may have5070 different feasible regions.5071 Return the appropriate CbcRangeCompare value (first argument being the5072 sub/superset if that's the case). In case of overlap (and if \c5073 replaceIfOverlap is true) replace the current branching object with one5074 whose feasible region is the overlap.5075 */5076 CbcRangeCompare5077 CbcOneGeneralBranchingObject::compareBranchingObject5078 (const CbcBranchingObject* /*brObj*/, const bool /*replaceIfOverlap*/)5079 {5080 throw("must implement");5081 }5082 // Default Constructor5083 CbcSubProblem::CbcSubProblem()5084 : objectiveValue_(0.0),5085 sumInfeasibilities_(0.0),5086 variables_(NULL),5087 newBounds_(NULL),5088 status_(NULL),5089 depth_(0),5090 numberChangedBounds_(0),5091 numberInfeasibilities_(0)5092 {5093 }5094 5095 // Useful constructor5096 CbcSubProblem::CbcSubProblem (const OsiSolverInterface * solver,5097 const double * lastLower,5098 const double * lastUpper,5099 const unsigned char * status,5100 int depth)5101 : objectiveValue_(0.0),5102 sumInfeasibilities_(0.0),5103 variables_(NULL),5104 newBounds_(NULL),5105 status_(NULL),5106 depth_(depth),5107 numberChangedBounds_(0),5108 numberInfeasibilities_(0)5109 {5110 const double * lower = solver>getColLower();5111 const double * upper = solver>getColUpper();5112 5113 numberChangedBounds_ = 0;5114 int numberColumns = solver>getNumCols();5115 int i;5116 for (i = 0; i < numberColumns; i++) {5117 if (lower[i] != lastLower[i])5118 numberChangedBounds_++;5119 if (upper[i] != lastUpper[i])5120 numberChangedBounds_++;5121 }5122 if (numberChangedBounds_) {5123 newBounds_ = new double [numberChangedBounds_] ;5124 variables_ = new int [numberChangedBounds_] ;5125 numberChangedBounds_ = 0;5126 for (i = 0; i < numberColumns; i++) {5127 if (lower[i] != lastLower[i]) {5128 variables_[numberChangedBounds_] = i;5129 newBounds_[numberChangedBounds_++] = lower[i];5130 }5131 if (upper[i] != lastUpper[i]) {5132 variables_[numberChangedBounds_] = i  0x80000000;5133 newBounds_[numberChangedBounds_++] = upper[i];5134 }5135 #ifdef CBC_DEBUG5136 if (lower[i] != lastLower[i]) {5137 std::cout5138 << "lower on " << i << " changed from "5139 << lastLower[i] << " to " << lower[i] << std::endl ;5140 }5141 if (upper[i] != lastUpper[i]) {5142 std::cout5143 << "upper on " << i << " changed from "5144 << lastUpper[i] << " to " << upper[i] << std::endl ;5145 }5146 #endif5147 }5148 #ifdef CBC_DEBUG5149 std::cout << numberChangedBounds_ << " changed bounds." << std::endl ;5150 #endif5151 }5152 const OsiClpSolverInterface * clpSolver5153 = dynamic_cast<const OsiClpSolverInterface *> (solver);5154 assert (clpSolver);5155 // Do difference5156 // Current basis5157 status_ = clpSolver>getBasis(status);5158 assert (status_>fullBasis());5159 //status_>print();5160 }5161 5162 // Copy constructor5163 CbcSubProblem::CbcSubProblem ( const CbcSubProblem & rhs)5164 : objectiveValue_(rhs.objectiveValue_),5165 sumInfeasibilities_(rhs.sumInfeasibilities_),5166 variables_(NULL),5167 newBounds_(NULL),5168 status_(NULL),5169 depth_(rhs.depth_),5170 numberChangedBounds_(rhs.numberChangedBounds_),5171 numberInfeasibilities_(rhs.numberInfeasibilities_)5172 {5173 if (numberChangedBounds_) {5174 variables_ = CoinCopyOfArray(rhs.variables_, numberChangedBounds_);5175 newBounds_ = CoinCopyOfArray(rhs.newBounds_, numberChangedBounds_);5176 }5177 if (rhs.status_) {5178 status_ = new CoinWarmStartBasis(*rhs.status_);5179 }5180 }5181 5182 // Assignment operator5183 CbcSubProblem &5184 CbcSubProblem::operator=( const CbcSubProblem & rhs)5185 {5186 if (this != &rhs) {5187 delete [] variables_;5188 delete [] newBounds_;5189 delete status_;5190 objectiveValue_ = rhs.objectiveValue_;5191 sumInfeasibilities_ = rhs.sumInfeasibilities_;5192 depth_ = rhs.depth_;5193 numberChangedBounds_ = rhs.numberChangedBounds_;5194 numberInfeasibilities_ = rhs.numberInfeasibilities_;5195 if (numberChangedBounds_) {5196 variables_ = CoinCopyOfArray(rhs.variables_, numberChangedBounds_);5197 newBounds_ = CoinCopyOfArray(rhs.newBounds_, numberChangedBounds_);5198 } else {5199 variables_ = NULL;5200 newBounds_ = NULL;5201 }5202 if (rhs.status_) {5203 status_ = new CoinWarmStartBasis(*rhs.status_);5204 } else {5205 status_ = NULL;5206 }5207 }5208 return *this;5209 }5210 5211 // Destructor5212 CbcSubProblem::~CbcSubProblem ()5213 {5214 delete [] variables_;5215 delete [] newBounds_;5216 delete status_;5217 }5218 // Apply subproblem5219 void5220 CbcSubProblem::apply(OsiSolverInterface * solver, int what) const5221 {5222 int i;5223 if ((what&1) != 0) {5224 int nSame = 0;5225 for (i = 0; i < numberChangedBounds_; i++) {5226 int variable = variables_[i];5227 int k = variable & 0x3fffffff;5228 if ((variable&0x80000000) == 0) {5229 // lower bound changing5230 //#define CBC_PRINT25231 #ifdef CBC_PRINT25232 if (solver>getColLower()[k] != newBounds_[i])5233 printf("lower change for column %d  from %g to %g\n", k, solver>getColLower()[k], newBounds_[i]);5234 #endif5235 #ifndef NDEBUG5236 if ((variable&0x40000000) == 0 && true) {5237 double oldValue = solver>getColLower()[k];5238 assert (newBounds_[i] > oldValue  1.0e8);5239 if (newBounds_[i] < oldValue + 1.0e8) {5240 #ifdef CBC_PRINT25241 printf("bad null lower change for column %d  bound %g\n", k, oldValue);5242 #endif5243 if (newBounds_[i] == oldValue)5244 nSame++;5245 }5246 }5247 #endif5248 solver>setColLower(k, newBounds_[i]);5249 } else {5250 // upper bound changing5251 #ifdef CBC_PRINT25252 if (solver>getColUpper()[k] != newBounds_[i])5253 printf("upper change for column %d  from %g to %g\n", k, solver>getColUpper()[k], newBounds_[i]);5254 #endif5255 #ifndef NDEBUG5256 if ((variable&0x40000000) == 0 && true) {5257 double oldValue = solver>getColUpper()[k];5258 assert (newBounds_[i] < oldValue + 1.0e8);5259 if (newBounds_[i] > oldValue  1.0e8) {5260 #ifdef CBC_PRINT25261 printf("bad null upper change for column %d  bound %g\n", k, oldValue);5262 #endif5263 if (newBounds_[i] == oldValue)5264 nSame++;5265 }5266 }5267 #endif5268 solver>setColUpper(k, newBounds_[i]);5269 }5270 }5271 if (nSame && (nSame < numberChangedBounds_  (what&3) != 3))5272 printf("%d changes out of %d redundant %d\n",5273 nSame, numberChangedBounds_, what);5274 else if (numberChangedBounds_ && what == 7 && !nSame)5275 printf("%d good changes %d\n",5276 numberChangedBounds_, what);5277 }5278 #if 05279 if ((what&2) != 0) {5280 OsiClpSolverInterface * clpSolver5281 = dynamic_cast<OsiClpSolverInterface *> (solver);5282 assert (clpSolver);5283 //assert (clpSolver>getNumRows()==numberRows_);5284 //clpSolver>setBasis(*status_);5285 // Current basis5286 CoinWarmStartBasis * basis = clpSolver>getPointerToWarmStart();5287 printf("BBBB\n");5288 basis>print();5289 assert (basis>fullBasis());5290 basis>applyDiff(status_);5291 printf("diff applied %x\n", status_);5292 printf("CCCC\n");5293 basis>print();5294 assert (basis>fullBasis());5295 #ifndef NDEBUG5296 if (!basis>fullBasis())5297 printf("Debug this basis!!\n");5298 #endif5299 clpSolver>setBasis(*basis);5300 }5301 #endif5302 if ((what&8) != 0) {5303 OsiClpSolverInterface * clpSolver5304 = dynamic_cast<OsiClpSolverInterface *> (solver);5305 assert (clpSolver);5306 clpSolver>setBasis(*status_);5307 delete status_;5308 status_ = NULL;5309 }5310 }5311 #endif5312 5313 /** Compare the original object of \c this with the original object of \c5314 brObj. Assumes that there is an ordering of the original objects.5315 This method should be invoked only if \c this and brObj are of the same5316 type.5317 Return negative/0/positive depending on whether \c this is5318 smaller/same/larger than the argument.5319 */5320 int5321 CbcDummyBranchingObject::compareOriginalObject5322 (const CbcBranchingObject* /*brObj*/) const5323 {5324 throw("must implement");5325 }5326 5327 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the5328 same type and must have the same original object, but they may have5329 different feasible regions.5330 Return the appropriate CbcRangeCompare value (first argument being the5331 sub/superset if that's the case). In case of overlap (and if \c5332 replaceIfOverlap is true) replace the current branching object with one5333 whose feasible region is the overlap.5334 */5335 CbcRangeCompare5336 CbcDummyBranchingObject::compareBranchingObject5337 (const CbcBranchingObject* /*brObj*/, const bool /*replaceIfOverlap*/)5338 {5339 throw("must implement");5340 }
Note: See TracChangeset
for help on using the changeset viewer.